FE-Project
All Classes Namespaces Files Functions Variables Pages
Functions/Subroutines
scale_atm_phy_tb_dgm_globalsmg Module Reference

module FElib / Atmosphere / Physics turbulence More...

Functions/Subroutines

subroutine, public atm_phy_tb_dgm_globalsmg_init (mesh, shallow_atm_approx)
 
subroutine, public atm_phy_tb_dgm_globalsmg_final ()
 
subroutine, public atm_phy_tb_dgm_globalsmg_cal_grad (t11, t12, t13, t21, t22, t23, t31, t32, t33, df1, df2, df3, tke, nu, kh, ddens_, momx_, momy_, momz_, drhot_, dens_hyd, pres_hyd, pres, pt, dx, dy, dz, sx, sy, sz, lift, lmesh, elem, lmesh2d, elem2d, is_bound)
 Calculate parameterized stress tensor and eddy heat flux with turbulent model.
 
subroutine, public atm_phy_tb_dgm_globalsmg_cal_grad_qtrc (dfq1, dfq2, dfq3, drdx, drdy, drdz, kh, qtrc, ddens, dens_hyd, dx, dy, dz, sx, sy, sz, lift, lmesh, elem, lmesh2d, elem2d, is_bound, cal_grad_dens)
 Calculate parameterized diffusive mass flux of tracer with turbulent model.
 
subroutine, public atm_phy_tb_dgm_globalsmg_cal_tend (momx_t, momy_t, momz_t, rhot_t, t11, t12, t13, t21, t22, t23, t31, t32, t33, df1, df2, df3, nu, kh, ddens_, momx_, momy_, momz_, drhot_, dens_hyd, pres_hyd, pres_, pt_, dx, dy, dz, sx, sy, sz, lift, lmesh, elem, lmesh2d, elem2d, is_bound)
 Calculate tendecies with turbulent model.
 
subroutine, public atm_phy_tb_dgm_globalsmg_cal_tend_qtrc (rhoq_t, dfq1, dfq2, dfq3, kh, ddens_, dens_hyd, dx, dy, dz, sx, sy, sz, lift, lmesh, elem, lmesh2d, elem2d, is_bound)
 Calculate tendecies of tracer density with turbulent model.
 

Detailed Description

module FElib / Atmosphere / Physics turbulence

Description
Sub-grid scale turbulnce process for global domain Smagorinsky-type
Author
Yuta Kawai, Team SCALE
Reference
  • Brown et al., 1994: Large-eddy simulaition of stable atmospheric boundary layers with a revised stochastic subgrid model. Roy. Meteor. Soc., 120, 1485-1512
  • Scotti et al., 1993: Generalized Smagorinsky model for anisotropic grids. Phys. Fluids A, 5, 2306-2308
  • Nishizawa et al., 2015: Influence of grid aspect ratio on planetary boundary layer turbulence in large-eddy simulations Geosci. Model Dev., 8, 3393–3419
NAMELIST
  • PARAM_ATMOS_PHY_TB_DGM_GLOBALSMG
    nametypedefault valuecomment
    CS real(RP) 0.13_RP Smagorinsky constant (Scotti et al. 1993)
    NU_MAX real(RP) 10000.0_RP
    FILTER_FAC real(RP) 2.0_RP
    CONSISTENT_TKE logical .true.

History Output
No history output

Function/Subroutine Documentation

◆ atm_phy_tb_dgm_globalsmg_init()

subroutine, public scale_atm_phy_tb_dgm_globalsmg::atm_phy_tb_dgm_globalsmg_init ( class(meshbase3d), intent(in) mesh,
logical, intent(in) shallow_atm_approx )

Definition at line 109 of file scale_atm_phy_tb_dgm_globalsmg.F90.

110 implicit none
111 class(MeshBase3D), intent(in) :: mesh
112 logical, intent(in) :: shallow_atm_approx
113
114 logical :: consistent_tke = .true.
115
116 namelist / param_atmos_phy_tb_dgm_globalsmg / &
117 cs, &
118 nu_max, &
119 filter_fac, &
120 consistent_tke
121
122 integer :: ierr
123 !--------------------------------------------------------------------
124
125 log_newline
126 log_info("ATMOS_PHY_TB_dgm_globalsmg_setup",*) 'Setup'
127 log_info("ATMOS_PHY_TB_dgm_globalsmg_setup",*) 'Smagorinsky-type Eddy Viscocity Model'
128
129 !--- read namelist
130 rewind(io_fid_conf)
131 read(io_fid_conf,nml=param_atmos_phy_tb_dgm_globalsmg,iostat=ierr)
132 if( ierr < 0 ) then !--- missing
133 log_info("ATMOS_PHY_TB_dgm_globalsmg_setup",*) 'Not found namelist. Default used.'
134 elseif( ierr > 0 ) then !--- fatal error
135 log_error("ATMOS_PHY_TB_dgm_globalsmg_setup",*) 'Not appropriate names in namelist PARAM_ATMOS_PHY_TB_DGM_SMG. Check!'
136 call prc_abort
137 endif
138 log_nml(param_atmos_phy_tb_dgm_globalsmg)
139
140 rprn = 1.0_rp / prn
141 rric = 1.0_rp / ric
142 onemprnovric = ( 1.0_rp - prn ) * rric
143
144 if ( consistent_tke ) then
145 tke_fac = 1.0_rp
146 else
147 tke_fac = 0.0_rp
148 end if
149
150 !--
151 shallow_atm_approx_flag = shallow_atm_approx
152 if ( shallow_atm_approx_flag ) then
153 shapro_coef = 0.0_rp
154 else
155 shapro_coef = 1.0_rp
156 end if
157
158 return

◆ atm_phy_tb_dgm_globalsmg_final()

subroutine, public scale_atm_phy_tb_dgm_globalsmg::atm_phy_tb_dgm_globalsmg_final

Definition at line 162 of file scale_atm_phy_tb_dgm_globalsmg.F90.

163 implicit none
164 !--------------------------------------------------------------------
165
166 return

◆ atm_phy_tb_dgm_globalsmg_cal_grad()

subroutine, public scale_atm_phy_tb_dgm_globalsmg::atm_phy_tb_dgm_globalsmg_cal_grad ( real(rp), dimension(elem%np,lmesh%nea), intent(out) t11,
real(rp), dimension(elem%np,lmesh%nea), intent(out) t12,
real(rp), dimension(elem%np,lmesh%nea), intent(out) t13,
real(rp), dimension(elem%np,lmesh%nea), intent(out) t21,
real(rp), dimension(elem%np,lmesh%nea), intent(out) t22,
real(rp), dimension(elem%np,lmesh%nea), intent(out) t23,
real(rp), dimension(elem%np,lmesh%nea), intent(out) t31,
real(rp), dimension(elem%np,lmesh%nea), intent(out) t32,
real(rp), dimension(elem%np,lmesh%nea), intent(out) t33,
real(rp), dimension(elem%np,lmesh%nea), intent(out) df1,
real(rp), dimension(elem%np,lmesh%nea), intent(out) df2,
real(rp), dimension(elem%np,lmesh%nea), intent(out) df3,
real(rp), dimension(elem%np,lmesh%nea), intent(out) tke,
real(rp), dimension(elem%np,lmesh%nea), intent(out) nu,
real(rp), dimension(elem%np,lmesh%nea), intent(out) kh,
real(rp), dimension(elem%np,lmesh%nea), intent(in) ddens_,
real(rp), dimension (elem%np,lmesh%nea), intent(in) momx_,
real(rp), dimension (elem%np,lmesh%nea), intent(in) momy_,
real(rp), dimension (elem%np,lmesh%nea), intent(in) momz_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) drhot_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) dens_hyd,
real(rp), dimension(elem%np,lmesh%nea), intent(in) pres_hyd,
real(rp), dimension(elem%np,lmesh%nea), intent(in) pres,
real(rp), dimension(elem%np,lmesh%nea), intent(in) pt,
type(sparsemat), intent(in) dx,
type(sparsemat), intent(in) dy,
type(sparsemat), intent(in) dz,
type(sparsemat), intent(in) sx,
type(sparsemat), intent(in) sy,
type(sparsemat), intent(in) sz,
type(sparsemat), intent(in) lift,
class(localmesh3d), intent(in) lmesh,
class(elementbase3d), intent(in) elem,
class(localmesh2d), intent(in) lmesh2d,
class(elementbase2d), intent(in) elem2d,
logical, dimension(elem%nfptot,lmesh%ne), intent(in) is_bound )

