Changeset 2997 for palm/trunk/SOURCE/pmc_interface_mod.f90
- Timestamp:
- Apr 19, 2018 1:35:17 PM (7 years ago)
- File:
-
- 1 edited
-
palm/trunk/SOURCE/pmc_interface_mod.f90 (modified) (18 diffs)
Legend:
- Unmodified
- Added
- Removed
-
palm/trunk/SOURCE/pmc_interface_mod.f90
r2984 r2997 25 25 ! ----------------- 26 26 ! $Id$ 27 ! Mask topography grid points in anterpolation 28 ! 29 ! 2984 2018-04-18 11:51:30Z hellstea 27 30 ! Bugfix in the log-law correction initialization. Pivot node cannot any more be 28 31 ! selected from outside the subdomain in order to prevent array under/overflows. … … 511 514 ! 512 515 !-- Child-array indices to be precomputed and stored for anterpolation. 513 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: iflu !< 514 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: ifuu !< 515 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: iflo !< 516 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: ifuo !< 517 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: jflv !< 518 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: jfuv !< 519 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: jflo !< 520 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: jfuo !< 521 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: kflw !< 522 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: kfuw !< 523 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: kflo !< 524 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: kfuo !< 516 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: iflu !< child index indicating left bound of parent grid box on u-grid 517 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: ifuu !< child index indicating right bound of parent grid box on u-grid 518 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: iflo !< child index indicating left bound of parent grid box on scalar-grid 519 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: ifuo !< child index indicating right bound of parent grid box on scalar-grid 520 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: jflv !< child index indicating south bound of parent grid box on v-grid 521 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: jfuv !< child index indicating north bound of parent grid box on v-grid 522 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: jflo !< child index indicating south bound of parent grid box on scalar-grid 523 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: jfuo !< child index indicating north bound of parent grid box on scalar-grid 524 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: kflw !< child index indicating lower bound of parent grid box on w-grid 525 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: kfuw !< child index indicating upper bound of parent grid box on w-grid 526 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: kflo !< child index indicating lower bound of parent grid box on scalar-grid 527 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: kfuo !< child index indicating upper bound of parent grid box on scalar-grid 525 528 ! 526 529 !-- Number of fine-grid nodes inside coarse-grid ij-faces 527 530 !-- to be precomputed for anterpolation. 528 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:,:) :: ijfc_u !< 529 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:,:) :: ijfc_v !< 530 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:,:) :: ijfc_s !< 531 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: kfc_w !< 532 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:) :: kfc_s !< 531 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: ijkfc_u !< number of child grid boxes contribution to a parent grid box, u-grid 532 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: ijkfc_v !< number of child grid boxes contribution to a parent grid box, v-grid 533 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: ijkfc_w !< number of child grid boxes contribution to a parent grid box, w-grid 534 INTEGER(iwp), SAVE, ALLOCATABLE, DIMENSION(:,:,:) :: ijkfc_s !< number of child grid boxes contribution to a parent grid box, scalar-grid 533 535 534 536 INTEGER(iwp), DIMENSION(3) :: parent_grid_info_int !< … … 2978 2980 2979 2981 INTEGER(iwp) :: i !< Fine-grid index 2980 INTEGER(iwp) :: ifc_o !<2981 INTEGER(iwp) :: ifc_u !<2982 2982 INTEGER(iwp) :: ii !< Coarse-grid index 2983 2983 INTEGER(iwp) :: istart !< … … 3039 3039 ALLOCATE( kfuo(0:kctu) ) 3040 3040 3041 ALLOCATE( ijfc_u(jcs:jcn,icl:icr) ) 3042 ALLOCATE( ijfc_v(jcs:jcn,icl:icr) ) 3043 ALLOCATE( ijfc_s(jcs:jcn,icl:icr) ) 3044 ALLOCATE( kfc_w(0:kctw) ) 3045 ALLOCATE( kfc_s(0:kctu) ) 3041 ALLOCATE( ijkfc_u(0:kctu,jcs:jcn,icl:icr) ) 3042 ALLOCATE( ijkfc_v(0:kctu,jcs:jcn,icl:icr) ) 3043 ALLOCATE( ijkfc_w(0:kctw,jcs:jcn,icl:icr) ) 3044 ALLOCATE( ijkfc_s(0:kctu,jcs:jcn,icl:icr) ) 3045 ijkfc_u = 0 3046 ijkfc_v = 0 3047 ijkfc_w = 0 3048 ijkfc_s = 0 3046 3049 ! 3047 3050 !-- i-indices of u for each ii-index value … … 3171 3174 !-- This information is needed in anterpolation. 3172 3175 DO ii = icl, icr 3173 ifc_u = ifuu(ii) - iflu(ii) + 13174 ifc_o = ifuo(ii) - iflo(ii) + 13175 3176 DO jj = jcs, jcn 3176 ijfc_u(jj,ii) = ifc_u * ( jfuo(jj) - jflo(jj) + 1 ) 3177 ijfc_v(jj,ii) = ifc_o * ( jfuv(jj) - jflv(jj) + 1 ) 3178 ijfc_s(jj,ii) = ifc_o * ( jfuo(jj) - jflo(jj) + 1 ) 3177 DO kk = 0, kctu 3178 ! 3179 !-- u-component 3180 DO i = iflu(ii), ifuu(ii) 3181 DO j = jflo(jj), jfuo(jj) 3182 DO k = kflo(kk), kfuo(kk) 3183 ijkfc_u(kk,jj,ii) = ijkfc_u(kk,jj,ii) + MERGE( 1, 0, & 3184 BTEST( wall_flags_0(k,j,i), 1 ) ) 3185 ENDDO 3186 ENDDO 3187 ENDDO 3188 ! 3189 !-- v-component 3190 DO i = iflo(ii), ifuo(ii) 3191 DO j = jflv(jj), jfuv(jj) 3192 DO k = kflo(kk), kfuo(kk) 3193 ijkfc_v(kk,jj,ii) = ijkfc_v(kk,jj,ii) + MERGE( 1, 0, & 3194 BTEST( wall_flags_0(k,j,i), 2 ) ) 3195 ENDDO 3196 ENDDO 3197 ENDDO 3198 ! 3199 !-- scalars 3200 DO i = iflo(ii), ifuo(ii) 3201 DO j = jflo(jj), jfuo(jj) 3202 DO k = kflo(kk), kfuo(kk) 3203 ijkfc_s(kk,jj,ii) = ijkfc_s(kk,jj,ii) + MERGE( 1, 0, & 3204 BTEST( wall_flags_0(k,j,i), 0 ) ) 3205 ENDDO 3206 ENDDO 3207 ENDDO 3208 ENDDO 3209 3210 DO kk = 0, kctw 3211 ! 3212 !-- w-component 3213 DO i = iflo(ii), ifuo(ii) 3214 DO j = jflo(jj), jfuo(jj) 3215 DO k = kflw(kk), kfuw(kk) 3216 ijkfc_w(kk,jj,ii) = ijkfc_w(kk,jj,ii) + MERGE( 1, 0, & 3217 BTEST( wall_flags_0(k,j,i), 3 ) ) 3218 ENDDO 3219 ENDDO 3220 ENDDO 3221 ENDDO 3222 3179 3223 ENDDO 3180 ENDDO3181 DO kk = 0, kctw3182 kfc_w(kk) = kfuw(kk) - kflw(kk) + 13183 ENDDO3184 DO kk = 0, kctu3185 kfc_s(kk) = kfuo(kk) - kflo(kk) + 13186 3224 ENDDO 3187 3225 ! … … 5012 5050 5013 5051 CALL pmci_anterp_tophat( u, uc, kctu, iflu, ifuu, jflo, jfuo, kflo, & 5014 kfuo, ij fc_u, kfc_s, 'u' )5052 kfuo, ijkfc_u, 'u' ) 5015 5053 CALL pmci_anterp_tophat( v, vc, kctu, iflo, ifuo, jflv, jfuv, kflo, & 5016 kfuo, ij fc_v, kfc_s, 'v' )5054 kfuo, ijkfc_v, 'v' ) 5017 5055 CALL pmci_anterp_tophat( w, wc, kctw, iflo, ifuo, jflo, jfuo, kflw, & 5018 kfuw, ij fc_s,kfc_w, 'w' )5056 kfuw, ijkfc_w, 'w' ) 5019 5057 ! 5020 5058 !-- Anterpolation of TKE and dissipation rate if parent and child are in … … 5022 5060 IF ( rans_mode_parent .AND. rans_mode ) THEN 5023 5061 CALL pmci_anterp_tophat( e, ec, kctu, iflo, ifuo, jflo, jfuo, kflo, & 5024 kfuo, ij fc_s,kfc_s, 'e' )5062 kfuo, ijkfc_s, 'e' ) 5025 5063 ! 5026 5064 !-- Anterpolation of dissipation rate only if TKE-e closure is applied. 5027 5065 IF ( rans_tke_e ) THEN 5028 5066 CALL pmci_anterp_tophat( diss, dissc, kctu, iflo, ifuo, jflo, jfuo,& 5029 kflo, kfuo, ij fc_s,kfc_s, 'diss' )5067 kflo, kfuo, ijkfc_s, 'diss' ) 5030 5068 ENDIF 5031 5069 … … 5034 5072 IF ( .NOT. neutral ) THEN 5035 5073 CALL pmci_anterp_tophat( pt, ptc, kctu, iflo, ifuo, jflo, jfuo, kflo, & 5036 kfuo, ij fc_s,kfc_s, 'pt' )5074 kfuo, ijkfc_s, 'pt' ) 5037 5075 ENDIF 5038 5076 … … 5040 5078 5041 5079 CALL pmci_anterp_tophat( q, q_c, kctu, iflo, ifuo, jflo, jfuo, kflo, & 5042 kfuo, ij fc_s,kfc_s, 'q' )5080 kfuo, ijkfc_s, 'q' ) 5043 5081 5044 5082 IF ( cloud_physics .AND. microphysics_morrison ) THEN 5045 5083 5046 5084 CALL pmci_anterp_tophat( qc, qcc, kctu, iflo, ifuo, jflo, jfuo, & 5047 kflo, kfuo, ij fc_s,kfc_s, 'qc' )5085 kflo, kfuo, ijkfc_s, 'qc' ) 5048 5086 5049 5087 CALL pmci_anterp_tophat( nc, ncc, kctu, iflo, ifuo, jflo, jfuo, & 5050 kflo, kfuo, ij fc_s, kfc_s,'nc' )5088 kflo, kfuo, ijkfc_s, 'nc' ) 5051 5089 5052 5090 ENDIF … … 5055 5093 5056 5094 CALL pmci_anterp_tophat( qr, qrc, kctu, iflo, ifuo, jflo, jfuo, & 5057 kflo, kfuo, ij fc_s,kfc_s, 'qr' )5095 kflo, kfuo, ijkfc_s, 'qr' ) 5058 5096 5059 5097 CALL pmci_anterp_tophat( nr, nrc, kctu, iflo, ifuo, jflo, jfuo, & 5060 kflo, kfuo, ij fc_s,kfc_s, 'nr' )5098 kflo, kfuo, ijkfc_s, 'nr' ) 5061 5099 5062 5100 ENDIF … … 5066 5104 IF ( passive_scalar ) THEN 5067 5105 CALL pmci_anterp_tophat( s, sc, kctu, iflo, ifuo, jflo, jfuo, kflo, & 5068 kfuo, ij fc_s,kfc_s, 's' )5106 kfuo, ijkfc_s, 's' ) 5069 5107 ENDIF 5070 5108 … … 5074 5112 chem_spec_c(:,:,:,n), & 5075 5113 kctu, iflo, ifuo, jflo, jfuo, kflo, & 5076 kfuo, ij fc_s,kfc_s, 's' )5114 kfuo, ijkfc_s, 's' ) 5077 5115 ENDDO 5078 5116 ENDIF … … 5625 5663 5626 5664 SUBROUTINE pmci_anterp_tophat( f, fc, kct, ifl, ifu, jfl, jfu, kfl, kfu, & 5627 ij fc,kfc, var )5665 ijkfc, var ) 5628 5666 ! 5629 5667 !-- Anterpolation of internal-node values to be used as the parent-domain … … 5635 5673 5636 5674 CHARACTER(LEN=*), INTENT(IN) :: var !< identifyer for treated variable 5675 5676 INTEGER(iwp), DIMENSION(0:kct,jcs:jcn,icl:icr), INTENT(IN) :: ijkfc !< number of child grid points contributing to a parent grid box 5637 5677 5638 5678 INTEGER(iwp) :: i !< Running index x-direction - fine-grid … … 5647 5687 INTEGER(iwp) :: kk !< Running index z-direction - coarse grid 5648 5688 INTEGER(iwp) :: kcb = 0 !< Bottom boundary index for anterpolation along z 5649 INTEGER(iwp) :: m !< Running index surface type 5650 INTEGER(iwp) :: nfc !< 5689 INTEGER(iwp) :: var_flag !< bit number used to flag topography on respective grid 5651 5690 5652 5691 INTEGER(iwp), INTENT(IN) :: kct !< Top boundary index for anterpolation along z 5653 5692 5654 INTEGER(iwp), DIMENSION(icl:icr), INTENT(IN) :: ifl !< Indicates start index of child cells belonging to certain parent cell - x direction 5655 INTEGER(iwp), DIMENSION(icl:icr), INTENT(IN) :: ifu !< Indicates end index of child cells belonging to certain parent cell - x direction 5656 INTEGER(iwp), DIMENSION(jcs:jcn,icl:icr), INTENT(IN) :: ijfc !< 5657 INTEGER(iwp), DIMENSION(jcs:jcn), INTENT(IN) :: jfl !< Indicates start index of child cells belonging to certain parent cell - y direction 5658 INTEGER(iwp), DIMENSION(jcs:jcn), INTENT(IN) :: jfu !< Indicates start index of child cells belonging to certain parent cell - y direction 5659 INTEGER(iwp), DIMENSION(0:kct), INTENT(IN) :: kfc !< 5660 INTEGER(iwp), DIMENSION(0:kct), INTENT(IN) :: kfl !< Indicates start index of child cells belonging to certain parent cell - z direction 5661 INTEGER(iwp), DIMENSION(0:kct), INTENT(IN) :: kfu !< Indicates start index of child cells belonging to certain parent cell - z direction 5693 INTEGER(iwp), DIMENSION(icl:icr), INTENT(IN) :: ifl !< Indicates start index of child cells belonging to certain parent cell - x direction 5694 INTEGER(iwp), DIMENSION(icl:icr), INTENT(IN) :: ifu !< Indicates end index of child cells belonging to certain parent cell - x direction 5695 INTEGER(iwp), DIMENSION(jcs:jcn), INTENT(IN) :: jfl !< Indicates start index of child cells belonging to certain parent cell - y direction 5696 INTEGER(iwp), DIMENSION(jcs:jcn), INTENT(IN) :: jfu !< Indicates start index of child cells belonging to certain parent cell - y direction 5697 INTEGER(iwp), DIMENSION(0:kct), INTENT(IN) :: kfl !< Indicates start index of child cells belonging to certain parent cell - z direction 5698 INTEGER(iwp), DIMENSION(0:kct), INTENT(IN) :: kfu !< Indicates start index of child cells belonging to certain parent cell - z direction 5662 5699 5663 5700 REAL(wp) :: cellsum !< sum of respective child cells belonging to parent cell … … 5719 5756 jcnm = jcn - nhln - 1 5720 5757 ENDIF 5721 ENDIF 5758 ENDIF 5759 ! 5760 !-- Set masking bit for topography flags 5761 IF ( var == 'u' ) THEN 5762 var_flag = 1 5763 ELSEIF ( var == 'v' ) THEN 5764 var_flag = 2 5765 ELSEIF ( var == 'w' ) THEN 5766 var_flag = 3 5767 ELSE 5768 var_flag = 0 5769 ENDIF 5770 5722 5771 ! 5723 5772 !-- Note that ii, jj, and kk are coarse-grid indices and i,j, and k … … 5729 5778 !-- terrain too 5730 5779 DO kk = kcb, kct - 1 5731 ! 5732 !-- ijfc and kfc are precomputed in pmci_init_anterp_tophat 5733 nfc = ijfc(jj,ii) * kfc(kk) 5780 5734 5781 cellsum = 0.0_wp 5735 5782 DO i = ifl(ii), ifu(ii) 5736 5783 DO j = jfl(jj), jfu(jj) 5737 5784 DO k = kfl(kk), kfu(kk) 5738 cellsum = cellsum + f(k,j,i) 5785 cellsum = cellsum + MERGE( f(k,j,i), 0.0_wp, & 5786 BTEST( wall_flags_0(k,j,i), var_flag ) ) 5739 5787 ENDDO 5740 5788 ENDDO … … 5743 5791 !-- Spatial under-relaxation. 5744 5792 fra = frax(ii) * fray(jj) * fraz(kk) 5745 fc(kk,jj,ii) = ( 1.0_wp - fra ) * fc(kk,jj,ii) + & 5746 fra * cellsum / REAL( nfc, KIND = wp ) 5793 ! 5794 !-- In case all child grid points are inside topography, i.e. 5795 !-- cellsum is zero, also parent solution would have zero values at 5796 !-- that point, which may cause problems in particular for the 5797 !-- temperature. Therefore, in case cellsum is zero, keep the 5798 !-- parent solution at this point. 5799 fc(kk,jj,ii) = MERGE( ( 1.0_wp - fra ) * fc(kk,jj,ii) + & 5800 fra * cellsum / & 5801 REAL( ijkfc(kk,jj,ii), KIND = wp ), & 5802 fc(kk,jj,ii), & 5803 cellsum /= 0.0_wp ) 5747 5804 5748 5805 ENDDO
Note: See TracChangeset
for help on using the changeset viewer.
