FE-Project
Loading...
Searching...
No Matches
scale_atm_dyn_dgm_nonhydro3d_etot_hevi_numflux.F90
Go to the documentation of this file.
1!-------------------------------------------------------------------------------
2!> module FElib / Fluid dyn solver / Atmosphere / Nonhydrostatic model / HEVI / Numflux
3!!
4!! @par Description
5!! HEVE DGM scheme for Atmospheric dynamical process.
6!!
7!! @author Yuta Kawai, Team SCALE
8!<
9!-------------------------------------------------------------------------------
10#include "scaleFElib.h"
12 !-----------------------------------------------------------------------------
13 !
14 !++ Used modules
15 !
16 use scale_precision
17 use scale_io
18 use scale_prc
19 use scale_prof
20 use scale_const, only: &
21 grav => const_grav, &
22 rdry => const_rdry, &
23 cpdry => const_cpdry, &
24 cvdry => const_cvdry, &
25 pres00 => const_pre00
26
28 use scale_element_base, only: &
36
38 dens_vid => prgvar_ddens_id, etot_vid => prgvar_etot_id, &
39 momx_vid => prgvar_momx_id, momy_vid => prgvar_momy_id, &
40 momz_vid => prgvar_momz_id, &
42
43 !-----------------------------------------------------------------------------
44 implicit none
45 private
46 !-----------------------------------------------------------------------------
47 !
48 !++ Public procedures
49 !
52
53 !-----------------------------------------------------------------------------
54 !
55 !++ Public parameters & variables
56 !
57
58 !-----------------------------------------------------------------------------
59 !
60 !++ Private procedures & variables
61 !
62 !-------------------
63
64contains
65
66!OCL SERIAL
68 del_flux, del_flux_hyd, & ! (out)
69 ddens_, momx_, momy_, momz_, etot_, dpres_, dens_hyd, pres_hyd, & ! (in)
70 rtot, cvtot, cptot, & ! (in)
71 gsqrt, g13, g23, zlev, nx, ny, nz, & ! (in)
72 vmapm, vmapp, lmesh, elem, lmesh2d, elem2d ) ! (in)
73
74 implicit none
75
76 class(localmesh3d), intent(in) :: lmesh
77 class(elementbase3d), intent(in) :: elem
78 class(localmesh2d), intent(in) :: lmesh2d
79 class(elementbase2d), intent(in) :: elem2d
80 real(rp), intent(out) :: del_flux(elem%nfptot,lmesh%ne,prgvar_num)
81 real(rp), intent(out) :: del_flux_hyd(elem%nfptot,lmesh%ne,2)
82 real(rp), intent(in) :: ddens_(elem%np*lmesh%nea)
83 real(rp), intent(in) :: momx_(elem%np*lmesh%nea)
84 real(rp), intent(in) :: momy_(elem%np*lmesh%nea)
85 real(rp), intent(in) :: momz_(elem%np*lmesh%nea)
86 real(rp), intent(in) :: etot_(elem%np*lmesh%nea)
87 real(rp), intent(in) :: dpres_(elem%np*lmesh%nea)
88 real(rp), intent(in) :: dens_hyd(elem%np*lmesh%nea)
89 real(rp), intent(in) :: pres_hyd(elem%np*lmesh%nea)
90 real(rp), intent(in) :: rtot (elem%np*lmesh%nea)
91 real(rp), intent(in) :: cvtot(elem%np*lmesh%nea)
92 real(rp), intent(in) :: cptot(elem%np*lmesh%nea)
93 real(rp), intent(in) :: gsqrt(elem%np*lmesh%nea)
94 real(rp), intent(in) :: g13(elem%np*lmesh%nea)
95 real(rp), intent(in) :: g23(elem%np*lmesh%nea)
96 real(rp), intent(in) :: zlev(elem%np*lmesh%ne)
97 real(rp), intent(in) :: nx(elem%nfptot,lmesh%ne)
98 real(rp), intent(in) :: ny(elem%nfptot,lmesh%ne)
99 real(rp), intent(in) :: nz(elem%nfptot,lmesh%ne)
100 integer, intent(in) :: vmapm(elem%nfptot,lmesh%ne)
101 integer, intent(in) :: vmapp(elem%nfptot,lmesh%ne)
102
103 integer :: ke, i, ip(elem%nfptot), im(elem%nfptot)
104 integer :: ke2d
105 real(rp) :: velp(elem%nfptot), velm(elem%nfptot), alpha(elem%nfptot)
106 real(rp) :: velhp(elem%nfptot), velhm(elem%nfptot)
107 real(rp) :: dpresp(elem%nfptot), dpresm(elem%nfptot)
108 real(rp) :: gsqrtdensm(elem%nfptot), gsqrtdensp(elem%nfptot)
109 real(rp) :: gsqrtenthalpym(elem%nfptot), gsqrtenthalpyp(elem%nfptot)
110 real(rp) :: gsqrtddens_p(elem%nfptot), gsqrtddens_m(elem%nfptot)
111 real(rp) :: gsqrtmomx_p(elem%nfptot), gsqrtmomx_m(elem%nfptot)
112 real(rp) :: gsqrtmomy_p(elem%nfptot), gsqrtmomy_m(elem%nfptot)
113 real(rp) :: gsqrtmomz_p(elem%nfptot), gsqrtmomz_m(elem%nfptot)
114 real(rp) :: gsqrtetot_p(elem%nfptot), gsqrtetot_m(elem%nfptot)
115 real(rp) :: phyd_p(elem%nfptot), phyd_m(elem%nfptot)
116 real(rp) :: gsqrt_p(elem%nfptot), gsqrt_m(elem%nfptot)
117 real(rp) :: gsqrtv_p(elem%nfptot), gsqrtv_m(elem%nfptot)
118 real(rp) :: g13_m(elem%nfptot), g13_p(elem%nfptot)
119 real(rp) :: g23_m(elem%nfptot), g23_p(elem%nfptot)
120 real(rp) :: swv(elem%nfptot)
121
122 real(rp) :: gamm, rgamm
123 real(rp) :: rp0
124 real(rp) :: rovp0, p0ovr
125 !------------------------------------------------------------------------
126
127 gamm = cpdry / cvdry
128 rgamm = cvdry / cpdry
129 rp0 = 1.0_rp / pres00
130 rovp0 = rdry * rp0
131 p0ovr = pres00 / rdry
132
133 !$omp parallel do private( &
134 !$omp ke, iM, iP, ke2d, &
135 !$omp alpha, VelM, VelP, VelhM, VelhP, &
136 !$omp dpresM, dpresP, GsqrtDensM, GsqrtDensP, &
137 !$omp GsqrtMOMX_M, GsqrtMOMX_P, GsqrtMOMY_M, GsqrtMOMY_P, GsqrtMOMZ_M, GsqrtMOMZ_P, &
138 !$omp GsqrtDDENS_M, GsqrtDDENS_P, GsqrtETOT_M, GsqrtETOT_P, &
139 !$omp GsqrtEnthalpyM, GsqrtEnthalpyP, &
140 !$omp Phyd_M, Phyd_P, &
141 !$omp Gsqrt_P, Gsqrt_M, GsqrtV_P, GsqrtV_M, G13_P, G13_M, G23_P, G23_M, &
142 !$omp swV )
143 do ke=lmesh%NeS, lmesh%NeE
144 im(:) = vmapm(:,ke); ip(:) = vmapp(:,ke)
145 ke2d = lmesh%EMap3Dto2D(ke)
146
147 gsqrt_m(:) = gsqrt(im)
148 gsqrt_p(:) = gsqrt(ip)
149 gsqrtv_m(:) = gsqrt_m(:)
150 gsqrtv_p(:) = gsqrt_p(:)
151
152 g13_m(:) = g13(im)
153 g13_p(:) = g13(ip)
154 g23_m(:) = g23(im)
155 g23_p(:) = g23(ip)
156
157 gsqrtddens_m(:) = gsqrt_m(:) * ddens_(im)
158 gsqrtddens_p(:) = gsqrt_p(:) * ddens_(ip)
159 gsqrtmomx_m(:) = gsqrt_m(:) * momx_(im)
160 gsqrtmomx_p(:) = gsqrt_p(:) * momx_(ip)
161 gsqrtmomy_m(:) = gsqrt_m(:) * momy_(im)
162 gsqrtmomy_p(:) = gsqrt_p(:) * momy_(ip)
163 gsqrtmomz_m(:) = gsqrt_m(:) * momz_(im)
164 gsqrtmomz_p(:) = gsqrt_p(:) * momz_(ip)
165 gsqrtetot_m(:) = gsqrt_m(:) * etot_(im)
166 gsqrtetot_p(:) = gsqrt_p(:) * etot_(ip)
167 phyd_m(:) = pres_hyd(im)
168 phyd_p(:) = pres_hyd(ip)
169 swv(:) = 1.0_rp - nz(:,ke)**2
170
171 gsqrtdensm(:) = gsqrtddens_m(:) + gsqrt_m(:) * dens_hyd(im)
172 gsqrtdensp(:) = gsqrtddens_p(:) + gsqrt_p(:) * dens_hyd(ip)
173
174 velhm(:) = ( gsqrtmomx_m(:) * nx(:,ke) + gsqrtmomy_m(:) * ny(:,ke) ) / gsqrtdensm(:)
175 velhp(:) = ( gsqrtmomx_p(:) * nx(:,ke) + gsqrtmomy_p(:) * ny(:,ke) ) / gsqrtdensp(:)
176
177 velm(:) = velhm(:) + ( &
178 gsqrtmomz_m(:) / gsqrtv_m(:) + g13_m(:) * gsqrtmomx_m(:) + g23_m(:) * gsqrtmomy_m(:) ) / gsqrtdensm(:) * nz(:,ke)
179 velp(:) = velhp(:) + ( &
180 gsqrtmomz_p(:) / gsqrtv_p(:) + g13_p(:) * gsqrtmomx_p(:) + g23_p(:) * gsqrtmomy_p(:) ) / gsqrtdensp(:) * nz(:,ke)
181
182
183 ! dpresM(:) = ( CPtot(iM) / CVtot(iM) - 1.0_RP ) &
184 ! * ( GsqrtETOT_M(:) - GsqrtDensM(:) * Grav * zlev(iM) &
185 ! - 0.5_RP * ( GsqrtMOMX_M(:)**2 + GsqrtMOMY_M(:)**2 + GsqrtMOMZ_M(:)**2 ) / GsqrtDENSM(:) ) / Gsqrt_M(:) &
186 ! - Phyd_M(:)
187 ! dpresP(:) = ( CPtot(iP) / CVtot(iP) - 1.0_RP ) &
188 ! * ( GsqrtETOT_P(:) - GsqrtDensP(:) * Grav * zlev(iM) &
189 ! - 0.5_RP * ( GsqrtMOMX_P(:)**2 + GsqrtMOMY_P(:)**2 + GsqrtMOMZ_P(:)**2 ) / GsqrtDENSP(:) ) / Gsqrt_P(:) &
190 ! - Phyd_P(:)
191 dpresm(:) = dpres_(im)
192 dpresp(:) = dpres_(ip)
193
194 gsqrtenthalpym(:) = gsqrtetot_m(:) + gsqrt_m(:) * ( phyd_m(:) + dpresm(:) )
195 gsqrtenthalpyp(:) = gsqrtetot_p(:) + gsqrt_p(:) * ( phyd_p(:) + dpresp(:) )
196
197 alpha(:) = swv(:) * max( sqrt( gamm * ( phyd_m(:) + dpresm(:) ) * gsqrt_m(:) / gsqrtdensm(:) ) + abs(velm(:)), &
198 sqrt( gamm * ( phyd_p(:) + dpresp(:) ) * gsqrt_p(:) / gsqrtdensp(:) ) + abs(velp(:)) )
199
200 del_flux(:,ke,dens_vid) = 0.5_rp * ( &
201 ( gsqrtdensp(:) * velhp(:) - gsqrtdensm(:) * velhm(:) ) &
202 - alpha(:) * ( gsqrtddens_p(:) - gsqrtddens_m(:) ) )
203
204 del_flux(:,ke,momx_vid ) = 0.5_rp * ( &
205 ( gsqrtmomx_p(:) * velp(:) - gsqrtmomx_m(:) * velm(:) ) &
206 + ( gsqrt_p(:) * ( nx(:,ke) + g13_p(:) * nz(:,ke)) * dpresp(:) &
207 - gsqrt_m(:) * ( nx(:,ke) + g13_m(:) * nz(:,ke)) * dpresm(:) ) &
208 - alpha(:) * ( gsqrtmomx_p(:) - gsqrtmomx_m(:) ) )
209
210 del_flux(:,ke,momy_vid ) = 0.5_rp * ( &
211 ( gsqrtmomy_p(:) * velp(:) - gsqrtmomy_m(:) * velm(:) ) &
212 + ( gsqrt_p(:) * ( ny(:,ke) + g23_p(:) * nz(:,ke)) * dpresp(:) &
213 - gsqrt_m(:) * ( ny(:,ke) + g23_m(:) * nz(:,ke)) * dpresm(:) ) &
214 - alpha(:) * ( gsqrtmomy_p(:) - gsqrtmomy_m(:) ) )
215
216 del_flux(:,ke,momz_vid ) = 0.5_rp * ( &
217 ( gsqrtmomz_p(:) * velp(:) - gsqrtmomz_m(:) * velm(:) ) &
218 - alpha(:) * ( gsqrtmomz_p(:) - gsqrtmomz_m(:) ) )
219
220 del_flux(:,ke,etot_vid) = 0.5_rp * ( &
221 ( gsqrtenthalpyp(:) * velhp(:) - gsqrtenthalpym(:) * velhm(:) ) &
222 - alpha(:) * ( gsqrtetot_p(:) - gsqrtetot_m(:) ) )
223
224 del_flux_hyd(:,ke,1) = 0.5_rp * ( &
225 gsqrtv_p(:) * ( nx(:,ke) + g13_p(:) * nz(:,ke) ) * phyd_p(:) &
226 - gsqrtv_m(:) * ( nx(:,ke) + g13_m(:) * nz(:,ke) ) * phyd_m(:) )
227
228 del_flux_hyd(:,ke,2) = 0.5_rp * ( &
229 gsqrtv_p(:) * ( ny(:,ke) + g23_p(:) * nz(:,ke) ) * phyd_p(:) &
230 - gsqrtv_m(:) * ( ny(:,ke) + g23_m(:) * nz(:,ke) ) * phyd_m(:) )
231 end do
232
233 return
235
236
237!OCL SERIAL
239 del_flux, del_flux_hyd, & ! (out)
240 ddens_, momx_, momy_, momz_, etot_, dpres_, dens_hyd, pres_hyd, & ! (in)
241 rtot, cvtot, cptot, & ! (in)
242 gsqrt, g11, g12, g22, g_11, g_12, g_22, gsqrth, g13, g23, zlev, & ! (in)
243 nx, ny, nz, vmapm, vmapp, im2dto3d, lmesh, elem, lmesh2d, elem2d ) ! (in)
244
245 implicit none
246
247 class(localmesh3d), intent(in) :: lmesh
248 class(elementbase3d), intent(in) :: elem
249 class(localmesh2d), intent(in) :: lmesh2d
250 class(elementbase2d), intent(in) :: elem2d
251 real(rp), intent(out) :: del_flux(elem%nfptot,lmesh%ne,prgvar_num)
252 real(rp), intent(out) :: del_flux_hyd(elem%nfptot,lmesh%ne,2)
253 real(rp), intent(in) :: ddens_(elem%np*lmesh%nea)
254 real(rp), intent(in) :: momx_(elem%np*lmesh%nea)
255 real(rp), intent(in) :: momy_(elem%np*lmesh%nea)
256 real(rp), intent(in) :: momz_(elem%np*lmesh%nea)
257 real(rp), intent(in) :: etot_(elem%np*lmesh%nea)
258 real(rp), intent(in) :: dpres_(elem%np*lmesh%nea)
259 real(rp), intent(in) :: dens_hyd(elem%np*lmesh%nea)
260 real(rp), intent(in) :: pres_hyd(elem%np*lmesh%nea)
261 real(rp), intent(in) :: rtot (elem%np*lmesh%nea)
262 real(rp), intent(in) :: cvtot(elem%np*lmesh%nea)
263 real(rp), intent(in) :: cptot(elem%np*lmesh%nea)
264 real(rp), intent(in) :: gsqrt(elem%np*lmesh%nea)
265 real(rp), intent(in) :: g11(elem2d%np,lmesh2d%ne)
266 real(rp), intent(in) :: g12(elem2d%np,lmesh2d%ne)
267 real(rp), intent(in) :: g22(elem2d%np,lmesh2d%ne)
268 real(rp), intent(in) :: g_11(elem2d%np,lmesh2d%ne)
269 real(rp), intent(in) :: g_12(elem2d%np,lmesh2d%ne)
270 real(rp), intent(in) :: g_22(elem2d%np,lmesh2d%ne)
271 real(rp), intent(in) :: gsqrth(elem2d%np,lmesh2d%ne)
272 real(rp), intent(in) :: g13(elem%np*lmesh%nea)
273 real(rp), intent(in) :: g23(elem%np*lmesh%nea)
274 real(rp), intent(in) :: zlev(elem%np*lmesh%ne)
275 real(rp), intent(in) :: nx(elem%nfptot,lmesh%ne)
276 real(rp), intent(in) :: ny(elem%nfptot,lmesh%ne)
277 real(rp), intent(in) :: nz(elem%nfptot,lmesh%ne)
278 integer, intent(in) :: vmapm(elem%nfptot,lmesh%ne)
279 integer, intent(in) :: vmapp(elem%nfptot,lmesh%ne)
280 integer, intent(in) :: im2dto3d(elem%nfptot)
281
282 integer :: ke, ip(elem%nfptot), im(elem%nfptot)
283 integer :: ke2d
284 real(rp) :: velp(elem%nfptot), velm(elem%nfptot), alpha(elem%nfptot)
285 real(rp) :: velhp(elem%nfptot), velhm(elem%nfptot)
286 real(rp) :: dpresp(elem%nfptot), dpresm(elem%nfptot)
287 real(rp) :: gsqrtdensm(elem%nfptot), gsqrtdensp(elem%nfptot)
288 real(rp) :: gsqrtenthalpym(elem%nfptot), gsqrtenthalpyp(elem%nfptot)
289 real(rp) :: gsqrtddens_p(elem%nfptot), gsqrtddens_m(elem%nfptot)
290 real(rp) :: gsqrtmomx_p(elem%nfptot), gsqrtmomx_m(elem%nfptot)
291 real(rp) :: gsqrtmomy_p(elem%nfptot), gsqrtmomy_m(elem%nfptot)
292 real(rp) :: gsqrtmomz_p(elem%nfptot), gsqrtmomz_m(elem%nfptot)
293 real(rp) :: gsqrtetot_p(elem%nfptot), gsqrtetot_m(elem%nfptot)
294 real(rp) :: phyd_p(elem%nfptot), phyd_m(elem%nfptot)
295 real(rp) :: gsqrt_p(elem%nfptot), gsqrt_m(elem%nfptot)
296 real(rp) :: gsqrtv_p(elem%nfptot), gsqrtv_m(elem%nfptot)
297 real(rp) :: g13_m(elem%nfptot), g13_p(elem%nfptot)
298 real(rp) :: g23_m(elem%nfptot), g23_p(elem%nfptot)
299 real(rp) :: g1n_m(elem%nfptot), g2n_m(elem%nfptot)
300 real(rp) :: gnn_m(elem%nfptot), gnn_p(elem%nfptot)
301 real(rp) :: gxz_m(elem%nfptot), gxz_p(elem%nfptot)
302 real(rp) :: gyz_m(elem%nfptot), gyz_p(elem%nfptot)
303 real(rp) :: swv(elem%nfptot)
304
305 real(rp) :: gsqrt_u1m(elem%nfptot), gsqrt_u2m(elem%nfptot)
306 real(rp) :: gsqrt_u1p(elem%nfptot), gsqrt_u2p(elem%nfptot)
307
308 real(rp) :: gamm, rgamm
309 real(rp) :: rp0
310 real(rp) :: rovp0, p0ovr
311
312 !------------------------------------------------------------------------
313
314 gamm = cpdry / cvdry
315 rgamm = cvdry / cpdry
316 rp0 = 1.0_rp / pres00
317 rovp0 = rdry * rp0
318 p0ovr = pres00 / rdry
319
320 !$omp parallel do private( &
321 !$omp ke, iM, iP, ke2d, &
322 !$omp alpha, VelM, VelP, VelhM, VelhP, &
323 !$omp dpresM, dpresP, GsqrtDensM, GsqrtDensP, &
324 !$omp GsqrtMOMX_M, GsqrtMOMX_P, GsqrtMOMY_M, GsqrtMOMY_P, GsqrtMOMZ_M, GsqrtMOMZ_P, &
325 !$omp GsqrtDDENS_M, GsqrtDDENS_P, GsqrtETOT_M, GsqrtETOT_P, &
326 !$omp GsqrtEnthalpyM, GsqrtEnthalpyP, &
327 !$omp Gsqrt_u1M, Gsqrt_u2M, Gsqrt_u1P, Gsqrt_u2P, &
328 !$omp Phyd_M, Phyd_P, &
329 !$omp Gsqrt_P, Gsqrt_M, GsqrtV_P, GsqrtV_M, G13_P, G13_M, G23_P, G23_M, &
330 !$omp Gxz_P, Gxz_M, Gyz_P, Gyz_M, G1n_M, G2n_M, Gnn_P, Gnn_M, swV )
331 do ke=lmesh%NeS, lmesh%NeE
332 im(:) = vmapm(:,ke); ip(:) = vmapp(:,ke)
333 ke2d = lmesh%EMap3Dto2D(ke)
334
335 gsqrt_m(:) = gsqrt(im)
336 gsqrt_p(:) = gsqrt(ip)
337 gsqrtv_m(:) = gsqrt_m(:) / gsqrth(im2dto3d(:),ke2d)
338 gsqrtv_p(:) = gsqrt_p(:) / gsqrth(im2dto3d(:),ke2d)
339
340 g13_m(:) = g13(im)
341 g13_p(:) = g13(ip)
342 g23_m(:) = g23(im)
343 g23_p(:) = g23(ip)
344
345 gsqrtddens_m(:) = gsqrt_m(:) * ddens_(im)
346 gsqrtddens_p(:) = gsqrt_p(:) * ddens_(ip)
347 gsqrtmomx_m(:) = gsqrt_m(:) * momx_(im)
348 gsqrtmomx_p(:) = gsqrt_p(:) * momx_(ip)
349 gsqrtmomy_m(:) = gsqrt_m(:) * momy_(im)
350 gsqrtmomy_p(:) = gsqrt_p(:) * momy_(ip)
351 gsqrtmomz_m(:) = gsqrt_m(:) * momz_(im)
352 gsqrtmomz_p(:) = gsqrt_p(:) * momz_(ip)
353 gsqrtetot_m(:) = gsqrt_m(:) * etot_(im)
354 gsqrtetot_p(:) = gsqrt_p(:) * etot_(ip)
355 phyd_m(:) = pres_hyd(im)
356 phyd_p(:) = pres_hyd(ip)
357 swv(:) = 1.0_rp - nz(:,ke)**2
358
359 gxz_m(:) = g11(im2dto3d(:),ke2d) * g13_m(:) + g12(im2dto3d(:),ke2d) * g23_m(:)
360 gxz_p(:) = g11(im2dto3d(:),ke2d) * g13_p(:) + g12(im2dto3d(:),ke2d) * g23_p(:)
361
362 gyz_m(:) = g12(im2dto3d(:),ke2d) * g13_m(:) + g22(im2dto3d(:),ke2d) * g23_m(:)
363 gyz_p(:) = g12(im2dto3d(:),ke2d) * g13_p(:) + g22(im2dto3d(:),ke2d) * g23_p(:)
364
365 g1n_m(:) = g11(im2dto3d(:),ke2d) * nx(:,ke) + g12(im2dto3d(:),ke2d) * ny(:,ke)
366 g2n_m(:) = g12(im2dto3d(:),ke2d) * nx(:,ke) + g22(im2dto3d(:),ke2d) * ny(:,ke)
367
368 gnn_m(:) = g11(im2dto3d(:),ke2d) * abs( nx(:,ke) ) + g22(im2dto3d(:),ke2d) * abs( ny(:,ke) ) &
369 + ( 1.0_rp / gsqrtv_m(:)**2 + g13_m(:) * gxz_m(:) + g23_m(:) * gyz_m(:) ) * abs( nz(:,ke) )
370 gnn_p(:) = g11(im2dto3d(:),ke2d) * abs( nx(:,ke) ) + g22(im2dto3d(:),ke2d) * abs( ny(:,ke) ) &
371 + ( 1.0_rp / gsqrtv_p(:)**2 + g13_p(:) * gxz_p(:) + g23_p(:) * gyz_p(:) ) * abs( nz(:,ke) )
372
373 gsqrtdensm(:) = gsqrtddens_m(:) + gsqrt_m(:) * dens_hyd(im)
374 gsqrtdensp(:) = gsqrtddens_p(:) + gsqrt_p(:) * dens_hyd(ip)
375
376 velhm(:) = ( gsqrtmomx_m(:) * nx(:,ke) + gsqrtmomy_m(:) * ny(:,ke) ) / gsqrtdensm(:)
377 velhp(:) = ( gsqrtmomx_p(:) * nx(:,ke) + gsqrtmomy_p(:) * ny(:,ke) ) / gsqrtdensp(:)
378
379 velm(:) = velhm(:) + ( &
380 gsqrtmomz_m(:) / gsqrtv_m(:) + g13_m(:) * gsqrtmomx_m(:) + g23_m(:) * gsqrtmomy_m(:) ) / gsqrtdensm(:) * nz(:,ke)
381 velp(:) = velhp(:) + ( &
382 gsqrtmomz_p(:) / gsqrtv_p(:) + g13_p(:) * gsqrtmomx_p(:) + g23_p(:) * gsqrtmomy_p(:) ) / gsqrtdensp(:) * nz(:,ke)
383
384 gsqrt_u1m(:) = ( g_11(im2dto3d(:),ke2d) * gsqrtmomx_m(:) + g_12(im2dto3d(:),ke2d) * gsqrtmomy_m(:) )
385 gsqrt_u2m(:) = ( g_12(im2dto3d(:),ke2d) * gsqrtmomx_m(:) + g_22(im2dto3d(:),ke2d) * gsqrtmomy_m(:) )
386 gsqrt_u1p(:) = ( g_11(im2dto3d(:),ke2d) * gsqrtmomx_p(:) + g_12(im2dto3d(:),ke2d) * gsqrtmomy_p(:) )
387 gsqrt_u2p(:) = ( g_12(im2dto3d(:),ke2d) * gsqrtmomx_p(:) + g_22(im2dto3d(:),ke2d) * gsqrtmomy_p(:) )
388
389 ! dpresM(:) = ( CPtot(iM) / CVtot(iM) - 1.0_RP ) &
390 ! * ( GsqrtETOT_M(:) - GsqrtDensM(:) * Grav * zlev(iM) &
391 ! - 0.5_RP * ( GsqrtMOMX_M(:) * Gsqrt_u1M(:) + GsqrtMOMY_M(:) * Gsqrt_u2M(:) + GsqrtMOMZ_M(:)**2 ) / GsqrtDENSM(:) ) / Gsqrt_M(:) &
392 ! - Phyd_M(:)
393 ! dpresP(:) = ( CPtot(iP) / CVtot(iP) - 1.0_RP ) &
394 ! * ( GsqrtETOT_P(:) - GsqrtDensP(:) * Grav * zlev(iM) &
395 ! - 0.5_RP * ( GsqrtMOMX_P(:) * Gsqrt_u1P(:) + GsqrtMOMY_P(:) * Gsqrt_u2P(:) + GsqrtMOMZ_P(:)**2 ) / GsqrtDENSP(:) ) / Gsqrt_P(:) &
396 ! - Phyd_P(:)
397 dpresm(:) = dpres_(im)
398 dpresp(:) = dpres_(ip)
399
400 gsqrtenthalpym(:) = gsqrtetot_m(:) + gsqrt_m(:) * ( phyd_m(:) + dpresm(:) )
401 gsqrtenthalpyp(:) = gsqrtetot_p(:) + gsqrt_p(:) * ( phyd_p(:) + dpresp(:) )
402
403 alpha(:) = swv(:) * max( sqrt( gnn_m(:) * gamm * ( phyd_m(:) + dpresm(:) ) * gsqrt_m(:) / gsqrtdensm(:) ) + abs(velm(:)), &
404 sqrt( gnn_p(:) * gamm * ( phyd_p(:) + dpresp(:) ) * gsqrt_p(:) / gsqrtdensp(:) ) + abs(velp(:)) )
405
406 del_flux(:,ke,dens_vid) = 0.5_rp * ( &
407 ( gsqrtdensp(:) * velhp(:) - gsqrtdensm(:) * velhm(:) ) &
408 - alpha(:) * ( gsqrtddens_p(:) - gsqrtddens_m(:) ) )
409
410 del_flux(:,ke,momx_vid ) = 0.5_rp * ( &
411 ( gsqrtmomx_p(:) * velp(:) - gsqrtmomx_m(:) * velm(:) ) &
412 + ( gsqrt_p(:) * ( g1n_m(:) + gxz_p(:) * nz(:,ke)) * dpresp(:) &
413 - gsqrt_m(:) * ( g1n_m(:) + gxz_m(:) * nz(:,ke)) * dpresm(:) ) &
414 - alpha(:) * ( gsqrtmomx_p(:) - gsqrtmomx_m(:) ) )
415
416 del_flux(:,ke,momy_vid ) = 0.5_rp * ( &
417 ( gsqrtmomy_p(:) * velp(:) - gsqrtmomy_m(:) * velm(:) ) &
418 + ( gsqrt_p(:) * ( g2n_m(:) + gyz_p(:) * nz(:,ke) ) * dpresp(:) &
419 - gsqrt_m(:) * ( g2n_m(:) + gyz_m(:) * nz(:,ke) ) * dpresm(:) ) &
420 - alpha(:) * ( gsqrtmomy_p(:) - gsqrtmomy_m(:) ) )
421
422 del_flux(:,ke,momz_vid ) = 0.5_rp * ( &
423 ( gsqrtmomz_p(:) * velp(:) - gsqrtmomz_m(:) * velm(:) ) &
424 - alpha(:) * ( gsqrtmomz_p(:) - gsqrtmomz_m(:) ) )
425
426 del_flux(:,ke,etot_vid) = 0.5_rp * ( &
427 ( gsqrtenthalpyp(:) * velhp(:) - gsqrtenthalpym(:) * velhm(:) ) &
428 - alpha(:) * ( gsqrtetot_p(:) - gsqrtetot_m(:) ) )
429
430 del_flux_hyd(:,ke,1) = 0.5_rp * ( &
431 gsqrtv_p(:) * ( nx(:,ke) + g13_p(:) * nz(:,ke) ) * phyd_p(:) &
432 - gsqrtv_m(:) * ( nx(:,ke) + g13_m(:) * nz(:,ke) ) * phyd_m(:) )
433
434 del_flux_hyd(:,ke,2) = 0.5_rp * ( &
435 gsqrtv_p(:) * ( ny(:,ke) + g23_p(:) * nz(:,ke) ) * phyd_p(:) &
436 - gsqrtv_m(:) * ( ny(:,ke) + g23_m(:) * nz(:,ke) ) * phyd_m(:) )
437 end do
438
439 return
441
module FElib / Fluid dyn solver / Atmosphere / Nonhydrostatic model / Common
module FElib / Fluid dyn solver / Atmosphere / Nonhydrostatic model / HEVI / Numflux
subroutine, public atm_dyn_dgm_nonhydro3d_etot_hevi_numflux_get_generalhvc(del_flux, del_flux_hyd, ddens_, momx_, momy_, momz_, etot_, dpres_, dens_hyd, pres_hyd, rtot, cvtot, cptot, gsqrt, g11, g12, g22, g_11, g_12, g_22, gsqrth, g13, g23, zlev, nx, ny, nz, vmapm, vmapp, im2dto3d, lmesh, elem, lmesh2d, elem2d)
subroutine, public atm_dyn_dgm_nonhydro3d_etot_hevi_numflux_get_generalvc(del_flux, del_flux_hyd, ddens_, momx_, momy_, momz_, etot_, dpres_, dens_hyd, pres_hyd, rtot, cvtot, cptot, gsqrt, g13, g23, zlev, nx, ny, nz, vmapm, vmapp, lmesh, elem, lmesh2d, elem2d)
module FElib / Element / Base
module FElib / Element / hexahedron
module FElib / Mesh / Local 2D
module FElib / Mesh / Local 3D
module FElib / Mesh / Base 3D
module FElib / Data / base
Module common / sparsemat.
Derived type representing a 2D reference element.
Derived type representing a 3D reference element.
Derived type representing a hexahedral element.
Derived type representing a local mesh for 2D domain.
Derived type to manage a local 3D computational domain.
Derived type representing a field with 3D local mesh.
Derived type to manage a computational mesh (base type for 3D domain)
Derived type representing a field with 3D mesh.