Calculate parameterized stress tensor and eddy heat flux with turbulent model.

Parameters
[out]t11(1,1) component of stress tensor
[out]t12(1,2) component of stress tensor
[out]t13(1,3) component of stress tensor
[out]t21(2,1) component of stress tensor
[out]t22(2,2) component of stress tensor
[out]t23(2,3) component of stress tensor
[out]t31(3,1) component of stress tensor
[out]t32(3,2) component of stress tensor
[out]t33(3,3) component of stress tensor
[out]df1Diffusive heat flux in x1 direction / density
[out]df2Diffusive heat flux in x2 direction / density
[out]df3Diffusive heat flux in x3 direction / density
[out]tkeParameterized turbulent kinetic energy
[out]nuEddy viscosity
[out]khEddy diffusivity
[in]ddens_Density perturbation
[in]momx_Momentum in x1 direction
[in]momy_Momentum in x2 direction
[in]momz_Momentum in x3 direction
[in]drhot_Density x potential temperature perturbation
[in]dens_hydReference pressure in hydrostatic balance
[in]pres_hydReference density in hydrostatic balance
[in]presPressure
[in]ptPotential temperature
[in]dzDifferential matrix managed by sparse matrix type
[in]szStiffness matrix managed by sparse matrix type
[in]liftLifting matrix managed by sparse matrix type
[in]is_boundFlag whether nodes are located at domain boundaries

Definition at line 172 of file scale_atm_phy_tb_dgm_globalsmg.F90.

