FE-Project
Loading...
Searching...
No Matches
scale_atm_phy_tb_dgm_globalsmg.F90
Go to the documentation of this file.
1!> module FElib / Atmosphere / Physics turbulence
2!!
3!! @par Description
4!! Sub-grid scale turbulence process for global domain
5!! Smagorinsky-type
6!!
7!! @author Yuta Kawai, Team SCALE
8!!
9!! @par Reference
10!! - Brown et al., 1994:
11!! Large-eddy simulation of stable atmospheric boundary layers with a revised stochastic subgrid model.
12!! Roy. Meteor. Soc., 120, 1485-1512
13!! - Scotti et al., 1993:
14!! Generalized Smagorinsky model for anisotropic grids.
15!! Phys. Fluids A, 5, 2306-2308
16!! - Nishizawa et al., 2015:
17!! Influence of grid aspect ratio on planetary boundary layer turbulence in large-eddy simulations
18!! Geosci. Model Dev., 8, 3393–3419
19!<
20!-------------------------------------------------------------------------------
21#include "scaleFElib.h"
22module scale_atm_phy_tb_dgm_globalsmg
23 !-----------------------------------------------------------------------------
24 !
25 !++ Used modules
26 !
27 use scale_precision
28 use scale_io
29 use scale_prc
30 use scale_prof
31 use scale_const, only: &
32 eps => const_eps, &
33 grav => const_grav, &
34 rdry => const_rdry, &
35 cpdry => const_cpdry, &
36 cvdry => const_cvdry, &
37 pres00 => const_pre00, &
38 karman => const_karman, &
39 rplanet => const_radius
40
41 use scale_sparsemat
42 use scale_element_base, only: &
43 elementbase2d, elementbase3d
44 use scale_element_hexahedral, only: hexahedralelement
45 use scale_localmesh_2d, only: localmesh2d
46 use scale_localmesh_3d, only: localmesh3d
47 use scale_mesh_base3d, only: meshbase3d
48 use scale_localmeshfield_base, only: localmeshfield3d
49 use scale_meshfield_base, only: meshfield3d
50
51 !-----------------------------------------------------------------------------
52 implicit none
53 private
54 !-----------------------------------------------------------------------------
55 !
56 !++ Public procedures
57 !
58 public :: atm_phy_tb_dgm_globalsmg_init
59 public :: atm_phy_tb_dgm_globalsmg_final
60 public :: atm_phy_tb_dgm_globalsmg_cal_grad
61 public :: atm_phy_tb_dgm_globalsmg_cal_grad_qtrc
62 public :: atm_phy_tb_dgm_globalsmg_cal_tend
63 public :: atm_phy_tb_dgm_globalsmg_cal_tend_qtrc
64
65 !-----------------------------------------------------------------------------
66 !
67 !++ Private procedure
68 !
69 private :: cal_del_flux_grad
70 private :: cal_del_flux_grad_qtrc
71 private :: cal_del_flux
72
73 !-----------------------------------------------------------------------------
74 !
75 !++ Private parameters & variables
76 !
77 real(RP), private, parameter :: OneOverThree = 1.0_rp / 3.0_rp
78 real(RP), private, parameter :: twoOverThree = 2.0_rp / 3.0_rp
79 real(RP), private, parameter :: FourOverThree = 4.0_rp / 3.0_rp
80
81 real(RP), private :: Cs = 0.13_rp ! Smagorinsky constant (Scotti et al. 1993)
82 real(RP), private, parameter :: PrN = 0.7_rp ! Prandtl number in neutral conditions
83 real(RP), private, parameter :: RiC = 0.25_rp ! critical Richardson number
84 real(RP), private, parameter :: FmC = 16.0_rp ! fum = sqrt(1 - c*Ri)
85 real(RP), private, parameter :: FhB = 40.0_rp ! fuh = sqrt(1 - b*Ri)/PrN
86 real(RP), private :: RPrN ! 1 / PrN
87 real(RP), private :: RRiC ! 1 / RiC
88 real(RP), private :: OnemPrNovRiC ! PrN / RiC
89
90 ! for backscatter
91 real(RP), private, parameter :: CB = 1.4_rp
92 real(RP), private, parameter :: CBt = 0.45_rp
93 real(RP), private, parameter :: aN = 0.47958315233127197_rp ! a_N = sqrt(0.23)
94 real(RP), private, parameter :: atN4 = 0.09_rp ! a_{\theta N}^4 = 0.3**2
95 real(RP), private, parameter :: C1o = an**3
96 real(RP), private, parameter :: D1o = prn * atn4 / an
97
98
99 real(RP), private :: filter_fac = 2.0_rp
100 real(RP), private :: NU_MAX = 10000.0_rp
101 real(RP), private :: tke_fac
102
103 !
104 logical, private :: shallow_atm_approx_flag
105 real(RP), private :: shapro_coef
106
107contains
108!OCL SERIAL
109 subroutine atm_phy_tb_dgm_globalsmg_init( mesh, shallow_atm_approx )
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
159 end subroutine atm_phy_tb_dgm_globalsmg_init
160
161!OCL SERIAL
162 subroutine atm_phy_tb_dgm_globalsmg_final()
163 implicit none
164 !--------------------------------------------------------------------
165
166 return
167 end subroutine atm_phy_tb_dgm_globalsmg_final
168
169!> Calculate parameterized stress tensor and eddy heat flux with turbulent model
170!!
171!OCL SERIAL
172 subroutine atm_phy_tb_dgm_globalsmg_cal_grad( &
173 T11, T12, T13, T21, T22, T23, T31, T32, T33, & ! (out)
174 df1, df2, df3, & ! (out)
175 tke, nu, kh, & ! (out)
176 ddens_, momx_, momy_, momz_, drhot_, dens_hyd, pres_hyd, & ! (in)
177 pres, pt, & ! (in)
178 dx, dy, dz, sx, sy, sz, lift, lmesh, elem, lmesh2d, elem2d, & ! (in)
179 is_bound ) ! (in)
180
181 use scale_atm_phy_tb_dgm_common, only: &
182 atm_phy_tb_dgm_common_calc_lambda
183 use scale_cubedsphere_coord_cnv, only: &
184 cubedspherecoordcnv_cs2lonlatvec
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) !< (1,1) component of stress tensor
193 real(rp), intent(out) :: t12(elem%np,lmesh%nea) !< (1,2) component of stress tensor
194 real(rp), intent(out) :: t13(elem%np,lmesh%nea) !< (1,3) component of stress tensor
195 real(rp), intent(out) :: t21(elem%np,lmesh%nea) !< (2,1) component of stress tensor
196 real(rp), intent(out) :: t22(elem%np,lmesh%nea) !< (2,2) component of stress tensor
197 real(rp), intent(out) :: t23(elem%np,lmesh%nea) !< (2,3) component of stress tensor
198 real(rp), intent(out) :: t31(elem%np,lmesh%nea) !< (3,1) component of stress tensor
199 real(rp), intent(out) :: t32(elem%np,lmesh%nea) !< (3,2) component of stress tensor
200 real(rp), intent(out) :: t33(elem%np,lmesh%nea) !< (3,3) component of stress tensor
201 real(rp), intent(out) :: df1(elem%np,lmesh%nea) !< Diffusive heat flux in x1 direction / density
202 real(rp), intent(out) :: df2(elem%np,lmesh%nea) !< Diffusive heat flux in x2 direction / density
203 real(rp), intent(out) :: df3(elem%np,lmesh%nea) !< Diffusive heat flux in x3 direction / density
204 real(rp), intent(out) :: tke(elem%np,lmesh%nea) !< Parameterized turbulent kinetic energy
205 real(rp), intent(out) :: nu(elem%np,lmesh%nea) !< Eddy viscosity
206 real(rp), intent(out) :: kh(elem%np,lmesh%nea) !< Eddy diffusivity
207 real(rp), intent(in) :: ddens_(elem%np,lmesh%nea) !< Density perturbation
208 real(rp), intent(in) :: momx_ (elem%np,lmesh%nea) !< Momentum in x1 direction
209 real(rp), intent(in) :: momy_ (elem%np,lmesh%nea) !< Momentum in x2 direction
210 real(rp), intent(in) :: momz_ (elem%np,lmesh%nea) !< Momentum in x3 direction
211 real(rp), intent(in) :: drhot_(elem%np,lmesh%nea) !< Density x potential temperature perturbation
212 real(rp), intent(in) :: dens_hyd(elem%np,lmesh%nea) !< Reference pressure in hydrostatic balance
213 real(rp), intent(in) :: pres_hyd(elem%np,lmesh%nea) !< Reference density in hydrostatic balance
214 real(rp), intent(in) :: pres(elem%np,lmesh%nea) !< Pressure
215 real(rp), intent(in) :: pt(elem%np,lmesh%nea) !< Potential temperature
216 type(sparsemat), intent(in) :: dx, dy, dz !< Differential matrix managed by sparse matrix type
217 type(sparsemat), intent(in) :: sx, sy, sz !< Stiffness matrix managed by sparse matrix type
218 type(sparsemat), intent(in) :: lift !< Lifting matrix managed by sparse matrix type
219 logical, intent(in) :: is_bound(elem%nfptot,lmesh%ne) !< Flag whether nodes are located at domain boundaries
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
525 end subroutine atm_phy_tb_dgm_globalsmg_cal_grad
526
527!> Calculate parameterized diffusive mass flux of tracer with turbulent model
528!!
529!OCL SERIAL
530 subroutine atm_phy_tb_dgm_globalsmg_cal_grad_qtrc( &
531 DFQ1, DFQ2, DFQ3, & ! (out)
532 drdx, drdy, drdz, & ! (inout)
533 kh, qtrc, ddens, dens_hyd, & ! (in)
534 dx, dy, dz, sx, sy, sz, lift, lmesh, elem, lmesh2d, elem2d, & ! (in)
535 is_bound, cal_grad_dens ) ! (in)
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) !< Diffusive mass flux in x1 direction / density (Kh dq/dx1)
544 real(rp), intent(out) :: dfq2(elem%np,lmesh%nea) !< Diffusive mass flux in x2 direction / density (Kh dq/dx2)
545 real(rp), intent(out) :: dfq3(elem%np,lmesh%nea) !< Diffusive mass flux in x3 direction / density (Kh dq/dx3)
546 real(rp), intent(inout) :: drdx(elem%np,lmesh%nea) !< Spatial gradient of density in x1 direction
547 real(rp), intent(inout) :: drdy(elem%np,lmesh%nea) !< Spatial gradient of density in x2 direction
548 real(rp), intent(inout) :: drdz(elem%np,lmesh%nea) !< Spatial gradient of density in x3 direction
549 real(rp), intent(in) :: kh(elem%np,lmesh%nea) !< Eddy diffusivity
550 real(rp), intent(in) :: qtrc(elem%np,lmesh%nea) !< Mass faction of tracer
551 real(rp), intent(in) :: ddens(elem%np,lmesh%nea) !< Density perturbation
552 real(rp), intent(in) :: dens_hyd(elem%np,lmesh%nea) !< Reference desity in hydrostatic state
553 type(sparsemat), intent(in) :: dx, dy, dz !< Differential matrix managed by sparse matrix type
554 type(sparsemat), intent(in) :: sx, sy, sz !< Stiffness matrix managed by sparse matrix type
555 type(sparsemat), intent(in) :: lift !< Lifting matrix managed by sparse matrix type
556 logical, intent(in) :: is_bound(elem%nfptot,lmesh%ne) !< Flag whether nodes are located at domain boundaries
557 logical, intent(in) :: cal_grad_dens !< Flag whether spatial gradients of density are calcuated
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
634 end subroutine atm_phy_tb_dgm_globalsmg_cal_grad_qtrc
635
636!OCL SERIAL
637 subroutine cal_del_flux_grad( del_flux_rho, del_flux_mom, del_flux_rhot, & ! (out)
638 ddens_, momx_, momy_, momz_, drhot_, dens_hyd, pres_hyd, pt_, & ! (in)
639 nx, ny, nz, vmapm, vmapp, lmesh, elem, is_bound ) ! (in)
640
641 implicit none
642
643 class(localmesh3d), intent(in) :: lmesh
644 class(elementbase3d), intent(in) :: elem
645 real(rp), intent(out) :: del_flux_rho(elem%nfptot*lmesh%ne,3)
646 real(rp), intent(out) :: del_flux_mom(elem%nfptot*lmesh%ne,3,3)
647 real(rp), intent(out) :: del_flux_rhot(elem%nfptot*lmesh%ne,3)
648 real(rp), intent(in) :: ddens_(elem%np*lmesh%nea)
649 real(rp), intent(in) :: momx_(elem%np*lmesh%nea)
650 real(rp), intent(in) :: momy_(elem%np*lmesh%nea)
651 real(rp), intent(in) :: momz_(elem%np*lmesh%nea)
652 real(rp), intent(in) :: drhot_(elem%np*lmesh%nea)
653 real(rp), intent(in) :: dens_hyd(elem%np*lmesh%nea)
654 real(rp), intent(in) :: pres_hyd(elem%np*lmesh%nea)
655 real(rp), intent(in) :: pt_(elem%np*lmesh%nea)
656 real(rp), intent(in) :: nx(elem%nfptot*lmesh%ne)
657 real(rp), intent(in) :: ny(elem%nfptot*lmesh%ne)
658 real(rp), intent(in) :: nz(elem%nfptot*lmesh%ne)
659 integer, intent(in) :: vmapm(elem%nfptot*lmesh%ne)
660 integer, intent(in) :: vmapp(elem%nfptot*lmesh%ne)
661 logical, intent(in) :: is_bound(elem%nfptot*lmesh%ne)
662
663 integer :: i, ip, im
664 real(rp) :: densm, densp
665 real(rp) :: del
666 real(rp) :: facx, facy, facz
667
668 real(rp) :: momz_p
669 !------------------------------------------------------------------------
670
671 !$omp parallel do private ( iM, iP, &
672 !$omp densM, densP, &
673 !$omp del, facx, facy, facz, MOMZ_P )
674 do i=1, elem%NfpTot * lmesh%Ne
675 im = vmapm(i); ip = vmapp(i)
676
677 densm = ddens_(im) + dens_hyd(im)
678 densp = ddens_(ip) + dens_hyd(ip)
679
680 if ( is_bound(i) ) then
681 facx = 1.0_rp
682 facy = 1.0_rp
683 facz = 1.0_rp
684 momz_p = - momz_(im)
685 else
686 ! facx = 1.0_RP - sign(1.0_RP,nx(i))
687 ! facy = 1.0_RP - sign(1.0_RP,ny(i))
688 ! facz = 1.0_RP - sign(1.0_RP,nz(i))
689 facx = 1.0_rp
690 facy = 1.0_rp
691 facz = 1.0_rp
692 momz_p = momz_(ip)
693 end if
694
695 del = 0.5_rp * ( densp - densm )
696 del_flux_rho(i,1) = facx * del * nx(i)
697 del_flux_rho(i,2) = facy * del * ny(i)
698 del_flux_rho(i,3) = facz * del * nz(i)
699
700 del = 0.5_rp * ( momx_(ip) - momx_(im) )
701 del_flux_mom(i,1,1) = facx * del * nx(i)
702 del_flux_mom(i,2,1) = facy * del * ny(i)
703 del_flux_mom(i,3,1) = facz * del * nz(i)
704
705 del = 0.5_rp * ( momy_(ip) - momy_(im) )
706 del_flux_mom(i,1,2) = facx * del * nx(i)
707 del_flux_mom(i,2,2) = facy * del * ny(i)
708 del_flux_mom(i,3,2) = facz * del * nz(i)
709
710 del = 0.5_rp * ( momz_p - momz_(im) )
711 del_flux_mom(i,1,3) = facx * del * nx(i)
712 del_flux_mom(i,2,3) = facy * del * ny(i)
713 del_flux_mom(i,3,3) = facz * del * nz(i)
714
715 del = 0.5_rp * ( densp * pt_(ip) - densm * pt_(im) )
716 del_flux_rhot(i,1) = facx * del * nx(i)
717 del_flux_rhot(i,2) = facy * del * ny(i)
718 del_flux_rhot(i,3) = facz * del * nz(i)
719 end do
720
721 return
722 end subroutine cal_del_flux_grad
723
724!OCL SERIAL
725 subroutine cal_del_flux_grad_qtrc( del_flux, del_flux_rho, & ! (out)
726 qtrc_, ddens_, dens_hyd_, & ! (in)
727 nx, ny, nz, vmapm, vmapp, lmesh, elem, is_bound, & ! (in)
728 cal_grad_dens ) ! (in)
729
730 implicit none
731
732 class(localmesh3d), intent(in) :: lmesh
733 class(elementbase3d), intent(in) :: elem
734 real(rp), intent(out) :: del_flux(elem%nfptot*lmesh%ne,3)
735 real(rp), intent(out) :: del_flux_rho(elem%nfptot*lmesh%ne,3)
736 real(rp), intent(in) :: qtrc_(elem%np*lmesh%nea)
737 real(rp), intent(in) :: ddens_(elem%np*lmesh%nea)
738 real(rp), intent(in) :: dens_hyd_(elem%np*lmesh%nea)
739 real(rp), intent(in) :: nx(elem%nfptot*lmesh%ne)
740 real(rp), intent(in) :: ny(elem%nfptot*lmesh%ne)
741 real(rp), intent(in) :: nz(elem%nfptot*lmesh%ne)
742 integer, intent(in) :: vmapm(elem%nfptot*lmesh%ne)
743 integer, intent(in) :: vmapp(elem%nfptot*lmesh%ne)
744 logical, intent(in) :: is_bound(elem%nfptot*lmesh%ne)
745 logical, intent(in) :: cal_grad_dens
746
747 integer :: i, ip, im
748 real(rp) :: del
749 real(rp) :: facx, facy, facz
750 real(rp) :: densm, densp
751
752 !------------------------------------------------------------------------
753
754 !$omp parallel do private ( iM, iP, &
755 !$omp del, facx, facy, facz, densM, densP )
756 do i=1, elem%NfpTot * lmesh%Ne
757 im = vmapm(i); ip = vmapp(i)
758
759 if ( is_bound(i) ) then
760 facx = 1.0_rp
761 facy = 1.0_rp
762 facz = 1.0_rp
763 else
764 ! facx = 1.0_RP - sign(1.0_RP,nx(i))
765 ! facy = 1.0_RP - sign(1.0_RP,ny(i))
766 ! facz = 1.0_RP - sign(1.0_RP,nz(i))
767 facx = 1.0_rp
768 facy = 1.0_rp
769 facz = 1.0_rp
770 end if
771
772 densm = ddens_(im) + dens_hyd_(im)
773 densp = ddens_(ip) + dens_hyd_(ip)
774
775 if ( cal_grad_dens ) then
776 del = 0.5_rp * ( densp - densm )
777 ! del = 0.5_RP * ( sqrt_DENS_Kh_P - sqrt_DENS_Kh_M )
778 del_flux_rho(i,1) = facx * del * nx(i)
779 del_flux_rho(i,2) = facy * del * ny(i)
780 del_flux_rho(i,3) = facz * del * nz(i)
781 end if
782
783 del = 0.5_rp * ( densp * qtrc_(ip) - densm * qtrc_(im) )
784 del_flux(i,1) = facx * del * nx(i)
785 del_flux(i,2) = facy * del * ny(i)
786 del_flux(i,3) = facz * del * nz(i)
787 end do
788
789 return
790 end subroutine cal_del_flux_grad_qtrc
791
792!> Calculate tendecies with turbulent model
793!!
794!OCL SERIAL
795 subroutine atm_phy_tb_dgm_globalsmg_cal_tend( &
796 MOMX_t, MOMY_t, MOMZ_t, RHOT_t, & ! (out)
797 t11, t12, t13, t21, t22, t23, t31, t32, t33, & ! (in)
798 df1, df2, df3, & ! (in)
799 nu, kh, & ! (in)
800 ddens_, momx_, momy_, momz_, drhot_, & ! (in)
801 dens_hyd, pres_hyd, pres_, pt_, & ! (in)
802 dx, dy, dz, sx, sy, sz, lift, lmesh, elem, lmesh2d, elem2d, & ! (in)
803 is_bound ) ! (in)
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) !< Tendency of momentum in x1 direction with turbulent model
812 real(rp), intent(out) :: momy_t(elem%np,lmesh%nea) !< Tendency of momentum in x2 direction with turbulent model
813 real(rp), intent(out) :: momz_t(elem%np,lmesh%nea) !< Tendency of momentum in x3 direction with turbulent model
814 real(rp), intent(out) :: rhot_t(elem%np,lmesh%nea) !< Tendency of density x potential temperature with turbulent model
815 real(rp), intent(in) :: t11(elem%np,lmesh%nea) !< (1,1) component of stress tensor
816 real(rp), intent(in) :: t12(elem%np,lmesh%nea) !< (1,2) component of stress tensor
817 real(rp), intent(in) :: t13(elem%np,lmesh%nea) !< (1,3) component of stress tensor
818 real(rp), intent(in) :: t21(elem%np,lmesh%nea) !< (2,1) component of stress tensor
819 real(rp), intent(in) :: t22(elem%np,lmesh%nea) !< (2,2) component of stress tensor
820 real(rp), intent(in) :: t23(elem%np,lmesh%nea) !< (2,3) component of stress tensor
821 real(rp), intent(in) :: t31(elem%np,lmesh%nea) !< (3,1) component of stress tensor
822 real(rp), intent(in) :: t32(elem%np,lmesh%nea) !< (3,2) component of stress tensor
823 real(rp), intent(in) :: t33(elem%np,lmesh%nea) !< (3,3) component of stress tensor
824 real(rp), intent(in) :: df1(elem%np,lmesh%nea) !< Diffusive heat flux in x1 direction / density
825 real(rp), intent(in) :: df2(elem%np,lmesh%nea) !< Diffusive heat flux in x2 direction / density
826 real(rp), intent(in) :: df3(elem%np,lmesh%nea) !< Diffusive heat flux in x3 direction / density
827 real(rp), intent(in) :: nu (elem%np,lmesh%nea) !< Eddy viscosity
828 real(rp), intent(in) :: kh (elem%np,lmesh%nea) !< Eddy diffusivity
829 real(rp), intent(in) :: ddens_(elem%np,lmesh%nea) !< Density perturbation
830 real(rp), intent(in) :: momx_ (elem%np,lmesh%nea) !< Momentum in x1 direction
831 real(rp), intent(in) :: momy_ (elem%np,lmesh%nea) !< Momentum in x2 direction
832 real(rp), intent(in) :: momz_ (elem%np,lmesh%nea) !< Momentum in x3 direction
833 real(rp), intent(in) :: drhot_(elem%np,lmesh%nea) !< Density x potential temperature perturbation
834 real(rp), intent(in) :: dens_hyd(elem%np,lmesh%nea) !< Reference pressure in hydrostatic balance
835 real(rp), intent(in) :: pres_hyd(elem%np,lmesh%nea) !< Reference density in hydrostatic balance
836 real(rp), intent(in) :: pres_(elem%np,lmesh%nea) !< Pressure
837 real(rp), intent(in) :: pt_ (elem%np,lmesh%nea) !< Potential temperature
838 type(sparsemat), intent(in) :: dx, dy, dz !< Differential matrix managed by sparse matrix type
839 type(sparsemat), intent(in) :: sx, sy, sz !< Stiffness matrix managed by sparse matrix type
840 type(sparsemat), intent(in) :: lift !< Lifting matrix managed by sparse matrix type
841 logical, intent(in) :: is_bound(elem%nfptot,lmesh%ne) !< Flag whether nodes are located at domain boundaries
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
943 end subroutine atm_phy_tb_dgm_globalsmg_cal_tend
944
945!> Calculate tendecies of tracer density with turbulent model
946!!
947!OCL SERIAL
948 subroutine atm_phy_tb_dgm_globalsmg_cal_tend_qtrc( &
949 RHOQ_t, & ! (out)
950 dfq1, dfq2, dfq3, kh, ddens_,dens_hyd, & ! (in)
951 dx, dy, dz, sx, sy, sz, lift, lmesh, elem, lmesh2d, elem2d, & ! (in)
952 is_bound ) ! (in)
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) !< Tendency of tracer mass fraction
961 real(rp), intent(in) :: dfq1(elem%np,lmesh%nea) !< Diffusive mass flux in x1 direction / density (Kh dq/dx1)
962 real(rp), intent(in) :: dfq2(elem%np,lmesh%nea) !< Diffusive mass flux in x2 direction / density (Kh dq/dx2)
963 real(rp), intent(in) :: dfq3(elem%np,lmesh%nea) !< Diffusive mass flux in x3 direction / density (Kh dq/dx3)
964 real(rp), intent(in) :: kh (elem%np,lmesh%nea) !< Eddy diffusivity
965 real(rp), intent(in) :: ddens_(elem%np,lmesh%nea) !< Density perturbation
966 real(rp), intent(in) :: dens_hyd(elem%np,lmesh%nea) !< Reference pressure in hydrostatic state
967 type(sparsemat), intent(in) :: dx, dy, dz !< Differential matrix managed by sparse matrix type
968 type(sparsemat), intent(in) :: sx, sy, sz !< Stiffness matrix managed by sparse matrix type
969 type(sparsemat), intent(in) :: lift !< Lifting matrix managed by sparse matrix type
970 logical, intent(in) :: is_bound(elem%nfptot,lmesh%ne) !< Flag whether nodes are located at domain boundaries
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
1003 end subroutine atm_phy_tb_dgm_globalsmg_cal_tend_qtrc
1004
1005!-- private --------------------------------------------------------
1006
1007!OCL SERIAL
1008 subroutine cal_del_flux( del_flux_mom, del_flux_rhot, & ! (out)
1009 t11, t12, t13, t21, t22, t23, t31, t32, t33, & ! (in)
1010 df1, df2, df3, & ! (in)
1011 nu, kh, & ! (in)
1012 ddens_, momx_, momy_, momz_, drhot_, dens_hyd, pres_hyd, & ! (in)
1013 gsqrt, nx, ny, nz, vmapm, vmapp, lmesh, elem, is_bound ) ! (in)
1014
1015 implicit none
1016
1017 class(localmesh3d), intent(in) :: lmesh
1018 class(elementbase3d), intent(in) :: elem
1019 real(rp), intent(out) :: del_flux_mom (elem%nfptot*lmesh%ne,3)
1020 real(rp), intent(out) :: del_flux_rhot(elem%nfptot*lmesh%ne)
1021 real(rp), intent(in) :: t11(elem%np*lmesh%nea), t12(elem%np*lmesh%nea), t13(elem%np*lmesh%nea)
1022 real(rp), intent(in) :: t21(elem%np*lmesh%nea), t22(elem%np*lmesh%nea), t23(elem%np*lmesh%nea)
1023 real(rp), intent(in) :: t31(elem%np*lmesh%nea), t32(elem%np*lmesh%nea), t33(elem%np*lmesh%nea)
1024 real(rp), intent(in) :: df1(elem%np*lmesh%nea)
1025 real(rp), intent(in) :: df2(elem%np*lmesh%nea)
1026 real(rp), intent(in) :: df3(elem%np*lmesh%nea)
1027 real(rp), intent(in) :: nu (elem%np*lmesh%nea) ! Eddy viscosity
1028 real(rp), intent(in) :: kh (elem%np*lmesh%nea) ! Eddy diffusivity
1029 real(rp), intent(in) :: ddens_(elem%np*lmesh%nea)
1030 real(rp), intent(in) :: momx_(elem%np*lmesh%nea)
1031 real(rp), intent(in) :: momy_(elem%np*lmesh%nea)
1032 real(rp), intent(in) :: momz_(elem%np*lmesh%nea)
1033 real(rp), intent(in) :: drhot_(elem%np*lmesh%nea)
1034 real(rp), intent(in) :: dens_hyd(elem%np*lmesh%nea)
1035 real(rp), intent(in) :: pres_hyd(elem%np*lmesh%nea)
1036 real(rp), intent(in) :: gsqrt(elem%np*lmesh%nea)
1037 real(rp), intent(in) :: nx(elem%nfptot*lmesh%ne)
1038 real(rp), intent(in) :: ny(elem%nfptot*lmesh%ne)
1039 real(rp), intent(in) :: nz(elem%nfptot*lmesh%ne)
1040 integer, intent(in) :: vmapm(elem%nfptot*lmesh%ne)
1041 integer, intent(in) :: vmapp(elem%nfptot*lmesh%ne)
1042 logical, intent(in) :: is_bound(elem%nfptot*lmesh%ne)
1043
1044 integer :: i, ip, im
1045 real(rp) :: gsqrtdensm, gsqrtdensp
1046 real(rp) :: taum_x, taup_x
1047 real(rp) :: taum_y, taup_y
1048 real(rp) :: taum_z, taup_z
1049 real(rp) :: nx_, ny_, nz_
1050 !------------------------------------------------------------------------
1051
1052 !$omp parallel do private( iM, iP, &
1053 !$omp GsqrtDensM, GsqrtDensP, &
1054 !$omp TauM_x, TauP_x, TauM_y, TauP_y, TauM_z, TauP_z, nx_, ny_, nz_ )
1055 do i=1, elem%NfpTot * lmesh%Ne
1056 im = vmapm(i); ip = vmapp(i)
1057
1058 gsqrtdensm = gsqrt(im) * ( ddens_(im) + dens_hyd(im) )
1059 gsqrtdensp = gsqrt(ip) * ( ddens_(ip) + dens_hyd(ip) )
1060
1061 if ( ip > elem%Np * lmesh%Ne .and. abs(nz(i)) > eps ) then ! Tentative implementation for the treatmnet of lower/upper boundary.
1062 nx_ = nx(i)
1063 ny_ = ny(i)
1064 nz_ = nz(i)
1065 else
1066 ! nx_ = ( 1.0_RP + sign(1.0_RP,nx(i)) ) * nx(i)
1067 ! ny_ = ( 1.0_RP + sign(1.0_RP,ny(i)) ) * ny(i)
1068 ! nz_ = ( 1.0_RP + sign(1.0_RP,nz(i)) ) * nz(i)
1069 nx_ = nx(i)
1070 ny_ = ny(i)
1071 nz_ = nz(i)
1072 end if
1073
1074 taum_x = gsqrt(im) * ( t11(im) * nx_ + t12(im) * ny_ + t13(im) * nz_ )
1075 taup_x = gsqrt(ip) * ( t11(ip) * nx_ + t12(ip) * ny_ + t13(ip) * nz_ )
1076
1077 taum_y = gsqrt(im) * ( t21(im) * nx_ + t22(im) * ny_ + t23(im) * nz_ )
1078 taup_y = gsqrt(ip) * ( t21(ip) * nx_ + t22(ip) * ny_ + t23(ip) * nz_ )
1079
1080 taum_z = gsqrt(im) * ( t31(im) * nx_ + t32(im) * ny_ + t33(im) * nz_ )
1081 taup_z = gsqrt(ip) * ( t31(ip) * nx_ + t32(ip) * ny_ + t33(ip) * nz_ )
1082
1083 if ( is_bound(i) ) then
1084 del_flux_mom(i,1) = - taum_x
1085 del_flux_mom(i,2) = - taum_y
1086 del_flux_mom(i,3) = 0.5_rp * ( taup_z - taum_z )
1087 del_flux_rhot(i) = - gsqrtdensm * ( df1(im) * nx(i) + df2(im) * ny(i) + df3(im) * nz(i) )
1088 else
1089 del_flux_mom(i,1) = 0.5_rp * ( taup_x - taum_x )
1090 del_flux_mom(i,2) = 0.5_rp * ( taup_y - taum_y )
1091 del_flux_mom(i,3) = 0.5_rp * ( taup_z - taum_z )
1092 del_flux_rhot(i) = 0.5_rp * ( gsqrtdensp * ( df1(ip) * nx_ + df2(ip) * ny_ + df3(ip) * nz_ ) &
1093 - gsqrtdensm * ( df1(im) * nx_ + df2(im) * ny_ + df3(im) * nz_ ) )
1094 end if
1095 end do
1096
1097 return
1098 end subroutine cal_del_flux
1099
1100!OCL SERIAL
1101 subroutine cal_del_flux_qtrc( del_flux, & ! (out)
1102 dfq1, dfq2, dfq3, & ! (in)
1103 kh, ddens_, dens_hyd, & ! (in)
1104 gsqrt, nx, ny, nz, vmapm, vmapp, lmesh, elem, is_bound ) ! (in)
1105
1106 implicit none
1107
1108 class(localmesh3d), intent(in) :: lmesh
1109 class(elementbase3d), intent(in) :: elem
1110 real(rp), intent(out) :: del_flux(elem%nfptot*lmesh%ne)
1111 real(rp), intent(in) :: dfq1(elem%np*lmesh%nea)
1112 real(rp), intent(in) :: dfq2(elem%np*lmesh%nea)
1113 real(rp), intent(in) :: dfq3(elem%np*lmesh%nea)
1114 real(rp), intent(in) :: kh (elem%np*lmesh%nea) ! Eddy diffusivity
1115 real(rp), intent(in) :: ddens_(elem%np*lmesh%nea)
1116 real(rp), intent(in) :: dens_hyd(elem%np*lmesh%nea)
1117 real(rp), intent(in) :: gsqrt(elem%np*lmesh%nea)
1118 real(rp), intent(in) :: nx(elem%nfptot*lmesh%ne)
1119 real(rp), intent(in) :: ny(elem%nfptot*lmesh%ne)
1120 real(rp), intent(in) :: nz(elem%nfptot*lmesh%ne)
1121 integer, intent(in) :: vmapm(elem%nfptot*lmesh%ne)
1122 integer, intent(in) :: vmapp(elem%nfptot*lmesh%ne)
1123 logical, intent(in) :: is_bound(elem%nfptot*lmesh%ne)
1124
1125 integer :: i, ip, im
1126 real(rp) :: gsqrtdensm, gsqrtdensp
1127 real(rp) :: nx_, ny_, nz_
1128 !------------------------------------------------------------------------
1129
1130 !$omp parallel do private( iM, iP, &
1131 !$omp GsqrtDensM, GsqrtDensP, nx_, ny_, nz_ )
1132 do i=1, elem%NfpTot * lmesh%Ne
1133 im = vmapm(i); ip = vmapp(i)
1134
1135 gsqrtdensm = gsqrt(im) * ( ddens_(im) + dens_hyd(im) )
1136 gsqrtdensp = gsqrt(ip) * ( ddens_(ip) + dens_hyd(ip) )
1137
1138 if ( ip > elem%Np * lmesh%Ne .and. abs(nz(i)) > eps ) then ! Tentative implementation for the treatmnet of lower/upper boundary.
1139 nx_ = nx(i)
1140 ny_ = ny(i)
1141 nz_ = nz(i)
1142 else
1143 ! nx_ = ( 1.0_RP + sign(1.0_RP,nx(i)) ) * nx(i)
1144 ! ny_ = ( 1.0_RP + sign(1.0_RP,ny(i)) ) * ny(i)
1145 ! nz_ = ( 1.0_RP + sign(1.0_RP,nz(i)) ) * nz(i)
1146 nx_ = nx(i)
1147 ny_ = ny(i)
1148 nz_ = nz(i)
1149 end if
1150
1151 if ( is_bound(i) ) then
1152 del_flux(i) = - gsqrtdensm * ( dfq1(im) * nx(i) + dfq2(im) * ny(i) + dfq3(im) * nz(i) )
1153 else
1154 del_flux(i) = 0.5_rp * ( gsqrtdensp * ( dfq1(ip) * nx_ + dfq2(ip) * ny_ + dfq3(ip) * nz_ ) &
1155 - gsqrtdensm * ( dfq1(im) * nx_ + dfq2(im) * ny_ + dfq3(im) * nz_ ) )
1156 end if
1157 end do
1158
1159 return
1160 end subroutine cal_del_flux_qtrc
1161
1162end module scale_atm_phy_tb_dgm_globalsmg
scale_sparsemat::sparsemat_matmul
Definition scale_sparsemat.F90:57
scale_atm_phy_tb_dgm_common
module FElib / Fluid dyn solver / Atmosphere / Physics turbulence / Common
Definition scale_atm_phy_tb_dgm_common.F90:11
scale_atm_phy_tb_dgm_common::atm_phy_tb_dgm_common_calc_lambda
subroutine, public atm_phy_tb_dgm_common_calc_lambda(lambda, cs, filter_fac, lmesh, elem, lmesh2d, elem2d)
Definition scale_atm_phy_tb_dgm_common.F90:226
scale_atm_phy_tb_dgm_globalsmg
module FElib / Atmosphere / Physics turbulence
Definition scale_atm_phy_tb_dgm_globalsmg.F90:22
scale_atm_phy_tb_dgm_globalsmg::atm_phy_tb_dgm_globalsmg_cal_tend
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.
Definition scale_atm_phy_tb_dgm_globalsmg.F90:804
scale_atm_phy_tb_dgm_globalsmg::atm_phy_tb_dgm_globalsmg_cal_tend_qtrc
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.
Definition scale_atm_phy_tb_dgm_globalsmg.F90:953
scale_atm_phy_tb_dgm_globalsmg::atm_phy_tb_dgm_globalsmg_cal_grad
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.
Definition scale_atm_phy_tb_dgm_globalsmg.F90:180
scale_atm_phy_tb_dgm_globalsmg::atm_phy_tb_dgm_globalsmg_cal_grad_qtrc
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.
Definition scale_atm_phy_tb_dgm_globalsmg.F90:536
scale_atm_phy_tb_dgm_globalsmg::atm_phy_tb_dgm_globalsmg_final
subroutine, public atm_phy_tb_dgm_globalsmg_final()
Definition scale_atm_phy_tb_dgm_globalsmg.F90:163
scale_atm_phy_tb_dgm_globalsmg::atm_phy_tb_dgm_globalsmg_init
subroutine, public atm_phy_tb_dgm_globalsmg_init(mesh, shallow_atm_approx)
Definition scale_atm_phy_tb_dgm_globalsmg.F90:110
scale_cubedsphere_coord_cnv
Module common / Coordinate conversion with cubed-sphere projection.
Definition scale_cubedsphere_coord_cnv.F90:21
scale_cubedsphere_coord_cnv::cubedspherecoordcnv_cs2lonlatvec
subroutine, public cubedspherecoordcnv_cs2lonlatvec(panelid, alpha, beta, gam, np, vecalpha, vecbeta, veclon, veclat, lat, gpu_async_id)
Convert the components of a vector in local coordinates with an equiangular gnomonic cubed-sphere pro...
Definition scale_cubedsphere_coord_cnv.F90:152
scale_element_base
module FElib / Element / Base
Definition scale_element_base.F90:11
scale_element_hexahedral
module FElib / Element / hexahedron
Definition scale_element_hexahedral.F90:10
scale_localmesh_2d
module FElib / Mesh / Local 2D
Definition scale_localmesh_2d.F90:10
scale_localmesh_3d
module FElib / Mesh / Local 3D
Definition scale_localmesh_3d.F90:10
scale_localmeshfield_base
module FElib / Data / base
Definition scale_localmeshfield_base.F90:10
scale_mesh_base3d
module FElib / Mesh / Base 3D
Definition scale_mesh_base3d.F90:10
scale_meshfield_base
module FElib / Data / base
Definition scale_meshfield_base.F90:10
scale_sparsemat
Module common / sparsemat.
Definition scale_sparsemat.F90:11
scale_element_base::elementbase2d
Derived type representing a 2D reference element.
Definition scale_element_base.F90:71
scale_element_base::elementbase3d
Derived type representing a 3D reference element.
Definition scale_element_base.F90:111
scale_element_hexahedral::hexahedralelement
Derived type representing a hexahedral element.
Definition scale_element_hexahedral.F90:31
scale_localmesh_2d::localmesh2d
Derived type representing a local mesh for 2D domain.
Definition scale_localmesh_2d.F90:30
scale_localmesh_3d::localmesh3d
Derived type to manage a local 3D computational domain.
Definition scale_localmesh_3d.F90:32
scale_localmeshfield_base::localmeshfield3d
Derived type representing a field with 3D local mesh.
Definition scale_localmeshfield_base.F90:69
scale_mesh_base3d::meshbase3d
Derived type to manage a computational mesh (base type for 3D domain)
Definition scale_mesh_base3d.F90:37
scale_meshfield_base::meshfield3d
Derived type representing a field with 3D mesh.
Definition scale_meshfield_base.F90:95
scale_sparsemat::sparsemat
Derived type to manage a sparse matrix.
Definition scale_sparsemat.F90:33
Generated by doxygen 1.14.0