Class | sltt_extarr |
In: |
sltt/sltt_extarr.f90
|
Subroutine : | |
y_ExtLatS(jexmins:jexmaxs) : | real(DP), intent(in ) |
y_ExtLatN(jexminn:jexmaxn) : | real(DP), intent(in ) |
x_SinLonS(0:imax-1) : | real(DP), intent(in ) |
x_CosLonS(0:imax-1) : | real(DP), intent(in ) |
x_SinLonN(0:imax-1) : | real(DP), intent(in ) |
x_CosLonN(0:imax-1) : | real(DP), intent(in ) |
xyzf_QMix(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) : | real(DP), intent(in ) |
xyz_U(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in ) |
xyz_V(0:imax-1, 1:jmax, 1:kmax) : | real(DP), intent(in ) |
xyzf_ExtDQMixDLatS(iexmin:iexmax, jexmins:jexmaxs, 1:kmax, 1:ncmax) : | real(DP), intent(in ) |
xyzf_ExtDQMixDLatN(iexmin:iexmax, jexminn:jexmaxn, 1:kmax, 1:ncmax) : | real(DP), intent(in ) |
xyzf_ExtQMixS(iexmin:iexmax, jexmins:jexmaxs, 1:kmax, 1:ncmax) : | real(DP), intent(out) |
xyzf_ExtQMixN(iexmin:iexmax, jexminn:jexmaxn, 1:kmax, 1:ncmax) : | real(DP), intent(out) |
xyz_ExtUS(iexmin:iexmax, jexmins:jexmaxs, 1:kmax) : | real(DP), intent(out) |
xyz_ExtUN(iexmin:iexmax, jexminn:jexmaxn, 1:kmax) : | real(DP), intent(out) |
xyz_ExtVS(iexmin:iexmax, jexmins:jexmaxs, 1:kmax) : | real(DP), intent(out) |
xyz_ExtVN(iexmin:iexmax, jexminn:jexmaxn, 1:kmax) : | real(DP), intent(out) |
subroutine SLTTExtArrExt( y_ExtLatS, y_ExtLatN, x_SinLonS, x_CosLonS, x_SinLonN, x_CosLonN, xyzf_QMix, xyz_U, xyz_V, xyzf_ExtDQMixDLatS, xyzf_ExtDQMixDLatN, xyzf_ExtQMixS, xyzf_ExtQMixN, xyz_ExtUS, xyz_ExtUN, xyz_ExtVS, xyz_ExtVN ) use mpi_wrapper, only : nprocs, myrank, MPIWrapperISend, MPIWrapperIRecv, MPIWrapperWait use sltt_const , only : dtjw, iexmin, iexmax, jexmins, jexmaxs, jexminn, jexmaxn use sltt_lagint, only : SLTTIrrHerIntQui1DNonUni real(DP), intent(in ) :: y_ExtLatS(jexmins:jexmaxs) real(DP), intent(in ) :: y_ExtLatN(jexminn:jexmaxn) real(DP), intent(in ) :: x_SinLonS(0:imax-1) real(DP), intent(in ) :: x_CosLonS(0:imax-1) real(DP), intent(in ) :: x_SinLonN(0:imax-1) real(DP), intent(in ) :: x_CosLonN(0:imax-1) real(DP), intent(in ) :: xyzf_QMix(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) real(DP), intent(in ) :: xyz_U (0:imax-1, 1:jmax, 1:kmax) real(DP), intent(in ) :: xyz_V (0:imax-1, 1:jmax, 1:kmax) real(DP), intent(in ) :: xyzf_ExtDQMixDLatS(iexmin:iexmax, jexmins:jexmaxs, 1:kmax, 1:ncmax) real(DP), intent(in ) :: xyzf_ExtDQMixDLatN(iexmin:iexmax, jexminn:jexmaxn, 1:kmax, 1:ncmax) real(DP), intent(out) :: xyzf_ExtQMixS(iexmin:iexmax, jexmins:jexmaxs, 1:kmax, 1:ncmax) real(DP), intent(out) :: xyzf_ExtQMixN(iexmin:iexmax, jexminn:jexmaxn, 1:kmax, 1:ncmax) real(DP), intent(out) :: xyz_ExtUS (iexmin:iexmax, jexmins:jexmaxs, 1:kmax) real(DP), intent(out) :: xyz_ExtUN (iexmin:iexmax, jexminn:jexmaxn, 1:kmax) real(DP), intent(out) :: xyz_ExtVS (iexmin:iexmax, jexmins:jexmaxs, 1:kmax) real(DP), intent(out) :: xyz_ExtVN (iexmin:iexmax, jexminn:jexmaxn, 1:kmax) ! ! local variables ! !!$ ! !!$ ! variables for data transfer using MPI !!$ ! !!$ real(DP) :: sendbuf( imax, dtjw, kmax, ntrc+2 ) !!$ real(DP) :: recvbuf( imax, dtjw, kmax, ntrc+2 ) !!$ integer :: irec_send, irec_recv !!$ integer :: istatus( mpi_status_size ) !!$ integer :: ierr ! ! variables for estimation of values at poles ! real(DP) :: Ave real(DP) :: SumC real(DP) :: SumS !!$ real(DP) :: & !!$ hcpif_q_dataz( imax*kmax, 6 ), & !!$ d1_1d ( imax*kmax ) , d2_1d( imax*kmax ) !!$ integer :: mm !!$ !!$ integer :: i, j, k, m, n, nt !!$ integer :: ii !!$ !!$ !!$ imaxh = imax / 2 !!$ mm = imax * kmax ! ! initialization for debug ! !!$ ex_gq( :, :, :, : ) = 1.0d100 integer :: idest integer :: idep integer :: a_ireq_send(4) integer :: a_ireq_recv(4) real(DP), allocatable :: xyza_SendBuf(:,:,:,:,:) real(DP), allocatable :: xyza_RecvBuf(:,:,:,:,:) real(DP), allocatable :: xyz_USPrev (:,:,:) real(DP), allocatable :: xyz_UNPrev (:,:,:) real(DP), allocatable :: xyz_USNext (:,:,:) real(DP), allocatable :: xyz_UNNext (:,:,:) real(DP), allocatable :: xyz_VSPrev (:,:,:) real(DP), allocatable :: xyz_VNPrev (:,:,:) real(DP), allocatable :: xyz_VSNext (:,:,:) real(DP), allocatable :: xyz_VNNext (:,:,:) real(DP), allocatable :: xyzf_QMixSPrev(:,:,:,:) real(DP), allocatable :: xyzf_QMixNPrev(:,:,:,:) real(DP), allocatable :: xyzf_QMixSNext(:,:,:,:) real(DP), allocatable :: xyzf_QMixNNext(:,:,:,:) real(DP) :: h, theta, thetasq integer :: i integer :: j integer :: k integer :: n integer :: ii !==================================================================== allocate( xyza_SendBuf(0:imax-1, 1:dtjw, 1:kmax, 2+ncmax, 4) ) allocate( xyza_RecvBuf(0:imax-1, 1:dtjw, 1:kmax, 2+ncmax, 4) ) if ( myrank > 0 ) then allocate( xyz_USPrev (0:imax-1, 1:dtjw, 1:kmax) ) allocate( xyz_UNPrev (0:imax-1, 1:dtjw, 1:kmax) ) allocate( xyz_VSPrev (0:imax-1, 1:dtjw, 1:kmax) ) allocate( xyz_VNPrev (0:imax-1, 1:dtjw, 1:kmax) ) allocate( xyzf_QMixSPrev(0:imax-1, 1:dtjw, 1:kmax, ncmax) ) allocate( xyzf_QMixNPrev(0:imax-1, 1:dtjw, 1:kmax, ncmax) ) end if if ( myrank < (nprocs-1) ) then allocate( xyz_USNext (0:imax-1, 1:dtjw, 1:kmax) ) allocate( xyz_UNNext (0:imax-1, 1:dtjw, 1:kmax) ) allocate( xyz_VSNext (0:imax-1, 1:dtjw, 1:kmax) ) allocate( xyz_VNNext (0:imax-1, 1:dtjw, 1:kmax) ) allocate( xyzf_QMixSNext(0:imax-1, 1:dtjw, 1:kmax, ncmax) ) allocate( xyzf_QMixNNext(0:imax-1, 1:dtjw, 1:kmax, ncmax) ) end if ! y_Array(1:dtjw) values are transfered (y_SPrev). ! if ( myrank < (nprocs-1) ) then xyza_SendBuf(:,:,:,1 ,1) = xyz_U (:,1:dtjw,:) xyza_SendBuf(:,:,:,2 ,1) = xyz_V (:,1:dtjw,:) xyza_SendBuf(:,:,:,3:2+ncmax,1) = xyzf_QMix(:,1:dtjw,:,:) idest = myrank + 1 call MPIWrapperISend( idest, imax, dtjw, kmax, 2+ncmax, xyza_SendBuf(:,:,:,:,1), a_ireq_send(1) ) end if if ( myrank > 0 ) then idep = myrank - 1 call MPIWrapperIRecv( idep, imax, dtjw, kmax, 2+ncmax, xyza_RecvBuf(:,:,:,:,1), a_ireq_recv(1) ) end if ! ! y_Array(jmax/2+1-dtjw:jmax/2) values are transfered (y_SNext). ! if ( myrank > 0 ) then xyza_SendBuf(:,:,:,1 ,2) = xyz_U (:,jmax/2+1-dtjw:jmax/2,:) xyza_SendBuf(:,:,:,2 ,2) = xyz_V (:,jmax/2+1-dtjw:jmax/2,:) xyza_SendBuf(:,:,:,3:2+ncmax,2) = xyzf_QMix(:,jmax/2+1-dtjw:jmax/2,:,:) idest = myrank - 1 call MPIWrapperISend( idest, imax, dtjw, kmax, 2+ncmax, xyza_SendBuf(:,:,:,:,2), a_ireq_send(2) ) end if if ( myrank < (nprocs-1) ) then idep = myrank + 1 call MPIWrapperIRecv( idep, imax, dtjw, kmax, 2+ncmax, xyza_RecvBuf(:,:,:,:,2), a_ireq_recv(2) ) end if ! ! y_Array(jmax/2+1:jmax/2+dtjw) values are transfered (y_NNext). ! if ( myrank > 0 ) then xyza_SendBuf(:,:,:,1 ,3) = xyz_U (:,jmax/2+1:jmax/2+dtjw,:) xyza_SendBuf(:,:,:,2 ,3) = xyz_V (:,jmax/2+1:jmax/2+dtjw,:) xyza_SendBuf(:,:,:,3:2+ncmax,3) = xyzf_QMix(:,jmax/2+1:jmax/2+dtjw,:,:) idest = myrank - 1 call MPIWrapperISend( idest, imax, dtjw, kmax, 2+ncmax, xyza_SendBuf(:,:,:,:,3), a_ireq_send(3) ) end if if ( myrank < (nprocs-1) ) then idep = myrank + 1 call MPIWrapperIRecv( idep, imax, dtjw, kmax, 2+ncmax, xyza_RecvBuf(:,:,:,:,3), a_ireq_recv(3) ) end if ! ! y_Array(jmax+1-dtjw:jmax) values are transfered (y_NPrev). ! if ( myrank < (nprocs-1) ) then xyza_SendBuf(:,:,:,1 ,4) = xyz_U (:,jmax+1-dtjw:jmax,:) xyza_SendBuf(:,:,:,2 ,4) = xyz_V (:,jmax+1-dtjw:jmax,:) xyza_SendBuf(:,:,:,3:2+ncmax,4) = xyzf_QMix(:,jmax+1-dtjw:jmax,:,:) idest = myrank + 1 call MPIWrapperISend( idest, imax, dtjw, kmax, 2+ncmax, xyza_SendBuf(:,:,:,:,4), a_ireq_send(4) ) end if if ( myrank > 0 ) then idep = myrank - 1 call MPIWrapperIRecv( idep, imax, dtjw, kmax, 2+ncmax, xyza_RecvBuf(:,:,:,:,4), a_ireq_recv(4) ) end if if ( myrank < (nprocs-1) ) call MPIWrapperWait( a_ireq_send(1) ) if ( myrank > 0 ) call MPIWrapperWait( a_ireq_recv(1) ) if ( myrank > 0 ) call MPIWrapperWait( a_ireq_send(2) ) if ( myrank < (nprocs-1) ) call MPIWrapperWait( a_ireq_recv(2) ) if ( myrank > 0 ) call MPIWrapperWait( a_ireq_send(3) ) if ( myrank < (nprocs-1) ) call MPIWrapperWait( a_ireq_recv(3) ) if ( myrank < (nprocs-1) ) call MPIWrapperWait( a_ireq_send(4) ) if ( myrank > 0 ) call MPIWrapperWait( a_ireq_recv(4) ) ! y_Array(1:dtjw) values are transfered (y_SPrev). if ( myrank > 0 ) then xyz_USPrev = xyza_RecvBuf(:,:,:,1 ,1) xyz_VSPrev = xyza_RecvBuf(:,:,:,2 ,1) xyzf_QMixSPrev = xyza_RecvBuf(:,:,:,3:2+ncmax,1) end if ! y_Array(jmax/2+1-dtjw:jmax/2) values are transfered (y_SNext). if ( myrank < (nprocs-1) ) then xyz_USNext = xyza_RecvBuf(:,:,:,1 ,2) xyz_VSNext = xyza_RecvBuf(:,:,:,2 ,2) xyzf_QMixSNext = xyza_RecvBuf(:,:,:,3:2+ncmax,2) end if ! y_Array(jmax/2+1:jmax/2+dtjw) values are transfered (y_NNext). if ( myrank < (nprocs-1) ) then xyz_UNNext = xyza_RecvBuf(:,:,:,1 ,3) xyz_VNNext = xyza_RecvBuf(:,:,:,2 ,3) xyzf_QMixNNext = xyza_RecvBuf(:,:,:,3:2+ncmax,3) end if ! y_Array(jmax+1-dtjw:jmax) values are transfered (y_NPrev). if ( myrank > 0 ) then xyz_UNPrev = xyza_RecvBuf(:,:,:,1 ,4) xyz_VNPrev = xyza_RecvBuf(:,:,:,2 ,4) xyzf_QMixNPrev = xyza_RecvBuf(:,:,:,3:2+ncmax,4) end if do k = 1, kmax do j = 1, jmax/2 do i = 0, imax-1 xyz_ExtUS(i,j,k) = xyz_U(i,j ,k) xyz_ExtUN(i,j,k) = xyz_U(i,j+jmax/2,k) xyz_ExtVS(i,j,k) = xyz_V(i,j ,k) xyz_ExtVN(i,j,k) = xyz_V(i,j+jmax/2,k) end do end do end do do n = 1, ncmax do k = 1, kmax do j = 1, jmax/2 do i = 0, imax-1 xyzf_ExtQMixS(i,j,k,n) = xyzf_QMix(i,j ,k,n) xyzf_ExtQMixN(i,j,k,n) = xyzf_QMix(i,j+jmax/2,k,n) end do end do end do end do ! southern edge of southern array if( myrank == (nprocs-1) ) then ! values at south pole do k = 1, kmax do j = 1, dtjw do i = 0, imax/2-1 ii = i + imax/2 xyz_ExtUS(ii,0-j,k) = - xyz_ExtUS(i,j,k) xyz_ExtVS(ii,0-j,k) = - xyz_ExtVS(i,j,k) end do do i = imax/2, imax-1 ii = i - imax/2 xyz_ExtUS(ii,0-j,k) = - xyz_ExtUS(i,j,k) xyz_ExtVS(ii,0-j,k) = - xyz_ExtVS(i,j,k) end do end do end do do k = 1, kmax do i = 0, imax-1 ! quadratic interpolation (old) !!$ xyz_ExtUS(i,0,k) = & !!$ & a_LQIFUVSP(-1) * xyz_ExtUS(i,-1,k) & !!$ & + a_LQIFUVSP( 1) * xyz_ExtUS(i, 1,k) & !!$ & + a_LQIFUVSP( 2) * xyz_ExtUS(i, 2,k) !!$ xyz_ExtVS(i,0,k) = & !!$ & a_LQIFUVSP(-1) * xyz_ExtVS(i,-1,k) & !!$ & + a_LQIFUVSP( 1) * xyz_ExtVS(i, 1,k) & !!$ & + a_LQIFUVSP( 2) * xyz_ExtVS(i, 2,k) ! cubic interpolation xyz_ExtUS(i,0,k) = a_LCIFUVSP(-2) * xyz_ExtUS(i,-2,k) + a_LCIFUVSP(-1) * xyz_ExtUS(i,-1,k) + a_LCIFUVSP( 1) * xyz_ExtUS(i, 1,k) + a_LCIFUVSP( 2) * xyz_ExtUS(i, 2,k) xyz_ExtVS(i,0,k) = a_LCIFUVSP(-2) * xyz_ExtVS(i,-2,k) + a_LCIFUVSP(-1) * xyz_ExtVS(i,-1,k) + a_LCIFUVSP( 1) * xyz_ExtVS(i, 1,k) + a_LCIFUVSP( 2) * xyz_ExtVS(i, 2,k) end do end do do n = 1, ncmax do k = 1, kmax do j = 1, dtjw do i = 0, imax/2-1 ii = i + imax/2 xyzf_ExtQMixS(ii,0-j,k,n) = xyzf_ExtQMixS(i,j,k,n) end do do i = imax/2, imax-1 ii = i - imax/2 xyzf_ExtQMixS(ii,0-j,k,n) = xyzf_ExtQMixS(i,j,k,n) end do end do end do do k = 1, kmax do i = 0, imax-1 ! quadratic interpolation (old) !!$ xyzf_ExtQMixS(i,0,k,n) = & !!$ & a_LQIFUVSP(-1) * xyzf_ExtQMixS(i,-1,k,n) & !!$ & + a_LQIFUVSP( 1) * xyzf_ExtQMixS(i, 1,k,n) & !!$ & + a_LQIFUVSP( 2) * xyzf_ExtQMixS(i, 2,k,n) !!$ ! cubic interpolation !!$ xyzf_ExtQMixS(i,0,k,n) = & !!$ & a_LCIFUVSP(-2) * xyzf_ExtQMixS(i,-2,k,n) & !!$ & + a_LCIFUVSP(-1) * xyzf_ExtQMixS(i,-1,k,n) & !!$ & + a_LCIFUVSP( 1) * xyzf_ExtQMixS(i, 1,k,n) & !!$ & + a_LCIFUVSP( 2) * xyzf_ExtQMixS(i, 2,k,n) !!$ ! Polar value estimated by 1D Hermite Quintic Interpolation xyzf_ExtQMixS(i,0,k,n) = SLTTIrrHerIntQui1DNonUni ( xyzf_ExtQMixS(i,-2,k,n), xyzf_ExtQMixS(i,-1,k,n), xyzf_ExtQMixS(i, 1,k,n), xyzf_ExtQMixS(i, 2,k,n), xyzf_ExtDQMixDLatS(i,-1,k,n), xyzf_ExtDQMixDLatS(i, 1,k,n), y_ExtLatS(-2)-y_ExtLatS(-1), y_ExtLatS( 1)-y_ExtLatS(-1), y_ExtLatS( 2)-y_ExtLatS(-1), y_ExtLatS( 0)-y_ExtLatS(-1) ) end do end do end do ! only wavenumber 1 component is retained for zonal and meridional ! wind velocities j = 0 do k = 1, kmax SumC = 0.0_DP SumS = 0.0_DP do i = 0, imax-1 SumC = SumC + xyz_ExtUS(i,j,k) * x_CosLonS(i) SumS = SumS + xyz_ExtUS(i,j,k) * x_SinLonS(i) end do SumC = SumC / SumSinSq SumS = SumS / SumSinSq do i = 0, imax-1 xyz_ExtUS(i,j,k) = SumC * x_CosLonS(i) + SumS * x_SinLonS(i) end do end do do k = 1, kmax SumC = 0.0_DP SumS = 0.0_DP do i = 0, imax-1 SumC = SumC + xyz_ExtVS(i,j,k) * x_CosLonS(i) SumS = SumS + xyz_ExtVS(i,j,k) * x_SinLonS(i) end do SumC = SumC / SumSinSq SumS = SumS / SumSinSq do i = 0, imax-1 xyz_ExtVS(i,j,k) = SumC * x_CosLonS(i) + SumS * x_SinLonS(i) end do end do ! zonal average is set for mixing ratio j = 0 do n = 1, ncmax do k = 1, kmax Ave = 0.0_DP do i = 0, imax-1 Ave = Ave + xyzf_ExtQMixS(i,j,k,n) end do Ave = Ave / dble( imax ) do i = 0, imax-1 xyzf_ExtQMixS(i,j,k,n) = Ave end do end do end do else do k = 1, kmax do j = 1, dtjw do i = 0, imax-1 xyz_ExtUS(i,1-j,k) = xyz_USNext(i,dtjw-(j-1),k) xyz_ExtVS(i,1-j,k) = xyz_VSNext(i,dtjw-(j-1),k) end do end do end do do n = 1, ncmax do k = 1, kmax do j = 1, dtjw do i = 0, imax-1 xyzf_ExtQMixS(i,1-j,k,n) = xyzf_QMixSNext(i,dtjw-(j-1),k,n) end do end do end do end do end if ! northern edge of southern array if ( myrank == 0 ) then do k = 1, kmax do j = 1, dtjw do i = 0, imax-1 xyz_ExtUS(i,jmax/2+j,k) = xyz_ExtUN(i,j,k) xyz_ExtVS(i,jmax/2+j,k) = xyz_ExtVN(i,j,k) end do end do end do do n = 1, ncmax do k = 1, kmax do j = 1, dtjw do i = 0, imax-1 xyzf_ExtQMixS(i,jmax/2+j,k,n) = xyzf_ExtQMixN(i,j,k,n) end do end do end do end do else do k = 1, kmax do j = 1, dtjw do i = 0, imax-1 xyz_ExtUS(i,jmax/2+j,k) = xyz_USPrev(i,j,k) xyz_ExtVS(i,jmax/2+j,k) = xyz_VSPrev(i,j,k) end do end do end do do n = 1, ncmax do k = 1, kmax do j = 1, dtjw do i = 0, imax-1 xyzf_ExtQMixS(i,jmax/2+j,k,n) = xyzf_QMixSPrev(i,j,k,n) end do end do end do end do end if ! ! southern edge of northern array if ( myrank == 0 ) then do k = 1, kmax do j = 1, dtjw do i = 0, imax-1 xyz_ExtUN(i,1-j,k) = xyz_ExtUS(i,jmax/2-(j-1),k) xyz_ExtVN(i,1-j,k) = xyz_ExtVS(i,jmax/2-(j-1),k) end do end do end do do n = 1, ncmax do k = 1, kmax do j = 1, dtjw do i = 0, imax-1 xyzf_ExtQMixN(i,1-j,k,n) = xyzf_ExtQMixS(i,jmax/2-(j-1),k,n) end do end do end do end do else do k = 1, kmax do j = 1, dtjw do i = 0, imax-1 xyz_ExtUN(i,1-j,k) = xyz_UNPrev(i,dtjw-(j-1),k) xyz_ExtVN(i,1-j,k) = xyz_VNPrev(i,dtjw-(j-1),k) end do end do end do do n = 1, ncmax do k = 1, kmax do j = 1, dtjw do i = 0, imax-1 xyzf_ExtQMixN(i,1-j,k,n) = xyzf_QMixNPrev(i,dtjw-(j-1),k,n) end do end do end do end do end if ! northern edge of northern array if( myrank == (nprocs-1) ) then ! values at north pole do k = 1, kmax do j = 1, dtjw do i = 0, imax/2-1 ii = i + imax/2 xyz_ExtUN(ii,jmax/2+1+j,k) = - xyz_ExtUN(i,jmax/2+1-j,k) xyz_ExtVN(ii,jmax/2+1+j,k) = - xyz_ExtVN(i,jmax/2+1-j,k) end do do i = imax/2, imax-1 ii = i - imax/2 xyz_ExtUN(ii,jmax/2+1+j,k) = - xyz_ExtUN(i,jmax/2+1-j,k) xyz_ExtVN(ii,jmax/2+1+j,k) = - xyz_ExtVN(i,jmax/2+1-j,k) end do end do end do do k = 1, kmax do i = 0, imax-1 ! quadratic interpolation (old) !!$ xyz_ExtUN(i,jmax/2+1,k) = & !!$ & a_LQIFUVNP(jmax/2-1) * xyz_ExtUN(i,jmax/2-1,k) & !!$ & + a_LQIFUVNP(jmax/2 ) * xyz_ExtUN(i,jmax/2 ,k) & !!$ & + a_LQIFUVNP(jmax/2+2) * xyz_ExtUN(i,jmax/2+2,k) !!$ xyz_ExtVN(i,jmax/2+1,k) = & !!$ & a_LQIFUVNP(jmax/2-1) * xyz_ExtVN(i,jmax/2-1,k) & !!$ & + a_LQIFUVNP(jmax/2 ) * xyz_ExtVN(i,jmax/2 ,k) & !!$ & + a_LQIFUVNP(jmax/2+2) * xyz_ExtVN(i,jmax/2+2,k) ! qcubic interpolation xyz_ExtUN(i,jmax/2+1,k) = a_LCIFUVNP(jmax/2-1) * xyz_ExtUN(i,jmax/2-1,k) + a_LCIFUVNP(jmax/2 ) * xyz_ExtUN(i,jmax/2 ,k) + a_LCIFUVNP(jmax/2+2) * xyz_ExtUN(i,jmax/2+2,k) + a_LCIFUVNP(jmax/2+3) * xyz_ExtUN(i,jmax/2+3,k) xyz_ExtVN(i,jmax/2+1,k) = a_LCIFUVNP(jmax/2-1) * xyz_ExtVN(i,jmax/2-1,k) + a_LCIFUVNP(jmax/2 ) * xyz_ExtVN(i,jmax/2 ,k) + a_LCIFUVNP(jmax/2+2) * xyz_ExtVN(i,jmax/2+2,k) + a_LCIFUVNP(jmax/2+3) * xyz_ExtVN(i,jmax/2+3,k) end do end do do n = 1, ncmax do k = 1, kmax do j = 1, dtjw do i = 0, imax/2-1 ii = i + imax/2 xyzf_ExtQMixN(ii,jmax/2+1+j,k,n) = xyzf_ExtQMixN(i,jmax/2+1-j,k,n) end do do i = imax/2, imax-1 ii = i - imax/2 xyzf_ExtQMixN(ii,jmax/2+1+j,k,n) = xyzf_ExtQMixN(i,jmax/2+1-j,k,n) end do end do end do do k = 1, kmax do i = 0, imax-1 !!$ ! quadratic interpolation (old) !!$ xyzf_ExtQMixN(i,jmax/2+1,k,n) = & !!$ & a_LQIFUVNP(jmax/2-1) * xyzf_ExtQMixN(i,jmax/2-1,k,n) & !!$ & + a_LQIFUVNP(jmax/2 ) * xyzf_ExtQMixN(i,jmax/2 ,k,n) & !!$ & + a_LQIFUVNP(jmax/2+2) * xyzf_ExtQMixN(i,jmax/2+2,k,n) !!$ ! cubic interpolation !!$ xyzf_ExtQMixN(i,jmax/2+1,k,n) = & !!$ & a_LCIFUVNP(jmax/2-1) * xyzf_ExtQMixN(i,jmax/2-1,k,n) & !!$ & + a_LCIFUVNP(jmax/2 ) * xyzf_ExtQMixN(i,jmax/2 ,k,n) & !!$ & + a_LCIFUVNP(jmax/2+2) * xyzf_ExtQMixN(i,jmax/2+2,k,n) & !!$ & + a_LCIFUVNP(jmax/2+3) * xyzf_ExtQMixN(i,jmax/2+3,k,n) !!$ ! Polar value estimated by 1D Hermite Quintic Interpolation xyzf_ExtQMixN(i,jmax/2+1,k,n) = SLTTIrrHerIntQui1DNonUni ( xyzf_ExtQMixN(i,jmax/2-1,k,n), xyzf_ExtQMixN(i,jmax/2 ,k,n), xyzf_ExtQMixN(i,jmax/2+2,k,n), xyzf_ExtQMixN(i,jmax/2+3,k,n), xyzf_ExtDQMixDLatN(i,jmax/2 ,k,n), xyzf_ExtDQMixDLatN(i,jmax/2+2,k,n), y_ExtLatN(jmax/2-1)-y_ExtLatN(jmax/2 ), y_ExtLatN(jmax/2+2)-y_ExtLatN(jmax/2 ), y_ExtLatN(jmax/2+3)-y_ExtLatN(jmax/2 ), y_ExtLatN(jmax/2+1)-y_ExtLatN(jmax/2 ) ) end do end do end do ! only wavenumber 1 component is retained for zonal and meridional ! wind velocities j = jmax/2+1 do k = 1, kmax SumC = 0.0_DP SumS = 0.0_DP do i = 0, imax-1 SumC = SumC + xyz_ExtUN(i,j,k) * x_CosLonN(i) SumS = SumS + xyz_ExtUN(i,j,k) * x_SinLonN(i) end do SumC = SumC / SumSinSq SumS = SumS / SumSinSq do i = 0, imax-1 xyz_ExtUN(i,j,k) = SumC * x_CosLonN(i) + SumS * x_SinLonN(i) end do end do do k = 1, kmax SumC = 0.0_DP SumS = 0.0_DP do i = 0, imax-1 SumC = SumC + xyz_ExtVN(i,j,k) * x_CosLonN(i) SumS = SumS + xyz_ExtVN(i,j,k) * x_SinLonN(i) end do SumC = SumC / SumSinSq SumS = SumS / SumSinSq do i = 0, imax-1 xyz_ExtVN(i,j,k) = SumC * x_CosLonN(i) + SumS * x_SinLonN(i) end do end do ! zonal average is set for mixing ratio j = jmax/2+1 do n = 1, ncmax do k = 1, kmax Ave = 0.0_DP do i = 0, imax-1 Ave = Ave + xyzf_ExtQMixN(i,j,k,n) end do Ave = Ave / dble( imax ) do i = 0, imax-1 xyzf_ExtQMixN(i,j,k,n) = Ave end do end do end do else do k = 1, kmax do j = 1, dtjw do i = 0, imax-1 xyz_ExtUN(i,jmax/2+j,k) = xyz_UNNext(i,j,k) xyz_ExtVN(i,jmax/2+j,k) = xyz_VNNext(i,j,k) end do end do end do do n = 1, ncmax do k = 1, kmax do j = 1, dtjw do i = 0, imax-1 xyzf_ExtQMixN(i,jmax/2+j,k,n) = xyzf_QMixNNext(i,j,k,n) end do end do end do end do end if deallocate( xyza_SendBuf ) deallocate( xyza_RecvBuf ) if ( myrank > 0 ) then deallocate( xyz_USPrev ) deallocate( xyz_UNPrev ) deallocate( xyz_VSPrev ) deallocate( xyz_VNPrev ) deallocate( xyzf_QMixSPrev ) deallocate( xyzf_QMixNPrev ) end if if ( myrank < (nprocs-1) ) then deallocate( xyz_USNext ) deallocate( xyz_UNNext ) deallocate( xyz_VSNext ) deallocate( xyz_VNNext ) deallocate( xyzf_QMixSNext ) deallocate( xyzf_QMixNNext ) end if !=========================================== ! set values at longitudinal edge !------------------------------------------- do k = 1, kmax do j = jexmins, jexmaxs do i = iexmin, 0-1 xyz_ExtUS(i,j,k) = xyz_ExtUS(imax+i,j,k) xyz_ExtVS(i,j,k) = xyz_ExtVS(imax+i,j,k) end do do i = imax-1+1, iexmax xyz_ExtUS(i,j,k) = xyz_ExtUS(i-imax,j,k) xyz_ExtVS(i,j,k) = xyz_ExtVS(i-imax,j,k) end do end do do j = jexminn, jexmaxn do i = iexmin, 0-1 xyz_ExtUN(i,j,k) = xyz_ExtUN(imax+i,j,k) xyz_ExtVN(i,j,k) = xyz_ExtVN(imax+i,j,k) end do do i = imax-1+1, iexmax xyz_ExtUN(i,j,k) = xyz_ExtUN(i-imax,j,k) xyz_ExtVN(i,j,k) = xyz_ExtVN(i-imax,j,k) end do end do end do do n = 1, ncmax do k = 1, kmax do j = jexmins, jexmaxs do i = iexmin, 0-1 xyzf_ExtQMixS(i,j,k,n) = xyzf_ExtQMixS(imax+i,j,k,n) end do do i = imax-1+1, iexmax xyzf_ExtQMixS(i,j,k,n) = xyzf_ExtQMixS(i-imax,j,k,n) end do end do do j = jexminn, jexmaxn do i = iexmin, 0-1 xyzf_ExtQMixN(i,j,k,n) = xyzf_ExtQMixN(imax+i,j,k,n) end do do i = imax-1+1, iexmax xyzf_ExtQMixN(i,j,k,n) = xyzf_ExtQMixN(i-imax,j,k,n) end do end do end do end do end subroutine SLTTExtArrExt
Subroutine : | |||
x_SinLonS(0:imax-1) : | real(DP), intent(in ) | ||
x_CosLonS(0:imax-1) : | real(DP), intent(in ) | ||
x_SinLonN(0:imax-1) : | real(DP), intent(in ) | ||
x_CosLonN(0:imax-1) : | real(DP), intent(in ) | ||
xyzf_QMix(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) : | real(DP), intent(in ) | ||
pm : | real(DP), intent(in )
| ||
xyzf_ExtQMixS(iexmin:iexmax, jexmins:jexmaxs, 1:kmax, 1:ncmax) : | real(DP), intent(out) | ||
xyzf_ExtQMixN(iexmin:iexmax, jexminn:jexmaxn, 1:kmax, 1:ncmax) : | real(DP), intent(out) | ||
PoleMethod : | character(*), intent(in )
|
subroutine SLTTExtArrExt2( x_SinLonS, x_CosLonS, x_SinLonN, x_CosLonN, xyzf_QMix, pm, xyzf_ExtQMixS, xyzf_ExtQMixN, PoleMethod ) ! メッセージ出力 ! Message output ! use dc_message, only: MessageNotify use mpi_wrapper, only : nprocs, myrank, MPIWrapperISend, MPIWrapperIRecv, MPIWrapperWait use sltt_const , only : dtjw, iexmin, iexmax, jexmins, jexmaxs, jexminn, jexmaxn real(DP), intent(in ) :: x_SinLonS(0:imax-1) real(DP), intent(in ) :: x_CosLonS(0:imax-1) real(DP), intent(in ) :: x_SinLonN(0:imax-1) real(DP), intent(in ) :: x_CosLonN(0:imax-1) real(DP), intent(in ) :: xyzf_QMix(0:imax-1, 1:jmax, 1:kmax, 1:ncmax) real(DP), intent(in ) :: pm ! 配列拡張する際、極ごえ後に符号が変わる場合は -1.0を与える。そうでない場合は1.0を与える。 real(DP), intent(out) :: xyzf_ExtQMixS(iexmin:iexmax, jexmins:jexmaxs, 1:kmax, 1:ncmax) real(DP), intent(out) :: xyzf_ExtQMixN(iexmin:iexmax, jexminn:jexmaxn, 1:kmax, 1:ncmax) character(*), intent(in ) :: PoleMethod ! "Mean" : Longitudinal mean ! "Wave1" : Only wave #1 is retained. ! ! local variables ! ! ! variables for estimation of values at poles ! real(DP) :: Ave real(DP) :: SumC real(DP) :: SumS integer :: idest integer :: idep integer :: a_ireq_send(4) integer :: a_ireq_recv(4) real(DP), allocatable :: xyza_SendBuf(:,:,:,:,:) real(DP), allocatable :: xyza_RecvBuf(:,:,:,:,:) real(DP), allocatable :: xyzf_QMixSPrev(:,:,:,:) real(DP), allocatable :: xyzf_QMixNPrev(:,:,:,:) real(DP), allocatable :: xyzf_QMixSNext(:,:,:,:) real(DP), allocatable :: xyzf_QMixNNext(:,:,:,:) real(DP) :: h, theta, thetasq integer :: i integer :: j integer :: k integer :: n integer :: ii !==================================================================== allocate( xyza_SendBuf(0:imax-1, 1:dtjw, 1:kmax, ncmax, 4) ) allocate( xyza_RecvBuf(0:imax-1, 1:dtjw, 1:kmax, ncmax, 4) ) if ( myrank > 0 ) then allocate( xyzf_QMixSPrev(0:imax-1, 1:dtjw, 1:kmax, ncmax) ) allocate( xyzf_QMixNPrev(0:imax-1, 1:dtjw, 1:kmax, ncmax) ) end if if ( myrank < (nprocs-1) ) then allocate( xyzf_QMixSNext(0:imax-1, 1:dtjw, 1:kmax, ncmax) ) allocate( xyzf_QMixNNext(0:imax-1, 1:dtjw, 1:kmax, ncmax) ) end if ! y_Array(1:dtjw) values are transfered (y_SPrev). ! if ( myrank < (nprocs-1) ) then xyza_SendBuf(:,:,:,:,1) = xyzf_QMix(:,1:dtjw,:,:) idest = myrank + 1 call MPIWrapperISend( idest, imax, dtjw, kmax, ncmax, xyza_SendBuf(:,:,:,:,1), a_ireq_send(1) ) end if if ( myrank > 0 ) then idep = myrank - 1 call MPIWrapperIRecv( idep, imax, dtjw, kmax, ncmax, xyza_RecvBuf(:,:,:,:,1), a_ireq_recv(1) ) end if ! ! y_Array(jmax/2+1-dtjw:jmax/2) values are transfered (y_SNext). ! if ( myrank > 0 ) then xyza_SendBuf(:,:,:,:,2) = xyzf_QMix(:,jmax/2+1-dtjw:jmax/2,:,:) idest = myrank - 1 call MPIWrapperISend( idest, imax, dtjw, kmax, ncmax, xyza_SendBuf(:,:,:,:,2), a_ireq_send(2) ) end if if ( myrank < (nprocs-1) ) then idep = myrank + 1 call MPIWrapperIRecv( idep, imax, dtjw, kmax, ncmax, xyza_RecvBuf(:,:,:,:,2), a_ireq_recv(2) ) end if ! ! y_Array(jmax/2+1:jmax/2+dtjw) values are transfered (y_NNext). ! if ( myrank > 0 ) then xyza_SendBuf(:,:,:,:,3) = xyzf_QMix(:,jmax/2+1:jmax/2+dtjw,:,:) idest = myrank - 1 call MPIWrapperISend( idest, imax, dtjw, kmax, ncmax, xyza_SendBuf(:,:,:,:,3), a_ireq_send(3) ) end if if ( myrank < (nprocs-1) ) then idep = myrank + 1 call MPIWrapperIRecv( idep, imax, dtjw, kmax, ncmax, xyza_RecvBuf(:,:,:,:,3), a_ireq_recv(3) ) end if ! ! y_Array(jmax+1-dtjw:jmax) values are transfered (y_NPrev). ! if ( myrank < (nprocs-1) ) then xyza_SendBuf(:,:,:,:,4) = xyzf_QMix(:,jmax+1-dtjw:jmax,:,:) idest = myrank + 1 call MPIWrapperISend( idest, imax, dtjw, kmax, ncmax, xyza_SendBuf(:,:,:,:,4), a_ireq_send(4) ) end if if ( myrank > 0 ) then idep = myrank - 1 call MPIWrapperIRecv( idep, imax, dtjw, kmax, ncmax, xyza_RecvBuf(:,:,:,:,4), a_ireq_recv(4) ) end if if ( myrank < (nprocs-1) ) call MPIWrapperWait( a_ireq_send(1) ) if ( myrank > 0 ) call MPIWrapperWait( a_ireq_recv(1) ) if ( myrank > 0 ) call MPIWrapperWait( a_ireq_send(2) ) if ( myrank < (nprocs-1) ) call MPIWrapperWait( a_ireq_recv(2) ) if ( myrank > 0 ) call MPIWrapperWait( a_ireq_send(3) ) if ( myrank < (nprocs-1) ) call MPIWrapperWait( a_ireq_recv(3) ) if ( myrank < (nprocs-1) ) call MPIWrapperWait( a_ireq_send(4) ) if ( myrank > 0 ) call MPIWrapperWait( a_ireq_recv(4) ) ! y_Array(1:dtjw) values are transfered (y_SPrev). if ( myrank > 0 ) then xyzf_QMixSPrev = xyza_RecvBuf(:,:,:,:,1) end if ! y_Array(jmax/2+1-dtjw:jmax/2) values are transfered (y_SNext). if ( myrank < (nprocs-1) ) then xyzf_QMixSNext = xyza_RecvBuf(:,:,:,:,2) end if ! y_Array(jmax/2+1:jmax/2+dtjw) values are transfered (y_NNext). if ( myrank < (nprocs-1) ) then xyzf_QMixNNext = xyza_RecvBuf(:,:,:,:,3) end if ! y_Array(jmax+1-dtjw:jmax) values are transfered (y_NPrev). if ( myrank > 0 ) then xyzf_QMixNPrev = xyza_RecvBuf(:,:,:,:,4) end if do n = 1, ncmax do k = 1, kmax do j = 1, jmax/2 do i = 0, imax-1 xyzf_ExtQMixS(i,j,k,n) = xyzf_QMix(i,j ,k,n) xyzf_ExtQMixN(i,j,k,n) = xyzf_QMix(i,j+jmax/2,k,n) end do end do end do end do ! southern edge of southern array if( myrank == (nprocs-1) ) then ! values at south pole do n = 1, ncmax do k = 1, kmax do j = 1, dtjw do i = 0, imax/2-1 ii = i + imax/2 xyzf_ExtQMixS(ii,0-j,k,n) = pm * xyzf_ExtQMixS(i,j,k,n) end do do i = imax/2, imax-1 ii = i - imax/2 xyzf_ExtQMixS(ii,0-j,k,n) = pm * xyzf_ExtQMixS(i,j,k,n) end do end do end do do k = 1, kmax do i = 0, imax-1 !!$ xyzf_ExtQMixS(i,0,k,n) = & !南極点の値 !!$ & a_LQIFUVSP(-1) * xyzf_ExtQMixS(i,-1,k,n) & !!$ & + a_LQIFUVSP( 1) * xyzf_ExtQMixS(i, 1,k,n) & !!$ & + a_LQIFUVSP( 2) * xyzf_ExtQMixS(i, 2,k,n) xyzf_ExtQMixS(i,0,k,n) = a_LCIFUVSP(-2) * xyzf_ExtQMixS(i,-2,k,n) + a_LCIFUVSP(-1) * xyzf_ExtQMixS(i,-1,k,n) + a_LCIFUVSP( 1) * xyzf_ExtQMixS(i, 1,k,n) + a_LCIFUVSP( 2) * xyzf_ExtQMixS(i, 2,k,n) end do end do end do select case ( PoleMethod ) case ( "Mean" ) ! Longitudinal mean j = 0 do n = 1, ncmax do k = 1, kmax Ave = 0.0_DP do i = 0, imax-1 Ave = Ave + xyzf_ExtQMixS(i,j,k,n) end do Ave = Ave / dble( imax ) do i = 0, imax-1 xyzf_ExtQMixS(i,j,k,n) = Ave !!南極点の値を各iで統一する。 end do end do end do case ( "Wave1" ) ! Only wave #1 is retained. j = 0 do n = 1, ncmax do k = 1, kmax SumC = 0.0_DP SumS = 0.0_DP do i = 0, imax-1 SumC = SumC + xyzf_ExtQMixS(i,j,k,n) * x_CosLonS(i) SumS = SumS + xyzf_ExtQMixS(i,j,k,n) * x_SinLonS(i) end do SumC = SumC / SumSinSq SumS = SumS / SumSinSq do i = 0, imax-1 xyzf_ExtQMixS(i,j,k,n) = SumC * x_CosLonS(i) + SumS * x_SinLonS(i) end do end do end do case default call MessageNotify( 'E', module_name, 'PoleMethod of %c is inappropriate.', c1 = trim(PoleMethod) ) end select else !!$ do j = 1, jew !!$ y_ExtLatS(1-j) = y_LatSNext(jew-(j-1)) !!$ end do do n = 1, ncmax do k = 1, kmax do j = 1, dtjw do i = 0, imax-1 xyzf_ExtQMixS(i,1-j,k,n) = xyzf_QMixSNext(i,dtjw-(j-1),k,n) end do end do end do end do end if ! northern edge of southern array if ( myrank == 0 ) then do n = 1, ncmax do k = 1, kmax do j = 1, dtjw do i = 0, imax-1 xyzf_ExtQMixS(i,jmax/2+j,k,n) = xyzf_ExtQMixN(i,j,k,n) end do end do end do end do else do n = 1, ncmax do k = 1, kmax do j = 1, dtjw do i = 0, imax-1 xyzf_ExtQMixS(i,jmax/2+j,k,n) = xyzf_QMixSPrev(i,j,k,n) end do end do end do end do end if ! ! southern edge of northern array if ( myrank == 0 ) then do n = 1, ncmax do k = 1, kmax do j = 1, dtjw do i = 0, imax-1 xyzf_ExtQMixN(i,-j+1,k,n) = xyzf_ExtQMixS(i,jmax/2-(j-1),k,n) end do end do end do end do else do n = 1, ncmax do k = 1, kmax do j = 1, dtjw do i = 0, imax-1 xyzf_ExtQMixN(i,-j+1,k,n) = xyzf_QMixNPrev(i,dtjw-(j-1),k,n) end do end do end do end do end if ! northern edge of northern array if( myrank == (nprocs-1) ) then ! values at north pole do n = 1, ncmax do k = 1, kmax do j = 1, dtjw do i = 0, imax/2-1 ii = i + imax/2 xyzf_ExtQMixN(ii,jmax/2+1+j,k,n) = pm * xyzf_ExtQMixN(i,jmax/2+1-j,k,n) end do do i = imax/2, imax-1 ii = i - imax/2 xyzf_ExtQMixN(ii,jmax/2+1+j,k,n) = pm * xyzf_ExtQMixN(i,jmax/2+1-j,k,n) end do end do end do do k = 1, kmax do i = 0, imax-1 !!$ xyzf_ExtQMixN(i,jmax/2+1,k,n) = & !北極点での値 !!$ & a_LQIFUVNP(jmax/2-1) * xyzf_ExtQMixN(i,jmax/2-1,k,n) & !!$ & + a_LQIFUVNP(jmax/2 ) * xyzf_ExtQMixN(i,jmax/2 ,k,n) & !!$ & + a_LQIFUVNP(jmax/2+2) * xyzf_ExtQMixN(i,jmax/2+2,k,n) xyzf_ExtQMixN(i,jmax/2+1,k,n) = a_LCIFUVNP(jmax/2-1) * xyzf_ExtQMixN(i,jmax/2-1,k,n) + a_LCIFUVNP(jmax/2 ) * xyzf_ExtQMixN(i,jmax/2 ,k,n) + a_LCIFUVNP(jmax/2+2) * xyzf_ExtQMixN(i,jmax/2+2,k,n) + a_LCIFUVNP(jmax/2+3) * xyzf_ExtQMixN(i,jmax/2+3,k,n) end do end do end do select case ( PoleMethod ) case ( "Mean" ) ! Longitudinal mean j = jmax/2+1 do n = 1, ncmax do k = 1, kmax Ave = 0.0_DP do i = 0, imax-1 Ave = Ave + xyzf_ExtQMixN(i,j,k,n) end do Ave = Ave / dble( imax ) do i = 0, imax-1 xyzf_ExtQMixN(i,j,k,n) = Ave !値を各iで統一する。 end do end do end do case ( "Wave1" ) ! Only wave #1 is retained. j = jmax/2+1 do n = 1, ncmax do k = 1, kmax SumC = 0.0_DP SumS = 0.0_DP do i = 0, imax-1 SumC = SumC + xyzf_ExtQMixN(i,j,k,n) * x_CosLonN(i) SumS = SumS + xyzf_ExtQMixN(i,j,k,n) * x_SinLonN(i) end do SumC = SumC / SumSinSq SumS = SumS / SumSinSq do i = 0, imax-1 xyzf_ExtQMixN(i,j,k,n) = SumC * x_CosLonN(i) + SumS * x_SinLonN(i) end do end do end do case default call MessageNotify( 'E', module_name, 'PoleMethod of %c is inappropriate.', c1 = trim(PoleMethod) ) end select else do n = 1, ncmax do k = 1, kmax do j = 1, dtjw do i = 0, imax-1 xyzf_ExtQMixN(i,jmax/2+j,k,n) = xyzf_QMixNNext(i,j,k,n) end do end do end do end do end if deallocate( xyza_SendBuf ) deallocate( xyza_RecvBuf ) if ( myrank > 0 ) then deallocate( xyzf_QMixSPrev ) deallocate( xyzf_QMixNPrev ) end if if ( myrank < (nprocs-1) ) then deallocate( xyzf_QMixSNext ) deallocate( xyzf_QMixNNext ) end if !=========================================== ! set values at longitudinal edge !------------------------------------------- do n = 1, ncmax do k = 1, kmax do j = jexmins, jexmaxs do i = iexmin, 0-1 xyzf_ExtQMixS(i,j,k,n) = xyzf_ExtQMixS(imax+i,j,k,n) end do do i = imax-1+1, iexmax xyzf_ExtQMixS(i,j,k,n) = xyzf_ExtQMixS(i-imax,j,k,n) end do end do do j = jexminn, jexmaxn do i = iexmin, 0-1 xyzf_ExtQMixN(i,j,k,n) = xyzf_ExtQMixN(imax+i,j,k,n) end do do i = imax-1+1, iexmax xyzf_ExtQMixN(i,j,k,n) = xyzf_ExtQMixN(i-imax,j,k,n) end do end do end do end do end subroutine SLTTExtArrExt2
Subroutine : | |
x_LonS( 0:imax-1 ) : | real(DP), intent(in ) |
y_LatS( 1:jmax/2 ) : | real(DP), intent(in ) |
x_LonN( 0:imax-1 ) : | real(DP), intent(in ) |
y_LatN( 1:jmax/2 ) : | real(DP), intent(in ) |
x_ExtLonS(iexmin :iexmax ) : | real(DP), intent(out) |
y_ExtLatS(jexmins:jexmaxs) : | real(DP), intent(out) |
x_ExtLonN(iexmin :iexmax ) : | real(DP), intent(out) |
y_ExtLatN(jexminn:jexmaxn) : | real(DP), intent(out) |
subroutine SLTTExtArrInit( x_LonS, y_LatS, x_LonN, y_LatN, x_ExtLonS, y_ExtLatS, x_ExtLonN, y_ExtLatN ) ! ! MPI ! use mpi_wrapper , only : myrank, nprocs, MPIWrapperISend, MPIWrapperIRecv, MPIWrapperWait use constants0, only : PI use axesset , only : y_Lat use sltt_const, only : PIx2, PIH, dtjw, iexmin, iexmax, jexmins, jexmaxs, jexminn, jexmaxn real(DP), intent(in ) :: x_LonS ( 0:imax-1 ) real(DP), intent(in ) :: y_LatS ( 1:jmax/2 ) real(DP), intent(in ) :: x_LonN ( 0:imax-1 ) real(DP), intent(in ) :: y_LatN ( 1:jmax/2 ) real(DP), intent(out) :: x_ExtLonS(iexmin :iexmax ) real(DP), intent(out) :: y_ExtLatS(jexmins:jexmaxs) real(DP), intent(out) :: x_ExtLonN(iexmin :iexmax ) real(DP), intent(out) :: y_ExtLatN(jexminn:jexmaxn) ! ! local variables ! integer :: idest integer :: idep integer :: a_ireq_send(4) integer :: a_ireq_recv(4) real(DP), allocatable :: y_LatSPrev(:) real(DP), allocatable :: y_LatNPrev(:) real(DP), allocatable :: y_LatSNext(:) real(DP), allocatable :: y_LatNNext(:) real(DP) :: h, theta, thetasq integer :: i, j, k, m logical, save :: sw_fs data sw_fs /.true./ if( .not. sw_fs ) return sw_fs = .false. x_ExtLonS(-2) = 0.0_DP - ( x_LonS(2) - x_LonS(0) ) x_ExtLonS(-1) = 0.0_DP - ( x_LonS(1) - x_LonS(0) ) do i = 0, imax-1 x_ExtLonS(i) = x_LonS(i) end do x_ExtLonS(imax-1+1) = PIx2 x_ExtLonS(imax-1+2) = PIx2 + ( x_LonS(1) - x_LonS(0) ) x_ExtLonS(imax-1+3) = PIx2 + ( x_LonS(2) - x_LonS(0) ) ! x_ExtLonN(-2) = 0.0_DP - ( x_LonN(2) - x_LonN(0) ) x_ExtLonN(-1) = 0.0_DP - ( x_LonN(1) - x_LonN(0) ) do i = 0, imax-1 x_ExtLonN(i) = x_LonN(i) end do x_ExtLonN(imax-1+1) = PIx2 x_ExtLonN(imax-1+2) = PIx2 + ( x_LonN(1) - x_LonN(0) ) x_ExtLonN(imax-1+3) = PIx2 + ( x_LonN(2) - x_LonN(0) ) !==================================================================== if ( myrank > 0 ) then allocate( y_LatSPrev(1:dtjw) ) allocate( y_LatNPrev(1:dtjw) ) end if if ( myrank < (nprocs-1) ) then allocate( y_LatSNext(1:dtjw) ) allocate( y_LatNNext(1:dtjw) ) end if ! y_Lat(1:dtjw) values are transfered (y_LatSPrev). ! if ( myrank < (nprocs-1) ) then idest = myrank + 1 call MPIWrapperISend( idest, dtjw, y_Lat(1:dtjw), a_ireq_send(1) ) end if if ( myrank > 0 ) then idep = myrank - 1 call MPIWrapperIRecv( idep, dtjw, y_LatSPrev, a_ireq_recv(1) ) end if ! ! y_Lat(jmax/2+1-dtjw:jmax/2) values are transfered (y_LatSNext). ! if ( myrank > 0 ) then idest = myrank - 1 call MPIWrapperISend( idest, dtjw, y_Lat(jmax/2+1-dtjw:jmax/2), a_ireq_send(2) ) end if if ( myrank < (nprocs-1) ) then idep = myrank + 1 call MPIWrapperIRecv( idep, dtjw, y_LatSNext, a_ireq_recv(2) ) end if ! ! y_Lat(jmax/2+1:jmax/2+dtjw) values are transfered (y_LatNNext). ! if ( myrank > 0 ) then idest = myrank - 1 call MPIWrapperISend( idest, dtjw, y_Lat(jmax/2+1:jmax/2+dtjw), a_ireq_send(3) ) end if if ( myrank < (nprocs-1) ) then idep = myrank + 1 call MPIWrapperIRecv( idep, dtjw, y_LatNNext, a_ireq_recv(3) ) end if ! ! y_Lat(jmax+1-dtjw:jmax) values are transfered (y_LatNPrev). ! if ( myrank < (nprocs-1) ) then idest = myrank + 1 call MPIWrapperISend( idest, dtjw, y_Lat(jmax+1-dtjw:jmax), a_ireq_send(4) ) end if if ( myrank > 0 ) then idep = myrank - 1 call MPIWrapperIRecv( idep, dtjw, y_LatNPrev, a_ireq_recv(4) ) end if if ( myrank < (nprocs-1) ) call MPIWrapperWait( a_ireq_send(1) ) if ( myrank > 0 ) call MPIWrapperWait( a_ireq_recv(1) ) if ( myrank > 0 ) call MPIWrapperWait( a_ireq_send(2) ) if ( myrank < (nprocs-1) ) call MPIWrapperWait( a_ireq_recv(2) ) if ( myrank > 0 ) call MPIWrapperWait( a_ireq_send(3) ) if ( myrank < (nprocs-1) ) call MPIWrapperWait( a_ireq_recv(3) ) if ( myrank < (nprocs-1) ) call MPIWrapperWait( a_ireq_send(4) ) if ( myrank > 0 ) call MPIWrapperWait( a_ireq_recv(4) ) do j = 1, jmax/2 y_ExtLatS(j) = y_LatS(j) y_ExtLatN(j) = y_LatN(j) end do ! southern edge of southern array if( myrank == (nprocs-1) ) then y_ExtLatS(0) = -PIH do j = 1, dtjw y_ExtLatS(0-j) = -PIH - ( y_LatS(j) - ( -PIH ) ) end do else do j = 1, dtjw y_ExtLatS(1-j) = y_LatSNext(dtjw-(j-1)) end do end if ! northern edge of southern array if ( myrank == 0 ) then do j = 1, dtjw y_ExtLatS(jmax/2+j) = y_LatN(j) end do else do j = 1, dtjw y_ExtLatS(jmax/2+j) = y_LatSPrev(j) end do end if ! ! southern edge of northern array if ( myrank == 0 ) then do j = 1, dtjw y_ExtLatN(-j+1) = y_LatS(jmax/2-(j-1)) end do else do j = 1, dtjw y_ExtLatN(-j+1) = y_LatNPrev(dtjw-(j-1)) end do end if ! northern edge of northern array if( myrank == (nprocs-1) ) then y_ExtLatN(jmax/2+1) = PIH do j = 1, dtjw y_ExtLatN(jmax/2+1+j) = PIH + ( PIH - y_LatN(jmax/2+1-j) ) end do else do j = 1, dtjw y_ExtLatN(jmax/2+j) = y_LatNNext(j) end do end if if ( myrank > 0 ) then deallocate( y_LatSPrev ) deallocate( y_LatNPrev ) end if if ( myrank < (nprocs-1) ) then deallocate( y_LatSNext ) deallocate( y_LatNNext ) end if !!$ allocate( a_LQIFUVNP(jmax/2-1:jmax/2+2) ) !!$ allocate( a_LQIFUVSP( -1:2 ) ) allocate( a_LCIFUVNP(jmax/2-1:jmax/2+3) ) allocate( a_LCIFUVSP( -2:2 ) ) ! ! calculation of factors for Lagrange cubic interpolation used ! to estimate mixing ratios at poles ! if( myrank == (nprocs-1) ) then !!$ ! !!$ ! calculation of factors for Lagrange quadratic interpolation used !!$ ! to estimate wind velocities at south pole (old) !!$ ! !!$ a_LQIFUVSP(-1) = & !!$ & ( y_ExtLatS( 0) - y_ExtLatS( 1) ) * ( y_ExtLatS( 0) - y_ExtLatS( 2) ) & !!$ & / ( ( y_ExtLatS(-1) - y_ExtLatS( 1) ) * ( y_ExtLatS(-1) - y_ExtLatS( 2) ) ) !!$ a_LQIFUVSP( 0) = 1.0d100 !!$ a_LQIFUVSP( 1) = & !!$ & ( y_ExtLatS( 0) - y_ExtLatS(-1) ) * ( y_ExtLatS( 0) - y_ExtLatS( 2) ) & !!$ & / ( ( y_ExtLatS( 1) - y_ExtLatS(-1) ) * ( y_ExtLatS( 1) - y_ExtLatS( 2) ) ) !!$ a_LQIFUVSP( 2) = & !!$ & ( y_ExtLatS( 0) - y_ExtLatS(-1) ) * ( y_ExtLatS( 0) - y_ExtLatS( 1) ) & !!$ & / ( ( y_ExtLatS( 2) - y_ExtLatS(-1) ) * ( y_ExtLatS( 2) - y_ExtLatS( 1) ) ) !!$ ! !!$ ! calculation of factors for Lagrange quadratic interpolation used !!$ ! to estimate wind velocities at north pole (old) !!$ ! !!$ a_LQIFUVNP(jmax/2-1) = & !!$ & ( y_ExtLatN(jmax/2+1) - y_ExtLatN(jmax/2 ) ) * ( y_ExtLatN(jmax/2+1) - y_ExtLatN(jmax/2+2) ) & !!$ & / ( ( y_ExtLatN(jmax/2-1) - y_ExtLatN(jmax/2 ) ) * ( y_ExtLatN(jmax/2-1) - y_ExtLatN(jmax/2+2) ) ) !!$ a_LQIFUVNP(jmax/2 ) = & !!$ & ( y_ExtLatN(jmax/2+1) - y_ExtLatN(jmax/2-1) ) * ( y_ExtLatN(jmax/2+1) - y_ExtLatN(jmax/2+2) ) & !!$ & / ( ( y_ExtLatN(jmax/2 ) - y_ExtLatN(jmax/2-1) ) * ( y_ExtLatN(jmax/2 ) - y_ExtLatN(jmax/2+2) ) ) !!$ a_LQIFUVNP(jmax/2+1) = 1.0d100 !!$ a_LQIFUVNP(jmax/2+2) = & !!$ & ( y_ExtLatN(jmax/2+1) - y_ExtLatN(jmax/2-1) ) * ( y_ExtLatN(jmax/2+1) - y_ExtLatN(jmax/2 ) ) & !!$ & / ( ( y_ExtLatN(jmax/2+2) - y_ExtLatN(jmax/2-1) ) * ( y_ExtLatN(jmax/2+2) - y_ExtLatN(jmax/2 ) ) ) ! ! calculation of factors for Lagrange cubic interpolation used ! to estimate wind velocities at south pole ! a_LCIFUVSP(-2) = ( y_ExtLatS( 0) - y_ExtLatS(-1) ) * ( y_ExtLatS( 0) - y_ExtLatS( 1) ) * ( y_ExtLatS( 0) - y_ExtLatS( 2) ) / ( ( y_ExtLatS(-2) - y_ExtLatS(-1) ) * ( y_ExtLatS(-2) - y_ExtLatS( 1) ) * ( y_ExtLatS(-2) - y_ExtLatS( 2) ) ) a_LCIFUVSP(-1) = ( y_ExtLatS( 0) - y_ExtLatS(-2) ) * ( y_ExtLatS( 0) - y_ExtLatS( 1) ) * ( y_ExtLatS( 0) - y_ExtLatS( 2) ) / ( ( y_ExtLatS(-1) - y_ExtLatS(-2) ) * ( y_ExtLatS(-1) - y_ExtLatS( 1) ) * ( y_ExtLatS(-1) - y_ExtLatS( 2) ) ) a_LCIFUVSP( 0) = 1.0d100 a_LCIFUVSP( 1) = ( y_ExtLatS( 0) - y_ExtLatS(-2) ) * ( y_ExtLatS( 0) - y_ExtLatS(-1) ) * ( y_ExtLatS( 0) - y_ExtLatS( 2) ) / ( ( y_ExtLatS( 1) - y_ExtLatS(-2) ) * ( y_ExtLatS( 1) - y_ExtLatS(-1) ) * ( y_ExtLatS( 1) - y_ExtLatS( 2) ) ) a_LCIFUVSP( 2) = ( y_ExtLatS( 0) - y_ExtLatS(-2) ) * ( y_ExtLatS( 0) - y_ExtLatS(-1) ) * ( y_ExtLatS( 0) - y_ExtLatS( 1) ) / ( ( y_ExtLatS( 2) - y_ExtLatS(-2) ) * ( y_ExtLatS( 2) - y_ExtLatS(-1) ) * ( y_ExtLatS( 2) - y_ExtLatS( 1) ) ) ! ! calculation of factors for Lagrange cubic interpolation used ! to estimate wind velocities at north pole ! a_LCIFUVNP(jmax/2-1) = ( y_ExtLatN(jmax/2+1) - y_ExtLatN(jmax/2 ) ) * ( y_ExtLatN(jmax/2+1) - y_ExtLatN(jmax/2+2) ) * ( y_ExtLatN(jmax/2+1) - y_ExtLatN(jmax/2+3) ) / ( ( y_ExtLatN(jmax/2-1) - y_ExtLatN(jmax/2 ) ) * ( y_ExtLatN(jmax/2-1) - y_ExtLatN(jmax/2+2) ) * ( y_ExtLatN(jmax/2-1) - y_ExtLatN(jmax/2+3) ) ) a_LCIFUVNP(jmax/2 ) = ( y_ExtLatN(jmax/2+1) - y_ExtLatN(jmax/2-1) ) * ( y_ExtLatN(jmax/2+1) - y_ExtLatN(jmax/2+2) ) * ( y_ExtLatN(jmax/2+1) - y_ExtLatN(jmax/2+3) ) / ( ( y_ExtLatN(jmax/2 ) - y_ExtLatN(jmax/2-1) ) * ( y_ExtLatN(jmax/2 ) - y_ExtLatN(jmax/2+2) ) * ( y_ExtLatN(jmax/2 ) - y_ExtLatN(jmax/2+3) ) ) a_LCIFUVNP(jmax/2+1) = 1.0d100 a_LCIFUVNP(jmax/2+2) = ( y_ExtLatN(jmax/2+1) - y_ExtLatN(jmax/2-1) ) * ( y_ExtLatN(jmax/2+1) - y_ExtLatN(jmax/2 ) ) * ( y_ExtLatN(jmax/2+1) - y_ExtLatN(jmax/2+3) ) / ( ( y_ExtLatN(jmax/2+2) - y_ExtLatN(jmax/2-1) ) * ( y_ExtLatN(jmax/2+2) - y_ExtLatN(jmax/2 ) ) * ( y_ExtLatN(jmax/2+2) - y_ExtLatN(jmax/2+3) ) ) a_LCIFUVNP(jmax/2+3) = ( y_ExtLatN(jmax/2+1) - y_ExtLatN(jmax/2-1) ) * ( y_ExtLatN(jmax/2+1) - y_ExtLatN(jmax/2 ) ) * ( y_ExtLatN(jmax/2+1) - y_ExtLatN(jmax/2+2) ) / ( ( y_ExtLatN(jmax/2+3) - y_ExtLatN(jmax/2-1) ) * ( y_ExtLatN(jmax/2+3) - y_ExtLatN(jmax/2 ) ) * ( y_ExtLatN(jmax/2+3) - y_ExtLatN(jmax/2+2) ) ) end if ! ! variable for estimating polar values ! if( myrank == ( nprocs-1 ) ) then SumSinSq = imax / 2 else SumSinSq = 1.0d100 end if end subroutine SLTTExtArrInit
Subroutine : |
subroutine SLTTExtArrMkJMAXTable ! メッセージ出力 ! Message output ! use dc_message, only: MessageNotify ! ! MPI ! use mpi_wrapper, only : nprocs, myrank, MPIWrapperISend, MPIWrapperIRecv, MPIWrapperWait ! ! local variables ! integer , save, allocatable :: a_TblJMAX(:) integer , save, allocatable :: a_NSSepTblRank(:) integer , save, allocatable :: a_NSSepTblNGrid(:) character(STRING), save, allocatable :: a_NSSepTblNS(:) integer :: a_SendBuf(1) integer, allocatable :: aa_RecvBuf(:,:) integer, allocatable :: a_iReqSend(:) integer, allocatable :: a_iReqRecv(:) integer :: n integer :: l allocate( a_TblJMAX(0:nprocs-1) ) ! Make a table containing jmax in all processes ! allocate( aa_RecvBuf(1,0:nprocs-1) ) allocate( a_iReqSend(0:nprocs-1) ) allocate( a_iReqRecv(0:nprocs-1) ) ! a_SendBuf = jmax do n = 0, nprocs-1 if ( n == myrank ) cycle call MPIWrapperISend( n, 1, a_SendBuf , a_iReqSend(n) ) call MPIWrapperIRecv( n, 1, aa_RecvBuf(:,n), a_iReqRecv(n) ) end do do n = 0, nprocs-1 if ( n == myrank ) cycle call MPIWrapperWait( a_iReqSend(n) ) call MPIWrapperWait( a_iReqRecv(n) ) end do ! aa_RecvBuf(:,myrank) = a_SendBuf do n = 0, nprocs-1 a_TblJMAX(n) = aa_RecvBuf(1,n) end do ! deallocate( aa_RecvBuf ) deallocate( a_iReqSend ) deallocate( a_iReqRecv ) ! Table is checked do n = 0, nprocs-1 if ( mod( a_TblJMAX(n), 2 ) /= 0 ) then call MessageNotify( 'E', module_name, 'Unexpected jmax in process %d, %d.', i = (/ n, a_TblJMAX(n) /) ) end if end do allocate( a_NSSepTblRank (1:nprocs*2) ) allocate( a_NSSepTblNGrid(1:nprocs*2) ) allocate( a_NSSepTblNS (1:nprocs*2) ) ! a_NSSepTblRank : rank included in a North-South separate table l = 1 do n = nprocs-1, 0, -1 a_NSSepTblRank(l) = n l = l + 1 end do do n = 0, nprocs-1 a_NSSepTblRank(l) = n l = l + 1 end do ! a_NSSepTblNGrid : number of grid included in a North-South separate table do l = 1, nprocs*2 a_NSSepTblNGrid(l) = a_TblJMAX(a_NSSepTblRank(l)) / 2 end do ! a_NSSepTblNS : Symbol representing south- or north-array ! : included in a North-South separate table do l = 1, nprocs*2 if ( l <= nprocs ) then a_NSSepTblNS(l) = 'S' else a_NSSepTblNS(l) = 'N' end if end do !!$ a_NSSepTblNToBeSent !!$ a_NSSepTblNToBeRecv deallocate( a_NSSepTblRank ) deallocate( a_NSSepTblNGrid ) deallocate( a_TblJMAX ) end subroutine SLTTExtArrMkJMAXTable
Constant : | |||
version = ’$Name: dcpam5-20130302 $’ // ’$Id: sltt_extarr.f90,v 1.2 2013-01-27 11:26:14 yot Exp $’ : | character(*), parameter
|