180
181 use scale_atm_phy_tb_dgm_common, only: &
183 use scale_cubedsphere_coord_cnv, only: &
185
186 implicit none
187
188 class(LocalMesh3D), intent(in) :: lmesh
189 class(ElementBase3D), intent(in) :: elem
190 class(LocalMesh2D), intent(in) :: lmesh2D
191 class(ElementBase2D), intent(in) :: elem2D
192 real(RP), intent(out) :: T11(elem%Np,lmesh%NeA)
193 real(RP), intent(out) :: T12(elem%Np,lmesh%NeA)
194 real(RP), intent(out) :: T13(elem%Np,lmesh%NeA)
195 real(RP), intent(out) :: T21(elem%Np,lmesh%NeA)
196 real(RP), intent(out) :: T22(elem%Np,lmesh%NeA)
197 real(RP), intent(out) :: T23(elem%Np,lmesh%NeA)
198 real(RP), intent(out) :: T31(elem%Np,lmesh%NeA)
199 real(RP), intent(out) :: T32(elem%Np,lmesh%NeA)
200 real(RP), intent(out) :: T33(elem%Np,lmesh%NeA)
201 real(RP), intent(out) :: DF1(elem%Np,lmesh%NeA)
202 real(RP), intent(out) :: DF2(elem%Np,lmesh%NeA)
203 real(RP), intent(out) :: DF3(elem%Np,lmesh%NeA)
204 real(RP), intent(out) :: TKE(elem%Np,lmesh%NeA)
205 real(RP), intent(out) :: Nu(elem%Np,lmesh%NeA)
206 real(RP), intent(out) :: Kh(elem%Np,lmesh%NeA)
207 real(RP), intent(in) :: DDENS_(elem%Np,lmesh%NeA)
208 real(RP), intent(in) :: MOMX_ (elem%Np,lmesh%NeA)
209 real(RP), intent(in) :: MOMY_ (elem%Np,lmesh%NeA)
210 real(RP), intent(in) :: MOMZ_ (elem%Np,lmesh%NeA)
211 real(RP), intent(in) :: DRHOT_(elem%Np,lmesh%NeA)
212 real(RP), intent(in) :: DENS_hyd(elem%Np,lmesh%NeA)
213 real(RP), intent(in) :: PRES_hyd(elem%Np,lmesh%NeA)
214 real(RP), intent(in) :: PRES(elem%Np,lmesh%NeA)
215 real(RP), intent(in) :: PT(elem%Np,lmesh%NeA)
216 type(SparseMat), intent(in) :: Dx, Dy, Dz
217 type(SparseMat), intent(in) :: Sx, Sy, Sz
218 type(SparseMat), intent(in) :: Lift
219 logical, intent(in) :: is_bound(elem%NfpTot,lmesh%Ne)
220
221 real(RP) :: Fx(elem%Np), Fy(elem%Np), Fz(elem%Np), LiftDelFlx(elem%Np)
222 real(RP) :: DENS(elem%Np), RDENS(elem%Np), RHOT(elem%Np), Q(elem%Np)
223 real(RP) :: DdensDxi(elem%Np,3)
224 real(RP) :: DqDxi_(elem%Np,2)
225 real(RP) :: DVelDxi(elem%Np,3,3)
226 real(RP) :: del_flux_rho (elem%NfpTot,lmesh%Ne,3)
227 real(RP) :: del_flux_mom (elem%NfpTot,lmesh%Ne,3,3)
228 real(RP) :: del_flux_rhot(elem%NfpTot,lmesh%Ne,3)
229
230 real(RP) :: Ri ! local gradient Richardson number
231 real(RP) :: S2 ! (2SijSij)^1/2
232 real(RP) :: fm ! factor in eddy viscosity which represents the stability dependence of
233 ! the Brown et al (1994)'s subgrid model
234 real(RP) :: Pr ! Parandtl number (=Nu/Kh= fm/fh)
235
236 real(RP) :: lambda (elem%Np,lmesh%Ne) ! basic mixing length
237 real(RP) :: lambda_r(elem%Np) ! characteristic subgrid length scale
238 real(RP) :: E(elem%Np) ! subgrid kinetic energy
239 real(RP) :: C1(elem%Np) ! factor in the relation with energy disspation rate, lambda_r, and E
240
241 real(RP) :: G11(elem%Np), G12(elem%Np), G22(elem%Np), G33(elem%Np)
242 real(RP) :: X2D(elem2D%Np,lmesh2D%Ne), Y2D(elem2D%Np,lmesh2D%Ne)
243 real(RP) :: X(elem%Np), Y(elem%Np), Rdel2(elem%Np)
244
245 real(RP) :: S11(elem%Np), S12(elem%Np), S22(elem%Np), S23(elem%Np), S31(elem%Np), S33(elem%Np)
246 real(RP) :: DivOvThree(elem%Np), TKEMulTwoOvThree(elem%Np)
247
248 real(RP) :: Sabs_tmp(elem%Np)
249 real(RP) :: SIJ(elem%Np,3,3)
250 real(RP) :: G_ij(elem%Np,3,3)
251
252 real(RP) :: rgam2(elem%Np)
253 real(RP) :: R(elem%Np)
254
255 integer :: ke, ke2D
256 integer :: p
257 integer :: i, j
258 !--------------------------------------------------------------------
259
260 call cal_del_flux_grad( del_flux_rho, del_flux_mom, del_flux_rhot, & ! (out)
261 ddens_, momx_, momy_, momz_, drhot_, dens_hyd, pres_hyd, pt, & ! (in)
262 lmesh%normal_fn(:,:,1), lmesh%normal_fn(:,:,2), lmesh%normal_fn(:,:,3), & ! (in)
263 lmesh%vmapM, lmesh%vmapP, & ! (in)
264 lmesh, elem, is_bound ) ! (in)
265
266 call atm_phy_tb_dgm_common_calc_lambda( lambda, & ! (out)
267 cs, filter_fac, lmesh, elem, lmesh2d, elem2d ) ! (in)
268
269 !$omp parallel private( ke, ke2D, &
270 !$omp Fx, Fy, Fz, LiftDelFlx, &
271 !$omp DENS, RHOT, RDENS, Q, DdensDxi, DqDxi_, DVelDxi, &
272 !$omp S11, S12, S22, S23, S31, S33, DivOvThree, TKEMulTwoOvThree, &
273 !$omp p, Ri, S2, fm, Pr, lambda_r, E, C1, &
274 !$omp Sabs_tmp, i, j, SIJ, G_ij, G11, G12, G22, G33, R, rgam2, &
275 !$omp X, Y, Rdel2 )
276
277 !$omp do
278 do ke2d = lmesh2d%NeS, lmesh2d%NeE
279 x2d(:,ke2d) = tan(lmesh2d%pos_en(:,ke2d,1))
280 y2d(:,ke2d) = tan(lmesh2d%pos_en(:,ke2d,2))
281 end do
282 !$omp end do
283
284 !$omp do
285 do ke=lmesh%NeS, lmesh%NeE
286 !---
287 ke2d = lmesh%EMap3Dto2D(ke)
288
289 r(:) = rplanet * lmesh%gam(:,ke)
290 rgam2(:) = 1.0_rp / lmesh%gam(:,ke)**2
291 g11(:) = lmesh%GIJ(elem%IndexH2Dto3D,ke2d,1,1) * rgam2(:)
292 g12(:) = lmesh%GIJ(elem%IndexH2Dto3D,ke2d,1,2) * rgam2(:)
293 g22(:) = lmesh%GIJ(elem%IndexH2Dto3D,ke2d,2,2) * rgam2(:)
294 g33(:) = 1.0_rp
295
296 g_ij(:,:,:) = 0.0_rp
297 g_ij(:,1,1) = lmesh%G_ij(elem%IndexH2Dto3D,ke2d,1,1) / rgam2(:)
298 g_ij(:,2,1) = lmesh%G_ij(elem%IndexH2Dto3D,ke2d,2,1) / rgam2(:)
299 g_ij(:,1,2) = g_ij(:,2,1)
300 g_ij(:,2,2) = lmesh%G_ij(elem%IndexH2Dto3D,ke2d,2,2) / rgam2(:)
301 g_ij(:,3,3) = 1.0_rp
302
303 x(:) = x2d(elem%IndexH2Dto3D,ke2d)
304 y(:) = y2d(elem%IndexH2Dto3D,ke2d)
305 rdel2(:) = 1.0_rp / ( 1.0_rp + x(:)**2 + y(:)**2 )
306
307 dens(:) = dens_hyd(:,ke) + ddens_(:,ke)
308 rdens(:) = 1.0_rp / dens(:)
309 rhot(:) = dens(:) * pt(:,ke)
310
311 ! gradient of density
312 call sparsemat_matmul( dx, dens, fx )
313 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_rho(:,ke,1), liftdelflx )
314 ddensdxi(:,1) = lmesh%Escale(:,ke,1,1) * fx(:) + liftdelflx(:)
315
316 call sparsemat_matmul( dy, dens, fy )
317 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_rho(:,ke,2), liftdelflx )
318 ddensdxi(:,2) = lmesh%Escale(:,ke,2,2) * fy(:) + liftdelflx(:)
319
320 call sparsemat_matmul( dz, dens, fz )
321 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_rho(:,ke,3), liftdelflx )
322 ddensdxi(:,3) = lmesh%Escale(:,ke,3,3) * fz(:) + liftdelflx(:)
323
324 ! gradient of u
325 q(:) = momx_(:,ke) * rdens(:)
326
327 call sparsemat_matmul( dx, momx_(:,ke), fx )
328 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_mom(:,ke,1,1), liftdelflx )
329 dqdxi_(:,1) = ( lmesh%Escale(:,ke,1,1) * fx(:) + liftdelflx(:) - q(:) * ddensdxi(:,1) &
330 + rdel2(:) * y(:) * &
331 ( 2.0_rp * x(:) * y(:) * momx_(:,ke) - ( 1.0_rp + y(:)**2 ) * momy_(:,ke) ) & !-> u^r Gam^1_1r
332 + shapro_coef * momz_(:,ke) / r(:) &
333 ) * rdens(:)
334
335 call sparsemat_matmul( dy, momx_(:,ke), fy )
336 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_mom(:,ke,2,1), liftdelflx )
337 dqdxi_(:,2) = ( lmesh%Escale(:,ke,2,2) * fy(:) + liftdelflx(:) - q(:) * ddensdxi(:,2) &
338 + rdel2(:) * y(:) * &
339 ( - ( 1.0_rp + y(:)**2 ) * momy_(:,ke) ) & ! u^r Gam^1_2r
340 ) * rdens(:)
341
342 call sparsemat_matmul( dz, momx_(:,ke), fz )
343 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_mom(:,ke,3,1), liftdelflx )
344 dveldxi(:,3,1) = ( lmesh%Escale(:,ke,3,3) * fz(:) + liftdelflx(:) - q(:) * ddensdxi(:,3) &
345 + shapro_coef * momx_(:,ke) / r(:) & ! u^r Gam^1_3r
346 ) * rdens(:)
347
348
349 divovthree(:) = dqdxi_(:,1)
350 dveldxi(:,1,1) = g11(:) * dqdxi_(:,1) + g12(:) * dqdxi_(:,2)
351 dveldxi(:,2,1) = g12(:) * dqdxi_(:,1) + g22(:) * dqdxi_(:,2)
352
353 ! gradient of v
354 q(:) = momy_(:,ke) * rdens(:)
355
356 call sparsemat_matmul( dx, momy_(:,ke), fx )
357 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_mom(:,ke,1,2), liftdelflx )
358 dqdxi_(:,1) = ( lmesh%Escale(:,ke,1,1) * fx(:) + liftdelflx(:) - q(:) * ddensdxi(:,1) &
359 + rdel2(:) * x(:) * &
360 ( - ( 1.0_rp + x(:)**2 ) * momy_(:,ke) ) & ! u^r Gam^2_1r
361 ) * rdens(:)
362
363 call sparsemat_matmul( dy, momy_(:,ke), fy )
364 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_mom(:,ke,2,2), liftdelflx )
365 dqdxi_(:,2) = ( lmesh%Escale(:,ke,2,2) * fy(:) + liftdelflx(:) - q(:) * ddensdxi(:,2) &
366 + rdel2(:) * x(:) * &
367 ( - ( 1.0_rp + x(:)**2 ) * momx_(:,ke) + 2.0_rp * x(:) * y(:) * momy_(:,ke) ) & ! u^r Gam^2_2r
368 + shapro_coef * momz_(:,ke) / r(:) &
369 ) * rdens(:)
370
371 call sparsemat_matmul( dz, momy_(:,ke), fz )
372 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_mom(:,ke,3,2), liftdelflx )
373 dveldxi(:,3,2) = ( lmesh%Escale(:,ke,3,3) * fz(:) + liftdelflx(:) - q(:) * ddensdxi(:,3) &
374 + shapro_coef * momy_(:,ke) / r(:) & ! u^r Gam^2_3r
375 ) * rdens(:)
376 divovthree(:) = divovthree(:) + dqdxi_(:,2)
377 dveldxi(:,1,2) = g11(:) * dqdxi_(:,1) + g12(:) * dqdxi_(:,2)
378 dveldxi(:,2,2) = g12(:) * dqdxi_(:,1) + g22(:) * dqdxi_(:,2)
379
380 ! gradient of w
381 q(:) = momz_(:,ke) * rdens(:)
382
383 call sparsemat_matmul( dx, momz_(:,ke), fx )
384 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) *del_flux_mom(:,ke,1,3), liftdelflx )
385 dqdxi_(:,1) = ( lmesh%Escale(:,ke,1,1) * fx(:) + liftdelflx(:) - q(:) * ddensdxi(:,1) &
386 + shapro_coef * r(:) * rdel2(:)**2 * ( 1.0_rp + x(:)**2 ) * ( 1.0_rp + y(:)**2 ) & ! u^r Gam^3_1r
387 * ( - ( 1.0_rp + x(:)**2 ) * momx_(:,ke) + x(:) * y(:) * momy_(:,ke) ) &
388 ) * rdens(:)
389
390
391 call sparsemat_matmul( dy, momz_(:,ke), fy )
392 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_mom(:,ke,2,3), liftdelflx )
393 dqdxi_(:,2) = ( lmesh%Escale(:,ke,2,2) * fy(:) + liftdelflx(:) - q(:) * ddensdxi(:,2) &
394 + shapro_coef * r(:) * rdel2(:)**2 * ( 1.0_rp + x(:)**2 ) * ( 1.0_rp + y(:)**2 ) & ! u^r Gam^3_2r
395 * ( x(:) * y(:) * momx_(:,ke) - ( 1.0_rp + y(:)**2 ) * momy_(:,ke) ) &
396 ) * rdens(:)
397
398 call sparsemat_matmul( dz, momz_(:,ke), fz )
399 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_mom(:,ke,3,3), liftdelflx )
400 dveldxi(:,3,3) = ( lmesh%Escale(:,ke,3,3) * fz(:) + liftdelflx(:) - q(:) * ddensdxi(:,3) ) * rdens(:)
401 ! u^r Gam^3_3r = 0
402
403 divovthree(:) = divovthree(:) + dveldxi(:,3,3)
404 dveldxi(:,1,3) = g11(:) * dqdxi_(:,1) + g12(:) * dqdxi_(:,2)
405 dveldxi(:,2,3) = g12(:) * dqdxi_(:,1) + g22(:) * dqdxi_(:,2)
406
407 ! gradient of pt
408 q(:) = rhot(:) * rdens(:)
409
410 call sparsemat_matmul( dx, rhot, fx )
411 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_rhot(:,ke,1), liftdelflx )
412 dqdxi_(:,1) = ( lmesh%Escale(:,ke,1,1) * fx(:) + liftdelflx(:) - q(:) * ddensdxi(:,1) ) * rdens(:)
413
414 call sparsemat_matmul( dy, rhot, fy )
415 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_rhot(:,ke,2), liftdelflx )
416 dqdxi_(:,2) = ( lmesh%Escale(:,ke,2,2) * fy(:) + liftdelflx(:) - q(:) * ddensdxi(:,2) ) * rdens(:)
417
418 call sparsemat_matmul( dz, rhot, fz )
419 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_rhot(:,ke,3), liftdelflx )
420 df3(:,ke) = ( lmesh%Escale(:,ke,3,3) * fz(:) + liftdelflx(:) - q(:) * ddensdxi(:,3) ) * rdens(:)
421
422 df1(:,ke) = g11(:) * dqdxi_(:,1) + g12(:) * dqdxi_(:,2)
423 df2(:,ke) = g12(:) * dqdxi_(:,1) + g22(:) * dqdxi_(:,2)
424
425 ! Calculate the component of strain velocity tensor
426 s11(:) = dveldxi(:,1,1)
427 s12(:) = 0.5_rp * ( dveldxi(:,1,2) + dveldxi(:,2,1) )
428 s22(:) = dveldxi(:,2,2)
429 s23(:) = 0.5_rp * ( dveldxi(:,2,3) + dveldxi(:,3,2) )
430 s31(:) = 0.5_rp * ( dveldxi(:,1,3) + dveldxi(:,3,1) )
431 s33(:) = dveldxi(:,3,3)
432
433 sij(:,1,1) = s11(:)
434 sij(:,1,2) = s12(:)
435 sij(:,1,3) = s31(:)
436 sij(:,2,1) = s12(:)
437 sij(:,2,2) = s22(:)
438 sij(:,2,3) = s23(:)
439 sij(:,3,1) = s31(:)
440 sij(:,3,2) = s23(:)
441 sij(:,3,3) = s33(:)
442
443 sabs_tmp(:) = 0.0_rp
444 do j=1, 3
445 do i=1, 3
446 sabs_tmp(:) = sabs_tmp(:) &
447 + sij(:,i,j) * ( &
448 g_ij(:,i,1) * ( g_ij(:,j,1) * sij(:,1,1) + g_ij(:,j,2) * sij(:,1,2) + g_ij(:,j,3) * sij(:,1,3) ) &
449 + g_ij(:,i,2) * ( g_ij(:,j,1) * sij(:,2,1) + g_ij(:,j,2) * sij(:,2,2) + g_ij(:,j,3) * sij(:,2,3) ) &
450 + g_ij(:,i,3) * ( g_ij(:,j,1) * sij(:,3,1) + g_ij(:,j,2) * sij(:,3,2) + g_ij(:,j,3) * sij(:,3,3) ) )
451 end do
452 end do
453
454 ! Caclulate eddy viscosity & eddy diffusivity
455
456 do p=1, elem%Np
457 s2 = 2.0_rp * sabs_tmp(p)
458
459 ri = grav / pt(p,ke) * df3(p,ke) / max( s2, eps )
460
461 ! The Stability functions fm and fh are given by the appendix A of Brown et al. (1994).
462 if (ri < 0.0_rp ) then ! unstable
463 fm = sqrt( 1.0_rp - fmc * ri )
464 nu(p,ke) = lambda(p,ke)**2 * sqrt( s2 ) * fm
465 pr = fm / sqrt( 1.0_rp - fhb * ri ) * prn
466 else if ( ri < ric ) then ! stable
467 fm = ( 1.0_rp - ri * rric )**4
468 nu(p,ke) = lambda(p,ke)**2 * sqrt( s2 ) * fm
469 pr = prn / ( 1.0_rp - onemprnovric * ri )
470 else ! strongly stable
471 fm = 0.0_rp
472 nu(p,ke) = 0.0_rp
473 kh(p,ke) = 0.0_rp
474 pr = 1.0_rp
475 end if
476
477 if ( ri < ric ) then
478 kh(p,ke) = max( min( nu(p,ke) / pr, nu_max ), eps )
479 nu(p,ke) = max( min( nu(p,ke), nu_max ), eps )
480 pr = nu(p,ke) / kh(p,ke)
481 lambda_r(p) = lambda(p,ke) * sqrt( fm / sqrt( 1.0_rp - ri/pr ) )
482 else
483 lambda_r(p) = 0.0_rp
484 end if
485 end do
486
487 ! if ( backscatter ) then
488 ! else
489 e(:) = nu(:,ke)**3 / ( lambda_r(:)**4 + eps )
490 c1(:) = c1o
491 ! end if
492
493 ! TKE
494 tke(:,ke) = ( e(:) * lambda_r(:) / c1(:) )**twooverthree
495
496
497 ! Calculate components of parameterized stress tensor
498 ! Tij = rho * { 2 Nu * [(G^im u^im + G^jm u^jm)/2 - G^ij/3 D] - G^ij 2/3 * TKE }
499 !
500 divovthree(:) = divovthree(:) * oneoverthree
501 tkemultwoovthree(:) = twooverthree * tke(:,ke) * tke_fac
502 do p=1, elem%Np
503 t11(p,ke) = dens(p) * ( 2.0_rp * nu(p,ke) * ( s11(p) - g11(p) * divovthree(p) ) - g11(p) * tkemultwoovthree(p) )
504 t12(p,ke) = dens(p) * ( 2.0_rp * nu(p,ke) * ( s12(p) - g12(p) * divovthree(p) ) - g12(p) * tkemultwoovthree(p) )
505 t13(p,ke) = dens(p) * 2.0_rp * nu(p,ke) * s31(p)
506 end do
507 do p=1, elem%Np
508 t21(p,ke) = dens(p) * ( 2.0_rp * nu(p,ke) * ( s12(p) - g12(p) * divovthree(p) ) - g12(p) * tkemultwoovthree(p) )
509 t22(p,ke) = dens(p) * ( 2.0_rp * nu(p,ke) * ( s22(p) - g22(p) * divovthree(p) ) - g22(p) * tkemultwoovthree(p) )
510 t23(p,ke) = dens(p) * 2.0_rp * nu(p,ke) * s23(p)
511 end do
512 do p=1, elem%Np
513 t31(p,ke) = dens(p) * 2.0_rp * nu(p,ke) * s31(p)
514 t32(p,ke) = dens(p) * 2.0_rp * nu(p,ke) * s23(p)
515 t33(p,ke) = dens(p) * ( 2.0_rp * nu(p,ke) * ( s33(p) - g33(p) * divovthree(p) ) - g33(p) * tkemultwoovthree(p) )
516 end do
517
518 df1(:,ke) = kh(:,ke) * df1(:,ke)
519 df2(:,ke) = kh(:,ke) * df2(:,ke)
520 df3(:,ke) = kh(:,ke) * df3(:,ke)
521 end do
522 !$omp end do
523 !$omp end parallel
524 return
module FElib / Fluid dyn solver / Atmosphere / Physics turbulence / Common
subroutine, public atm_phy_tb_dgm_common_calc_lambda(lambda, cs, filter_fac, lmesh, elem, lmesh2d, elem2d)
module common / Coordinate conversion with a cubed-sphere
subroutine, public cubedspherecoordcnv_cs2lonlatvec(panelid, alpha, beta, gam, np, vecalpha, vecbeta, veclon, veclat, lat)

