Changeset 3708 for palm/trunk/SOURCE
- Timestamp:
- Jan 30, 2019 12:58:13 PM (6 years ago)
- File:
-
- 1 edited
-
palm/trunk/SOURCE/pmc_interface_mod.f90 (modified) (14 diffs)
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/SOURCE/pmc_interface_mod.f90
r3697 r3708 25 25 ! ----------------- 26 26 ! $Id$ 27 ! Checks for parent / child grid line matching introduced. 28 ! Interpolation of nest-boundary-tangential velocity components revised. 29 ! 30 ! 3697 2019-01-24 17:16:13Z hellstea 27 31 ! Bugfix: upper k-bound in the child initialization interpolation 28 32 ! pmci_interp_1sto_all corrected. … … 1475 1479 !-- included in the interpolation algorithms. 1476 1480 CALL pmci_init_anterp_tophat 1481 ! 1482 !-- Check that the child and parent grid lines match 1483 CALL pmci_check_grid_matching 1477 1484 1478 1485 ENDIF … … 1928 1935 CEILING( parentdzmax / dzmin ) 1929 1936 ALLOCATE( celltmpd(1:acsize) ) 1930 1937 1931 1938 END SUBROUTINE pmci_init_interp_tril 1932 1939 … … 1941 1948 igsr = NINT( cg%dx / dx, iwp ) 1942 1949 jgsr = NINT( cg%dy / dy, iwp ) 1943 kgsr = NINT( cg%dz u(cg%nz+1) / dzu(nzt+1), iwp )1950 kgsr = NINT( cg%dzw(1) / dzw(1), iwp ) 1944 1951 write(9,"('igsr, jgsr, kgsr: ',3(i3,2x))") igsr, jgsr, kgsr 1945 1952 flush(9) … … 3899 3906 END SUBROUTINE pmci_init_anterp_tophat 3900 3907 3901 #endif 3908 3909 3910 SUBROUTINE pmci_check_grid_matching 3911 ! 3912 !-- Check that the grid lines of child and parent do match 3913 IMPLICIT NONE 3914 REAL(wp), PARAMETER :: tolefac = 1.0E-6_wp 3915 REAL(wp) :: tolex 3916 REAL(wp) :: toley 3917 REAL(wp) :: tolez 3918 INTEGER(iwp) :: non_matching_lower_left_corner 3919 INTEGER(iwp) :: non_int_gsr_x 3920 INTEGER(iwp) :: non_int_gsr_y 3921 INTEGER(iwp) :: non_int_gsr_z 3922 3923 3924 non_matching_lower_left_corner = 0 3925 non_int_gsr_x = 0 3926 non_int_gsr_y = 0 3927 non_int_gsr_z = 0 3928 3929 IF ( myid == 0 ) THEN 3930 3931 tolex = tolefac * dx 3932 toley = tolefac * dy 3933 tolez = tolefac * MINVAL( dzw ) 3934 ! 3935 !-- First check that the child grid lower left corner matches the paren grid lines. 3936 IF ( MOD( lower_left_coord_x, cg%dx ) > tolex ) non_matching_lower_left_corner = 1 3937 IF ( MOD( lower_left_coord_y, cg%dy ) > toley ) non_matching_lower_left_corner = 1 3938 ! 3939 !-- Then check that the grid-spacing ratios in each direction are integer valued. 3940 IF ( MOD( cg%dx, dx ) > tolex ) non_int_gsr_x = 1 3941 IF ( MOD( cg%dy, dy ) > toley ) non_int_gsr_y = 1 3942 DO n = 0, kctw+1 3943 IF ( ABS( cg%zw(n) - zw(kflw(n)) ) > tolez ) non_int_gsr_z = 1 3944 ENDDO 3945 3946 ENDIF 3947 3948 CALL MPI_BCAST( non_matching_lower_left_corner, 1, MPI_INTEGER, 0, & 3949 comm2d, ierr ) 3950 IF ( non_matching_lower_left_corner > 0 ) THEN 3951 WRITE ( message_string, * ) 'nested child domain lower left ', & 3952 'corner must match its parent grid lines' 3953 CALL message( 'pmci_check_grid_matching', 'PA0414', 3, 2, 0, 6, 0 ) 3954 ENDIF 3955 3956 CALL MPI_BCAST( non_int_gsr_x, 1, MPI_INTEGER, 0, comm2d, ierr ) 3957 IF ( non_int_gsr_x > 0 ) THEN 3958 WRITE ( message_string, * ) 'nesting grid-spacing ratio ', & 3959 '( parent dx / child dx ) must have an integer value' 3960 CALL message( 'pmci_check_grid_matching', 'PA0416', 3, 2, 0, 6, 0 ) 3961 ENDIF 3962 3963 CALL MPI_BCAST( non_int_gsr_y, 1, MPI_INTEGER, 0, comm2d, ierr ) 3964 IF ( non_int_gsr_y > 0 ) THEN 3965 WRITE ( message_string, * ) 'nesting grid-spacing ratio ', & 3966 '( parent dy / child dy ) must have an integer value' 3967 CALL message( 'pmci_check_grid_matching', 'PA0416', 3, 2, 0, 6, 0 ) 3968 ENDIF 3969 3970 CALL MPI_BCAST( non_int_gsr_z, 1, MPI_INTEGER, 0, comm2d, ierr ) 3971 IF ( non_int_gsr_z > 0 ) THEN 3972 WRITE ( message_string, * ) 'nesting grid-spacing ratio ', & 3973 '( parent dz / child dz ) must have an integer value for each z-level' 3974 CALL message( 'pmci_check_grid_matching', 'PA0416', 3, 2, 0, 6, 0 ) 3975 ENDIF 3976 3977 END SUBROUTINE pmci_check_grid_matching 3978 3979 3980 #endif 3902 3981 END SUBROUTINE pmci_setup_child 3903 3982 … … 4557 4636 DO m = mb, me 4558 4637 DO n = 0, kct + 1 4559 4638 4560 4639 DO i = ifl(l), ifu(l) 4561 4640 DO j = jfl(m), jfl(m+1)-1 … … 4594 4673 DO m = mb, me 4595 4674 DO n = 0, kct + 1 4596 4675 4597 4676 DO i = ifl(l), ifu(l) 4598 DO j = jfl(m), jfu(m) 4677 DO j = jfl(m), jfu(m) 4599 4678 DO k = kfl(n), MIN( kfu(n), nzt+1 ) 4600 4679 f(k,j,i) = fc(n,m,l) … … 5785 5864 REAL(wp) :: rcorr !< Average reversibility correction for the whole anterpolation cell 5786 5865 REAL(wp) :: rcorr_ijk !< Reversibility correction distributed to the individual child-grid nodes 5866 5867 real(wp), parameter :: c1 = 2.0_wp / 6.0_wp 5868 real(wp), parameter :: c2 = 5.0_wp / 6.0_wp 5869 real(wp), parameter :: c3 = -1.0_wp / 6.0_wp 5787 5870 5788 5871 ! … … 5858 5941 ENDIF 5859 5942 5860 IF ( var == 'u' ) THEN 5861 !AH! 5862 !AH!-- Substitute the necessary parent-grid data to the work array workarrc_lr. 5863 !AH workarrc_lr = 0.0_wp 5864 !AH IF ( pdims(2) > 1 ) THEN 5865 !AH#if defined( __parallel ) 5866 !AH IF ( nys == 0 ) THEN 5867 !AH workarrc_lr(0:cg%nz+1,jcsw:jcnw-1,0:2) & 5868 !AH = fc(0:cg%nz+1,jcsw:jcnw-1,lbeg:lbeg+2) 5869 !AH ELSE IF ( nyn == ny ) THEN 5870 !AH workarrc_lr(0:cg%nz+1,jcsw+1:jcnw,0:2) & 5871 !AH = fc(0:cg%nz+1,jcsw+1:jcnw,lbeg:lbeg+2) 5872 !AH ELSE 5873 !AH workarrc_lr(0:cg%nz+1,jcsw+1:jcnw-1,0:2) & 5874 !AH = fc(0:cg%nz+1,jcsw+1:jcnw-1,lbeg:lbeg+2) 5875 !AH ENDIF 5876 !AH! 5877 !AH!-- South-north exchange if more than one PE subdomain in the y-direction. 5878 !AH!-- Note that in case of 3-D nesting the south (psouth == MPI_PROC_NULL) 5879 !AH!-- and north (pnorth == MPI_PROC_NULL) boundaries are not exchanged 5880 !AH!-- because the nest domain is not cyclic. 5881 !AH!-- From south to north 5882 !AH CALL MPI_SENDRECV( workarrc_lr(0,jcsw+1,0), 1, & 5883 !AH workarrc_lr_exchange_type, psouth, 0, & 5884 !AH workarrc_lr(0,jcnw,0), 1, & 5885 !AH workarrc_lr_exchange_type, pnorth, 0, & 5886 !AH comm2d, status, ierr ) 5887 !AH! 5888 !AH!-- From north to south 5889 !AH CALL MPI_SENDRECV( workarrc_lr(0,jcnw-1,0), 1, & 5890 !AH workarrc_lr_exchange_type, pnorth, 1, & 5891 !AH workarrc_lr(0,jcsw,0), 1, & 5892 !AH workarrc_lr_exchange_type, psouth, 1, & 5893 !AH comm2d, status, ierr ) 5894 !AH#endif 5895 !AH ELSE 5896 !AH workarrc_lr(0:cg%nz+1,jcsw:jcnw,0:2) & 5897 !AH = fc(0:cg%nz+1,jcsw:jcnw,lbeg:lbeg+2) 5898 !AH ENDIF 5899 5900 DO m = jcsw+1, jcnw-1 5901 DO n = 0, kct 5902 5903 DO j = jfl(m), jfu(m) 5904 DO k = kfl(n), kfu(n) 5905 f(k,j,ibc) = fc(n,m,l) 5943 IF ( var == 'v' .OR. var == 'w' ) THEN 5944 ! 5945 !-- Substitute the necessary parent-grid data to the work array workarrc_lr. 5946 workarrc_lr = 0.0_wp 5947 IF ( pdims(2) > 1 ) THEN 5948 #if defined( __parallel ) 5949 IF ( nys == 0 ) THEN 5950 workarrc_lr(0:cg%nz+1,jcsw:jcnw-1,0:2) & 5951 = fc(0:cg%nz+1,jcsw:jcnw-1,lbeg:lbeg+2) 5952 ELSE IF ( nyn == ny ) THEN 5953 workarrc_lr(0:cg%nz+1,jcsw+1:jcnw,0:2) & 5954 = fc(0:cg%nz+1,jcsw+1:jcnw,lbeg:lbeg+2) 5955 ELSE 5956 workarrc_lr(0:cg%nz+1,jcsw+1:jcnw-1,0:2) & 5957 = fc(0:cg%nz+1,jcsw+1:jcnw-1,lbeg:lbeg+2) 5958 ENDIF 5959 ! 5960 !-- South-north exchange if more than one PE subdomain in the y-direction. 5961 !-- Note that in case of 3-D nesting the south (psouth == MPI_PROC_NULL) 5962 !-- and north (pnorth == MPI_PROC_NULL) boundaries are not exchanged 5963 !-- because the nest domain is not cyclic. 5964 !-- From south to north 5965 CALL MPI_SENDRECV( workarrc_lr(0,jcsw+1,0), 1, & 5966 workarrc_lr_exchange_type, psouth, 0, & 5967 workarrc_lr(0,jcnw,0), 1, & 5968 workarrc_lr_exchange_type, pnorth, 0, & 5969 comm2d, status, ierr ) 5970 ! 5971 !-- From north to south 5972 CALL MPI_SENDRECV( workarrc_lr(0,jcnw-1,0), 1, & 5973 workarrc_lr_exchange_type, pnorth, 1, & 5974 workarrc_lr(0,jcsw,0), 1, & 5975 workarrc_lr_exchange_type, psouth, 1, & 5976 comm2d, status, ierr ) 5977 #endif 5978 ELSE 5979 workarrc_lr(0:cg%nz+1,jcsw:jcnw,0:2) & 5980 = fc(0:cg%nz+1,jcsw:jcnw,lbeg:lbeg+2) 5981 ENDIF 5982 5983 IF ( var == 'v' ) THEN 5984 5985 DO m = jcsw+1, jcnw 5986 DO n = 0, kct 5987 ! 5988 !-- First substitute only the matching-node values 5989 f(kfl(n):kfu(n),jfl(m),ibc) = workarrc_lr(n,m,lw) 5990 5991 ENDDO 5992 ENDDO 5993 5994 DO m = jcsw+1, jcnw-1 5995 DO n = 0, kct 5996 ! 5997 !-- Then fill up the nodes in between with the averages 5998 DO j = jfl(m)+1, jfl(m+1)-1 5999 f(kfl(n):kfu(n),j,ibc) = 0.5_wp * ( f(kfl(n):kfu(n),jfl(m),ibc) & 6000 + f(kfl(n):kfu(n),jfl(m+1),ibc) ) 5906 6001 ENDDO 6002 5907 6003 ENDDO 5908 5909 6004 ENDDO 5910 ENDDO 5911 5912 ELSE IF ( var == 'v' ) THEN 5913 5914 DO m = jcsw+1, jcnw-1 5915 DO n = 0, kct 5916 5917 DO j = jfl(m), jfl(m+1)-1 5918 DO k = kfl(n), kfu(n) 5919 f(k,j,ibc) = fc(n,m,l) 6005 6006 ELSE IF ( var == 'w' ) THEN 6007 6008 DO m = jcsw+1, jcnw-1 6009 DO n = 0, kct + 1 ! It is important to go up to kct+1 6010 ! 6011 !-- First substitute only the matching-node values 6012 f(kfu(n),jfl(m):jfu(m),ibc) = workarrc_lr(n,m,lw) 6013 6014 ENDDO 6015 ENDDO 6016 6017 DO m = jcsw+1, jcnw-1 6018 DO n = 1, kct + 1 ! It is important to go up to kct+1 6019 ! 6020 !-- Then fill up the nodes in between with the averages 6021 DO k = kfu(n-1)+1, kfu(n)-1 6022 f(k,jfl(m):jfu(m),ibc) = 0.5_wp * ( f(kfu(n-1),jfl(m):jfu(m),ibc) & 6023 + f(kfu(n),jfl(m):jfu(m),ibc) ) 5920 6024 ENDDO 6025 5921 6026 ENDDO 5922 5923 6027 ENDDO 5924 ENDDO 5925 5926 ELSE IF ( var == 'w' ) THEN 5927 5928 DO m = jcsw+1, jcnw-1 5929 DO n = 1, kct + 1 ! It is important to go up to kct+1 5930 5931 DO j = jfl(m), jfu(m) 5932 f(nzb,j,ibc) = 0.0_wp ! Because the n-loop starts from n=1 instead of 0 5933 DO k = kfu(n-1)+1, kfu(n) 5934 f(k,j,ibc) = fc(n,m,l) 5935 ENDDO 5936 ENDDO 5937 5938 ENDDO 5939 ENDDO 5940 5941 ELSE ! scalars 6028 6029 ENDIF 6030 6031 ELSE ! u or scalars 5942 6032 5943 6033 DO m = jcsw+1, jcnw-1 … … 6081 6171 REAL(wp) :: rcorr !< Average reversibility correction for the whole anterpolation cell 6082 6172 REAL(wp) :: rcorr_ijk !< Reversibility correction distributed to the individual child-grid nodes 6173 6174 real(wp), parameter :: c1 = 2.0_wp / 6.0_wp 6175 real(wp), parameter :: c2 = 5.0_wp / 6.0_wp 6176 real(wp), parameter :: c3 = -1.0_wp / 6.0_wp 6083 6177 6084 6178 ! … … 6154 6248 ENDIF 6155 6249 6156 IF ( var == 'v' ) THEN 6157 !AH! 6158 !AH!-- Substitute the necessary parent-grid data to the work array workarrc_sn. 6159 !AH workarrc_sn = 0.0_wp 6160 !AH IF ( pdims(1) > 1 ) THEN 6161 !AH#if defined( __parallel ) 6162 !AH IF ( nxl == 0 ) THEN ! if ( bc_dirichlet_l ) 6163 !AH workarrc_sn(0:cg%nz+1,0:2,iclw:icrw-1) & 6164 !AH = fc(0:cg%nz+1,mbeg:mbeg+2,iclw:icrw-1) 6165 !AH ELSE IF ( nxr == nx ) THEN ! if ( bc_dirichlet_r ) 6166 !AH workarrc_sn(0:cg%nz+1,0:2,iclw+1:icrw) & 6167 !AH = fc(0:cg%nz+1,mbeg:mbeg+2,iclw+1:icrw) 6168 !AH ELSE 6169 !AH workarrc_sn(0:cg%nz+1,0:2,iclw+1:icrw-1) & 6170 !AH = fc(0:cg%nz+1,mbeg:mbeg+2,iclw+1:icrw-1) 6171 !AH ENDIF 6172 !AH! 6173 !AH!-- Left-right exchange if more than one PE subdomain in the x-direction. 6174 !AH!-- Note that in case of 3-D nesting the left (pleft == MPI_PROC_NULL) and 6175 !AH!-- right (pright == MPI_PROC_NULL) boundaries are not exchanged because 6176 !AH!-- the nest domain is not cyclic. 6177 !AH!-- From left to right 6178 !AH CALL MPI_SENDRECV( workarrc_sn(0,0,iclw+1), 1, & 6179 !AH workarrc_sn_exchange_type, pleft, 0, & 6180 !AH workarrc_sn(0,0,icrw), 1, & 6181 !AH workarrc_sn_exchange_type, pright, 0, & 6182 !AH comm2d, status, ierr ) 6183 !AH! 6184 !AH!-- From right to left 6185 !AH CALL MPI_SENDRECV( workarrc_sn(0,0,icrw-1), 1, & 6186 !AH workarrc_sn_exchange_type, pright, 1, & 6187 !AH workarrc_sn(0,0,iclw), 1, & 6188 !AH workarrc_sn_exchange_type, pleft, 1, & 6189 !AH comm2d, status, ierr ) 6190 !AH#endif 6191 !AH ELSE 6192 !AH workarrc_sn(0:cg%nz+1,0:2,iclw+1:icrw-1) & 6193 !AH = fc(0:cg%nz+1,mbeg:mbeg+2,iclw+1:icrw-1) 6194 !AH ENDIF 6195 6196 DO l = iclw+1, icrw-1 6197 DO n = 0, kct 6198 6199 DO i = ifl(l), ifu(l) 6200 DO k = kfl(n), kfu(n) 6201 f(k,jbc,i) = fc(n,m,l) 6250 IF ( var == 'u' .OR. var == 'w' ) THEN 6251 ! 6252 !-- Substitute the necessary parent-grid data to the work array workarrc_sn. 6253 workarrc_sn = 0.0_wp 6254 IF ( pdims(1) > 1 ) THEN 6255 #if defined( __parallel ) 6256 IF ( nxl == 0 ) THEN ! if ( bc_dirichlet_l ) 6257 workarrc_sn(0:cg%nz+1,0:2,iclw:icrw-1) & 6258 = fc(0:cg%nz+1,mbeg:mbeg+2,iclw:icrw-1) 6259 ELSE IF ( nxr == nx ) THEN ! if ( bc_dirichlet_r ) 6260 workarrc_sn(0:cg%nz+1,0:2,iclw+1:icrw) & 6261 = fc(0:cg%nz+1,mbeg:mbeg+2,iclw+1:icrw) 6262 ELSE 6263 workarrc_sn(0:cg%nz+1,0:2,iclw+1:icrw-1) & 6264 = fc(0:cg%nz+1,mbeg:mbeg+2,iclw+1:icrw-1) 6265 ENDIF 6266 ! 6267 !-- Left-right exchange if more than one PE subdomain in the x-direction. 6268 !-- Note that in case of 3-D nesting the left (pleft == MPI_PROC_NULL) and 6269 !-- right (pright == MPI_PROC_NULL) boundaries are not exchanged because 6270 !-- the nest domain is not cyclic. 6271 !-- From left to right 6272 CALL MPI_SENDRECV( workarrc_sn(0,0,iclw+1), 1, & 6273 workarrc_sn_exchange_type, pleft, 0, & 6274 workarrc_sn(0,0,icrw), 1, & 6275 workarrc_sn_exchange_type, pright, 0, & 6276 comm2d, status, ierr ) 6277 ! 6278 !-- From right to left 6279 CALL MPI_SENDRECV( workarrc_sn(0,0,icrw-1), 1, & 6280 workarrc_sn_exchange_type, pright, 1, & 6281 workarrc_sn(0,0,iclw), 1, & 6282 workarrc_sn_exchange_type, pleft, 1, & 6283 comm2d, status, ierr ) 6284 #endif 6285 ELSE 6286 workarrc_sn(0:cg%nz+1,0:2,iclw+1:icrw-1) & 6287 = fc(0:cg%nz+1,mbeg:mbeg+2,iclw+1:icrw-1) 6288 ENDIF 6289 6290 IF ( var == 'u' ) THEN 6291 6292 DO l = iclw+1, icrw 6293 DO n = 0, kct 6294 ! 6295 !-- First substitute only the matching-node values 6296 f(kfl(n):kfu(n),jbc,ifl(l)) = workarrc_sn(n,mw,l) 6297 6298 ENDDO 6299 ENDDO 6300 6301 DO l = iclw+1, icrw-1 6302 DO n = 0, kct 6303 ! 6304 !-- Then fill up the nodes in between with the averages 6305 DO i = ifl(l)+1, ifl(l+1)-1 6306 f(kfl(n):kfu(n),jbc,i) = 0.5_wp * ( f(kfl(n):kfu(n),jbc,ifl(l)) & 6307 + f(kfl(n):kfu(n),jbc,ifl(l+1)) ) 6202 6308 ENDDO 6309 6203 6310 ENDDO 6204 6205 6311 ENDDO 6206 ENDDO 6207 6208 ELSE IF ( var == 'u' ) THEN 6209 6210 DO l = iclw+1, icrw-1 6211 DO n = 0, kct 6212 6213 DO i = ifl(l), ifl(l+1)-1 6214 DO k = kfl(n), kfu(n) 6215 f(k,jbc,i) = fc(n,m,l) 6216 !AH 6217 ! write(9,"('sn u: ',6(i3,2x),2(e12.5,2x))") n, m, l, k, jbc, i, fc(n,m,l), f(k,jbc,i) 6218 ! flush(9) 6219 !AH 6312 6313 ELSE IF ( var == 'w' ) THEN 6314 6315 DO l = iclw+1, icrw-1 6316 DO n = 0, kct + 1 ! It is important to go up to kct+1 6317 ! 6318 !-- First substitute only the matching-node values 6319 f(kfu(n),jbc,ifl(l):ifu(l)) = workarrc_sn(n,mw,l) 6320 6321 ENDDO 6322 ENDDO 6323 6324 DO l = iclw+1, icrw-1 6325 DO n = 1, kct + 1 ! It is important to go up to kct+1 6326 ! 6327 !-- Then fill up the nodes in between with the averages 6328 DO k = kfu(n-1)+1, kfu(n)-1 6329 f(k,jbc,ifl(l):ifu(l)) = 0.5_wp * ( f(kfu(n-1),jbc,ifl(l):ifu(l)) & 6330 + f(kfu(n),jbc,ifl(l):ifu(l)) ) 6220 6331 ENDDO 6332 6221 6333 ENDDO 6222 6223 6334 ENDDO 6224 ENDDO 6225 6226 ELSE IF ( var == 'w' ) THEN 6227 6228 DO l = iclw+1, icrw-1 6229 DO n = 1, kct + 1 ! It is important to go up to kct+1 6230 6231 DO i = ifl(l), ifu(l) 6232 f(nzb,jbc,i) = 0.0_wp ! Because the n-loop starts from n=1 instead of 0 6233 DO k = kfu(n-1)+1, kfu(n) 6234 f(k,jbc,i) = fc(n,m,l) 6235 ENDDO 6236 ENDDO 6237 6238 ENDDO 6239 ENDDO 6240 6241 ELSE ! scalars 6335 6336 ENDIF 6337 6338 ELSE ! v or scalars 6242 6339 6243 6340 DO l = iclw+1, icrw-1 … … 6358 6455 REAL(wp) :: fkpjp !< 6359 6456 REAL(wp) :: fkp !< 6457 REAL(wp) :: f_interp !< 6360 6458 REAL(wp) :: rcorr !< 6361 6459 REAL(wp) :: rcorr_ijk !< 6460 6461 real(wp), parameter :: c1 = 2.0_wp / 6.0_wp 6462 real(wp), parameter :: c2 = 5.0_wp / 6.0_wp 6463 real(wp), parameter :: c3 = -1.0_wp / 6.0_wp 6362 6464 6363 6465 … … 6380 6482 ENDIF 6381 6483 n = kc(k) + noff 6382 nw = 1 6383 6384 IF ( var == 'w' ) THEN 6385 !AH! 6386 !AH!-- Substitute the necessary parent-grid data to the work array. 6387 !AH!-- Note that the dimension of workarrc_t is (0:2,jcsw:jcnw,iclw:icrw), 6388 !AH!-- And the jc?w and ic?w-index bounds depend on the location of the PE- 6389 !AH!-- subdomain relative to the side boundaries. 6390 !AH iclc = iclw + 1 6391 !AH icrc = icrw - 1 6392 !AH jcsc = jcsw + 1 6393 !AH jcnc = jcnw - 1 6394 !AH IF ( bc_dirichlet_l ) THEN 6395 !AH iclc = iclw 6396 !AH ENDIF 6397 !AH IF ( bc_dirichlet_r ) THEN 6398 !AH icrc = icrw 6399 !AH ENDIF 6400 !AH IF ( bc_dirichlet_s ) THEN 6401 !AH jcsc = jcsw 6402 !AH ENDIF 6403 !AH IF ( bc_dirichlet_n ) THEN 6404 !AH jcnc = jcnw 6405 !AH ENDIF 6406 !AH workarrc_t = 0.0_wp 6407 !AH workarrc_t(0:2,jcsc:jcnc,iclc:icrc) & 6408 !AH = fc(kc(k):kc(k)+2,jcsc:jcnc,iclc:icrc) 6409 !AH! 6410 !AH!-- Left-right exchange if more than one PE subdomain in the x-direction. 6411 !AH!-- Note that in case of 3-D nesting the left and right boundaries are 6412 !AH!-- not exchanged because the nest domain is not cyclic. 6413 !AH#if defined( __parallel ) 6414 !AH IF ( pdims(1) > 1 ) THEN 6415 !AH! 6416 !AH!-- From left to right 6417 !AH CALL MPI_SENDRECV( workarrc_t(0,jcsw,iclw+1), 1, & 6418 !AH workarrc_t_exchange_type_y, pleft, 0, & 6419 !AH workarrc_t(0,jcsw,icrw), 1, & 6420 !AH workarrc_t_exchange_type_y, pright, 0, & 6421 !AH comm2d, status, ierr ) 6422 !AH! 6423 !AH!-- From right to left 6424 !AH CALL MPI_SENDRECV( workarrc_t(0,jcsw,icrw-1), 1, & 6425 !AH workarrc_t_exchange_type_y, pright, 1, & 6426 !AH workarrc_t(0,jcsw,iclw), 1, & 6427 !AH workarrc_t_exchange_type_y, pleft, 1, & 6428 !AH comm2d, status, ierr ) 6429 !AH ENDIF 6430 !AH! 6431 !AH!-- South-north exchange if more than one PE subdomain in the y-direction. 6432 !AH!-- Note that in case of 3-D nesting the south and north boundaries are 6433 !AH!-- not exchanged because the nest domain is not cyclic. 6434 !AH IF ( pdims(2) > 1 ) THEN 6435 !AH! 6436 !AH!-- From south to north 6437 !AH CALL MPI_SENDRECV( workarrc_t(0,jcsw+1,iclw), 1, & 6438 !AH workarrc_t_exchange_type_x, psouth, 2, & 6439 !AH workarrc_t(0,jcnw,iclw), 1, & 6440 !AH workarrc_t_exchange_type_x, pnorth, 2, & 6441 !AH comm2d, status, ierr ) 6442 !AH! 6443 !AH!-- From north to south 6444 !AH CALL MPI_SENDRECV( workarrc_t(0,jcnw-1,iclw), 1, & 6445 !AH workarrc_t_exchange_type_x, pnorth, 3, & 6446 !AH workarrc_t(0,jcsw,iclw), 1, & 6447 !AH workarrc_t_exchange_type_x, psouth, 3, & 6448 !AH comm2d, status, ierr ) 6449 !AH ENDIF 6450 !AH#endif 6451 !AH 6452 6484 nw = noff 6485 ! write(9,"(a1,2x,5(i3,2x))") var, k, kc(k), noff, n, nw 6486 ! flush(9) 6487 ! 6488 !-- Substitute the necessary parent-grid data to the work array. 6489 !-- Note that the dimension of workarrc_t is (0:2,jcsw:jcnw,iclw:icrw), 6490 !-- And the jc?w and ic?w-index bounds depend on the location of the PE- 6491 !-- subdomain relative to the side boundaries. 6492 iclc = iclw + 1 6493 icrc = icrw - 1 6494 jcsc = jcsw + 1 6495 jcnc = jcnw - 1 6496 IF ( bc_dirichlet_l ) THEN 6497 iclc = iclw 6498 ENDIF 6499 IF ( bc_dirichlet_r ) THEN 6500 icrc = icrw 6501 ENDIF 6502 IF ( bc_dirichlet_s ) THEN 6503 jcsc = jcsw 6504 ENDIF 6505 IF ( bc_dirichlet_n ) THEN 6506 jcnc = jcnw 6507 ENDIF 6508 workarrc_t = 0.0_wp 6509 workarrc_t(0:2,jcsc:jcnc,iclc:icrc) = fc(kc(k):kc(k)+2,jcsc:jcnc,iclc:icrc) 6510 ! 6511 !-- Left-right exchange if more than one PE subdomain in the x-direction. 6512 !-- Note that in case of 3-D nesting the left and right boundaries are 6513 !-- not exchanged because the nest domain is not cyclic. 6514 #if defined( __parallel ) 6515 IF ( pdims(1) > 1 ) THEN 6516 ! 6517 !-- From left to right 6518 CALL MPI_SENDRECV( workarrc_t(0,jcsw,iclw+1), 1, & 6519 workarrc_t_exchange_type_y, pleft, 0, & 6520 workarrc_t(0,jcsw,icrw), 1, & 6521 workarrc_t_exchange_type_y, pright, 0, & 6522 comm2d, status, ierr ) 6523 ! 6524 !-- From right to left 6525 CALL MPI_SENDRECV( workarrc_t(0,jcsw,icrw-1), 1, & 6526 workarrc_t_exchange_type_y, pright, 1, & 6527 workarrc_t(0,jcsw,iclw), 1, & 6528 workarrc_t_exchange_type_y, pleft, 1, & 6529 comm2d, status, ierr ) 6530 ENDIF 6531 ! 6532 !-- South-north exchange if more than one PE subdomain in the y-direction. 6533 !-- Note that in case of 3-D nesting the south and north boundaries are 6534 !-- not exchanged because the nest domain is not cyclic. 6535 IF ( pdims(2) > 1 ) THEN 6536 ! 6537 !-- From south to north 6538 CALL MPI_SENDRECV( workarrc_t(0,jcsw+1,iclw), 1, & 6539 workarrc_t_exchange_type_x, psouth, 2, & 6540 workarrc_t(0,jcnw,iclw), 1, & 6541 workarrc_t_exchange_type_x, pnorth, 2, & 6542 comm2d, status, ierr ) 6543 ! 6544 !-- From north to south 6545 CALL MPI_SENDRECV( workarrc_t(0,jcnw-1,iclw), 1, & 6546 workarrc_t_exchange_type_x, pnorth, 3, & 6547 workarrc_t(0,jcsw,iclw), 1, & 6548 workarrc_t_exchange_type_x, psouth, 3, & 6549 comm2d, status, ierr ) 6550 ENDIF 6551 #endif 6552 6553 IF ( var == 'w' ) THEN 6554 ! DO l = iclw+1, icrw-1 6555 ! DO m = jcsw+1, jcnw-1 6556 DO l = iclw, icrw ! Also fill up beyond edges 6557 DO m = jcsw, jcnw ! Also fill up beyond edges 6558 6559 DO i = ifl(l), ifu(l) 6560 DO j = jfl(m), jfu(m) 6561 f(k,j,i) = workarrc_t(nw,m,l) 6562 ENDDO 6563 ENDDO 6564 6565 ENDDO 6566 ENDDO 6567 ENDIF 6568 6569 IF ( var == 'u' ) THEN 6570 6571 DO l = iclw+1, icrw ! Note that we have to go up to icrw here. 6572 DO m = jcsw+1, jcnw-1 6573 ! 6574 !-- 1st-order upwind interpolation 6575 f_interp = workarrc_t(nw,m,l) 6576 ! 6577 !-- 3rd-order upwind biased interpolation. Note that w is available 6578 !-- here since it is already interpolated. 6579 ! IF ( w(k,jfl(m),ifl(l)) + w(k,jfl(m),ifl(l)-1) < 0.0_wp ) THEN 6580 ! f_interp = c1 * workarrc_t(0,m,l) + c2 * workarrc_t(1,m,l) & 6581 ! + c3 * workarrc_t(2,m,l) 6582 ! ELSE 6583 ! f_interp = workarrc_t(nw,m,l) 6584 ! ENDIF 6585 6586 ! DO i = ifl(l), ifl(l+1)-1 6587 ! DO j = jfl(m), jfu(m) 6588 ! f(k,j,i) = fc(n,m,l) 6589 ! ENDDO 6590 ! ENDDO 6591 ! 6592 !-- First substitute only the matching-node values 6593 f(k,jfl(m):jfu(m),ifl(l)) = f_interp 6594 6595 ENDDO 6596 ENDDO 6597 ! 6598 !-- Then fill up the nodes in between with the averages 6453 6599 DO l = iclw+1, icrw-1 6454 6600 DO m = jcsw+1, jcnw-1 6455 6601 6602 DO i = ifl(l)+1, ifl(l+1)-1 6603 f(k,jfl(m):jfu(m),i) = 0.5_wp * ( f(k,jfl(m):jfu(m),ifl(l)) & 6604 + f(k,jfl(m):jfu(m),ifl(l+1)) ) 6605 ENDDO 6606 6607 ENDDO 6608 ENDDO 6609 6610 ENDIF 6611 6612 IF ( var == 'v' ) THEN 6613 6614 DO l = iclw+1, icrw-1 6615 DO m = jcsw+1, jcnw ! Note that we have to go up to jcnw here. 6616 ! 6617 !-- 1st-order upwind interpolation 6618 f_interp = workarrc_t(nw,m,l) 6619 ! 6620 !-- 3rd-order upwind biased interpolation. Note that w is available 6621 !-- here since it is already interpolated. 6622 ! IF ( w(k,jfl(m),ifl(l)) + w(k,jfl(m)-1,ifl(l)) < 0.0_wp ) THEN 6623 ! f_interp = c1 * workarrc_t(0,m,l) + c2 * workarrc_t(1,m,l) + c3 * workarrc_t(2,m,l) 6624 ! ELSE 6625 ! f_interp = workarrc_t(nw,m,l) 6626 ! ENDIF 6627 6628 ! DO i = ifl(l), ifu(l) 6629 ! DO j = jfl(m), jfl(m+1)-1 6630 ! f(k,j,i) = fc(n,m,l) 6631 ! ENDDO 6632 ! ENDDO 6633 ! 6634 !-- First substitute only the matching-node values 6635 f(k,jfl(m),ifl(l):ifu(l)) = f_interp 6636 6637 ENDDO 6638 ENDDO 6639 ! 6640 !-- Then fill up the nodes in between with the averages 6641 DO l = iclw+1, icrw-1 6642 DO m = jcsw+1, jcnw-1 6643 6644 DO j = jfl(m)+1, jfl(m+1)-1 6645 f(k,j,ifl(l):ifu(l)) = 0.5_wp * ( f(k,jfl(m),ifl(l):ifu(l)) & 6646 + f(k,jfl(m+1),ifl(l):ifu(l)) ) 6647 ENDDO 6648 6649 ENDDO 6650 ENDDO 6651 6652 ENDIF 6653 6654 IF ( var == 's' ) THEN 6655 DO l = iclw+1, icrw-1 6656 DO m = jcsw+1, jcnw-1 6657 ! 6658 !-- 1st-order upwind interpolation 6659 f_interp = workarrc_t(nw,m,l) 6660 ! 6661 !-- 3rd-order upwind biased interpolation. Note that w is available 6662 !-- here since it is already interpolated. 6663 ! IF ( w(k,jfl(m),ifl(l)) < 0.0_wp ) THEN 6664 ! f_interp = c1 * workarrc_t(0,m,l) + c2 * workarrc_t(1,m,l) + c3 * workarrc_t(2,m,l) 6665 ! ELSE 6666 ! f_interp = workarrc_t(nw,m,l) 6667 ! ENDIF 6668 6456 6669 DO i = ifl(l), ifu(l) 6457 6670 DO j = jfl(m), jfu(m) 6458 f(k,j,i) = f c(n,m,l)6671 f(k,j,i) = f_interp 6459 6672 ENDDO 6460 6673 ENDDO … … 6462 6675 ENDDO 6463 6676 ENDDO 6464 6465 ELSE IF ( var == 'u' ) THEN 6466 6467 DO l = iclw+1, icrw-1 6468 DO m = jcsw+1, jcnw-1 6469 6470 DO i = ifl(l), ifl(l+1)-1 6471 DO j = jfl(m), jfu(m) 6472 f(k,j,i) = fc(n,m,l) 6473 ENDDO 6474 ENDDO 6475 6476 ENDDO 6477 ENDDO 6478 6479 ELSE IF ( var == 'v' ) THEN 6480 6481 DO l = iclw+1, icrw-1 6482 DO m = jcsw+1, jcnw-1 6483 6484 DO i = ifl(l), ifu(l) 6485 DO j = jfl(m), jfl(m+1)-1 6486 f(k,j,i) = fc(n,m,l) 6487 ENDDO 6488 ENDDO 6489 6490 ENDDO 6491 ENDDO 6492 6493 ELSE ! scalars 6494 6495 DO l = iclw+1, icrw-1 6496 DO m = jcsw+1, jcnw-1 6497 6498 DO i = ifl(l), ifu(l) 6499 DO j = jfl(m), jfu(m) 6500 f(k,j,i) = fc(n,m,l) 6501 ENDDO 6502 ENDDO 6503 6504 ENDDO 6505 ENDDO 6506 6507 ENDIF ! var 6677 ENDIF 6508 6678 ! 6509 6679 !-- Just fill up the redundant second ghost-node layer in case of var == w.
Note: See TracChangeset
for help on using the changeset viewer.
