SCALE-RM
Functions/Subroutines
scale_atmos_phy_mp_common Module Reference

module ATMOSPHERE / Physics Cloud Microphysics - Common More...

Functions/Subroutines

subroutine, public atmos_phy_mp_negative_fixer (DENS, RHOT, QTRC, I_QV, limit_negative)
 Negative fixer. More...
 
subroutine, public atmos_phy_mp_saturation_adjustment (RHOE_t, QTRC_t, RHOE0, QTRC0, DENS0, I_QV, I_QC, I_QI, flag_liquid)
 Saturation adjustment. More...
 
subroutine, public atmos_phy_mp_precipitation (QA_MP, QS_MP, qflx, vterm, DENS, MOMZ, MOMX, MOMY, RHOE, QTRC, temp, CVq, dt, vt_fixed)
 precipitation transport More...
 

Detailed Description

module ATMOSPHERE / Physics Cloud Microphysics - Common

Description
Common module for Cloud Microphysics Sedimentation/Precipitation and Saturation adjustment
Author
Team SCALE
History
  • 2012-12-23 (H.Yashiro) [new]

Function/Subroutine Documentation

◆ atmos_phy_mp_negative_fixer()

subroutine, public scale_atmos_phy_mp_common::atmos_phy_mp_negative_fixer ( real(rp), dimension(ka,ia,ja), intent(inout)  DENS,
real(rp), dimension(ka,ia,ja), intent(inout)  RHOT,
real(rp), dimension(ka,ia,ja,qa), intent(inout)  QTRC,
integer, intent(in)  I_QV,
real(rp), intent(in)  limit_negative 
)

Negative fixer.

Definition at line 64 of file scale_atmos_phy_mp_common.F90.

References scale_grid_index::ieb, scale_stdio::io_fid_log, scale_stdio::io_l, scale_grid_index::isb, scale_grid_index::jeb, scale_grid_index::jsb, scale_grid_index::ke, scale_grid_index::ks, scale_process::prc_mpistop(), scale_process::prc_myrank, scale_prof::prof_rapend(), scale_prof::prof_rapstart(), scale_atmos_hydrometeor::qhe, and scale_atmos_hydrometeor::qhs.

Referenced by scale_atmos_phy_mp_kessler::atmos_phy_mp_kessler(), scale_atmos_phy_mp_suzuki10::atmos_phy_mp_suzuki10(), and scale_atmos_phy_mp_tomita08::atmos_phy_mp_tomita08().