References scale_atm_phy_tb_dgm_common::atm_phy_tb_dgm_common_calc_lambda(), and scale_cubedsphere_coord_cnv::cubedspherecoordcnv_cs2lonlatvec().

◆ atm_phy_tb_dgm_globalsmg_cal_grad_qtrc()

subroutine, public scale_atm_phy_tb_dgm_globalsmg::atm_phy_tb_dgm_globalsmg_cal_grad_qtrc ( real(rp), dimension(elem%np,lmesh%nea), intent(out) dfq1,
real(rp), dimension(elem%np,lmesh%nea), intent(out) dfq2,
real(rp), dimension(elem%np,lmesh%nea), intent(out) dfq3,
real(rp), dimension(elem%np,lmesh%nea), intent(inout) drdx,
real(rp), dimension(elem%np,lmesh%nea), intent(inout) drdy,
real(rp), dimension(elem%np,lmesh%nea), intent(inout) drdz,
real(rp), dimension(elem%np,lmesh%nea), intent(in) kh,
real(rp), dimension(elem%np,lmesh%nea), intent(in) qtrc,
real(rp), dimension(elem%np,lmesh%nea), intent(in) ddens,
real(rp), dimension(elem%np,lmesh%nea), intent(in) dens_hyd,
type(sparsemat), intent(in) dx,
type(sparsemat), intent(in) dy,
type(sparsemat), intent(in) dz,
type(sparsemat), intent(in) sx,
type(sparsemat), intent(in) sy,
type(sparsemat), intent(in) sz,
type(sparsemat), intent(in) lift,
class(localmesh3d), intent(in) lmesh,
class(elementbase3d), intent(in) elem,
class(localmesh2d), intent(in) lmesh2d,
class(elementbase2d), intent(in) elem2d,
logical, dimension(elem%nfptot,lmesh%ne), intent(in) is_bound,
logical, intent(in) cal_grad_dens )

