Calculate the correlation function over the extended domain
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
complex, | intent(out), | dimension (:,:,:) | :: | cfnc | ||
real, | intent(out), | dimension (:,:,:) | :: | phi2extend |
subroutine correlation_extend (cfnc, phi2extend)
! use constants, only: pi
use fields_arrays, only: phinew
use theta_grid, only: ntgrid, jacob, delthet, ntgrid
! use theta_grid, only: theta
use kt_grids, only: ntheta0, naky, jtwist_out
implicit none
real, dimension (:,:,:), intent (out) :: phi2extend
complex, dimension (:,:,:), intent (out) :: cfnc
integer :: ig, it, ik, im, igmod
integer :: itshift, nconnect, offset
! real :: fac
real, dimension (:), allocatable :: dl_over_b
complex, dimension (:,:,:), allocatable :: phiextend
complex, dimension (:,:), allocatable :: phir
allocate (dl_over_b(2*ntgrid+1))
dl_over_b = delthet*jacob
allocate (phiextend(ntg_extend,ntheta0,naky)) ; phiextend = 0.0
allocate (phir(-ntgrid:ntgrid,ntheta0))
cfnc = 0.0 ; phiextend = 0.0
offset = ((ntheta0-1)/jtwist_out)*(2*ntgrid+1)/2
do ig = -ntgrid, ntgrid
call reorder_kx (phinew(ig,:,1), phir(ig,:))
end do
phiextend(offset+1:offset+(2*ntgrid+1),:,1) = phir
do it = 1, ntheta0
do im = 1, 2*ntgrid+1
do ig = im+offset, offset+(2*ntgrid+1)
igmod = mod(ig-offset-1,2*ntgrid+1)+1
cfnc(im,it,1) = cfnc(im,it,1) &
+ phiextend(ig,it,1)*conjg(phiextend(ig-im+1,it,1)) &
* dl_over_b(igmod)
end do
end do
if (abs(cfnc(1, it, 1)) > epsilon(0.0)) then
cfnc(:, it, 1) = cfnc(:, it, 1) / cfnc(1, it, 1)
end if
end do
do ik = 2, naky
do ig = -ntgrid, ntgrid
call reorder_kx (phinew(ig,:,ik), phir(ig,:))
end do
! shift in kx due to parallel boundary condition
! also the number of independent theta0s
itshift = jtwist_out*(ik-1)
do it = 1, min(itshift,ntheta0)
! number of connections between kx's
nconnect = (ntheta0-it)/itshift
! shift of theta index to account for fact that not all ky's
! have same length in extended theta
offset = (2*ntgrid+1)*((ntheta0-1)/jtwist_out-nconnect)/2
do ig = 1, nconnect+1
phiextend(offset+(ig-1)*(2*ntgrid+1)+1:offset+ig*(2*ntgrid+1),it,ik) &
= phir(:,ntheta0-it-(ig-1)*itshift+1)
end do
do im = 1, (2*ntgrid+1)*(nconnect+1)
do ig = im+offset, offset+(2*ntgrid+1)*(nconnect+1)
igmod = mod(ig-offset-1,2*ntgrid+1)+1
cfnc(im,it,ik) = cfnc(im,it,ik) &
+ phiextend(ig,it,ik)*conjg(phiextend(ig-im+1,it,ik)) &
* dl_over_b(igmod)
end do
end do
cfnc(:,it,ik) = cfnc(:,it,ik) / cfnc(1,it,ik)
end do
end do
phi2extend = real(phiextend*conjg(phiextend))
deallocate (dl_over_b, phir, phiextend)
end subroutine correlation_extend