Homec4science
Diffusion miniapp_gbs_gradients (master) Browse

R11986/gradients_mod.F90
77d2bb5ebd48master

Last Change
Raw File

gradients_mod.F90
View Options

!****************************************************
!************** FIRST DERIVATIVE ********************
!****************************************************
!-----------------------------------------------------------------------------
module gradients
use space_grid
use prec_const
real(dp), DIMENSION(1:4), public, PROTECTED :: coef_int, coef_der1_stag, coef_der1_n2n, coef_der2_stag
real(dp), DIMENSION(1:5), public, PROTECTED :: coef_der2
!real(dp), DIMENSION(:,:), pointer :: gradpar_x_n, gradpar_y_n
real(dp), dimension(:,:), allocatable :: gradpar_y_n,gradpar_x_n
real(dp), PUBLIC :: gradpar_z=0.5 ! It contains the sign of the toroidal magnetic field
contains
subroutine intialize_coefficients
coef_int = (/ -1.0_dp/16.0_dp, 9.0_dp/16.0_dp, 9.0_dp/16.0_dp, -1.0_dp/16.0_dp /)
coef_der1_stag = (/ +1.0_dp/24.0_dp, -9.0_dp/8.0_dp, 9.0_dp/8.0_dp, -1.0_dp/24.0_dp /)
coef_der2 = (/-1.0_dp/12.0_dp, 4.0_dp/3.0_dp, -5.0_dp/2.0_dp, +4.0_dp/3.0_dp, -1.0_dp/12.0_dp /)
coef_der1_n2n = (/ 1.0_dp/12.0_dp, -2.0_dp/3.0_dp, 2.0_dp/3.0_dp, -1.0_dp/12.0_dp /)
coef_der2_stag = (/ 1.0_dp/2._dp, -1.0_dp/2.0_dp, -1.0_dp/2._dp, 1.0_dp/2.0_dp /)
end subroutine intialize_coefficients
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
subroutine gradz_n2n_fd4(f , f_z)
use prec_const
IMPLICIT none
real(dp), DIMENSION(iysg:iyeg,ixsg:ixeg,izsg:izeg), INTENT(in) :: f
real(dp), DIMENSION(iysg:iyeg,ixsg:ixeg,izsg:izeg), intent(out) :: f_z
integer :: iz
real(dp), DIMENSION(1:4) :: coef_der
coef_der(:) = deltazi*coef_der1_n2n(:)
! f_z(:,:,:) = nan_
do iz=izs,ize
f_z(:, :, iz) = coef_der(1)*f(:, :, iz-2) &
+ coef_der(2)*f(:, :, iz-1) &
+ coef_der(3)*f(:, :, iz+1) &
+ coef_der(4)*f(:, :, iz+2)
end do
end subroutine gradz_n2n_fd4
subroutine gradz_n2n_fd4_omp(f , f_z)
use prec_const
implicit none
real(dp), dimension(iysg:iyeg,ixsg:ixeg,izsg:izeg), intent(in) :: f
real(dp), dimension(iysg:iyeg,ixsg:ixeg,izsg:izeg), intent(out) :: f_z
integer :: iz,iy,ix
real(dp) :: toto
real(dp), dimension(1:4) :: coef_der
coef_der(:) = deltazi*coef_der1_n2n(:)
! f_z(:,:,:) = nan_
!$omp target is_device_ptr(f_z)
!$omp target teams distribute parallel do simd collapse(3) map(from: f_z)
!$acc parallel loop collapse(3) copyout(f_z)
do iz=izs,ize
do ix=ixsg,ixeg
do iy=iysg,iyeg
toto = f_z(iy, ix, iz) * 1.2_dp
f_z(iy, ix, iz) = 1.*f_z(iy, ix, iz) ! coef_der(1)*f(iy, ix, iz-2) &
! + coef_der(2)*f(iy, ix, iz-1) &
! + coef_der(3)*f(iy, ix, iz+1) &
! + coef_der(4)*f(iy, ix, iz+2)
end do
end do
end do
!$omp end target teams distribute parallel do simd
!$omp end target
end subroutine gradz_n2n_fd4_omp
subroutine gradz_v2n_fd4(f,f_z_v2n)
use prec_const
real(dp), DIMENSION(iysg:iyeg,ixsg:ixeg,izsg:izeg), INTENT(in) :: f
real(dp), DIMENSION(iysg:iyeg,ixsg:ixeg,izsg:izeg), intent(out) :: f_z_v2n
integer :: iy, iz
real(dp), DIMENSION(1:4) :: coef_der
coef_der(:) = deltazi*coef_der1_stag(:)
! f_z_v2n(:,:,:) = nan_
do iy=iys,iye
do iz=izs,ize
f_z_v2n(iy,:,iz) = coef_int(1)*( coef_der(1)*f(iy-1, :, iz-1) &
+ coef_der(2)*f(iy-1, :, iz) &
+ coef_der(3)*f(iy-1, :, iz +1) &
+ coef_der(4)*f(iy-1, :, iz+2) )&
+ coef_int(2)*( coef_der(1)*f(iy, :, iz-1) &
+ coef_der(2)*f(iy, :, iz) &
+ coef_der(3)*f(iy, :, iz+1) &
+ coef_der(4)*f(iy, :, iz+2) )&
+ coef_int(3)*( coef_der(1)*f(iy+1, :, iz-1) &
+ coef_der(2)*f(iy+1, :, iz) &
+ coef_der(3)*f(iy+1, :, iz+1) &
+ coef_der(4)*f(iy+1, :, iz+2) )&
+ coef_int(4)*( coef_der(1)*f(iy+2, :, iz-1) &
+ coef_der(2)*f(iy+2, :, iz) &
+ coef_der(3)*f(iy+2, :, iz+1) &
+ coef_der(4)*f(iy+2, :, iz+2) )
end do
end do
end subroutine gradz_v2n_fd4
subroutine gradz_v2n_fd4_omp(f,f_z_v2n)
use prec_const
real(dp), dimension(iysg:iyeg,ixsg:ixeg,izsg:izeg), intent(in) :: f
real(dp), dimension(iysg:iyeg,ixsg:ixeg,izsg:izeg), intent(out) :: f_z_v2n
integer :: iy, iz, ix
real(dp), dimension(1:4) :: coef_der
coef_der(:) = deltazi*coef_der1_stag(:)
! f_z_v2n(:,:,:) = nan_
!$omp target teams distribute parallel do simd collapse(3) map(from: f_z_v2n)
!$acc parallel loop collapse(3) copyout(f_z_v2n)
do iz=izs,ize
do ix=ixsg,ixeg
do iy=iys,iye
f_z_v2n(iy,ix,iz) = coef_int(1)*( coef_der(1)*f(iy-1, ix, iz-1) &
+ coef_der(2)*f(iy-1, ix, iz) &
+ coef_der(3)*f(iy-1, ix, iz +1) &
+ coef_der(4)*f(iy-1, ix, iz+2) )&
+ coef_int(2)*( coef_der(1)*f(iy, ix, iz-1) &
+ coef_der(2)*f(iy, ix, iz) &
+ coef_der(3)*f(iy, ix, iz+1) &
+ coef_der(4)*f(iy, ix, iz+2) )&
+ coef_int(3)*( coef_der(1)*f(iy+1, ix, iz-1) &
+ coef_der(2)*f(iy+1, ix, iz) &
+ coef_der(3)*f(iy+1, ix, iz+1) &
+ coef_der(4)*f(iy+1, ix, iz+2) )&
+ coef_int(4)*( coef_der(1)*f(iy+2, ix, iz-1) &
+ coef_der(2)*f(iy+2, ix, iz) &
+ coef_der(3)*f(iy+2, ix, iz+1) &
+ coef_der(4)*f(iy+2, ix, iz+2) )
end do
end do
end do
!$omp end target teams distribute parallel do simd
end subroutine gradz_v2n_fd4_omp
subroutine grady_v2n_fd4(f,f_y_v2n)
use prec_const
real(dp), DIMENSION(iysg:iyeg,ixsg:ixeg,izsg:izeg), INTENT(in) :: f
real(dp), DIMENSION(iysg:iyeg,ixsg:ixeg,izsg:izeg), intent(out) :: f_y_v2n
integer :: iy, iz
real(dp), DIMENSION(1:4) :: coef_der
coef_der(:) = deltayi*coef_der1_stag(:)
! f_y_v2n(:,:,:) = nan_
do iy=iys,iye
do iz=izs,ize
f_y_v2n(iy,:,iz) = coef_int(1)*( coef_der(1)*f(iy-1, :, iz-1) &
+ coef_der(2)*f(iy, :, iz-1) &
+ coef_der(3)*f(iy+1, :, iz-1) &
+ coef_der(4)*f(iy+2, :, iz-1) )&
+ coef_int(2)*( coef_der(1)*f(iy-1, :, iz) &
+ coef_der(2)*f(iy, :, iz) &
+ coef_der(3)*f(iy+1, :, iz) &
+ coef_der(4)*f(iy+2, :, iz) )&
+ coef_int(3)*( coef_der(1)*f(iy-1, :, iz+1) &
+ coef_der(2)*f(iy, :, iz+1) &
+ coef_der(3)*f(iy+1, :, iz+1) &
+ coef_der(4)*f(iy+2, :, iz+1) )&
+ coef_int(4)*( coef_der(1)*f(iy-1, :, iz+2) &
+ coef_der(2)*f(iy, :, iz+2) &
+ coef_der(3)*f(iy+1, :, iz+2) &
+ coef_der(4)*f(iy+2, :, iz+2) )
end do
end do
end subroutine grady_v2n_fd4
subroutine grady_v2n_fd4_omp(f,f_y_v2n)
use prec_const
real(dp), dimension(iysg:iyeg,ixsg:ixeg,izsg:izeg), intent(in) :: f
real(dp), dimension(iysg:iyeg,ixsg:ixeg,izsg:izeg), intent(out) :: f_y_v2n
integer :: iy, iz, ix
real(dp), dimension(1:4) :: coef_der
coef_der(:) = deltayi*coef_der1_stag(:)
! f_y_v2n(:,:,:) = nan_
!$omp target teams distribute parallel do simd collapse(3) map(from: f_y_v2n)
!$acc parallel loop collapse(3) copyout(f_y_v2n)
do iz=izs,ize
do ix=ixsg,ixeg
do iy=iys,iye
f_y_v2n(iy,ix,iz) = coef_int(1)*( coef_der(1)*f(iy-1, ix, iz-1) &
+ coef_der(2)*f(iy, ix, iz-1) &
+ coef_der(3)*f(iy+1, ix, iz-1) &
+ coef_der(4)*f(iy+2, ix, iz-1) )&
+ coef_int(2)*( coef_der(1)*f(iy-1, ix, iz) &
+ coef_der(2)*f(iy, ix, iz) &
+ coef_der(3)*f(iy+1, ix, iz) &
+ coef_der(4)*f(iy+2, ix, iz) )&
+ coef_int(3)*( coef_der(1)*f(iy-1, ix, iz+1) &
+ coef_der(2)*f(iy, ix, iz+1) &
+ coef_der(3)*f(iy+1, ix, iz+1) &
+ coef_der(4)*f(iy+2, ix, iz+1) )&
+ coef_int(4)*( coef_der(1)*f(iy-1, ix, iz+2) &
+ coef_der(2)*f(iy, ix, iz+2) &
+ coef_der(3)*f(iy+1, ix, iz+2) &
+ coef_der(4)*f(iy+2, ix, iz+2) )
enddo
end do
end do
!$omp end target teams distribute parallel do simd
end subroutine grady_v2n_fd4_omp
subroutine interp_v2n_fd4 ( f, f_int_v2n )
use prec_const
implicit none
real(dp), dimension(iysg:iyeg,ixsg:ixeg,izsg:izeg),intent(in) :: f
real(dp), dimension(iysg:iyeg,ixsg:ixeg,izsg:izeg), intent(out) :: f_int_v2n
integer :: iy, iz
! f_int_v2n(:,:,:) = nan_
do iy=iys,iye
do iz=izs,ize
f_int_v2n(iy,:,iz) = coef_int(1)*(coef_int(1)*f(iy-1, :, iz-1) &
+ coef_int(2)*f(iy-1, : , iz) &
+ coef_int(3)*f(iy-1, :, iz +1) &
+ coef_int(4)*f(iy-1, :, iz+2) )&
+ coef_int(2)*( coef_int(1)*f(iy, :, iz-1) &
+ coef_int(2)*f(iy, :, iz) &
+ coef_int(3)*f(iy, :, iz+1) &
+ coef_int(4)*f(iy, :, iz+2) )&
+ coef_int(3)*( coef_int(1)*f(iy+1, :, iz-1) &
+ coef_int(2)*f(iy+1, :, iz) &
+ coef_int(3)*f(iy+1, :, iz+1) &
+ coef_int(4)*f(iy+1, :, iz+2) )&
+ coef_int(4)*( coef_int(1)*f(iy+2, :, iz-1) &
+ coef_int(2)*f(iy+2, :, iz) &
+ coef_int(3)*f(iy+2, :, iz+1) &
+ coef_int(4)*f(iy+2, :, iz+2) )
end do
end do
end subroutine interp_v2n_fd4
subroutine interp_v2n_fd4_omp ( f , f_int_v2n )
use prec_const
implicit none
real(dp), dimension(iysg:iyeg,ixsg:ixeg,izsg:izeg),intent(in) :: f
real(dp), dimension(iysg:iyeg,ixsg:ixeg,izsg:izeg), intent(out) :: f_int_v2n
integer :: iy, iz, ix
! f_int_v2n(:,:,:) = nan_
!$omp target teams distribute parallel do simd collapse(3) map(from: f_int_v2n)
!$acc parallel loop collapse(3) copyout(f_int_v2n)
do iz=izs,ize
do ix=ixsg,ixeg
do iy=iys,iye
f_int_v2n(iy,ix,iz) = coef_int(1)*(coef_int(1)*f(iy-1, ix, iz-1) &
+ coef_int(2)*f(iy-1, ix , iz) &
+ coef_int(3)*f(iy-1, ix, iz +1) &
+ coef_int(4)*f(iy-1, ix, iz+2) )&
+ coef_int(2)*( coef_int(1)*f(iy, ix, iz-1) &
+ coef_int(2)*f(iy, ix, iz) &
+ coef_int(3)*f(iy, ix, iz+1) &
+ coef_int(4)*f(iy, ix, iz+2) )&
+ coef_int(3)*( coef_int(1)*f(iy+1, ix, iz-1) &
+ coef_int(2)*f(iy+1, ix, iz) &
+ coef_int(3)*f(iy+1, ix, iz+1) &
+ coef_int(4)*f(iy+1, ix, iz+2) )&
+ coef_int(4)*( coef_int(1)*f(iy+2, ix, iz-1) &
+ coef_int(2)*f(iy+2, ix, iz) &
+ coef_int(3)*f(iy+2, ix, iz+1) &
+ coef_int(4)*f(iy+2, ix, iz+2) )
end do
end do
end do
!$omp end target teams distribute parallel do simd
end subroutine interp_v2n_fd4_omp
subroutine gradx_n2n_fd4(f,f_x)
use prec_const
IMPLICIT none
real(dp), DIMENSION(iysg:iyeg,ixsg:ixeg,izsg:izeg), INTENT(in) :: f
real(dp), DIMENSION(iysg:iyeg,ixsg:ixeg,izsg:izeg), intent(out) :: f_x
integer :: ix
real(dp), DIMENSION(1:4) :: coef_der
coef_der(:)= deltaxi*coef_der1_n2n(:)
! f_x (:,:,: ) = nan_
do ix=ixs,ixe
f_x(:, ix, :) = coef_der(1)*f(:, ix-2, :) &
+ coef_der(2)*f(:, ix-1, :) &
+ coef_der(3)*f(:, ix+1, :) &
+ coef_der(4)*f(:, ix+2, :)
end do
end subroutine gradx_n2n_fd4
subroutine gradx_n2n_fd4_omp(f,f_x)
use prec_const
implicit none
real(dp), dimension(iysg:iyeg,ixsg:ixeg,izsg:izeg), intent(in) :: f
real(dp), dimension(iysg:iyeg,ixsg:ixeg,izsg:izeg), intent(out) :: f_x
integer :: ix,iy,iz
real(dp), dimension(1:4) :: coef_der
coef_der(:)= deltaxi*coef_der1_n2n(:)
! f_x (:,:,: ) = nan_
!$omp target teams distribute parallel do simd collapse(3) map(from: f_x)
!$acc parallel loop collapse(3) copyout(f_x)
do iz=izsg,izeg
do ix=ixs,ixe
do iy=iysg,iyeg
f_x(iy, ix, iz) = coef_der(1)*f(iy, ix-2, iz) &
+ coef_der(2)*f(iy, ix-1, iz) &
+ coef_der(3)*f(iy, ix+1, iz) &
+ coef_der(4)*f(iy, ix+2, iz)
end do
end do
end do
!$omp end target teams distribute parallel do simd
end subroutine gradx_n2n_fd4_omp
subroutine gradpar_v2n_fd4 ( f, f_grad)
!flu and frd are useless here for periodic BC
use prec_const
IMPLICIT none
real(dp), DIMENSION(iysg:iyeg,ixsg:ixeg,izsg:izeg), intent(in) :: f
real(dp), DIMENSION(iysg:iyeg,ixsg:ixeg,izsg:izeg), intent(out) :: f_grad
real(dp), DIMENSION(iysg:iyeg,ixsg:ixeg,izsg:izeg) :: f_z, f_y, f_n, f_x_n !f_x_v
integer :: ix, iy
! f_grad(:,:,:) = nan_
call gradz_v2n_fd4(f, f_z)
call grady_v2n_fd4(f, f_y)
! its important to do the int first
call interp_v2n_fd4(f, f_n )
call gradx_n2n_fd4(f_n, f_x_n)
!call gradx_n2n_fd4(f, f_x_n)
do ix = ixs,ixe
do iy = iys, iye
f_grad(iy,ix,izs:ize) = gradpar_z*f_z(iy,ix,izs:ize) + gradpar_y_n(iy,ix)*f_y(iy,ix,izs:ize) &
+ gradpar_x_n(iy,ix)*f_x_n(iy,ix,izs:ize)
end do
end do
end subroutine gradpar_v2n_fd4
! parallel gradient for finite differences 4rth order from v grid to n grid
subroutine gradpar_v2n_fd4_omp ( f, f_grad)
!flu and frd are useless here for periodic BC
use prec_const
implicit none
real(dp), dimension(iysg:iyeg,ixsg:ixeg,izsg:izeg), intent(in) :: f
real(dp), dimension(iysg:iyeg,ixsg:ixeg,izsg:izeg), intent(out) :: f_grad
real(dp), dimension(iysg:iyeg,ixsg:ixeg,izsg:izeg) :: f_z, f_y, f_n, f_x_n !f_x_v
integer :: ix, iy,iz
!$acc enter data create(f_z,f_y,f_n,f_x_n)
!$omp target enter data map(alloc:f_z,f_y,f_n,f_x_n)
! f_grad(:,:,:) = nan_
call gradz_v2n_fd4_omp(f, f_z)
call grady_v2n_fd4_omp(f, f_y)
! its important to do the int first
call interp_v2n_fd4_omp(f, f_n )
call gradx_n2n_fd4_omp(f_n, f_x_n)
!$omp target teams distribute parallel do simd collapse(3) map(from: f_grad)
!$acc parallel loop collapse(3) copyout(f_grad)
do iz = izs,ize
do ix = ixs,ixe
do iy = iys, iye
f_grad(iy,ix,iz) = gradpar_z*f_z(iy,ix,iz) + gradpar_y_n(iy,ix)*f_y(iy,ix,iz) &
+ gradpar_x_n(iy,ix)*f_x_n(iy,ix,iz)
end do
end do
end do
!$omp end target teams distribute parallel do simd
!$omp target exit data map(delete:f_z,f_y,f_n,f_x_n)
!$acc exit data delete(f_z,f_y,f_n,f_x_n)
end subroutine gradpar_v2n_fd4_omp
end module gradients
c4science · Help