Calculate parameterized diffusive mass flux of tracer with turbulent model.

Parameters
[out]dfq1Diffusive mass flux in x1 direction / density (Kh dq/dx1)
[out]dfq2Diffusive mass flux in x2 direction / density (Kh dq/dx2)
[out]dfq3Diffusive mass flux in x3 direction / density (Kh dq/dx3)
[in,out]drdxSpatial gradient of density in x1 direction
[in,out]drdySpatial gradient of density in x2 direction
[in,out]drdzSpatial gradient of density in x3 direction
[in]khEddy diffusivity
[in]qtrcMass faction of tracer
[in]ddensDensity perturbation
[in]dens_hydReference desity in hydrostatic state
[in]dzDifferential matrix managed by sparse matrix type
[in]szStiffness matrix managed by sparse matrix type
[in]liftLifting matrix managed by sparse matrix type
[in]is_boundFlag whether nodes are located at domain boundaries
[in]cal_grad_densFlag whether spatial gradients of density are calcuated

Definition at line 530 of file scale_atm_phy_tb_dgm_globalsmg.F90.

536
537 implicit none
538
539 class(LocalMesh3D), intent(in) :: lmesh
540 class(ElementBase3D), intent(in) :: elem
541 class(LocalMesh2D), intent(in) :: lmesh2D
542 class(ElementBase2D), intent(in) :: elem2D
543 real(RP), intent(out) :: DFQ1(elem%Np,lmesh%NeA)
544 real(RP), intent(out) :: DFQ2(elem%Np,lmesh%NeA)
545 real(RP), intent(out) :: DFQ3(elem%Np,lmesh%NeA)
546 real(RP), intent(inout) :: dRdx(elem%Np,lmesh%NeA)
547 real(RP), intent(inout) :: dRdy(elem%Np,lmesh%NeA)
548 real(RP), intent(inout) :: dRdz(elem%Np,lmesh%NeA)
549 real(RP), intent(in) :: Kh(elem%Np,lmesh%NeA)
550 real(RP), intent(in) :: QTRC(elem%Np,lmesh%NeA)
551 real(RP), intent(in) :: DDENS(elem%Np,lmesh%NeA)
552 real(RP), intent(in) :: DENS_hyd(elem%Np,lmesh%NeA)
553 type(SparseMat), intent(in) :: Dx, Dy, Dz
554 type(SparseMat), intent(in) :: Sx, Sy, Sz
555 type(SparseMat), intent(in) :: Lift
556 logical, intent(in) :: is_bound(elem%NfpTot,lmesh%Ne)
557 logical, intent(in) :: cal_grad_dens
558
559 real(RP) :: G11(elem%Np), G12(elem%Np), G22(elem%Np)
560 real(RP) :: Fx(elem%Np), Fy(elem%Np), Fz(elem%Np), LiftDelFlx(elem%Np)
561 real(RP) :: DqDxi_(elem%Np,2)
562 real(RP) :: del_flux(elem%NfpTot,lmesh%Ne,3)
563 real(RP) :: del_flux_rho (elem%NfpTot,lmesh%Ne,3)
564
565 real(RP) :: DENS(elem%Np), RDENS(elem%Np)
566 real(RP) :: RHOxQTRC(elem%Np)
567
568 integer :: ke, ke2D
569 integer :: p
570 !--------------------------------------------------------------------
571
572 call cal_del_flux_grad_qtrc( del_flux, del_flux_rho, & ! (out)
573 qtrc, ddens, dens_hyd, & ! (in)
574 lmesh%normal_fn(:,:,1), lmesh%normal_fn(:,:,2), lmesh%normal_fn(:,:,3), & ! (in)
575 lmesh%vmapM, lmesh%vmapP, & ! (in)
576 lmesh, elem, is_bound, cal_grad_dens ) ! (in)
577
578 !$omp parallel private( ke, ke2D, &
579 !$omp Fx, Fy, Fz, LiftDelFlx, &
580 !$omp DqDxi_, RHOxQTRC, DENS, RDENS )
581
582 ! Calculate gradient of density
583 if ( cal_grad_dens ) then
584 !$omp do
585 do ke=lmesh%NeS, lmesh%NeE
586 dens(:) = dens_hyd(:,ke) + ddens(:,ke)
587
588 call sparsemat_matmul( dx, dens, fx )
589 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_rho(:,ke,1), liftdelflx )
590 drdx(:,ke) = lmesh%Escale(:,ke,1,1) * fx(:) + liftdelflx(:)
591
592 call sparsemat_matmul( dy, dens, fy )
593 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_rho(:,ke,2), liftdelflx )
594 drdy(:,ke) = lmesh%Escale(:,ke,2,2) * fy(:) + liftdelflx(:)
595
596 call sparsemat_matmul( dz, dens, fz )
597 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_rho(:,ke,3), liftdelflx )
598 drdz(:,ke) = lmesh%Escale(:,ke,3,3) * fz(:) + liftdelflx(:)
599 end do
600 end if
601
602 ! Calculate gradient of tracer
603 !$omp do
604 do ke=lmesh%NeS, lmesh%NeE
605 ke2d = lmesh%EMap3Dto2D(ke)
606 g11(:) = lmesh%GIJ(elem%IndexH2Dto3D,ke2d,1,1)
607 g12(:) = lmesh%GIJ(elem%IndexH2Dto3D,ke2d,1,2)
608 g22(:) = lmesh%GIJ(elem%IndexH2Dto3D,ke2d,2,2)
609
610 dens(:) = dens_hyd(:,ke) + ddens(:,ke)
611 rdens(:) = 1.0_rp / dens(:)
612 rhoxqtrc(:) = dens(:) * qtrc(:,ke)
613
614 !---
615 call sparsemat_matmul( dx, rhoxqtrc(:), fx )
616 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux(:,ke,1), liftdelflx )
617 dqdxi_(:,1) = ( lmesh%Escale(:,ke,1,1) * fx(:) + liftdelflx(:) - qtrc(:,ke) * drdx(:,ke) ) * rdens(:)
618
619 call sparsemat_matmul( dy, rhoxqtrc(:), fy )
620 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux(:,ke,2), liftdelflx )
621 dqdxi_(:,2) = ( lmesh%Escale(:,ke,2,2) * fy(:) + liftdelflx(:) - qtrc(:,ke) * drdy(:,ke) ) * rdens(:)
622
623 call sparsemat_matmul( dz, rhoxqtrc(:), fz )
624 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux(:,ke,3), liftdelflx )
625 dfq3(:,ke) = kh(:,ke) * ( lmesh%Escale(:,ke,3,3) * fz(:) + liftdelflx(:) - qtrc(:,ke) * drdz(:,ke) ) * rdens(:)
626
627 dfq1(:,ke) = kh(:,ke) * ( g11(:) * dqdxi_(:,1) + g12(:) * dqdxi_(:,2) )
628 dfq2(:,ke) = kh(:,ke) * ( g12(:) * dqdxi_(:,1) + g22(:) * dqdxi_(:,2) )
629 end do
630
631 !$omp end parallel
632
633 return

