FE-Project
Loading...
Searching...
No Matches
scale_atm_dyn_dgm_nonhydro3d_rhot_heve.F90
Go to the documentation of this file.
1!-------------------------------------------------------------------------------
2!> module FElib / Fluid dyn solver / Atmosphere / Regional nonhydrostatic model / HEVE
3!!
4!! @par Description
5!! HEVE DGM scheme for Atmospheric dynamical process.
6!! The governing equations is a fully compressible nonhydrostatic equations,
7!! which consist of mass, momentum, and thermodynamics (density * potential temperature conservation) equations.
8!!
9!! @author Yuta Kawai, Xuanzhengbo Ren, and Team SCALE
10!<
11!-------------------------------------------------------------------------------
12#include "scaleFElib.h"
13module scale_atm_dyn_dgm_nonhydro3d_rhot_heve
14 !-----------------------------------------------------------------------------
15 !
16 !++ Used modules
17 !
18 use scale_precision
19 use scale_io
20 use scale_prc
21 use scale_prof
22 use scale_const, only: &
23 grav => const_grav, &
24 rdry => const_rdry, &
25 cpdry => const_cpdry, &
26 cvdry => const_cvdry, &
27 pres00 => const_pre00
28
29 use scale_sparsemat
30 use scale_element_base, only: &
31 elementbase2d, elementbase3d
32 use scale_element_hexahedral, only: hexahedralelement
33 use scale_element_operation_base, only: elementoperationbase3d
34 use scale_localmesh_2d, only: localmesh2d
35 use scale_localmesh_3d, only: localmesh3d
36 use scale_mesh_base3d, only: meshbase3d
37 use scale_localmeshfield_base, only: localmeshfield3d
38 use scale_meshfield_base, only: meshfield3d
39
40 use scale_atm_dyn_dgm_nonhydro3d_common, only: &
41 atm_dyn_dgm_nonhydro3d_common_init, &
42 atm_dyn_dgm_nonhydro3d_common_final, &
43 dens_vid => prgvar_ddens_id, rhot_vid => prgvar_drhot_id, &
44 momx_vid => prgvar_momx_id, momy_vid => prgvar_momy_id, &
45 momz_vid => prgvar_momz_id, &
46 prgvar_num, intrpmat_vpordm1
47
48
49 !-----------------------------------------------------------------------------
50 implicit none
51 private
52 !-----------------------------------------------------------------------------
53 !
54 !++ Public procedures
55 !
56 public :: atm_dyn_dgm_nonhydro3d_rhot_heve_init
57 public :: atm_dyn_dgm_nonhydro3d_rhot_heve_final
58 public :: atm_dyn_dgm_nonhydro3d_rhot_heve_cal_tend_asis
59 public :: atm_dyn_dgm_nonhydro3d_rhot_heve_cal_tend
60 public :: atm_dyn_dgm_nonhydro3d_rhot_heve_cal_tend_cco
61
62 !-----------------------------------------------------------------------------
63 !
64 !++ Public parameters & variables
65 !
66
67 !-----------------------------------------------------------------------------
68 !
69 !++ Private procedures & variables
70 !
71 !-------------------
72
73contains
74 subroutine atm_dyn_dgm_nonhydro3d_rhot_heve_init( mesh )
75 implicit none
76 class(meshbase3d), intent(in) :: mesh
77 !--------------------------------------------
78
79 call atm_dyn_dgm_nonhydro3d_common_init( mesh )
80 return
81 end subroutine atm_dyn_dgm_nonhydro3d_rhot_heve_init
82
83
84 subroutine atm_dyn_dgm_nonhydro3d_rhot_heve_final()
85 implicit none
86 !--------------------------------------------
87
88 call atm_dyn_dgm_nonhydro3d_common_final()
89 return
90 end subroutine atm_dyn_dgm_nonhydro3d_rhot_heve_final
91
92 !-------------------------------
93
94!OCL SERIAL
95 subroutine atm_dyn_dgm_nonhydro3d_rhot_heve_cal_tend_asis( &
96 DENS_dt, MOMX_dt, MOMY_dt, MOMZ_dt, RHOT_dt, & ! (out)
97 ddens_, momx_, momy_, momz_, drhot_, dpres_, & ! (in)
98 dens_hyd, pres_hyd, pres_hyd_ref, therm_hyd, & ! (in)
99 coriolis, rtot, cvtot, cptot, dphyddx, dphyddy, & ! (in)
100 element3d_operation, dx, dy, dz, sx, sy, sz, lift, & ! (in)
101 lmesh, elem, lmesh2d, elem2d ) ! (in)
102
103 use scale_atm_dyn_dgm_nonhydro3d_rhot_heve_numflux, only: &
104 get_ebnd_flux => atm_dyn_dgm_nonhydro3d_rhot_heve_numflux_get_generalvc_asis
105
106 implicit none
107
108 class(localmesh3d), intent(in) :: lmesh
109 class(elementbase3d), intent(in) :: elem
110 class(localmesh2d), intent(in) :: lmesh2d
111 class(elementbase2d), intent(in) :: elem2d
112 class(elementoperationbase3d), intent(in) :: element3d_operation
113 type(sparsemat), intent(in) :: dx, dy, dz, sx, sy, sz, lift
114 real(rp), intent(out) :: dens_dt(elem%np,lmesh%nea)
115 real(rp), intent(out) :: momx_dt(elem%np,lmesh%nea)
116 real(rp), intent(out) :: momy_dt(elem%np,lmesh%nea)
117 real(rp), intent(out) :: momz_dt(elem%np,lmesh%nea)
118 real(rp), intent(out) :: rhot_dt(elem%np,lmesh%nea)
119 real(rp), intent(in) :: ddens_(elem%np,lmesh%nea)
120 real(rp), intent(in) :: momx_(elem%np,lmesh%nea)
121 real(rp), intent(in) :: momy_(elem%np,lmesh%nea)
122 real(rp), intent(in) :: momz_(elem%np,lmesh%nea)
123 real(rp), intent(in) :: drhot_(elem%np,lmesh%nea)
124 real(rp), intent(in) :: dpres_(elem%np,lmesh%nea)
125 real(rp), intent(in) :: dens_hyd(elem%np,lmesh%nea)
126 real(rp), intent(in) :: pres_hyd(elem%np,lmesh%nea)
127 real(rp), intent(in) :: pres_hyd_ref(elem%np,lmesh%nea)
128 real(rp), intent(in) :: therm_hyd(elem%np,lmesh%nea)
129 real(rp), intent(in) :: coriolis(elem2d%np,lmesh2d%nea)
130 real(rp), intent(in) :: rtot(elem%np,lmesh%nea)
131 real(rp), intent(in) :: cvtot(elem%np,lmesh%nea)
132 real(rp), intent(in) :: cptot(elem%np,lmesh%nea)
133 real(rp), intent(in) :: dphyddx(elem%np,lmesh%nea)
134 real(rp), intent(in) :: dphyddy(elem%np,lmesh%nea)
135
136 real(rp) :: fx(elem%np), fy(elem%np), fz(elem%np), liftdelflx(elem%np)
137 real(rp) :: dpres_hyd(elem%np), gradphyd_x(elem%np), gradphyd_y(elem%np)
138 real(rp) :: del_flux(elem%nfptot,lmesh%ne,prgvar_num)
139 real(rp) :: del_flux_hyd(elem%nfptot,lmesh%ne,2)
140 real(rp) :: rhot_(elem%np)
141 real(rp) :: rdens_(elem%np), u_(elem%np), v_(elem%np), w_(elem%np), wt_(elem%np)
142 real(rp) :: cori(elem%np)
143 real(rp) :: drho(elem%np)
144 real(rp) :: gsqrtv(elem%np), rgsqrtv(elem%np)
145
146 integer :: ke, ke2d
147
148 real(rp) :: gamm, rgamm
149 real(rp) :: rp0
150 real(rp) :: rovp0, p0ovr
151 !------------------------------------------------------------------------
152
153 call prof_rapstart('cal_dyn_tend_bndflux', 3)
154 call get_ebnd_flux( &
155 del_flux, del_flux_hyd, & ! (out)
156 ddens_, momx_, momy_, momz_, drhot_, dpres_, dens_hyd, pres_hyd, & ! (in)
157 rtot, cvtot, cptot, & ! (in)
158 lmesh%Gsqrt, lmesh%GI3(:,:,1), lmesh%GI3(:,:,2), & ! (in)
159 lmesh%normal_fn(:,:,1), lmesh%normal_fn(:,:,2), lmesh%normal_fn(:,:,3), & ! (in)
160 lmesh%vmapM, lmesh%vmapP, & ! (in)
161 lmesh, elem, lmesh2d, elem2d ) ! (in)
162 call prof_rapend('cal_dyn_tend_bndflux', 3)
163
164 !-----
165 call prof_rapstart('cal_dyn_tend_interior', 3)
166 gamm = cpdry / cvdry
167 rgamm = cvdry / cpdry
168 rp0 = 1.0_rp / pres00
169 rovp0 = rdry * rp0
170 p0ovr = pres00 / rdry
171
172 !$omp parallel do private( ke, ke2d, Cori, &
173 !$omp RHOT_, rdens_, u_, v_, w_, wt_, &
174 !$omp drho, DPRES_hyd, GradPhyd_x, GradPhyd_y, &
175 !$omp GsqrtV, RGsqrtV, &
176 !$omp Fx, Fy, Fz, LiftDelFlx )
177 do ke = lmesh%NeS, lmesh%NeE
178 !--
179 ke2d = lmesh%EMap3Dto2D(ke)
180 cori(:) = coriolis(elem%IndexH2Dto3D(:),ke2d)
181
182 gsqrtv(:) = lmesh%Gsqrt(:,ke) / lmesh%GsqrtH(elem%IndexH2Dto3D,ke2d)
183 rgsqrtv(:) = 1.0_rp / gsqrtv(:)
184
185 !--
186 rhot_(:) = p0ovr * (pres_hyd(:,ke) * rp0)**rgamm + drhot_(:,ke)
187 ! DPRES_(:) = PRES00 * ( Rtot(:,ke) * rP0 * RHOT_(:) )**( CPtot(:,ke) / CVtot(:,ke) ) &
188 ! - PRES_hyd(:,ke)
189
190 rdens_(:) = 1.0_rp / (ddens_(:,ke) + dens_hyd(:,ke))
191 u_(:) = momx_(:,ke) * rdens_(:)
192 v_(:) = momy_(:,ke) * rdens_(:)
193 w_(:) = momz_(:,ke) * rdens_(:)
194 wt_(:) = w_(:) * rgsqrtv(:) + lmesh%GI3(:,ke,1) * u_(:) + lmesh%GI3(:,ke,2) * v_(:)
195
196 drho(:) = matmul(intrpmat_vpordm1, ddens_(:,ke))
197
198 !-- Gradient hydrostatic pressure
199
200 dpres_hyd(:) = pres_hyd(:,ke) - pres_hyd_ref(:,ke)
201
202 call sparsemat_matmul(dx, gsqrtv(:) * dpres_hyd(:), fx)
203 call sparsemat_matmul(dz, gsqrtv(:) * lmesh%GI3(:,ke,1) * dpres_hyd(:), fz)
204 call sparsemat_matmul(lift, lmesh%Fscale(:,ke) * del_flux_hyd(:,ke,1), liftdelflx)
205 gradphyd_x(:) = lmesh%Escale(:,ke,1,1) * fx(:) &
206 + lmesh%Escale(:,ke,3,3) * fz(:) &
207 + liftdelflx(:)
208
209 call sparsemat_matmul(dy, gsqrtv(:) * dpres_hyd(:), fy)
210 call sparsemat_matmul(dz, gsqrtv(:) * lmesh%GI3(:,ke,2) * dpres_hyd(:), fz)
211 call sparsemat_matmul(lift, lmesh%Fscale(:,ke) * del_flux_hyd(:,ke,2), liftdelflx)
212 gradphyd_y(:) = lmesh%Escale(:,ke,2,2) * fy(:) &
213 + lmesh%Escale(:,ke,3,3) * fz(:) &
214 + liftdelflx(:)
215
216 !-- DENS
217 call sparsemat_matmul(dx, lmesh%Gsqrt(:,ke) * momx_(:,ke), fx)
218 call sparsemat_matmul(dy, lmesh%Gsqrt(:,ke) * momy_(:,ke), fy)
219 call sparsemat_matmul(dz, lmesh%Gsqrt(:,ke) * ( ddens_(:,ke) + dens_hyd(:,ke) ) * wt_(:), fz)
220 call sparsemat_matmul(lift, lmesh%Fscale(:,ke) * del_flux(:,ke,dens_vid), liftdelflx)
221
222 dens_dt(:,ke) = - ( &
223 lmesh%Escale(:,ke,1,1) * fx(:) &
224 + lmesh%Escale(:,ke,2,2) * fy(:) &
225 + lmesh%Escale(:,ke,3,3) * fz(:) &
226 + liftdelflx(:) ) / lmesh%Gsqrt(:,ke)
227
228 !-- MOMX
229 call sparsemat_matmul(dx, lmesh%Gsqrt(:,ke) * ( u_(:) * momx_(:,ke) + dpres_(:,ke) ), fx)
230 call sparsemat_matmul(dy, lmesh%Gsqrt(:,ke) * v_(:) * momx_(:,ke) , fy)
231 call sparsemat_matmul(dz, lmesh%Gsqrt(:,ke) * ( wt_(:) * momx_(:,ke) &
232 + lmesh%GI3(:,ke,1) * dpres_(:,ke) ), fz)
233 call sparsemat_matmul(lift, lmesh%Fscale(:,ke) * del_flux(:,ke,momx_vid), liftdelflx)
234
235 momx_dt(:,ke) = &
236 - ( lmesh%Escale(:,ke,1,1) * fx(:) &
237 + lmesh%Escale(:,ke,2,2) * fy(:) &
238 + lmesh%Escale(:,ke,3,3) * fz(:) &
239 + liftdelflx(:) ) / lmesh%Gsqrt(:,ke) &
240 - gradphyd_x(:) * rgsqrtv(:) &
241 + cori(:) * momy_(:,ke)
242
243 !-- MOMY
244 call sparsemat_matmul(dx, lmesh%Gsqrt(:,ke) * u_(:) * momy_(:,ke) , fx)
245 call sparsemat_matmul(dy, lmesh%Gsqrt(:,ke) * ( v_(:) * momy_(:,ke) + dpres_(:,ke) ), fy)
246 call sparsemat_matmul(dz, lmesh%Gsqrt(:,ke) * ( wt_(:) * momy_(:,ke) &
247 + lmesh%GI3(:,ke,2) * dpres_(:,ke) ), fz)
248 call sparsemat_matmul(lift, lmesh%Fscale(:,ke) * del_flux(:,ke,momy_vid), liftdelflx)
249
250 momy_dt(:,ke) = &
251 - ( lmesh%Escale(:,ke,1,1) * fx(:) &
252 + lmesh%Escale(:,ke,2,2) * fy(:) &
253 + lmesh%Escale(:,ke,3,3) * fz(:) &
254 + liftdelflx(:) ) / lmesh%Gsqrt(:,ke) &
255 - gradphyd_y(:) * rgsqrtv(:) &
256 - cori(:) * momx_(:,ke)
257
258 !-- MOMZ
259 call sparsemat_matmul(dx, lmesh%Gsqrt(:,ke) * u_(:) * momz_(:,ke), fx)
260 call sparsemat_matmul(dy, lmesh%Gsqrt(:,ke) * v_(:) * momz_(:,ke), fy)
261 call sparsemat_matmul(dz, lmesh%Gsqrt(:,ke) * ( wt_(:) * momz_(:,ke) &
262 + rgsqrtv(:) * dpres_(:,ke) ), fz)
263 call sparsemat_matmul(lift, lmesh%Fscale(:,ke) * del_flux(:,ke,momz_vid), liftdelflx)
264
265 momz_dt(:,ke) = &
266 - ( lmesh%Escale(:,ke,1,1) * fx(:) &
267 + lmesh%Escale(:,ke,2,2) * fy(:) &
268 + lmesh%Escale(:,ke,3,3) * fz(:) &
269 + liftdelflx(:) ) / lmesh%Gsqrt(:,ke) &
270 - grav * drho(:)
271
272 !-- RHOT
273 call sparsemat_matmul(dx, lmesh%Gsqrt(:,ke) * u_(:) * rhot_(:), fx)
274 call sparsemat_matmul(dy, lmesh%Gsqrt(:,ke) * v_(:) * rhot_(:), fy)
275 call sparsemat_matmul(dz, lmesh%Gsqrt(:,ke) * wt_(:) * rhot_(:), fz)
276 call sparsemat_matmul(lift, lmesh%Fscale(:,ke) * del_flux(:,ke,rhot_vid), liftdelflx)
277
278 rhot_dt(:,ke) = &
279 - ( lmesh%Escale(:,ke,1,1) * fx(:) &
280 + lmesh%Escale(:,ke,2,2) * fy(:) &
281 + lmesh%Escale(:,ke,3,3) * fz(:) &
282 + liftdelflx(:) ) / lmesh%Gsqrt(:,ke)
283 end do
284
285 call prof_rapend('cal_dyn_tend_interior', 3)
286
287 return
288 end subroutine atm_dyn_dgm_nonhydro3d_rhot_heve_cal_tend_asis
289
290
291!OCL SERIAL
292 subroutine atm_dyn_dgm_nonhydro3d_rhot_heve_cal_tend( &
293 DENS_dt, MOMX_dt, MOMY_dt, MOMZ_dt, RHOT_dt, & ! (out)
294 ddens_, momx_, momy_, momz_, drhot_, dpres_, & ! (in)
295 dens_hyd, pres_hyd, pres_hyd_ref, therm_hyd, & ! (in)
296 coriolis, rtot, cvtot, cptot, dphyddx, dphyddy, & ! (in)
297 element3d_operation, dx, dy, dz, sx, sy, sz, lift, & ! (in)
298 lmesh, elem, lmesh2d, elem2d ) ! (in)
299
300 use scale_atm_dyn_dgm_nonhydro3d_rhot_heve_numflux, only: &
301 get_ebnd_flux => atm_dyn_dgm_nonhydro3d_rhot_heve_numflux_get_generalvc
302
303 implicit none
304
305 class(localmesh3d), intent(in) :: lmesh
306 class(elementbase3d), intent(in) :: elem
307 class(localmesh2d), intent(in) :: lmesh2d
308 class(elementbase2d), intent(in) :: elem2d
309 class(elementoperationbase3d), intent(in) :: element3d_operation
310 type(sparsemat), intent(in) :: dx, dy, dz, sx, sy, sz, lift
311 real(rp), intent(out) :: dens_dt(elem%np,lmesh%nea)
312 real(rp), intent(out) :: momx_dt(elem%np,lmesh%nea)
313 real(rp), intent(out) :: momy_dt(elem%np,lmesh%nea)
314 real(rp), intent(out) :: momz_dt(elem%np,lmesh%nea)
315 real(rp), intent(out) :: rhot_dt(elem%np,lmesh%nea)
316 real(rp), intent(in) :: ddens_(elem%np,lmesh%nea)
317 real(rp), intent(in) :: momx_(elem%np,lmesh%nea)
318 real(rp), intent(in) :: momy_(elem%np,lmesh%nea)
319 real(rp), intent(in) :: momz_(elem%np,lmesh%nea)
320 real(rp), intent(in) :: drhot_(elem%np,lmesh%nea)
321 real(rp), intent(in) :: dpres_(elem%np,lmesh%nea)
322 real(rp), intent(in) :: dens_hyd(elem%np,lmesh%nea)
323 real(rp), intent(in) :: pres_hyd(elem%np,lmesh%nea)
324 real(rp), intent(in) :: pres_hyd_ref(elem%np,lmesh%nea)
325 real(rp), intent(in) :: therm_hyd(elem%np,lmesh%nea)
326 real(rp), intent(in) :: coriolis(elem2d%np,lmesh2d%nea)
327 real(rp), intent(in) :: rtot(elem%np,lmesh%nea)
328 real(rp), intent(in) :: cvtot(elem%np,lmesh%nea)
329 real(rp), intent(in) :: cptot(elem%np,lmesh%nea)
330 real(rp), intent(in) :: dphyddx(elem%np,lmesh%nea)
331 real(rp), intent(in) :: dphyddy(elem%np,lmesh%nea)
332
333 real(rp) :: flux(elem%np,3,5), dflux(elem%np,4,5)
334 real(rp) :: del_flux(elem%nfptot,prgvar_num,lmesh%ne)
335 real(rp) :: u_, v_, w_, pt_
336 real(rp) :: cor
337 real(rp) :: drho(elem%np)
338 real(rp) :: rdens_(elem%np), gsqrtv(elem%np), rgsqrtv(elem%np), rgsqrt(elem%np)
339 real(rp) :: gsqrt_, gsqrtdpres_, e11, e22, e33
340
341 integer :: ke, ke2d
342 integer :: p
343
344 real(rp) :: gamm, rgamm
345 real(rp) :: rp0
346 real(rp) :: rovp0, p0ovr
347 !------------------------------------------------------------------------
348
349 call prof_rapstart('cal_dyn_tend_bndflux', 3)
350 call get_ebnd_flux( &
351 del_flux, & ! (out)
352 ddens_, momx_, momy_, momz_, drhot_, dpres_, dens_hyd, pres_hyd, therm_hyd, & ! (in)
353 rtot, cvtot, cptot, & ! (in)
354 lmesh%Gsqrt, lmesh%GI3(:,:,1), lmesh%GI3(:,:,2), & ! (in)
355 lmesh%normal_fn(:,:,1), lmesh%normal_fn(:,:,2), lmesh%normal_fn(:,:,3), & ! (in)
356 lmesh%vmapM, lmesh%vmapP, & ! (in)
357 lmesh, elem, lmesh2d, elem2d ) ! (in)
358 call prof_rapend('cal_dyn_tend_bndflux', 3)
359
360 call prof_rapstart('cal_dyn_tend_interior', 3)
361 gamm = cpdry / cvdry
362 rgamm = cvdry / cpdry
363 rp0 = 1.0_rp / pres00
364 rovp0 = rdry * rp0
365 p0ovr = pres00 / rdry
366
367 !$omp parallel do private( ke, ke2d, p, cor, &
368 !$omp u_, v_, w_, pt_, &
369 !$omp drho, &
370 !$omp RDENS_, RGsqrt, GsqrtV, RGsqrtV, Gsqrt_, GsqrtDPRES_, &
371 !$omp E11, E22, E33, DFlux, Flux )
372 do ke = lmesh%NeS, lmesh%NeE
373 !--
374 ke2d = lmesh%EMap3Dto2D(ke)
375
376 do p=1, elem%Np
377 gsqrtv(p) = lmesh%Gsqrt(p,ke) / lmesh%GsqrtH(elem%IndexH2Dto3D(p),ke2d)
378 rgsqrtv(p) = 1.0_rp / gsqrtv(p)
379 rgsqrt(p) = 1.0_rp / lmesh%Gsqrt(p,ke)
380 rdens_(p) = 1.0_rp / (ddens_(p,ke) + dens_hyd(p,ke))
381 end do
382
383 !--
384
385 do p=1, elem%Np
386 gsqrt_ = lmesh%Gsqrt(p,ke)
387 flux(p,1,dens_vid) = gsqrt_ * momx_(p,ke)
388 flux(p,2,dens_vid) = gsqrt_ * momy_(p,ke)
389 flux(p,3,dens_vid) = gsqrt_ * ( &
390 momz_(p,ke) * rgsqrtv(p) &
391 + lmesh%GI3(p,ke,1) * momx_(p,ke) &
392 + lmesh%GI3(p,ke,2) * momy_(p,ke) )
393 end do
394 do p=1, elem%Np
395 pt_ = ( therm_hyd(p,ke) + drhot_(p,ke) ) * rdens_(p)
396
397 flux(p,1,rhot_vid) = flux(p,1,dens_vid) * pt_
398 flux(p,2,rhot_vid) = flux(p,2,dens_vid) * pt_
399 flux(p,3,rhot_vid) = flux(p,3,dens_vid) * pt_
400
401 w_ = momz_(p,ke) * rdens_(p)
402 flux(p,1,momz_vid) = flux(p,1,dens_vid) * w_
403 flux(p,2,momz_vid) = flux(p,2,dens_vid) * w_
404 flux(p,3,momz_vid) = flux(p,3,dens_vid) * w_ + lmesh%Gsqrt(p,ke) * dpres_(p,ke) * rgsqrtv(p)
405 end do
406 do p=1, elem%Np
407 gsqrtdpres_ = lmesh%Gsqrt(p,ke) * dpres_(p,ke)
408
409 u_ = momx_(p,ke) * rdens_(p)
410 flux(p,1,momx_vid) = flux(p,1,dens_vid) * u_ + gsqrtdpres_
411 flux(p,2,momx_vid) = flux(p,2,dens_vid) * u_
412 flux(p,3,momx_vid) = flux(p,3,dens_vid) * u_ + gsqrtdpres_ * lmesh%GI3(p,ke,1)
413
414 v_ = momy_(p,ke) * rdens_(p)
415 flux(p,1,momy_vid) = flux(p,1,dens_vid) * v_
416 flux(p,2,momy_vid) = flux(p,2,dens_vid) * v_ + gsqrtdpres_
417 flux(p,3,momy_vid) = flux(p,3,dens_vid) * v_ + gsqrtdpres_ * lmesh%GI3(p,ke,2)
418 end do
419
420 call element3d_operation%Div_var5( &
421 flux, del_flux(:,:,ke), & ! (in)
422 dflux ) ! (out)
423
424 do p=1, elem%Np
425 e11 = lmesh%Escale(p,ke,1,1)
426 e22 = lmesh%Escale(p,ke,2,2)
427 e33 = lmesh%Escale(p,ke,3,3)
428
429 dens_dt(p,ke) = - ( &
430 e11 * dflux(p,1,dens_vid) &
431 + e22 * dflux(p,2,dens_vid) &
432 + e33 * dflux(p,3,dens_vid) &
433 + dflux(p,4,dens_vid) ) * rgsqrt(p)
434
435 rhot_dt(p,ke) = - ( &
436 e11 * dflux(p,1,rhot_vid) &
437 + e22 * dflux(p,2,rhot_vid) &
438 + e33 * dflux(p,3,rhot_vid) &
439 + dflux(p,4,rhot_vid) ) * rgsqrt(p)
440 end do
441
442 call element3d_operation%VFilterPM1( ddens_(:,ke), & ! (in)
443 drho ) ! (out)
444
445 do p=1, elem%Np
446 e11 = lmesh%Escale(p,ke,1,1)
447 e22 = lmesh%Escale(p,ke,2,2)
448 e33 = lmesh%Escale(p,ke,3,3)
449
450 momz_dt(p,ke) = - ( &
451 e11 * dflux(p,1,momz_vid) &
452 + e22 * dflux(p,2,momz_vid) &
453 + e33 * dflux(p,3,momz_vid) &
454 + dflux(p,4,momz_vid) ) * rgsqrt(p) &
455 - grav * drho(p)
456 end do
457
458 !--
459 do p=1, elem%Np
460 cor = coriolis(elem%IndexH2Dto3D(p),ke2d)
461 momx_dt(p,ke) = - dphyddx(p,ke) &
462 + cor * momy_(p,ke)
463 momy_dt(p,ke) = - dphyddy(p,ke) &
464 - cor * momx_(p,ke)
465 end do
466
467 do p=1, elem%Np
468 e11 = lmesh%Escale(p,ke,1,1)
469 e22 = lmesh%Escale(p,ke,2,2)
470 e33 = lmesh%Escale(p,ke,3,3)
471
472 momx_dt(p,ke) = momx_dt(p,ke) - ( &
473 e11 * dflux(p,1,momx_vid) &
474 + e22 * dflux(p,2,momx_vid) &
475 + e33 * dflux(p,3,momx_vid) &
476 + dflux(p,4,momx_vid) ) * rgsqrt(p)
477
478 momy_dt(p,ke) = momy_dt(p,ke) - ( &
479 e11 * dflux(p,1,momy_vid) &
480 + e22 * dflux(p,2,momy_vid) &
481 + e33 * dflux(p,3,momy_vid) &
482 + dflux(p,4,momy_vid) ) * rgsqrt(p)
483 end do
484 end do
485
486 call prof_rapend('cal_dyn_tend_interior', 3)
487
488 return
489 end subroutine atm_dyn_dgm_nonhydro3d_rhot_heve_cal_tend
490
491
492!OCL SERIAL
493 subroutine atm_dyn_dgm_nonhydro3d_rhot_heve_cal_tend_cco( &
494 DENS_dt, MOMX_dt, MOMY_dt, MOMZ_dt, RHOT_dt, & ! (out)
495 ddens_, momx_, momy_, momz_, drhot_, dpres_, & ! (in)
496 dens_hyd, pres_hyd, pres_hyd_ref, therm_hyd, & ! (in)
497 coriolis, rtot, cvtot, cptot, dphyddx, dphyddy, & ! (in)
498 element3d_operation, dx, dy, dz, sx, sy, sz, lift, & ! (in)
499 lmesh, elem, lmesh2d, elem2d ) ! (in)
500
501 use scale_atm_dyn_dgm_nonhydro3d_rhot_heve_numflux, only: &
502 get_ebnd_flux => atm_dyn_dgm_nonhydro3d_rhot_heve_numflux_get_generalvc
503
504 implicit none
505
506 class(localmesh3d), intent(in) :: lmesh
507 class(elementbase3d), intent(in) :: elem
508 class(localmesh2d), intent(in) :: lmesh2d
509 class(elementbase2d), intent(in) :: elem2d
510 class(elementoperationbase3d), intent(in) :: element3d_operation
511 type(sparsemat), intent(in) :: dx, dy, dz, sx, sy, sz, lift
512 real(rp), intent(out) :: dens_dt(elem%np,lmesh%nea)
513 real(rp), intent(out) :: momx_dt(elem%np,lmesh%nea)
514 real(rp), intent(out) :: momy_dt(elem%np,lmesh%nea)
515 real(rp), intent(out) :: momz_dt(elem%np,lmesh%nea)
516 real(rp), intent(out) :: rhot_dt(elem%np,lmesh%nea)
517 real(rp), intent(in) :: ddens_(elem%np,lmesh%nea)
518 real(rp), intent(in) :: momx_(elem%np,lmesh%nea)
519 real(rp), intent(in) :: momy_(elem%np,lmesh%nea)
520 real(rp), intent(in) :: momz_(elem%np,lmesh%nea)
521 real(rp), intent(in) :: drhot_(elem%np,lmesh%nea)
522 real(rp), intent(in) :: dpres_(elem%np,lmesh%nea)
523 real(rp), intent(in) :: dens_hyd(elem%np,lmesh%nea)
524 real(rp), intent(in) :: pres_hyd(elem%np,lmesh%nea)
525 real(rp), intent(in) :: pres_hyd_ref(elem%np,lmesh%nea)
526 real(rp), intent(in) :: therm_hyd(elem%np,lmesh%nea)
527 real(rp), intent(in) :: coriolis(elem2d%np,lmesh2d%nea)
528 real(rp), intent(in) :: rtot(elem%np,lmesh%nea)
529 real(rp), intent(in) :: cvtot(elem%np,lmesh%nea)
530 real(rp), intent(in) :: cptot(elem%np,lmesh%nea)
531 real(rp), intent(in) :: dphyddx(elem%np,lmesh%nea)
532 real(rp), intent(in) :: dphyddy(elem%np,lmesh%nea)
533
534 real(rp) :: flux(elem%np,3,5), dflux(elem%np,3,5)
535 real(rp) :: u_, v_, w_, pt_
536 real(rp) :: cor
537 real(rp) :: drho(elem%np)
538 real(rp) :: rdens_(elem%np), gsqrtv(elem%np), rgsqrtv(elem%np), rgsqrt(elem%np)
539 real(rp) :: gsqrt_, gsqrtdpres_, e11, e22, e33
540
541 integer :: ke, ke2d
542 integer :: p
543
544 real(rp) :: gamm, rgamm
545 real(rp) :: rp0
546 real(rp) :: rovp0, p0ovr
547 !------------------------------------------------------------------------
548
549 call prof_rapstart('cal_dyn_tend_interior', 3)
550 gamm = cpdry / cvdry
551 rgamm = cvdry / cpdry
552 rp0 = 1.0_rp / pres00
553 rovp0 = rdry * rp0
554 p0ovr = pres00 / rdry
555
556 !$omp parallel do private( ke, ke2d, p, cor, &
557 !$omp u_, v_, w_, pt_, &
558 !$omp drho, &
559 !$omp RDENS_, RGsqrt, GsqrtV, RGsqrtV, Gsqrt_, GsqrtDPRES_, &
560 !$omp E11, E22, E33, DFlux, Flux )
561 do ke = lmesh%NeS, lmesh%NeE
562 !--
563 ke2d = lmesh%EMap3Dto2D(ke)
564
565 do p=1, elem%Np
566 gsqrtv(p) = lmesh%Gsqrt(p,ke) / lmesh%GsqrtH(elem%IndexH2Dto3D(p),ke2d)
567 rgsqrtv(p) = 1.0_rp / gsqrtv(p)
568 rgsqrt(p) = 1.0_rp / lmesh%Gsqrt(p,ke)
569 rdens_(p) = 1.0_rp / (ddens_(p,ke) + dens_hyd(p,ke))
570 end do
571
572 !--
573
574 do p=1, elem%Np
575 gsqrt_ = lmesh%Gsqrt(p,ke)
576 flux(p,1,dens_vid) = gsqrt_ * momx_(p,ke)
577 flux(p,2,dens_vid) = gsqrt_ * momy_(p,ke)
578 flux(p,3,dens_vid) = gsqrt_ * ( &
579 momz_(p,ke) * rgsqrtv(p) &
580 + lmesh%GI3(p,ke,1) * momx_(p,ke) &
581 + lmesh%GI3(p,ke,2) * momy_(p,ke) )
582 end do
583 do p=1, elem%Np
584 pt_ = ( therm_hyd(p,ke) + drhot_(p,ke) ) * rdens_(p)
585
586 flux(p,1,rhot_vid) = flux(p,1,dens_vid) * pt_
587 flux(p,2,rhot_vid) = flux(p,2,dens_vid) * pt_
588 flux(p,3,rhot_vid) = flux(p,3,dens_vid) * pt_
589
590 w_ = momz_(p,ke) * rdens_(p)
591 flux(p,1,momz_vid) = flux(p,1,dens_vid) * w_
592 flux(p,2,momz_vid) = flux(p,2,dens_vid) * w_
593 flux(p,3,momz_vid) = flux(p,3,dens_vid) * w_ + lmesh%Gsqrt(p,ke) * dpres_(p,ke) * rgsqrtv(p)
594 end do
595 do p=1, elem%Np
596 gsqrtdpres_ = lmesh%Gsqrt(p,ke) * dpres_(p,ke)
597
598 u_ = momx_(p,ke) * rdens_(p)
599 flux(p,1,momx_vid) = flux(p,1,dens_vid) * u_ + gsqrtdpres_
600 flux(p,2,momx_vid) = flux(p,2,dens_vid) * u_
601 flux(p,3,momx_vid) = flux(p,3,dens_vid) * u_ + gsqrtdpres_ * lmesh%GI3(p,ke,1)
602
603 v_ = momy_(p,ke) * rdens_(p)
604 flux(p,1,momy_vid) = flux(p,1,dens_vid) * v_
605 flux(p,2,momy_vid) = flux(p,2,dens_vid) * v_ + gsqrtdpres_
606 flux(p,3,momy_vid) = flux(p,3,dens_vid) * v_ + gsqrtdpres_ * lmesh%GI3(p,ke,2)
607 end do
608
609 call element3d_operation%Div_var5_2( &
610 flux, & ! (in)
611 dflux ) ! (out)
612
613 do p=1, elem%Np
614 e11 = lmesh%Escale(p,ke,1,1)
615 e22 = lmesh%Escale(p,ke,2,2)
616 e33 = lmesh%Escale(p,ke,3,3)
617
618 dens_dt(p,ke) = - ( &
619 e11 * dflux(p,1,dens_vid) &
620 + e22 * dflux(p,2,dens_vid) &
621 + e33 * dflux(p,3,dens_vid) &
622 ) * rgsqrt(p)
623
624 rhot_dt(p,ke) = - ( &
625 e11 * dflux(p,1,rhot_vid) &
626 + e22 * dflux(p,2,rhot_vid) &
627 + e33 * dflux(p,3,rhot_vid) &
628 ) * rgsqrt(p)
629 end do
630
631 call element3d_operation%VFilterPM1( ddens_(:,ke), & ! (in)
632 drho ) ! (out)
633
634 do p=1, elem%Np
635 e11 = lmesh%Escale(p,ke,1,1)
636 e22 = lmesh%Escale(p,ke,2,2)
637 e33 = lmesh%Escale(p,ke,3,3)
638
639 momz_dt(p,ke) = - ( &
640 e11 * dflux(p,1,momz_vid) &
641 + e22 * dflux(p,2,momz_vid) &
642 + e33 * dflux(p,3,momz_vid) &
643 ) * rgsqrt(p) &
644 - grav * drho(p)
645 end do
646
647 !--
648 do p=1, elem%Np
649 cor = coriolis(elem%IndexH2Dto3D(p),ke2d)
650 momx_dt(p,ke) = - dphyddx(p,ke) &
651 + cor * momy_(p,ke)
652 momy_dt(p,ke) = - dphyddy(p,ke) &
653 - cor * momx_(p,ke)
654 end do
655
656 do p=1, elem%Np
657 e11 = lmesh%Escale(p,ke,1,1)
658 e22 = lmesh%Escale(p,ke,2,2)
659 e33 = lmesh%Escale(p,ke,3,3)
660
661 momx_dt(p,ke) = momx_dt(p,ke) - ( &
662 e11 * dflux(p,1,momx_vid) &
663 + e22 * dflux(p,2,momx_vid) &
664 + e33 * dflux(p,3,momx_vid) &
665 ) * rgsqrt(p)
666
667 momy_dt(p,ke) = momy_dt(p,ke) - ( &
668 e11 * dflux(p,1,momy_vid) &
669 + e22 * dflux(p,2,momy_vid) &
670 + e33 * dflux(p,3,momy_vid) &
671 ) * rgsqrt(p)
672 end do
673 end do
674
675 call prof_rapend('cal_dyn_tend_interior', 3)
676
677 return
678 end subroutine atm_dyn_dgm_nonhydro3d_rhot_heve_cal_tend_cco
679
680end module scale_atm_dyn_dgm_nonhydro3d_rhot_heve
scale_sparsemat::sparsemat_matmul
Definition scale_sparsemat.F90:57
scale_atm_dyn_dgm_nonhydro3d_common
module FElib / Fluid dyn solver / Atmosphere / Nonhydrostatic model / Common
Definition scale_atm_dyn_dgm_nonhydro3d_common.F90:11
scale_atm_dyn_dgm_nonhydro3d_common::prgvar_momy_id
integer, parameter, public prgvar_momy_id
Definition scale_atm_dyn_dgm_nonhydro3d_common.F90:74
scale_atm_dyn_dgm_nonhydro3d_common::prgvar_ddens_id
integer, parameter, public prgvar_ddens_id
Definition scale_atm_dyn_dgm_nonhydro3d_common.F90:68
scale_atm_dyn_dgm_nonhydro3d_common::prgvar_momz_id
integer, parameter, public prgvar_momz_id
Definition scale_atm_dyn_dgm_nonhydro3d_common.F90:72
scale_atm_dyn_dgm_nonhydro3d_common::atm_dyn_dgm_nonhydro3d_common_init
subroutine, public atm_dyn_dgm_nonhydro3d_common_init(mesh)
Initialize a common module for atmospheric nonhydrostatic dynamical core.
Definition scale_atm_dyn_dgm_nonhydro3d_common.F90:113
scale_atm_dyn_dgm_nonhydro3d_common::prgvar_momx_id
integer, parameter, public prgvar_momx_id
Definition scale_atm_dyn_dgm_nonhydro3d_common.F90:73
scale_atm_dyn_dgm_nonhydro3d_common::prgvar_drhot_id
integer, parameter, public prgvar_drhot_id
Definition scale_atm_dyn_dgm_nonhydro3d_common.F90:70
scale_atm_dyn_dgm_nonhydro3d_common::prgvar_num
integer, parameter, public prgvar_num
Definition scale_atm_dyn_dgm_nonhydro3d_common.F90:77
scale_atm_dyn_dgm_nonhydro3d_common::intrpmat_vpordm1
real(rp), dimension(:,:), allocatable, public intrpmat_vpordm1
Definition scale_atm_dyn_dgm_nonhydro3d_common.F90:102
scale_atm_dyn_dgm_nonhydro3d_common::atm_dyn_dgm_nonhydro3d_common_final
subroutine, public atm_dyn_dgm_nonhydro3d_common_final()
Finalize a common module for atmospheric nonhydrostatic dynamical core.
Definition scale_atm_dyn_dgm_nonhydro3d_common.F90:129
scale_atm_dyn_dgm_nonhydro3d_rhot_heve_numflux
module FElib / Fluid dyn solver / Atmosphere / Nonhydrostatic model / HEVE / Numflux
Definition scale_atm_dyn_dgm_nonhydro3d_rhot_heve_numflux.F90:11
scale_atm_dyn_dgm_nonhydro3d_rhot_heve_numflux::atm_dyn_dgm_nonhydro3d_rhot_heve_numflux_get_generalvc
subroutine, public atm_dyn_dgm_nonhydro3d_rhot_heve_numflux_get_generalvc(del_flux, ddens_, momx_, momy_, momz_, drhot_, dpres, dens_hyd, pres_hyd, therm_hyd, rtot, cvtot, cptot, gsqrt, g13, g23, nx, ny, nz, vmapm, vmapp, lmesh, elem, lmesh2d, elem2d)
Definition scale_atm_dyn_dgm_nonhydro3d_rhot_heve_numflux.F90:952
scale_atm_dyn_dgm_nonhydro3d_rhot_heve_numflux::atm_dyn_dgm_nonhydro3d_rhot_heve_numflux_get_generalvc_asis
subroutine, public atm_dyn_dgm_nonhydro3d_rhot_heve_numflux_get_generalvc_asis(del_flux, del_flux_hyd, ddens_, momx_, momy_, momz_, drhot_, dpres_, dens_hyd, pres_hyd, rtot, cvtot, cptot, gsqrt, g13, g23, nx, ny, nz, vmapm, vmapp, lmesh, elem, lmesh2d, elem2d)
Definition scale_atm_dyn_dgm_nonhydro3d_rhot_heve_numflux.F90:78
scale_atm_dyn_dgm_nonhydro3d_rhot_heve
module FElib / Fluid dyn solver / Atmosphere / Regional nonhydrostatic model / HEVE
Definition scale_atm_dyn_dgm_nonhydro3d_rhot_heve.F90:13
scale_atm_dyn_dgm_nonhydro3d_rhot_heve::atm_dyn_dgm_nonhydro3d_rhot_heve_init
subroutine, public atm_dyn_dgm_nonhydro3d_rhot_heve_init(mesh)
Definition scale_atm_dyn_dgm_nonhydro3d_rhot_heve.F90:75
scale_atm_dyn_dgm_nonhydro3d_rhot_heve::atm_dyn_dgm_nonhydro3d_rhot_heve_final
subroutine, public atm_dyn_dgm_nonhydro3d_rhot_heve_final()
Definition scale_atm_dyn_dgm_nonhydro3d_rhot_heve.F90:85
scale_atm_dyn_dgm_nonhydro3d_rhot_heve::atm_dyn_dgm_nonhydro3d_rhot_heve_cal_tend_cco
subroutine, public atm_dyn_dgm_nonhydro3d_rhot_heve_cal_tend_cco(dens_dt, momx_dt, momy_dt, momz_dt, rhot_dt, ddens_, momx_, momy_, momz_, drhot_, dpres_, dens_hyd, pres_hyd, pres_hyd_ref, therm_hyd, coriolis, rtot, cvtot, cptot, dphyddx, dphyddy, element3d_operation, dx, dy, dz, sx, sy, sz, lift, lmesh, elem, lmesh2d, elem2d)
Definition scale_atm_dyn_dgm_nonhydro3d_rhot_heve.F90:500
scale_atm_dyn_dgm_nonhydro3d_rhot_heve::atm_dyn_dgm_nonhydro3d_rhot_heve_cal_tend
subroutine, public atm_dyn_dgm_nonhydro3d_rhot_heve_cal_tend(dens_dt, momx_dt, momy_dt, momz_dt, rhot_dt, ddens_, momx_, momy_, momz_, drhot_, dpres_, dens_hyd, pres_hyd, pres_hyd_ref, therm_hyd, coriolis, rtot, cvtot, cptot, dphyddx, dphyddy, element3d_operation, dx, dy, dz, sx, sy, sz, lift, lmesh, elem, lmesh2d, elem2d)
Definition scale_atm_dyn_dgm_nonhydro3d_rhot_heve.F90:299
scale_atm_dyn_dgm_nonhydro3d_rhot_heve::atm_dyn_dgm_nonhydro3d_rhot_heve_cal_tend_asis
subroutine, public atm_dyn_dgm_nonhydro3d_rhot_heve_cal_tend_asis(dens_dt, momx_dt, momy_dt, momz_dt, rhot_dt, ddens_, momx_, momy_, momz_, drhot_, dpres_, dens_hyd, pres_hyd, pres_hyd_ref, therm_hyd, coriolis, rtot, cvtot, cptot, dphyddx, dphyddy, element3d_operation, dx, dy, dz, sx, sy, sz, lift, lmesh, elem, lmesh2d, elem2d)
Definition scale_atm_dyn_dgm_nonhydro3d_rhot_heve.F90:102
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_element_operation_base
module FElib / Element / Operation / Base
Definition scale_element_operation_base.F90:12
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_element_operation_base::elementoperationbase3d
Base type for elementwise operations.
Definition scale_element_operation_base.F90:33
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