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), PUBLIC :: gradpar_z 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_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 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 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 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 gradpar_v2n_fd4 ( f, f_grad) 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) call interp_v2n_fd4(f, f_n ) call gradx_n2n_fd4(f_n, 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 end module gradients