◆ atm_phy_tb_dgm_globalsmg_cal_tend()

subroutine, public scale_atm_phy_tb_dgm_globalsmg::atm_phy_tb_dgm_globalsmg_cal_tend ( real(rp), dimension(elem%np,lmesh%nea), intent(out) momx_t,
real(rp), dimension(elem%np,lmesh%nea), intent(out) momy_t,
real(rp), dimension(elem%np,lmesh%nea), intent(out) momz_t,
real(rp), dimension(elem%np,lmesh%nea), intent(out) rhot_t,
real(rp), dimension(elem%np,lmesh%nea), intent(in) t11,
real(rp), dimension(elem%np,lmesh%nea), intent(in) t12,
real(rp), dimension(elem%np,lmesh%nea), intent(in) t13,
real(rp), dimension(elem%np,lmesh%nea), intent(in) t21,
real(rp), dimension(elem%np,lmesh%nea), intent(in) t22,
real(rp), dimension(elem%np,lmesh%nea), intent(in) t23,
real(rp), dimension(elem%np,lmesh%nea), intent(in) t31,
real(rp), dimension(elem%np,lmesh%nea), intent(in) t32,
real(rp), dimension(elem%np,lmesh%nea), intent(in) t33,
real(rp), dimension(elem%np,lmesh%nea), intent(in) df1,
real(rp), dimension(elem%np,lmesh%nea), intent(in) df2,
real(rp), dimension(elem%np,lmesh%nea), intent(in) df3,
real(rp), dimension (elem%np,lmesh%nea), intent(in) nu,
real(rp), dimension (elem%np,lmesh%nea), intent(in) kh,
real(rp), dimension(elem%np,lmesh%nea), intent(in) ddens_,
real(rp), dimension (elem%np,lmesh%nea), intent(in) momx_,
real(rp), dimension (elem%np,lmesh%nea), intent(in) momy_,
real(rp), dimension (elem%np,lmesh%nea), intent(in) momz_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) drhot_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) dens_hyd,
real(rp), dimension(elem%np,lmesh%nea), intent(in) pres_hyd,
real(rp), dimension(elem%np,lmesh%nea), intent(in) pres_,
real(rp), dimension (elem%np,lmesh%nea), intent(in) pt_,
type(sparsemat), intent(in) dx,
type(sparsemat), intent(in) dy,
type(sparsemat), intent(in) dz,
type(sparsemat), intent(in) sx,
type(sparsemat), intent(in) sy,
type(sparsemat), intent(in) sz,
type(sparsemat), intent(in) lift,
class(localmesh3d), intent(in) lmesh,
class(elementbase3d), intent(in) elem,
class(localmesh2d), intent(in) lmesh2d,
class(elementbase2d), intent(in) elem2d,
logical, dimension(elem%nfptot,lmesh%ne), intent(in) is_bound )

Calculate tendecies with turbulent model.

Parameters
[out]momx_tTendency of momentum in x1 direction with turbulent model
[out]momy_tTendency of momentum in x2 direction with turbulent model
[out]momz_tTendency of momentum in x3 direction with turbulent model
[out]rhot_tTendency of density x potential temperature with turbulent model
[in]t11(1,1) component of stress tensor
[in]t12(1,2) component of stress tensor
[in]t13(1,3) component of stress tensor
[in]t21(2,1) component of stress tensor
[in]t22(2,2) component of stress tensor
[in]t23(2,3) component of stress tensor
[in]t31(3,1) component of stress tensor
[in]t32(3,2) component of stress tensor
[in]t33(3,3) component of stress tensor
[in]df1Diffusive heat flux in x1 direction / density
[in]df2Diffusive heat flux in x2 direction / density
[in]df3Diffusive heat flux in x3 direction / density
[in]nuEddy viscosity
[in]khEddy diffusivity
[in]ddens_Density perturbation
[in]momx_Momentum in x1 direction
[in]momy_Momentum in x2 direction
[in]momz_Momentum in x3 direction
[in]drhot_Density x potential temperature perturbation
[in]dens_hydReference pressure in hydrostatic balance
[in]pres_hydReference density in hydrostatic balance
[in]pres_Pressure
[in]pt_Potential temperature
[in]dzDifferential matrix managed by sparse matrix type
[in]szStiffness matrix managed by sparse matrix type
[in]liftLifting matrix managed by sparse matrix type
[in]is_boundFlag whether nodes are located at domain boundaries