64  use scale_process, only: &
65  prc_myrank, &
67  use scale_atmos_hydrometeor, only: &
68  qhs, &
69  qhe
70  implicit none
71 
72  real(RP), intent(inout) :: DENS(KA,IA,JA)
73  real(RP), intent(inout) :: RHOT(KA,IA,JA)
74  real(RP), intent(inout) :: QTRC(KA,IA,JA,QA)
75  integer, intent(in) :: I_QV
76  real(RP), intent(in) :: limit_negative
77 
78  real(RP) :: diffq
79  real(RP) :: diffq_check(KA,IA,JA)
80  real(RP) :: diffq_min
81 
82  integer :: k, i, j, iq
83  !---------------------------------------------------------------------------
84 
85  call prof_rapstart('MP_filter', 3)
86 
87  !$omp parallel do default(none) private(i,j,k,iq,diffq) OMP_SCHEDULE_ collapse(2) &
88  !$omp shared(JSB,JEB,ISB,IEB,KS,KE,QHS,QHE,QTRC,DENS,RHOT,I_QV,diffq_check)
89  do j = jsb, jeb
90  do i = isb, ieb
91  do k = ks, ke
92 
93  diffq = 0.0_rp
94  do iq = qhs, qhe
95  ! total hydrometeor (before correction)
96  diffq = diffq + qtrc(k,i,j,iq)
97  ! remove negative value of hydrometeors (mass)
98  qtrc(k,i,j,iq) = max( qtrc(k,i,j,iq), 0.0_rp )
99  enddo
100 
101  do iq = qhs, qhe
102  ! difference between before and after correction
103  diffq = diffq - qtrc(k,i,j,iq)
104  enddo
105 
106  ! Compensate for the lack of hydrometeors by the water vapor
107  qtrc(k,i,j,i_qv) = qtrc(k,i,j,i_qv) + diffq
108  diffq_check(k,i,j) = diffq
109 
110  ! TODO: We have to consider energy conservation (but very small)
111 
112  ! remove negative value of water vapor (mass)
113  diffq = qtrc(k,i,j,i_qv)
114  qtrc(k,i,j,i_qv) = max( qtrc(k,i,j,i_qv), 0.0_rp )
115  diffq = diffq - qtrc(k,i,j,i_qv)
116 
117  ! Apply correction to total density
118  ! TODO: We have to consider energy conservation (but very small)
119  dens(k,i,j) = dens(k,i,j) * ( 1.0_rp - diffq ) ! diffq is negative
120  rhot(k,i,j) = rhot(k,i,j) * ( 1.0_rp - diffq )
121 
122  enddo
123  enddo
124  enddo
125 
126  diffq_min = minval( diffq_check(ks:ke,isb:ieb,jsb:jeb) )
127 
128  if ( abs(limit_negative) > 0.0_rp &
129  .AND. abs(limit_negative) < abs(diffq_min) ) then
130  if( io_l ) write(io_fid_log,*) 'xxx [MP_negative_fixer] large negative is found.'
131  write(*,*) 'xxx [MP_negative_fixer] large negative is found. rank = ', prc_myrank
132 
133  do j = jsb, jeb
134  do i = isb, ieb
135  do k = ks, ke
136  if ( abs(limit_negative) < abs(diffq_check(k,i,j)) &
137  .OR. abs(qtrc(k,i,j,i_qv)) < abs(diffq_check(k,i,j)) ) then
138  if( io_l ) write(io_fid_log,*) &
139  'xxx k,i,j,value(QHYD,QV) = ', k, i, j, diffq_check(k,i,j), qtrc(k,i,j,i_qv)
140  endif
141  enddo
142  enddo
143  enddo
144  if( io_l ) write(io_fid_log,*) 'xxx criteria: total negative hydrometeor < ', abs(limit_negative)
145 
146  call prc_mpistop
147  endif
148 
149  call prof_rapend('MP_filter', 3)
150 
151  return
subroutine, public prc_mpistop
Abort MPI.
integer, public jeb
integer, public ke
end point of inner domain: z, local
integer, public ieb
module PROCESS
integer, public ks
start point of inner domain: z, local
integer, public prc_myrank
process num in local communicator
integer, public isb
integer, public jsb
Here is the call graph for this function:
Here is the caller graph for this function:

◆ atmos_phy_mp_saturation_adjustment()

subroutine, public scale_atmos_phy_mp_common::atmos_phy_mp_saturation_adjustment ( real(rp), dimension(ka,ia,ja), intent(inout)  RHOE_t,
real(rp), dimension(ka,ia,ja,qa), intent(inout)  QTRC_t,
real(rp), dimension (ka,ia,ja), intent(inout)  RHOE0,
real(rp), dimension (ka,ia,ja,qa), intent(inout)  QTRC0,
real(rp), dimension (ka,ia,ja), intent(in)  DENS0,
integer, intent(in)  I_QV,
integer, intent(in)  I_QC,
integer, intent(in)  I_QI,
logical, intent(in)  flag_liquid 
)

Saturation adjustment.

Definition at line 167 of file scale_atmos_phy_mp_common.F90.

References scale_const::const_undef, scale_atmos_saturation::cvovr_ice, scale_atmos_saturation::cvovr_liq, scale_grid_index::ie, scale_grid_index::ieb, scale_grid_index::is, scale_grid_index::isb, scale_grid_index::je, scale_grid_index::jeb, scale_grid_index::js, scale_grid_index::jsb, scale_grid_index::ke, scale_grid_index::ks, scale_atmos_hydrometeor::lhf, scale_atmos_hydrometeor::lhv, scale_atmos_saturation::lovr_ice, scale_atmos_saturation::lovr_liq, scale_process::prc_mpistop(), scale_prof::prof_rapend(), scale_prof::prof_rapstart(), scale_tracer::qa, scale_precision::rp, scale_time::time_dtsec_atmos_phy_mp, scale_tracer::tracer_cv, scale_tracer::tracer_mass, and scale_tracer::tracer_r.

Referenced by scale_atmos_phy_mp_kessler::atmos_phy_mp_kessler(), and scale_atmos_phy_mp_tomita08::atmos_phy_mp_tomita08().