Definition at line 795 of file scale_atm_phy_tb_dgm_globalsmg.F90.

804
805 implicit none
806
807 class(LocalMesh3D), intent(in) :: lmesh
808 class(ElementBase3D), intent(in) :: elem
809 class(LocalMesh2D), intent(in) :: lmesh2D
810 class(ElementBase2D), intent(in) :: elem2D
811 real(RP), intent(out) :: MOMX_t(elem%Np,lmesh%NeA)
812 real(RP), intent(out) :: MOMY_t(elem%Np,lmesh%NeA)
813 real(RP), intent(out) :: MOMZ_t(elem%Np,lmesh%NeA)
814 real(RP), intent(out) :: RHOT_t(elem%Np,lmesh%NeA)
815 real(RP), intent(in) :: T11(elem%Np,lmesh%NeA)
816 real(RP), intent(in) :: T12(elem%Np,lmesh%NeA)
817 real(RP), intent(in) :: T13(elem%Np,lmesh%NeA)
818 real(RP), intent(in) :: T21(elem%Np,lmesh%NeA)
819 real(RP), intent(in) :: T22(elem%Np,lmesh%NeA)
820 real(RP), intent(in) :: T23(elem%Np,lmesh%NeA)
821 real(RP), intent(in) :: T31(elem%Np,lmesh%NeA)
822 real(RP), intent(in) :: T32(elem%Np,lmesh%NeA)
823 real(RP), intent(in) :: T33(elem%Np,lmesh%NeA)
824 real(RP), intent(in) :: DF1(elem%Np,lmesh%NeA)
825 real(RP), intent(in) :: DF2(elem%Np,lmesh%NeA)
826 real(RP), intent(in) :: DF3(elem%Np,lmesh%NeA)
827 real(RP), intent(in) :: Nu (elem%Np,lmesh%NeA)
828 real(RP), intent(in) :: Kh (elem%Np,lmesh%NeA)
829 real(RP), intent(in) :: DDENS_(elem%Np,lmesh%NeA)
830 real(RP), intent(in) :: MOMX_ (elem%Np,lmesh%NeA)
831 real(RP), intent(in) :: MOMY_ (elem%Np,lmesh%NeA)
832 real(RP), intent(in) :: MOMZ_ (elem%Np,lmesh%NeA)
833 real(RP), intent(in) :: DRHOT_(elem%Np,lmesh%NeA)
834 real(RP), intent(in) :: DENS_hyd(elem%Np,lmesh%NeA)
835 real(RP), intent(in) :: PRES_hyd(elem%Np,lmesh%NeA)
836 real(RP), intent(in) :: PRES_(elem%Np,lmesh%NeA)
837 real(RP), intent(in) :: PT_ (elem%Np,lmesh%NeA)
838 type(SparseMat), intent(in) :: Dx, Dy, Dz
839 type(SparseMat), intent(in) :: Sx, Sy, Sz
840 type(SparseMat), intent(in) :: Lift
841 logical, intent(in) :: is_bound(elem%NfpTot,lmesh%Ne)
842
843 integer :: ke, ke2d
844
845 real(RP) :: X2D(elem2D%Np,lmesh2D%Ne), Y2D(elem2D%Np,lmesh2D%Ne)
846 real(RP) :: X(elem%Np), Y(elem%Np), twoOVdel2(elem%Np)
847 real(RP) :: R(elem%Np)
848
849 real(RP) :: Fx(elem%Np), Fy(elem%Np), Fz(elem%Np), LiftDelFlx(elem%Np)
850 real(RP) :: GsqrtDENS(elem%Np), RHOT(elem%Np)
851 real(RP) :: del_flux_mom(elem%NfpTot,lmesh%Ne,3)
852 real(RP) :: del_flux_rhot(elem%NfpTot,lmesh%Ne)
853 !--------------------------------------------------------------------
854
855 call cal_del_flux( del_flux_mom, del_flux_rhot, & ! (out)
856 t11, t12, t13, t21, t22, t23, t31, t32, t33, & ! (in)
857 df1, df2, df3, & ! (in)
858 nu, kh, & ! (in)
859 ddens_, momx_, momy_, momz_, drhot_, dens_hyd, pres_hyd, & ! (in)
860 lmesh%Gsqrt, & ! (in)
861 lmesh%normal_fn(:,:,1), lmesh%normal_fn(:,:,2), lmesh%normal_fn(:,:,3), & ! (in)
862 lmesh%vmapM, lmesh%vmapP, & ! (in)
863 lmesh, elem, is_bound ) ! (in)
864
865 !$omp parallel private( ke, ke2D, &
866 !$omp Fx, Fy, Fz, LiftDelFlx, &
867 !$omp GsqrtDENS, X, Y, twoOVdel2, R )
868
869 !$omp do
870 do ke2d = lmesh2d%NeS, lmesh2d%NeE
871 x2d(:,ke2d) = tan(lmesh2d%pos_en(:,ke2d,1))
872 y2d(:,ke2d) = tan(lmesh2d%pos_en(:,ke2d,2))
873 end do
874
875 !$omp do
876 do ke=lmesh%NeS, lmesh%NeE
877 ke2d = lmesh%EMap3Dto2D(ke)
878 x(:) = x2d(elem%IndexH2Dto3D,ke2d)
879 y(:) = y2d(elem%IndexH2Dto3D,ke2d)
880 twoovdel2(:) = 2.0_rp / ( 1.0_rp + x(:)**2 + y(:)**2 )
881 r(:) = rplanet * lmesh%gam(:,ke)
882
883 gsqrtdens(:) = lmesh%Gsqrt(:,ke) * ( dens_hyd(:,ke) + ddens_(:,ke) )
884
885 ! MOMX
886 call sparsemat_matmul( dx, lmesh%Gsqrt(:,ke) * t11(:,ke), fx )
887 call sparsemat_matmul( dy, lmesh%Gsqrt(:,ke) * t12(:,ke), fy )
888 call sparsemat_matmul( dz, lmesh%Gsqrt(:,ke) * t13(:,ke), fz )
889 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_mom(:,ke,1), liftdelflx )
890
891 momx_t(:,ke) = ( lmesh%Escale(:,ke,1,1) * fx(:) &
892 + lmesh%Escale(:,ke,2,2) * fy(:) &
893 + lmesh%Escale(:,ke,3,3) * fz(:) &
894 + liftdelflx(:) ) / lmesh%Gsqrt(:,ke) &
895 + twoovdel2(:) * y(:) * &
896 ( x(:) * y(:) * t11(:,ke) - ( 1.0_rp + y(:)**2 ) * t12(:,ke) ) &
897 + shapro_coef * 2.0_rp * t13(:,ke) / r(:)
898
899 ! MOMY
900 call sparsemat_matmul( dx, lmesh%Gsqrt(:,ke) * t21(:,ke), fx )
901 call sparsemat_matmul( dy, lmesh%Gsqrt(:,ke) * t22(:,ke), fy )
902 call sparsemat_matmul( dz, lmesh%Gsqrt(:,ke) * t23(:,ke), fz )
903 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_mom(:,ke,2), liftdelflx )
904
905 momy_t(:,ke) = ( lmesh%Escale(:,ke,1,1) * fx(:) &
906 + lmesh%Escale(:,ke,2,2) * fy(:) &
907 + lmesh%Escale(:,ke,3,3) * fz(:) &
908 + liftdelflx(:) ) / lmesh%Gsqrt(:,ke) &
909 + twoovdel2(:) * x(:) * &
910 ( - ( 1.0_rp + x(:)**2 ) * t21(:,ke) + x(:) * y(:) * t22(:,ke) ) &
911 + shapro_coef * 2.0_rp * t23(:,ke) / r(:)
912
913 ! MOMZ
914 call sparsemat_matmul( dx, lmesh%Gsqrt(:,ke) * t31(:,ke), fx )
915 call sparsemat_matmul( dy, lmesh%Gsqrt(:,ke) * t32(:,ke), fy )
916 call sparsemat_matmul( dz, lmesh%Gsqrt(:,ke) * t33(:,ke), fz )
917 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_mom(:,ke,3), liftdelflx )
918
919 momz_t(:,ke) = ( lmesh%Escale(:,ke,1,1) * fx(:) &
920 + lmesh%Escale(:,ke,2,2) * fy(:) &
921 + lmesh%Escale(:,ke,3,3) * fz(:) &
922 + liftdelflx(:) ) / lmesh%Gsqrt(:,ke) &
923 + shapro_coef * 0.25_rp * r(:) * twoovdel2(:)**2 * ( 1.0_rp * x(:)**2 ) * ( 1.0_rp * y(:)**2 ) & !-> metric terms
924 * ( - ( 1.0_rp + x(:)**2 ) * t11(:,ke) + 2.0_rp * x(:) * y(:) * t12(:,ke) & !
925 - ( 1.0_rp + y(:)**2 ) * t22(:,ke) )
926
927 ! RHOT
928 call sparsemat_matmul( dx, gsqrtdens(:) * df1(:,ke), fx )
929 call sparsemat_matmul( dy, gsqrtdens(:) * df2(:,ke), fy )
930 call sparsemat_matmul( dz, gsqrtdens(:) * df3(:,ke), fz )
931 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_rhot(:,ke), liftdelflx )
932
933 rhot_t(:,ke) = ( lmesh%Escale(:,ke,1,1) * fx(:) &
934 + lmesh%Escale(:,ke,2,2) * fy(:) &
935 + lmesh%Escale(:,ke,3,3) * fz(:) &
936 + liftdelflx(:) ) / lmesh%Gsqrt(:,ke)
937
938 end do
939 !$omp end do
940 !$omp end parallel
941
942 return

◆ atm_phy_tb_dgm_globalsmg_cal_tend_qtrc()

subroutine, public scale_atm_phy_tb_dgm_globalsmg::atm_phy_tb_dgm_globalsmg_cal_tend_qtrc ( real(rp), dimension(elem%np,lmesh%nea), intent(out) rhoq_t,
real(rp), dimension(elem%np,lmesh%nea), intent(in) dfq1,
real(rp), dimension(elem%np,lmesh%nea), intent(in) dfq2,
real(rp), dimension(elem%np,lmesh%nea), intent(in) dfq3,
real(rp), dimension (elem%np,lmesh%nea), intent(in) kh,
real(rp), dimension(elem%np,lmesh%nea), intent(in) ddens_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) dens_hyd,
type(sparsemat), intent(in) dx,
type(sparsemat), intent(in) dy,
type(sparsemat), intent(in) dz,
type(sparsemat), intent(in) sx,
type(sparsemat), intent(in) sy,
type(sparsemat), intent(in) sz,
type(sparsemat), intent(in) lift,
class(localmesh3d), intent(in) lmesh,
class(elementbase3d), intent(in) elem,
class(localmesh2d), intent(in) lmesh2d,
class(elementbase2d), intent(in) elem2d,
logical, dimension(elem%nfptot,lmesh%ne), intent(in) is_bound )

Calculate tendecies of tracer density with turbulent model.

Parameters
[out]rhoq_tTendency of tracer mass fraction
[in]dfq1Diffusive mass flux in x1 direction / density (Kh dq/dx1)
[in]dfq2Diffusive mass flux in x2 direction / density (Kh dq/dx2)
[in]dfq3Diffusive mass flux in x3 direction / density (Kh dq/dx3)
[in]khEddy diffusivity
[in]ddens_Density perturbation
[in]dens_hydReference pressure in hydrostatic state
[in]dzDifferential matrix managed by sparse matrix type
[in]szStiffness matrix managed by sparse matrix type
[in]liftLifting matrix managed by sparse matrix type
[in]is_boundFlag whether nodes are located at domain boundaries

Definition at line 948 of file scale_atm_phy_tb_dgm_globalsmg.F90.

953
954 implicit none
955
956 class(LocalMesh3D), intent(in) :: lmesh
957 class(ElementBase3D), intent(in) :: elem
958 class(LocalMesh2D), intent(in) :: lmesh2D
959 class(ElementBase2D), intent(in) :: elem2D
960 real(RP), intent(out) :: RHOQ_t(elem%Np,lmesh%NeA)
961 real(RP), intent(in) :: DFQ1(elem%Np,lmesh%NeA)
962 real(RP), intent(in) :: DFQ2(elem%Np,lmesh%NeA)
963 real(RP), intent(in) :: DFQ3(elem%Np,lmesh%NeA)
964 real(RP), intent(in) :: Kh (elem%Np,lmesh%NeA)
965 real(RP), intent(in) :: DDENS_(elem%Np,lmesh%NeA)
966 real(RP), intent(in) :: DENS_hyd(elem%Np,lmesh%NeA)
967 type(SparseMat), intent(in) :: Dx, Dy, Dz
968 type(SparseMat), intent(in) :: Sx, Sy, Sz
969 type(SparseMat), intent(in) :: Lift
970 logical, intent(in) :: is_bound(elem%NfpTot,lmesh%Ne)
971
972 integer :: ke
973
974 real(RP) :: Fx(elem%Np), Fy(elem%Np), Fz(elem%Np), LiftDelFlx(elem%Np)
975 real(RP) :: GsqrtDENS(elem%Np), RHOT(elem%Np)
976 real(RP) :: del_flux(elem%NfpTot,lmesh%Ne)
977 !--------------------------------------------------------------------
978
979 call cal_del_flux_qtrc( del_flux, & ! (out)
980 dfq1, dfq2, dfq3, kh, ddens_, dens_hyd, & ! (in)
981 lmesh%Gsqrt, lmesh%normal_fn(:,:,1), lmesh%normal_fn(:,:,2), lmesh%normal_fn(:,:,3), & ! (in)
982 lmesh%vmapM, lmesh%vmapP, lmesh, elem, is_bound ) ! (in)
983
984 !$omp parallel do private( &
985 !$omp Fx, Fy, Fz, LiftDelFlx, &
986 !$omp GsqrtDENS )
987 do ke=lmesh%NeS, lmesh%NeE
988 gsqrtdens(:) = lmesh%Gsqrt(:,ke) * ( dens_hyd(:,ke) + ddens_(:,ke) )
989
990 ! RHOQ
991 call sparsemat_matmul( dx, gsqrtdens(:) * dfq1(:,ke), fx )
992 call sparsemat_matmul( dy, gsqrtdens(:) * dfq2(:,ke), fy )
993 call sparsemat_matmul( dz, gsqrtdens(:) * dfq3(:,ke), fz )
994 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux(:,ke), liftdelflx )
995
996 rhoq_t(:,ke) = ( lmesh%Escale(:,ke,1,1) * fx(:) &
997 + lmesh%Escale(:,ke,2,2) * fy(:) &
998 + lmesh%Escale(:,ke,3,3) * fz(:) &
999 + liftdelflx(:) ) / lmesh%Gsqrt(:,ke)
1000 end do
1001
1002 return