167 #ifdef DRY
168  use scale_const, only: &
169  undef => const_undef
170 #endif
171  use scale_time, only: &
173  use scale_atmos_thermodyn, only: &
174  thermodyn_qd => atmos_thermodyn_qd, &
175  thermodyn_cv => atmos_thermodyn_cv, &
176  thermodyn_temp_pres_e => atmos_thermodyn_temp_pres_e
177  use scale_atmos_hydrometeor, only: &
178  lhv, &
179  lhf
180  use scale_atmos_saturation, only: &
181  saturation_dens2qsat_liq => atmos_saturation_dens2qsat_liq, &
182  saturation_dens2qsat_all => atmos_saturation_dens2qsat_all
183  implicit none
184 
185  real(RP), intent(inout) :: RHOE_t(KA,IA,JA) ! tendency rhoe [J/m3/s]
186  real(RP), intent(inout) :: QTRC_t(KA,IA,JA,QA) ! tendency tracer [kg/kg/s]
187  real(RP), intent(inout) :: RHOE0 (KA,IA,JA) ! density * internal energy [J/m3]
188  real(RP), intent(inout) :: QTRC0 (KA,IA,JA,QA) ! mass concentration [kg/kg]
189  real(RP), intent(in) :: DENS0 (KA,IA,JA) ! density [kg/m3]
190  integer, intent(in) :: I_QV ! index for water vapor
191  integer, intent(in) :: I_QC ! index for water cloud
192  integer, intent(in) :: I_QI ! index for ice cloud
193  logical, intent(in) :: flag_liquid ! use scheme only for the liquid water?
194 
195  ! working
196  real(RP) :: TEMP0 (KA,IA,JA)
197  real(RP) :: PRES0 (KA,IA,JA)
198  real(RP) :: QDRY0 (KA,IA,JA)
199  real(RP) :: CVtot (KA,IA,JA)
200 
201  real(RP) :: Emoist(KA,IA,JA) ! moist internal energy
202  real(RP) :: QSUM1 (KA,IA,JA) ! QV+QC+QI
203  real(RP) :: TEMP1 (KA,IA,JA)
204 
205  real(RP) :: RHOE1 (KA,IA,JA)
206  real(RP) :: QTRC1 (KA,IA,JA,QA)
207  real(RP) :: rdt
208 
209  integer :: k, i, j, iq
210  !---------------------------------------------------------------------------
211 
212 #ifndef DRY
213 
214  call prof_rapstart('MP_Saturation_adjustment', 2)
215 
216  rdt = 1.0_rp / dt
217 
218  !$omp parallel do private(i,j,k,iq) OMP_SCHEDULE_ collapse(4)
219  do iq = 1, qa
220  do j = jsb, jeb
221  do i = isb, ieb
222  do k = ks, ke
223  qtrc1(k,i,j,iq) = qtrc0(k,i,j,iq)
224  enddo
225  enddo
226  enddo
227  enddo
228 
229  call thermodyn_temp_pres_e( temp0(:,:,:), & ! [OUT]
230  pres0(:,:,:), & ! [OUT]
231  dens0(:,:,:), & ! [IN]
232  rhoe0(:,:,:), & ! [IN]
233  qtrc0(:,:,:,:), & ! [IN]
234  tracer_cv(:), & ! [IN]
235  tracer_r(:), & ! [IN]
236  tracer_mass(:) ) ! [IN]
237 
238  ! qdry dont change through the process
239  call thermodyn_qd( qdry0(:,:,:), & ! [OUT]
240  qtrc0(:,:,:,:), & ! [IN]
241  tracer_mass(:) ) ! [IN]
242 
243  call thermodyn_cv( cvtot(:,:,:), & ! [OUT]
244  qtrc0(:,:,:,:), & ! [IN]
245  tracer_cv(:) , & ! [IN]
246  qdry0(:,:,:) ) ! [IN]
247 
248  if ( i_qi <= 0 .OR. flag_liquid ) then ! warm rain
249 
250  ! Turn QC into QV with consistency of moist internal energy
251  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
252  do j = jsb, jeb
253  do i = isb, ieb
254  do k = ks, ke
255  emoist(k,i,j) = temp0(k,i,j) * cvtot(k,i,j) &
256  + qtrc1(k,i,j,i_qv) * lhv
257 
258  qsum1(k,i,j) = qtrc1(k,i,j,i_qv) &
259  + qtrc1(k,i,j,i_qc)
260 
261  qtrc1(k,i,j,i_qv) = qsum1(k,i,j)
262  qtrc1(k,i,j,i_qc) = 0.0_rp
263  enddo
264  enddo
265  enddo
266 
267  call thermodyn_cv( cvtot(:,:,:), & ! [OUT]
268  qtrc1(:,:,:,:), & ! [IN]
269  tracer_cv(:), & ! [IN]
270  qdry0(:,:,:) ) ! [IN]
271 
272  ! new temperature (after QC evaporation)
273  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
274  do j = jsb, jeb
275  do i = isb, ieb
276  do k = ks, ke
277  temp1(k,i,j) = ( emoist(k,i,j) - qtrc1(k,i,j,i_qv) * lhv ) / cvtot(k,i,j)
278  enddo
279  enddo
280  enddo
281 
282  call moist_conversion_liq( temp1(:,:,:), & ! [INOUT]
283  qtrc1(:,:,:,:), & ! [INOUT]
284  dens0(:,:,:), & ! [IN]
285  qsum1(:,:,:), & ! [IN]
286  qdry0(:,:,:), & ! [IN]
287  emoist(:,:,:), & ! [IN]
288  i_qv, & ! [IN]
289  i_qc ) ! [IN]
290 
291  else ! cold rain
292 
293  ! Turn QC & QI into QV with consistency of moist internal energy
294  !$omp parallel do default(none) private(i,j,k) OMP_SCHEDULE_ collapse(2) &
295  !$omp shared(JSB,JEB,ISB,IEB,KS,KE,Emoist,TEMP0,CVtot,QTRC1,LHV,LHF,QSUM1,I_QV,I_QC,I_QI)
296  do j = jsb, jeb
297  do i = isb, ieb
298  do k = ks, ke
299  emoist(k,i,j) = temp0(k,i,j) * cvtot(k,i,j) &
300  + qtrc1(k,i,j,i_qv) * lhv &
301  - qtrc1(k,i,j,i_qi) * lhf
302 
303  qsum1(k,i,j) = qtrc1(k,i,j,i_qv) &
304  + qtrc1(k,i,j,i_qc) &
305  + qtrc1(k,i,j,i_qi)
306 
307  qtrc1(k,i,j,i_qv) = qsum1(k,i,j)
308  qtrc1(k,i,j,i_qc) = 0.0_rp
309  qtrc1(k,i,j,i_qi) = 0.0_rp
310  enddo
311  enddo
312  enddo
313 
314  call thermodyn_cv( cvtot(:,:,:), & ! [OUT]
315  qtrc1(:,:,:,:), & ! [IN]
316  tracer_cv(:), & ! [IN]
317  qdry0(:,:,:) ) ! [IN]
318 
319  ! new temperature (after QC & QI evaporation)
320  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
321  do j = jsb, jeb
322  do i = isb, ieb
323  do k = ks, ke
324  temp1(k,i,j) = ( emoist(k,i,j) - qtrc1(k,i,j,i_qv) * lhv ) / cvtot(k,i,j)
325  enddo
326  enddo
327  enddo
328 
329  call moist_conversion_all( temp1(:,:,:), & ! [INOUT]
330  qtrc1(:,:,:,:), & ! [INOUT]
331  dens0(:,:,:), & ! [IN]
332  qsum1(:,:,:), & ! [IN]
333  qdry0(:,:,:), & ! [IN]
334  emoist(:,:,:), & ! [IN]
335  i_qv, & ! [IN]
336  i_qc, & ! [IN]
337  i_qi )
338 
339  endif
340 
341  call thermodyn_cv( cvtot(:,:,:), & ! [OUT]
342  qtrc1(:,:,:,:), & ! [IN]
343  tracer_cv(:), & ! [IN]
344  qdry0(:,:,:) ) ! [IN]
345 
346 
347  ! mass & energy update
348 
349  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(3)
350  do j = js, je
351  do i = is, ie
352  do k = ks, ke
353  qtrc_t(k,i,j,i_qv) = qtrc_t(k,i,j,i_qv) + ( qtrc1(k,i,j,i_qv) - qtrc0(k,i,j,i_qv) ) * rdt
354 
355  qtrc0(k,i,j,i_qv) = qtrc1(k,i,j,i_qv)
356  enddo
357  enddo
358  enddo
359 
360  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(3)
361  do j = js, je
362  do i = is, ie
363  do k = ks, ke
364  qtrc_t(k,i,j,i_qc) = qtrc_t(k,i,j,i_qc) + ( qtrc1(k,i,j,i_qc) - qtrc0(k,i,j,i_qc) ) * rdt
365 
366  qtrc0(k,i,j,i_qc) = qtrc1(k,i,j,i_qc)
367  enddo
368  enddo
369  enddo
370 
371  ! mass & energy update
372  if ( i_qi > 0 .AND. (.not. flag_liquid) ) then
373  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(3)
374  do j = js, je
375  do i = is, ie
376  do k = ks, ke
377  qtrc_t(k,i,j,i_qi) = qtrc_t(k,i,j,i_qi) + ( qtrc1(k,i,j,i_qi) - qtrc0(k,i,j,i_qi) ) * rdt
378 
379  qtrc0(k,i,j,i_qi) = qtrc1(k,i,j,i_qi)
380  enddo
381  enddo
382  enddo
383  end if
384 
385 
386  !$omp parallel do default(none) private(i,j,k) OMP_SCHEDULE_ collapse(2) &
387  !$omp shared(JS,JE,IS,IE,KS,KE,RHOE1,DENS0,TEMP1,CVtot,RHOE_t,RHOE0,rdt)
388  do j = js, je
389  do i = is, ie
390  do k = ks, ke
391  rhoe1(k,i,j) = dens0(k,i,j) * temp1(k,i,j) * cvtot(k,i,j)
392 
393  rhoe_t(k,i,j) = rhoe_t(k,i,j) + ( rhoe1(k,i,j) - rhoe0(k,i,j) ) * rdt
394 
395  rhoe0(k,i,j) = rhoe1(k,i,j)
396  enddo
397  enddo
398  enddo
399 
400  call prof_rapend ('MP_Saturation_adjustment', 2)
401 
402 #else
403  rhoe_t = undef
404  qtrc_t = undef
405  rhoe0 = undef
406  qtrc0 = undef
407 #endif
408  return
integer, public is
start point of inner domain: x, local
integer, public je
end point of inner domain: y, local
module ATMOSPHERE / Saturation adjustment
integer, public jeb
real(rp), dimension(qa_max), public tracer_r
real(dp), public time_dtsec_atmos_phy_mp
time interval of physics(microphysics) [sec]
Definition: scale_time.F90:41
integer, public ke
end point of inner domain: z, local
real(rp), dimension(qa_max), public tracer_cv
real(rp), public const_undef
Definition: scale_const.F90:43
integer, public ieb
integer, public js
start point of inner domain: y, local
module TIME
Definition: scale_time.F90:15
module CONSTANT
Definition: scale_const.F90:14
integer, public ks
start point of inner domain: z, local
integer, public ie
end point of inner domain: x, local
module ATMOSPHERE / Thermodynamics
integer, public isb
integer, public jsb
real(rp), dimension(qa_max), public tracer_mass
Here is the call graph for this function:
Here is the caller graph for this function:

◆ atmos_phy_mp_precipitation()

subroutine, public scale_atmos_phy_mp_common::atmos_phy_mp_precipitation ( integer, intent(in)  QA_MP,
integer, intent(in)  QS_MP,
real(rp), dimension (ka,ia,ja,qa_mp-1), intent(out)  qflx,
real(rp), dimension(ka,ia,ja,qa_mp-1), intent(inout)  vterm,
real(rp), dimension (ka,ia,ja), intent(inout)  DENS,
real(rp), dimension (ka,ia,ja), intent(inout)  MOMZ,
real(rp), dimension (ka,ia,ja), intent(inout)  MOMX,
real(rp), dimension (ka,ia,ja), intent(inout)  MOMY,
real(rp), dimension (ka,ia,ja), intent(inout)  RHOE,
real(rp), dimension (ka,ia,ja,qa), intent(inout)  QTRC,
real(rp), dimension (ka,ia,ja), intent(in)  temp,
real(rp), dimension (qa), intent(in)  CVq,
real(dp), intent(in)  dt,
logical, intent(in), optional  vt_fixed 
)

precipitation transport

Definition at line 763 of file scale_atmos_phy_mp_common.F90.

References scale_const::const_grav, scale_gridtrans::gtrans_j33g, scale_grid_index::ie, scale_grid_index::is, scale_grid_index::je, scale_grid_index::js, scale_grid_index::ke, scale_grid_index::ks, scale_prof::prof_rapend(), scale_prof::prof_rapstart(), scale_tracer::qa, scale_grid_real::real_cz, scale_grid_real::real_fz, and scale_tracer::tracer_mass.

Referenced by scale_atmos_phy_mp_kessler::atmos_phy_mp_kessler(), scale_atmos_phy_mp_sn14::atmos_phy_mp_sn14(), scale_atmos_phy_mp_suzuki10::atmos_phy_mp_suzuki10(), and scale_atmos_phy_mp_tomita08::atmos_phy_mp_tomita08().

763  use scale_const, only: &
764  grav => const_grav
765  use scale_grid_real, only: &
766  real_cz, &
767  real_fz
768  use scale_gridtrans, only: &
769  j33g => gtrans_j33g
770  use scale_comm, only: &
771  comm_vars8, &
772  comm_wait
773  implicit none
774 
775  integer, intent(in) :: QA_MP
776  integer, intent(in) :: QS_MP
777  real(RP), intent(out) :: qflx (KA,IA,JA,QA_MP-1)
778  real(RP), intent(inout) :: vterm(KA,IA,JA,QA_MP-1) ! terminal velocity of cloud mass
779  real(RP), intent(inout) :: DENS (KA,IA,JA)
780  real(RP), intent(inout) :: MOMZ (KA,IA,JA)
781  real(RP), intent(inout) :: MOMX (KA,IA,JA)
782  real(RP), intent(inout) :: MOMY (KA,IA,JA)
783  real(RP), intent(inout) :: RHOE (KA,IA,JA)
784  real(RP), intent(inout) :: QTRC (KA,IA,JA,QA)
785  real(RP), intent(in) :: temp (KA,IA,JA)
786  real(RP), intent(in) :: CVq (QA)
787  real(DP), intent(in) :: dt
788  logical, intent(in), optional :: vt_fixed
789 
790  real(RP) :: rhoq (KA,IA,JA,QA) ! rho * q before precipitation
791  real(RP) :: eflx (KA,IA,JA)
792  real(RP) :: rfdz (KA,IA,JA)
793  real(RP) :: rcdz (KA,IA,JA)
794  real(RP) :: rcdz_u(KA,IA,JA)
795  real(RP) :: rcdz_v(KA,IA,JA)
796 
797  integer :: k, i, j
798  integer :: iq, iqa
799  logical :: vt_fixed_
800  !---------------------------------------------------------------------------
801 
802  call prof_rapstart('MP_Precipitation', 2)
803 
804  if ( present(vt_fixed) ) then
805  vt_fixed_ = vt_fixed
806  else
807  vt_fixed_ = .false.
808  end if
809 
810  do iq = 1, qa_mp-1
811  iqa = qs_mp + iq
812 
813  if( tracer_mass(iqa) == 0.0_rp ) cycle
814 
815  if( .NOT. vt_fixed_ ) call comm_vars8( vterm(:,:,:,iq), iq )
816 
817  call comm_vars8( qtrc(:,:,:,iqa), qa_mp+iq )
818  enddo
819 
820  ! tracer/energy transport by falldown
821  ! 1st order upwind, forward euler, velocity is always negative
822 
823  !$omp parallel do private(i,j,k) OMP_SCHEDULE_ collapse(2)
824  do j = js, je
825  do i = is, ie
826  rfdz(ks-1,i,j) = 1.0_rp / ( real_cz(ks,i,j) - real_fz(ks-1,i,j) )
827  do k = ks, ke
828  rfdz(k,i,j) = 1.0_rp / ( real_cz(k+1,i,j) - real_cz(k ,i,j) )
829  rcdz(k,i,j) = 1.0_rp / ( real_fz(k ,i,j) - real_fz(k-1,i,j) )
830 
831  rcdz_u(k,i,j) = 2.0_rp / ( ( real_fz(k,i+1,j) - real_fz(k-1,i+1,j) ) &
832  + ( real_fz(k,i ,j) - real_fz(k-1,i ,j) ) )
833  rcdz_v(k,i,j) = 2.0_rp / ( ( real_fz(k,i,j+1) - real_fz(k-1,i,j+1) ) &
834  + ( real_fz(k,i,j ) - real_fz(k-1,i,j ) ) )
835  enddo
836  enddo
837  enddo
838 
839  do iq = 1, qa_mp-1
840  iqa = qs_mp + iq
841  if ( tracer_mass(iqa) == 0.0_rp ) cycle
842 
843  if ( .not. vt_fixed_ ) then
844  call comm_wait( vterm(:,:,:,iq), iq )
845  endif
846  call comm_wait( qtrc(:,:,:,iqa), qa_mp+iq )
847 
848  !$omp parallel do default(none) &
849  !$omp shared(JS,JE,IS,IE,KS,KE,qflx,iq,vterm,DENS,QTRC,iqa,J33G,eflx,temp,CVq,RHOE,dt) &
850  !$omp shared(rcdz,GRAV,rfdz,MOMZ,MOMX,rcdz_u,MOMY,rcdz_v) &
851  !$omp private(i,j,k) OMP_SCHEDULE_ collapse(2)
852  do j = js, je
853  do i = is, ie
854 
855  !--- mass flux for each mass tracer, upwind with vel < 0
856  do k = ks-1, ke-1
857  qflx(k,i,j,iq) = vterm(k+1,i,j,iq) * dens(k+1,i,j) * qtrc(k+1,i,j,iqa) * j33g
858  enddo
859  qflx(ke,i,j,iq) = 0.0_rp
860 
861  !--- internal energy
862  eflx(ks-1,i,j) = qflx(ks-1,i,j,iq) * temp(ks,i,j) * cvq(iqa)
863  do k = ks, ke-1
864  eflx(k,i,j) = qflx(k,i,j,iq) * temp(k+1,i,j) * cvq(iqa)
865  rhoe(k,i,j) = rhoe(k,i,j) - dt * ( eflx(k,i,j) - eflx(k-1,i,j) ) * rcdz(k,i,j)
866  enddo
867  eflx(ke,i,j) = 0.0_rp
868  rhoe(ke,i,j) = rhoe(ke,i,j) - dt * ( - eflx(ke-1,i,j) ) * rcdz(ke,i,j)
869 
870  !--- potential energy
871  eflx(ks-1,i,j) = qflx(ks-1,i,j,iq) * grav / rfdz(ks-1,i,j)
872  do k = ks, ke-1
873  eflx(k,i,j) = qflx(k,i,j,iq) * grav / rfdz(k,i,j)
874  rhoe(k,i,j) = rhoe(k,i,j) - dt * ( eflx(k,i,j) - eflx(k-1,i,j) ) * rcdz(k,i,j)
875  enddo
876  rhoe(ke,i,j) = rhoe(ke,i,j) - dt * ( - eflx(ke-1,i,j) ) * rcdz(ke,i,j)
877 
878  !--- momentum z (half level)
879  do k = ks-1, ke-1
880  eflx(k,i,j) = 0.25_rp * ( vterm(k+1,i,j,iq ) + vterm(k,i,j,iq ) ) &
881  * ( qtrc(k+1,i,j,iqa) + qtrc(k,i,j,iqa) ) &
882  * momz(k,i,j)
883  enddo
884  do k = ks, ke-1
885  momz(k,i,j) = momz(k,i,j) - dt * ( eflx(k+1,i,j) - eflx(k,i,j) ) * rfdz(k,i,j)
886  enddo
887 
888  !--- momentum x
889  eflx(ks-1,i,j) = 0.25_rp * ( vterm(ks,i,j,iq ) + vterm(ks,i+1,j,iq ) ) &
890  * ( qtrc(ks,i,j,iqa) + qtrc(ks,i+1,j,iqa) ) &
891  * momx(ks,i,j)
892  do k = ks, ke-1
893  eflx(k,i,j) = 0.25_rp * ( vterm(k+1,i,j,iq ) + vterm(k+1,i+1,j,iq ) ) &
894  * ( qtrc(k+1,i,j,iqa) + qtrc(k+1,i+1,j,iqa) ) &
895  * momx(k+1,i,j)
896  momx(k,i,j) = momx(k,i,j) - dt * ( eflx(k,i,j) - eflx(k-1,i,j) ) * rcdz_u(k,i,j)
897  enddo
898  momx(ke,i,j) = momx(ke,i,j) - dt * ( - eflx(ke-1,i,j) ) * rcdz_u(ke,i,j)
899 
900  !--- momentum y
901  eflx(ks-1,i,j) = 0.25_rp * ( vterm(ks,i,j,iq ) + vterm(ks,i,j+1,iq ) ) &
902  * ( qtrc(ks,i,j,iqa) + qtrc(ks,i,j+1,iqa) ) &
903  * momy(ks,i,j)
904  do k = ks, ke-1
905  eflx(k,i,j) = 0.25_rp * ( vterm(k+1,i,j,iq ) + vterm(k+1,i,j+1,iq ) ) &
906  * ( qtrc(k+1,i,j,iqa) + qtrc(k+1,i,j+1,iqa) ) &
907  * momy(k+1,i,j)
908  momy(k,i,j) = momy(k,i,j) - dt * ( eflx(k,i,j) - eflx(k-1,i,j) ) * rcdz_v(k,i,j)
909  enddo
910  momy(ke,i,j) = momy(ke,i,j) - dt * ( - eflx(ke-1,i,j) ) * rcdz_v(ke,i,j)
911 
912  enddo
913  enddo
914 
915  enddo ! falling (water mass & number) tracer
916 
917  !--- save previous value
918  do iqa = 1, qa
919  do j = js, je
920  do i = is, ie
921  do k = ks, ke
922  rhoq(k,i,j,iqa) = qtrc(k,i,j,iqa) * dens(k,i,j)
923  enddo
924  enddo
925  enddo
926  enddo ! all tracer
927 
928  do iq = 1, qa_mp-1
929  iqa = qs_mp + iq
930 
931  !--- mass flux for each tracer, upwind with vel < 0
932  do j = js, je
933  do i = is, ie
934  do k = ks-1, ke-1
935  qflx(k,i,j,iq) = vterm(k+1,i,j,iq) * rhoq(k+1,i,j,iqa)
936  enddo
937  qflx(ke,i,j,iq) = 0.0_rp
938  enddo
939  enddo
940 
941  if ( tracer_mass(iqa) == 0.0_rp ) cycle
942 
943  !--- update total density
944  do j = js, je
945  do i = is, ie
946  do k = ks, ke
947  dens(k,i,j) = dens(k,i,j) - dt * ( qflx(k,i,j,iq) - qflx(k-1,i,j,iq) ) * rcdz(k,i,j)
948  enddo
949  enddo
950  enddo
951 
952  enddo ! falling (water mass & number) tracer
953 
954  !--- update falling tracer
955  do iq = 1, qa_mp-1
956  iqa = qs_mp + iq
957  do j = js, je
958  do i = is, ie
959  do k = ks, ke
960  rhoq(k,i,j,iqa) = rhoq(k,i,j,iqa) - dt * ( qflx(k,i,j,iq) - qflx(k-1,i,j,iq) ) * rcdz(k,i,j)
961  enddo
962  enddo
963  enddo
964  enddo ! falling (water mass & number) tracer
965 
966  !--- update tracer ratio with updated total density)
967  do iqa = 1, qa
968  do j = js, je
969  do i = is, ie
970  do k = ks, ke
971  qtrc(k,i,j,iqa) = rhoq(k,i,j,iqa) / dens(k,i,j)
972  enddo
973  enddo
974  enddo
975  enddo ! all tracer
976 
977  call prof_rapend ('MP_Precipitation', 2)
978 
979  return
integer, public is
start point of inner domain: x, local
integer, public je
end point of inner domain: y, local
integer, public ke
end point of inner domain: z, local
real(rp), public gtrans_j33g
(3,3) element of Jacobian matrix * {G}^1/2
real(rp), dimension(:,:,:), allocatable, public real_fz
geopotential height [m] (cell face )
real(rp), dimension(:,:,:), allocatable, public real_cz
geopotential height [m] (cell center)
module GRIDTRANS
module GRID (real space)
module COMMUNICATION
Definition: scale_comm.F90:23
real(rp), public const_grav
standard acceleration of gravity [m/s2]
Definition: scale_const.F90:48
integer, public js
start point of inner domain: y, local
module CONSTANT
Definition: scale_const.F90:14
integer, public ks
start point of inner domain: z, local
integer, public ie
end point of inner domain: x, local
real(rp), dimension(qa_max), public tracer_mass
Here is the call graph for this function:
Here is the caller graph for this function: