FE-Project
Loading...
Searching...
No Matches
scale_atm_dyn_dgm_globalnonhydro3d_rhot_heve Module Reference

module FElib / Fluid dyn solver / Atmosphere / Global nonhydrostatic model / HEVE More...

Functions/Subroutines

subroutine, public atm_dyn_dgm_globalnonhydro3d_rhot_heve_init (mesh)
subroutine, public atm_dyn_dgm_globalnonhydro3d_rhot_heve_final ()
subroutine, public atm_dyn_dgm_globalnh3d_rhot_heve_cal_tend_shallow_atm_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)
subroutine, public atm_dyn_dgm_globalnonhydro3d_rhot_heve_cal_tend_shallow_atm (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)
subroutine, public atm_dyn_dgm_globalnonhydro3d_rhot_heve_cal_tend_deep_atm (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)

Detailed Description

module FElib / Fluid dyn solver / Atmosphere / Global nonhydrostatic model / HEVE

Description
HEVE DGM scheme for Global Atmospheric Dynamical process. The governing equations is a fully compressible nonhydrostatic equations, which consist of mass, momentum, and thermodynamics (density * potential temperature conservation) equations.
Author
Yuta Kawai, Xuanzhengbo Ren, and Team SCALE

Function/Subroutine Documentation

◆ atm_dyn_dgm_globalnonhydro3d_rhot_heve_init()

subroutine, public scale_atm_dyn_dgm_globalnonhydro3d_rhot_heve::atm_dyn_dgm_globalnonhydro3d_rhot_heve_init ( class(meshbase3d), intent(in) mesh)

Definition at line 76 of file scale_atm_dyn_dgm_globalnonhydro3d_rhot_heve.F90.

77 implicit none
78 class(MeshBase3D), intent(in) :: mesh
79 !--------------------------------------------
80
81 call atm_dyn_dgm_nonhydro3d_common_init( mesh )
82
83 return

References scale_atm_dyn_dgm_nonhydro3d_common::atm_dyn_dgm_nonhydro3d_common_init().

◆ atm_dyn_dgm_globalnonhydro3d_rhot_heve_final()

subroutine, public scale_atm_dyn_dgm_globalnonhydro3d_rhot_heve::atm_dyn_dgm_globalnonhydro3d_rhot_heve_final

Definition at line 87 of file scale_atm_dyn_dgm_globalnonhydro3d_rhot_heve.F90.

88 implicit none
89 !--------------------------------------------
90
91 call atm_dyn_dgm_nonhydro3d_common_final()
92 return

References scale_atm_dyn_dgm_nonhydro3d_common::atm_dyn_dgm_nonhydro3d_common_final().

◆ atm_dyn_dgm_globalnh3d_rhot_heve_cal_tend_shallow_atm_asis()

subroutine, public scale_atm_dyn_dgm_globalnonhydro3d_rhot_heve::atm_dyn_dgm_globalnh3d_rhot_heve_cal_tend_shallow_atm_asis ( real(rp), dimension(elem%np,lmesh%nea), intent(out) dens_dt,
real(rp), dimension(elem%np,lmesh%nea), intent(out) momx_dt,
real(rp), dimension(elem%np,lmesh%nea), intent(out) momy_dt,
real(rp), dimension(elem%np,lmesh%nea), intent(out) momz_dt,
real(rp), dimension(elem%np,lmesh%nea), intent(out) rhot_dt,
real(rp), dimension(elem%np,lmesh%nea), intent(in) ddens_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) momx_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) momy_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) momz_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) drhot_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) dpres_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) dens_hyd,
real(rp), dimension(elem%np,lmesh%nea), intent(in) pres_hyd,
real(rp), dimension(elem%np,lmesh%nea), intent(in) pres_hyd_ref,
real(rp), dimension(elem%np,lmesh%nea), intent(in) therm_hyd,
real(rp), dimension(elem2d%np,lmesh2d%nea), intent(in) coriolis,
real(rp), dimension (elem%np,lmesh%nea), intent(in) rtot,
real(rp), dimension(elem%np,lmesh%nea), intent(in) cvtot,
real(rp), dimension(elem%np,lmesh%nea), intent(in) cptot,
real(rp), dimension(elem%np,lmesh%nea), intent(in) dphyddx,
real(rp), dimension(elem%np,lmesh%nea), intent(in) dphyddy,
class(elementoperationbase3d), intent(in) element3d_operation,
type(sparsemat), intent(in) dx,
type(sparsemat), intent(in) dy,
type(sparsemat), intent(in) dz,
type(sparsemat), intent(in) sx,
type(sparsemat), intent(in) sy,
type(sparsemat), intent(in) sz,
type(sparsemat), intent(in) lift,
class(localmesh3d), intent(in) lmesh,
class(elementbase3d), intent(in) elem,
class(localmesh2d), intent(in) lmesh2d,
class(elementbase2d), intent(in) elem2d )

Definition at line 98 of file scale_atm_dyn_dgm_globalnonhydro3d_rhot_heve.F90.

105
108 use scale_const, only: &
109 ohm => const_ohm
110 implicit none
111
112 class(LocalMesh3D), intent(in) :: lmesh
113 class(ElementBase3D), intent(in) :: elem
114 class(LocalMesh2D), intent(in) :: lmesh2D
115 class(ElementBase2D), intent(in) :: elem2D
116 class(ElementOperationBase3D), intent(in) :: element3D_operation
117 type(SparseMat), intent(in) :: Dx, Dy, Dz, Sx, Sy, Sz, Lift
118 real(RP), intent(out) :: DENS_dt(elem%Np,lmesh%NeA)
119 real(RP), intent(out) :: MOMX_dt(elem%Np,lmesh%NeA)
120 real(RP), intent(out) :: MOMY_dt(elem%Np,lmesh%NeA)
121 real(RP), intent(out) :: MOMZ_dt(elem%Np,lmesh%NeA)
122 real(RP), intent(out) :: RHOT_dt(elem%Np,lmesh%NeA)
123 real(RP), intent(in) :: DDENS_(elem%Np,lmesh%NeA)
124 real(RP), intent(in) :: MOMX_(elem%Np,lmesh%NeA)
125 real(RP), intent(in) :: MOMY_(elem%Np,lmesh%NeA)
126 real(RP), intent(in) :: MOMZ_(elem%Np,lmesh%NeA)
127 real(RP), intent(in) :: DRHOT_(elem%Np,lmesh%NeA)
128 real(RP), intent(in) :: DPRES_(elem%Np,lmesh%NeA)
129 real(RP), intent(in) :: DENS_hyd(elem%Np,lmesh%NeA)
130 real(RP), intent(in) :: PRES_hyd(elem%Np,lmesh%NeA)
131 real(RP), intent(in) :: PRES_hyd_ref(elem%Np,lmesh%NeA)
132 real(RP), intent(in) :: THERM_hyd(elem%Np,lmesh%NeA)
133 real(RP), intent(in) :: CORIOLIS(elem2D%Np,lmesh2D%NeA)
134 real(RP), intent(in) :: Rtot (elem%Np,lmesh%NeA)
135 real(RP), intent(in) :: CVtot(elem%Np,lmesh%NeA)
136 real(RP), intent(in) :: CPtot(elem%Np,lmesh%NeA)
137 real(RP), intent(in) :: DPhydDx(elem%Np,lmesh%NeA)
138 real(RP), intent(in) :: DPhydDy(elem%Np,lmesh%NeA)
139
140 real(RP) :: Fx(elem%Np), Fy(elem%Np), Fz(elem%Np), LiftDelFlx(elem%Np)
141 real(RP) :: DPRES_hyd(elem%Np), GradPhyd_x(elem%Np), GradPhyd_y(elem%Np)
142 real(RP) :: del_flux(elem%NfpTot,lmesh%Ne,PRGVAR_NUM)
143 real(RP) :: del_flux_hyd(elem%NfpTot,lmesh%Ne,2)
144 real(RP) :: RHOT_(elem%Np)
145 real(RP) :: rdens_(elem%Np), u_(elem%Np), v_(elem%Np), w_(elem%Np), wt_(elem%np), drho(elem%Np)
146
147 real(RP) :: G11(elem%Np), G12(elem%Np), G22(elem%Np)
148 real(RP) :: GsqrtV(elem%Np), RGsqrtV(elem%Np)
149 real(RP) :: X2D(elem2D%Np,lmesh2D%Ne), Y2D(elem2D%Np,lmesh2D%Ne)
150 real(RP) :: X(elem%Np), Y(elem%Np), twoOVdel2(elem%Np)
151 real(RP) :: CORI(elem%Np,2)
152 logical :: is_panel1to4
153 real(RP) :: s
154
155 integer :: ke, ke2d
156 integer :: p, p12, p3
157
158 real(RP) :: gamm, rgamm
159 real(RP) :: rP0
160 real(RP) :: RovP0, P0ovR
161 !------------------------------------------------------------------------
162
163 call prof_rapstart('cal_dyn_tend_bndflux', 3)
164 call get_ebnd_flux( &
165 del_flux, del_flux_hyd, & ! (out)
166 ddens_, momx_, momy_, momz_, drhot_, dpres_, dens_hyd, pres_hyd, & ! (in)
167 rtot, cvtot, cptot, & ! (in)
168 lmesh%Gsqrt, lmesh%GIJ(:,:,1,1), lmesh%GIJ(:,:,1,2), lmesh%GIJ(:,:,2,2), & ! (in)
169 lmesh%GsqrtH, lmesh%gam, lmesh%GI3(:,:,1), lmesh%GI3(:,:,2), & ! (in)
170 lmesh%normal_fn(:,:,1), lmesh%normal_fn(:,:,2), lmesh%normal_fn(:,:,3), & ! (in)
171 lmesh%vmapM, lmesh%vmapP, elem%IndexH2Dto3D_bnd, & ! (in)
172 lmesh, elem, lmesh2d, elem2d ) ! (in)
173 call prof_rapend('cal_dyn_tend_bndflux', 3)
174
175 !-----
176 call prof_rapstart('cal_dyn_tend_interior', 3)
177 gamm = cpdry / cvdry
178 rgamm = cvdry / cpdry
179 rp0 = 1.0_rp / pres00
180 rovp0 = rdry * rp0
181 p0ovr = pres00 / rdry
182
183 s = 1.0_rp
184 is_panel1to4 = .true.
185 if ( lmesh%panelID == 5 ) then
186 is_panel1to4 = .false.
187 else if ( lmesh%panelID == 6 ) then
188 is_panel1to4 = .false.
189 s = - 1.0_rp
190 end if
191
192 !$omp parallel private( &
193 !$omp RHOT_, rdens_, u_, v_, w_, wt_, &
194 !$omp Fx, Fy, Fz, LiftDelFlx, &
195 !$omp drho, DPRES_hyd, GradPhyd_x, GradPhyd_y, &
196 !$omp G11, G12, G22, GsqrtV, RGsqrtV, &
197 !$omp X, Y, twoOVdel2, &
198 !$omp CORI, ke, ke2D )
199
200 !$omp do
201 do ke2d = lmesh2d%NeS, lmesh2d%NeE
202 x2d(:,ke2d) = tan(lmesh2d%pos_en(:,ke2d,1))
203 y2d(:,ke2d) = tan(lmesh2d%pos_en(:,ke2d,2))
204 end do
205
206 !$omp do
207 do ke = lmesh%NeS, lmesh%NeE
208 !--
209 ke2d = lmesh%EMap3Dto2D(ke)
210 g11(:) = lmesh%GIJ(elem%IndexH2Dto3D,ke2d,1,1)
211 g12(:) = lmesh%GIJ(elem%IndexH2Dto3D,ke2d,1,2)
212 g22(:) = lmesh%GIJ(elem%IndexH2Dto3D,ke2d,2,2)
213 gsqrtv(:) = lmesh%Gsqrt(:,ke) / lmesh%GsqrtH(elem%IndexH2Dto3D,ke2d)
214 rgsqrtv(:) = 1.0_rp / gsqrtv(:)
215
216 !--
217 rhot_(:) = p0ovr * ( pres_hyd(:,ke) * rp0 )**rgamm + drhot_(:,ke)
218 ! DPRES_(:) = PRES00 * ( Rtot(:,ke) * rP0 * RHOT_(:) )**( CPtot(:,ke) / CVtot(:,ke) ) &
219 ! - PRES_hyd(:,ke)
220
221 rdens_(:) = 1.0_rp / ( ddens_(:,ke) + dens_hyd(:,ke) )
222 u_(:) = momx_(:,ke) * rdens_(:)
223 v_(:) = momy_(:,ke) * rdens_(:)
224 w_(:) = momz_(:,ke) * rdens_(:)
225 wt_(:) = w_(:) * rgsqrtv(:) + lmesh%GI3(:,ke,1) * u_(:) + lmesh%GI3(:,ke,2) * v_(:)
226
227 x(:) = x2d(elem%IndexH2Dto3D,ke2d)
228 y(:) = y2d(elem%IndexH2Dto3D,ke2d)
229 twoovdel2(:) = 2.0_rp / ( 1.0_rp + x(:)**2 + y(:)**2 )
230
231 cori(:,1) = s * ohm * twoovdel2(:) * ( - x(:) * y(:) * momx_(:,ke) + ( 1.0_rp + y(:)**2 ) * momy_(:,ke) )
232 cori(:,2) = s * ohm * twoovdel2(:) * ( - ( 1.0_rp + x(:)**2 ) * momx_(:,ke) + x(:) * y(:) * momy_(:,ke) )
233 if ( is_panel1to4 ) then
234 cori(:,1) = s * y(:) * cori(:,1)
235 cori(:,2) = s * y(:) * cori(:,2)
236 end if
237
238 drho(:) = matmul(intrpmat_vpordm1, ddens_(:,ke))
239
240 !-- Gradient hydrostatic pressure
241
242 dpres_hyd(:) = pres_hyd(:,ke) - pres_hyd_ref(:,ke)
243
244 call sparsemat_matmul(dx, gsqrtv(:) * dpres_hyd(:), fx)
245 call sparsemat_matmul(dz, gsqrtv(:) * lmesh%GI3(:,ke,1) * dpres_hyd(:), fz)
246 call sparsemat_matmul(lift, lmesh%Fscale(:,ke) * del_flux_hyd(:,ke,1), liftdelflx)
247 gradphyd_x(:) = lmesh%Escale(:,ke,1,1) * fx(:) &
248 + lmesh%Escale(:,ke,3,3) * fz(:) &
249 + liftdelflx(:)
250
251 call sparsemat_matmul(dy, gsqrtv(:) * dpres_hyd(:), fy)
252 call sparsemat_matmul(dz, gsqrtv(:) * lmesh%GI3(:,ke,2) * dpres_hyd(:), fz)
253 call sparsemat_matmul(lift, lmesh%Fscale(:,ke) * del_flux_hyd(:,ke,2), liftdelflx)
254 gradphyd_y(:) = lmesh%Escale(:,ke,2,2) * fy(:) &
255 + lmesh%Escale(:,ke,3,3) * fz(:) &
256 + liftdelflx(:)
257
258 !-- DENS
259 call sparsemat_matmul(dx, lmesh%Gsqrt(:,ke) * momx_(:,ke), fx)
260 call sparsemat_matmul(dy, lmesh%Gsqrt(:,ke) * momy_(:,ke), fy)
261 call sparsemat_matmul(dz, lmesh%Gsqrt(:,ke) * ( ddens_(:,ke) + dens_hyd(:,ke) ) * wt_(:), fz)
262 call sparsemat_matmul(lift, lmesh%Fscale(:,ke) * del_flux(:,ke,dens_vid), liftdelflx)
263
264 dens_dt(:,ke) = - ( &
265 lmesh%Escale(:,ke,1,1) * fx(:) &
266 + lmesh%Escale(:,ke,2,2) * fy(:) &
267 + lmesh%Escale(:,ke,3,3) * fz(:) &
268 + liftdelflx(:) ) / lmesh%Gsqrt(:,ke)
269
270 !-- MOMX
271 call sparsemat_matmul(dx, lmesh%Gsqrt(:,ke) * ( u_(:) * momx_(:,ke) + g11(:) * dpres_(:,ke) ), fx)
272 call sparsemat_matmul(dy, lmesh%Gsqrt(:,ke) * ( v_(:) * momx_(:,ke) + g12(:) * dpres_(:,ke) ), fy)
273 call sparsemat_matmul(dz, lmesh%Gsqrt(:,ke) * ( wt_(:) * momx_(:,ke) &
274 + ( lmesh%GI3(:,ke,1) * g11(:) + lmesh%GI3(:,ke,2) * g12(:) ) * dpres_(:,ke) ), fz)
275 call sparsemat_matmul(lift, lmesh%Fscale(:,ke) * del_flux(:,ke,momx_vid), liftdelflx)
276
277 momx_dt(:,ke) = &
278 - ( lmesh%Escale(:,ke,1,1) * fx(:) &
279 + lmesh%Escale(:,ke,2,2) * fy(:) &
280 + lmesh%Escale(:,ke,3,3) * fz(:) &
281 + liftdelflx(:) ) / lmesh%Gsqrt(:,ke) &
282 - twoovdel2(:) * y(:) * &
283 ( x(:) * y(:) * u_(:) - (1.0_rp + y(:)**2) * v_(:) ) * momx_(:,ke) &
284 - ( g11(:) * gradphyd_x(:) + g12(:) * gradphyd_y(:) ) * rgsqrtv(:) &
285 + cori(:,1)
286
287 !-- MOMY
288 call sparsemat_matmul(dx, lmesh%Gsqrt(:,ke) * ( u_(:) * momy_(:,ke) + g12(:) * dpres_(:,ke) ), fx)
289 call sparsemat_matmul(dy, lmesh%Gsqrt(:,ke) * ( v_(:) * momy_(:,ke) + g22(:) * dpres_(:,ke) ), fy)
290 call sparsemat_matmul(dz, lmesh%Gsqrt(:,ke) * ( wt_(:) * momy_(:,ke) &
291 + ( lmesh%GI3(:,ke,1) * g12(:) + lmesh%GI3(:,ke,2) * g22(:) ) * dpres_(:,ke) ), fz)
292 call sparsemat_matmul(lift, lmesh%Fscale(:,ke) * del_flux(:,ke,momy_vid), liftdelflx)
293
294 momy_dt(:,ke) = &
295 - ( lmesh%Escale(:,ke,1,1) * fx(:) &
296 + lmesh%Escale(:,ke,2,2) * fy(:) &
297 + lmesh%Escale(:,ke,3,3) * fz(:) &
298 + liftdelflx(:) ) / lmesh%Gsqrt(:,ke) &
299 - twoovdel2(:) * x(:) * &
300 ( - (1.0_rp + x(:)**2) * u_(:) + x(:) * y(:) * v_(:) ) * momy_(:,ke) &
301 - ( g12(:) * gradphyd_x(:) + g22(:) * gradphyd_y(:) ) * rgsqrtv(:) &
302 + cori(:,2)
303
304 !-- MOMZ
305 call sparsemat_matmul(dx, lmesh%Gsqrt(:,ke) * u_(:) * momz_(:,ke), fx)
306 call sparsemat_matmul(dy, lmesh%Gsqrt(:,ke) * v_(:) * momz_(:,ke), fy)
307 call sparsemat_matmul(dz, lmesh%Gsqrt(:,ke) * ( wt_(:) * momz_(:,ke) + rgsqrtv(:) * dpres_(:,ke) ), fz)
308 call sparsemat_matmul(lift, lmesh%Fscale(:,ke) * del_flux(:,ke,momz_vid), liftdelflx)
309
310 momz_dt(:,ke) = &
311 - ( lmesh%Escale(:,ke,1,1) * fx(:) &
312 + lmesh%Escale(:,ke,2,2) * fy(:) &
313 + lmesh%Escale(:,ke,3,3) * fz(:) &
314 + liftdelflx(:) ) / lmesh%Gsqrt(:,ke) &
315 - grav * drho(:)
316
317 !-- RHOT
318 call sparsemat_matmul(dx, lmesh%Gsqrt(:,ke) * u_(:) * rhot_(:), fx)
319 call sparsemat_matmul(dy, lmesh%Gsqrt(:,ke) * v_(:) * rhot_(:), fy)
320 call sparsemat_matmul(dz, lmesh%Gsqrt(:,ke) * wt_(:) * rhot_(:), fz)
321 call sparsemat_matmul(lift, lmesh%Fscale(:,ke) * del_flux(:,ke,rhot_vid), liftdelflx)
322
323 rhot_dt(:,ke) = &
324 - ( lmesh%Escale(:,ke,1,1) * fx(:) &
325 + lmesh%Escale(:,ke,2,2) * fy(:) &
326 + lmesh%Escale(:,ke,3,3) * fz(:) &
327 + liftdelflx(:) ) / lmesh%Gsqrt(:,ke)
328
329 end do
330 !$omp end do
331 !$omp end parallel
332 call prof_rapend('cal_dyn_tend_interior', 3)
333
334 return
module FElib / Fluid dyn solver / Atmosphere / Nonhydrostatic model / HEVE / Numflux
subroutine, public atm_dyn_dgm_nonhydro3d_rhot_heve_numflux_get_generalhvc_asis(del_flux, del_flux_hyd, ddens_, momx_, momy_, momz_, drhot_, dpres_, dens_hyd, pres_hyd, rtot, cvtot, cptot, gsqrt, g11, g12, g22, gsqrth, gam, g13, g23, nx, ny, nz, vmapm, vmapp, im2dto3d, lmesh, elem, lmesh2d, elem2d)

References scale_atm_dyn_dgm_nonhydro3d_rhot_heve_numflux::atm_dyn_dgm_nonhydro3d_rhot_heve_numflux_get_generalhvc_asis(), and scale_atm_dyn_dgm_nonhydro3d_common::intrpmat_vpordm1.

◆ atm_dyn_dgm_globalnonhydro3d_rhot_heve_cal_tend_shallow_atm()

subroutine, public scale_atm_dyn_dgm_globalnonhydro3d_rhot_heve::atm_dyn_dgm_globalnonhydro3d_rhot_heve_cal_tend_shallow_atm ( real(rp), dimension(elem%np,lmesh%nea), intent(out) dens_dt,
real(rp), dimension(elem%np,lmesh%nea), intent(out) momx_dt,
real(rp), dimension(elem%np,lmesh%nea), intent(out) momy_dt,
real(rp), dimension(elem%np,lmesh%nea), intent(out) momz_dt,
real(rp), dimension(elem%np,lmesh%nea), intent(out) rhot_dt,
real(rp), dimension(elem%np,lmesh%nea), intent(in) ddens_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) momx_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) momy_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) momz_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) drhot_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) dpres_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) dens_hyd,
real(rp), dimension(elem%np,lmesh%nea), intent(in) pres_hyd,
real(rp), dimension(elem%np,lmesh%nea), intent(in) pres_hyd_ref,
real(rp), dimension(elem%np,lmesh%nea), intent(in) therm_hyd,
real(rp), dimension(elem2d%np,lmesh2d%nea), intent(in) coriolis,
real(rp), dimension (elem%np,lmesh%nea), intent(in) rtot,
real(rp), dimension(elem%np,lmesh%nea), intent(in) cvtot,
real(rp), dimension(elem%np,lmesh%nea), intent(in) cptot,
real(rp), dimension(elem%np,lmesh%nea), intent(in) dphyddx,
real(rp), dimension(elem%np,lmesh%nea), intent(in) dphyddy,
class(elementoperationbase3d), intent(in) element3d_operation,
type(sparsemat), intent(in) dx,
type(sparsemat), intent(in) dy,
type(sparsemat), intent(in) dz,
type(sparsemat), intent(in) sx,
type(sparsemat), intent(in) sy,
type(sparsemat), intent(in) sz,
type(sparsemat), intent(in) lift,
class(localmesh3d), intent(in) lmesh,
class(elementbase3d), intent(in) elem,
class(localmesh2d), intent(in) lmesh2d,
class(elementbase2d), intent(in) elem2d )

Definition at line 338 of file scale_atm_dyn_dgm_globalnonhydro3d_rhot_heve.F90.

345
348 use scale_const, only: &
349 ohm => const_ohm
350 implicit none
351
352 class(LocalMesh3D), intent(in) :: lmesh
353 class(ElementBase3D), intent(in) :: elem
354 class(LocalMesh2D), intent(in) :: lmesh2D
355 class(ElementBase2D), intent(in) :: elem2D
356 class(ElementOperationBase3D), intent(in) :: element3D_operation
357 type(SparseMat), intent(in) :: Dx, Dy, Dz, Sx, Sy, Sz, Lift
358 real(RP), intent(out) :: DENS_dt(elem%Np,lmesh%NeA)
359 real(RP), intent(out) :: MOMX_dt(elem%Np,lmesh%NeA)
360 real(RP), intent(out) :: MOMY_dt(elem%Np,lmesh%NeA)
361 real(RP), intent(out) :: MOMZ_dt(elem%Np,lmesh%NeA)
362 real(RP), intent(out) :: RHOT_dt(elem%Np,lmesh%NeA)
363 real(RP), intent(in) :: DDENS_(elem%Np,lmesh%NeA)
364 real(RP), intent(in) :: MOMX_(elem%Np,lmesh%NeA)
365 real(RP), intent(in) :: MOMY_(elem%Np,lmesh%NeA)
366 real(RP), intent(in) :: MOMZ_(elem%Np,lmesh%NeA)
367 real(RP), intent(in) :: DRHOT_(elem%Np,lmesh%NeA)
368 real(RP), intent(in) :: DPRES_(elem%Np,lmesh%NeA)
369 real(RP), intent(in) :: DENS_hyd(elem%Np,lmesh%NeA)
370 real(RP), intent(in) :: PRES_hyd(elem%Np,lmesh%NeA)
371 real(RP), intent(in) :: PRES_hyd_ref(elem%Np,lmesh%NeA)
372 real(RP), intent(in) :: THERM_hyd(elem%Np,lmesh%NeA)
373 real(RP), intent(in) :: CORIOLIS(elem2D%Np,lmesh2D%NeA)
374 real(RP), intent(in) :: Rtot (elem%Np,lmesh%NeA)
375 real(RP), intent(in) :: CVtot(elem%Np,lmesh%NeA)
376 real(RP), intent(in) :: CPtot(elem%Np,lmesh%NeA)
377 real(RP), intent(in) :: DPhydDx(elem%Np,lmesh%NeA)
378 real(RP), intent(in) :: DPhydDy(elem%Np,lmesh%NeA)
379
380 real(RP) :: Flux(elem%Np,3,5), DFlux(elem%Np,4,5)
381 real(RP) :: del_flux(elem%NfpTot,PRGVAR_NUM,lmesh%Ne)
382 real(RP) :: u_, v_, w_, pt_
383 real(RP) :: drho(elem%Np)
384 real(RP) :: RDENS_(elem%Np), GsqrtV(elem%Np), RGsqrtV(elem%Np), RGsqrt(elem%Np)
385 real(RP) :: Gsqrt_, GsqrtDPRES_, E11, E22, E33
386
387 real(RP) :: G11(elem%Np), G12(elem%Np), G22(elem%Np)
388 real(RP) :: X2D(elem2D%Np,lmesh2D%Ne), Y2D(elem2D%Np,lmesh2D%Ne)
389 real(RP) :: X(elem%Np), Y(elem%Np), twoOVdel2(elem%Np)
390 real(RP) :: CORI(elem%Np,2)
391 logical :: is_panel1to4
392 real(RP) :: s
393
394 integer :: ke, ke2d
395 integer :: p
396
397 real(RP) :: gamm, rgamm
398 real(RP) :: rP0
399 real(RP) :: RovP0, P0ovR
400 !------------------------------------------------------------------------
401
402 call prof_rapstart('cal_dyn_tend_bndflux', 3)
403 call get_ebnd_flux( &
404 del_flux, & ! (out)
405 ddens_, momx_, momy_, momz_, drhot_, dpres_, & ! (in)
406 dens_hyd, pres_hyd, therm_hyd, rtot, cvtot, cptot, & ! (in)
407 lmesh%Gsqrt, lmesh%GIJ(:,:,1,1), lmesh%GIJ(:,:,1,2), lmesh%GIJ(:,:,2,2), & ! (in)
408 lmesh%GsqrtH, lmesh%gam, lmesh%GI3(:,:,1), lmesh%GI3(:,:,2), & ! (in)
409 lmesh%normal_fn(:,:,1), lmesh%normal_fn(:,:,2), lmesh%normal_fn(:,:,3), & ! (in)
410 lmesh%vmapM, lmesh%vmapP, elem%IndexH2Dto3D_bnd, & ! (in)
411 lmesh, elem, lmesh2d, elem2d ) ! (in)
412 call prof_rapend('cal_dyn_tend_bndflux', 3)
413
414 !-----
415 call prof_rapstart('cal_dyn_tend_interior', 3)
416 gamm = cpdry / cvdry
417 rgamm = cvdry / cpdry
418 rp0 = 1.0_rp / pres00
419 rovp0 = rdry * rp0
420 p0ovr = pres00 / rdry
421
422 s = 1.0_rp
423 is_panel1to4 = .true.
424 if ( lmesh%panelID == 5 ) then
425 is_panel1to4 = .false.
426 else if ( lmesh%panelID == 6 ) then
427 is_panel1to4 = .false.
428 s = - 1.0_rp
429 end if
430
431 !$omp parallel private( ke, ke2d, p, &
432 !$omp u_, v_, w_, pt_, &
433 !$omp drho, CORI, &
434 !$omp RDENS_, RGsqrt, GsqrtV, RGsqrtV, Gsqrt_, GsqrtDPRES_, &
435 !$omp G11, G12, G22, X, Y, twoOVdel2, &
436 !$omp E11, E22, E33, DFlux, Flux )
437
438 !$omp do
439 do ke2d = lmesh2d%NeS, lmesh2d%NeE
440 x2d(:,ke2d) = tan(lmesh2d%pos_en(:,ke2d,1))
441 y2d(:,ke2d) = tan(lmesh2d%pos_en(:,ke2d,2))
442 end do
443
444 !$omp do
445 do ke = lmesh%NeS, lmesh%NeE
446 !--
447 ke2d = lmesh%EMap3Dto2D(ke)
448
449 do p=1, elem%Np
450 g11(p) = lmesh%GIJ(elem%IndexH2Dto3D(p),ke2d,1,1)
451 g12(p) = lmesh%GIJ(elem%IndexH2Dto3D(p),ke2d,1,2)
452 g22(p) = lmesh%GIJ(elem%IndexH2Dto3D(p),ke2d,2,2)
453 end do
454 do p=1, elem%Np
455 gsqrtv(p) = lmesh%Gsqrt(p,ke) / lmesh%GsqrtH(elem%IndexH2Dto3D(p),ke2d)
456 rgsqrtv(p) = 1.0_rp / gsqrtv(p)
457 rgsqrt(p) = 1.0_rp / lmesh%Gsqrt(p,ke)
458 rdens_(p) = 1.0_rp / ( ddens_(p,ke) + dens_hyd(p,ke) )
459 end do
460
461 !--
462
463 do p=1, elem%Np
464 gsqrt_ = lmesh%Gsqrt(p,ke)
465 flux(p,1,dens_vid) = gsqrt_ * momx_(p,ke)
466 flux(p,2,dens_vid) = gsqrt_ * momy_(p,ke)
467 flux(p,3,dens_vid) = gsqrt_ * ( &
468 momz_(p,ke) * rgsqrtv(p) &
469 + lmesh%GI3(p,ke,1) * momx_(p,ke) &
470 + lmesh%GI3(p,ke,2) * momy_(p,ke) )
471 end do
472 do p=1, elem%Np
473 pt_ = ( therm_hyd(p,ke) + drhot_(p,ke) ) * rdens_(p)
474
475 flux(p,1,rhot_vid) = flux(p,1,dens_vid) * pt_
476 flux(p,2,rhot_vid) = flux(p,2,dens_vid) * pt_
477 flux(p,3,rhot_vid) = flux(p,3,dens_vid) * pt_
478
479 w_ = momz_(p,ke) * rdens_(p)
480 flux(p,1,momz_vid) = flux(p,1,dens_vid) * w_
481 flux(p,2,momz_vid) = flux(p,2,dens_vid) * w_
482 flux(p,3,momz_vid) = flux(p,3,dens_vid) * w_ + lmesh%Gsqrt(p,ke) * rgsqrtv(p) * dpres_(p,ke)
483 end do
484 do p=1, elem%Np
485 gsqrtdpres_ = lmesh%Gsqrt(p,ke) * dpres_(p,ke)
486
487 u_ = momx_(p,ke) * rdens_(p)
488 v_ = momy_(p,ke) * rdens_(p)
489
490 flux(p,1,momx_vid) = flux(p,1,dens_vid) * u_ + g11(p) * gsqrtdpres_
491 flux(p,2,momx_vid) = flux(p,2,dens_vid) * u_ + g12(p) * gsqrtdpres_
492 flux(p,3,momx_vid) = flux(p,3,dens_vid) * u_ + gsqrtdpres_ * ( g11(p) * lmesh%GI3(p,ke,1) + g12(p) * lmesh%GI3(p,ke,2) )
493
494 flux(p,1,momy_vid) = flux(p,1,dens_vid) * v_ + g12(p) * gsqrtdpres_
495 flux(p,2,momy_vid) = flux(p,2,dens_vid) * v_ + g22(p) * gsqrtdpres_
496 flux(p,3,momy_vid) = flux(p,3,dens_vid) * v_ + gsqrtdpres_ * ( g12(p) * lmesh%GI3(p,ke,1) + g22(p) * lmesh%GI3(p,ke,2) )
497 end do
498
499 call element3d_operation%Div_var5( &
500 flux, del_flux(:,:,ke), & ! (in)
501 dflux ) ! (out)
502
503
504 do p=1, elem%Np
505 e11 = lmesh%Escale(p,ke,1,1)
506 e22 = lmesh%Escale(p,ke,2,2)
507 e33 = lmesh%Escale(p,ke,3,3)
508
509 dens_dt(p,ke) = - ( &
510 e11 * dflux(p,1,dens_vid) &
511 + e22 * dflux(p,2,dens_vid) &
512 + e33 * dflux(p,3,dens_vid) &
513 + dflux(p,4,dens_vid) ) * rgsqrt(p)
514
515 rhot_dt(p,ke) = - ( &
516 e11 * dflux(p,1,rhot_vid) &
517 + e22 * dflux(p,2,rhot_vid) &
518 + e33 * dflux(p,3,rhot_vid) &
519 + dflux(p,4,rhot_vid) ) * rgsqrt(p)
520 end do
521
522 call element3d_operation%VFilterPM1( ddens_(:,ke), & ! (in)
523 drho ) ! (out)
524
525 do p=1, elem%Np
526 e11 = lmesh%Escale(p,ke,1,1)
527 e22 = lmesh%Escale(p,ke,2,2)
528 e33 = lmesh%Escale(p,ke,3,3)
529
530 momz_dt(p,ke) = - ( &
531 e11 * dflux(p,1,momz_vid) &
532 + e22 * dflux(p,2,momz_vid) &
533 + e33 * dflux(p,3,momz_vid) &
534 + dflux(p,4,momz_vid) ) * rgsqrt(p) &
535 - grav * drho(p)
536 end do
537
538 !--
539 do p=1, elem%Np
540 x(p) = x2d(elem%IndexH2Dto3D(p),ke2d)
541 y(p) = y2d(elem%IndexH2Dto3D(p),ke2d)
542 twoovdel2(p) = 2.0_rp / ( 1.0_rp + x(p)**2 + y(p)**2 )
543 end do
544
545 do p=1, elem%Np
546 cori(p,1) = s * ohm * twoovdel2(p) * ( - x(p) * y(p) * momx_(p,ke) + ( 1.0_rp + y(p)**2 ) * momy_(p,ke) )
547 cori(p,2) = s * ohm * twoovdel2(p) * ( - ( 1.0_rp + x(p)**2 ) * momx_(p,ke) + x(p) * y(p) * momy_(p,ke) )
548 end do
549 if ( is_panel1to4 ) then
550 do p=1, elem%Np
551 cori(p,1) = s * y(p) * cori(p,1)
552 cori(p,2) = s * y(p) * cori(p,2)
553 end do
554 end if
555
556 do p=1, elem%Np
557 u_ = momx_(p,ke) * rdens_(p)
558 v_ = momy_(p,ke) * rdens_(p)
559
560 momx_dt(p,ke) = - ( g11(p) * dphyddx(p,ke) + g12(p) * dphyddy(p,ke) ) &
561 - twoovdel2(p) * y(p) * &
562 ( x(p) * y(p) * u_ - (1.0_rp + y(p)**2) * v_ ) * momx_(p,ke) &
563 + cori(p,1)
564
565 momy_dt(p,ke) = - ( g12(p) * dphyddx(p,ke) + g22(p) * dphyddy(p,ke) ) &
566 - twoovdel2(p) * x(p) * &
567 ( - (1.0_rp + x(p)**2) * u_ + x(p) * y(p) * v_ ) * momy_(p,ke) &
568 + cori(p,2)
569 end do
570
571 do p=1, elem%Np
572 e11 = lmesh%Escale(p,ke,1,1)
573 e22 = lmesh%Escale(p,ke,2,2)
574 e33 = lmesh%Escale(p,ke,3,3)
575
576 momx_dt(p,ke) = momx_dt(p,ke) - ( &
577 e11 * dflux(p,1,momx_vid) &
578 + e22 * dflux(p,2,momx_vid) &
579 + e33 * dflux(p,3,momx_vid) &
580 + dflux(p,4,momx_vid) ) * rgsqrt(p)
581
582 momy_dt(p,ke) = momy_dt(p,ke) - ( &
583 e11 * dflux(p,1,momy_vid) &
584 + e22 * dflux(p,2,momy_vid) &
585 + e33 * dflux(p,3,momy_vid) &
586 + dflux(p,4,momy_vid) ) * rgsqrt(p)
587 end do
588 end do
589 !$omp end do
590 !$omp end parallel
591 call prof_rapend('cal_dyn_tend_interior', 3)
592
593 return
subroutine, public atm_dyn_dgm_nonhydro3d_rhot_heve_numflux_get_generalhvc(del_flux, ddens_, momx_, momy_, momz_, drhot_, dpres, dens_hyd, pres_hyd, therm_hyd, rtot, cvtot, cptot, gsqrt, g11, g12, g22, gsqrth, gam, g13, g23, nx, ny, nz, vmapm, vmapp, im2dto3d, lmesh, elem, lmesh2d, elem2d)

References scale_atm_dyn_dgm_nonhydro3d_rhot_heve_numflux::atm_dyn_dgm_nonhydro3d_rhot_heve_numflux_get_generalhvc().

◆ atm_dyn_dgm_globalnonhydro3d_rhot_heve_cal_tend_deep_atm()

subroutine, public scale_atm_dyn_dgm_globalnonhydro3d_rhot_heve::atm_dyn_dgm_globalnonhydro3d_rhot_heve_cal_tend_deep_atm ( real(rp), dimension(elem%np,lmesh%nea), intent(out) dens_dt,
real(rp), dimension(elem%np,lmesh%nea), intent(out) momx_dt,
real(rp), dimension(elem%np,lmesh%nea), intent(out) momy_dt,
real(rp), dimension(elem%np,lmesh%nea), intent(out) momz_dt,
real(rp), dimension(elem%np,lmesh%nea), intent(out) rhot_dt,
real(rp), dimension(elem%np,lmesh%nea), intent(in) ddens_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) momx_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) momy_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) momz_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) drhot_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) dpres_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) dens_hyd,
real(rp), dimension(elem%np,lmesh%nea), intent(in) pres_hyd,
real(rp), dimension(elem%np,lmesh%nea), intent(in) pres_hyd_ref,
real(rp), dimension(elem%np,lmesh%nea), intent(in) therm_hyd,
real(rp), dimension(elem2d%np,lmesh2d%nea), intent(in) coriolis,
real(rp), dimension (elem%np,lmesh%nea), intent(in) rtot,
real(rp), dimension(elem%np,lmesh%nea), intent(in) cvtot,
real(rp), dimension(elem%np,lmesh%nea), intent(in) cptot,
real(rp), dimension(elem%np,lmesh%nea), intent(in) dphyddx,
real(rp), dimension(elem%np,lmesh%nea), intent(in) dphyddy,
class(elementoperationbase3d), intent(in) element3d_operation,
type(sparsemat), intent(in) dx,
type(sparsemat), intent(in) dy,
type(sparsemat), intent(in) dz,
type(sparsemat), intent(in) sx,
type(sparsemat), intent(in) sy,
type(sparsemat), intent(in) sz,
type(sparsemat), intent(in) lift,
class(localmesh3d), intent(in) lmesh,
class(elementbase3d), intent(in) elem,
class(localmesh2d), intent(in) lmesh2d,
class(elementbase2d), intent(in) elem2d )

Definition at line 597 of file scale_atm_dyn_dgm_globalnonhydro3d_rhot_heve.F90.

604
607 use scale_const, only: &
608 ohm => const_ohm
609 implicit none
610
611 class(LocalMesh3D), intent(in) :: lmesh
612 class(ElementBase3D), intent(in) :: elem
613 class(LocalMesh2D), intent(in) :: lmesh2D
614 class(ElementBase2D), intent(in) :: elem2D
615 class(ElementOperationBase3D), intent(in) :: element3D_operation
616 type(SparseMat), intent(in) :: Dx, Dy, Dz, Sx, Sy, Sz, Lift
617 real(RP), intent(out) :: DENS_dt(elem%Np,lmesh%NeA)
618 real(RP), intent(out) :: MOMX_dt(elem%Np,lmesh%NeA)
619 real(RP), intent(out) :: MOMY_dt(elem%Np,lmesh%NeA)
620 real(RP), intent(out) :: MOMZ_dt(elem%Np,lmesh%NeA)
621 real(RP), intent(out) :: RHOT_dt(elem%Np,lmesh%NeA)
622 real(RP), intent(in) :: DDENS_(elem%Np,lmesh%NeA)
623 real(RP), intent(in) :: MOMX_(elem%Np,lmesh%NeA)
624 real(RP), intent(in) :: MOMY_(elem%Np,lmesh%NeA)
625 real(RP), intent(in) :: MOMZ_(elem%Np,lmesh%NeA)
626 real(RP), intent(in) :: DRHOT_(elem%Np,lmesh%NeA)
627 real(RP), intent(in) :: DPRES_(elem%Np,lmesh%NeA)
628 real(RP), intent(in) :: DENS_hyd(elem%Np,lmesh%NeA)
629 real(RP), intent(in) :: PRES_hyd(elem%Np,lmesh%NeA)
630 real(RP), intent(in) :: PRES_hyd_ref(elem%Np,lmesh%NeA)
631 real(RP), intent(in) :: THERM_hyd(elem%Np,lmesh%NeA)
632 real(RP), intent(in) :: CORIOLIS(elem2D%Np,lmesh2D%NeA)
633 real(RP), intent(in) :: Rtot (elem%Np,lmesh%NeA)
634 real(RP), intent(in) :: CVtot(elem%Np,lmesh%NeA)
635 real(RP), intent(in) :: CPtot(elem%Np,lmesh%NeA)
636 real(RP), intent(in) :: DPhydDx(elem%Np,lmesh%NeA)
637 real(RP), intent(in) :: DPhydDy(elem%Np,lmesh%NeA)
638
639 real(RP) :: Fx(elem%Np), Fy(elem%Np), Fz(elem%Np), LiftDelFlx(elem%Np)
640 real(RP) :: DPRES_hyd(elem%Np), GradPhyd_x(elem%Np), GradPhyd_y(elem%Np)
641 real(RP) :: del_flux(elem%NfpTot,lmesh%Ne,PRGVAR_NUM)
642 real(RP) :: del_flux_hyd(elem%NfpTot,lmesh%Ne,2)
643 real(RP) :: RHOT_(elem%Np)
644 real(RP) :: rdens_(elem%Np), u_(elem%Np), v_(elem%Np), w_(elem%Np), wt_(elem%np), drho(elem%Np)
645
646 real(RP) :: G11(elem%Np), G12(elem%Np), G22(elem%Np)
647 real(RP) :: GsqrtV(elem%Np), RGsqrtV(elem%Np), Rgam2(elem%Np)
648 real(RP) :: X2D(elem2D%Np,lmesh2D%Ne), Y2D(elem2D%Np,lmesh2D%Ne)
649 real(RP) :: X(elem%Np), Y(elem%Np), twoOVdel2(elem%Np)
650 real(RP) :: OM1(elem%Np), OM2(elem%Np), OM3(elem%Np), DEL(elem%Np), R(elem%Np)
651 logical :: is_panel1to4
652 real(RP) :: s
653
654 integer :: ke, ke2d
655 integer :: p
656
657 real(RP) :: rgamm
658 real(RP) :: rP0
659 real(RP) :: P0ovR
660 !------------------------------------------------------------------------
661
662 call prof_rapstart('cal_dyn_tend_bndflux', 3)
663 call get_ebnd_flux( &
664 del_flux, del_flux_hyd, & ! (out)
665 ddens_, momx_, momy_, momz_, drhot_, dpres_, dens_hyd, pres_hyd, & ! (in)
666 rtot, cvtot, cptot, & ! (in)
667 lmesh%Gsqrt, lmesh%GIJ(:,:,1,1), lmesh%GIJ(:,:,1,2), lmesh%GIJ(:,:,2,2), & ! (in)
668 lmesh%GsqrtH, lmesh%gam, lmesh%GI3(:,:,1), lmesh%GI3(:,:,2), & ! (in)
669 lmesh%normal_fn(:,:,1), lmesh%normal_fn(:,:,2), lmesh%normal_fn(:,:,3), & ! (in)
670 lmesh%vmapM, lmesh%vmapP, elem%IndexH2Dto3D_bnd, & ! (in)
671 lmesh, elem, lmesh2d, elem2d ) ! (in)
672 call prof_rapend('cal_dyn_tend_bndflux', 3)
673
674 !-----
675 call prof_rapstart('cal_dyn_tend_interior', 3)
676 rgamm = cvdry / cpdry
677 rp0 = 1.0_rp / pres00
678 p0ovr = pres00 / rdry
679
680 s = 2.0_rp * ohm
681 is_panel1to4 = .true.
682 if ( lmesh%panelID == 5 ) then
683 is_panel1to4 = .false.
684 else if ( lmesh%panelID == 6 ) then
685 is_panel1to4 = .false.
686 s = - s
687 end if
688
689 !$omp parallel private( &
690 !$omp RHOT_, rdens_, u_, v_, w_, wt_, &
691 !$omp Fx, Fy, Fz, LiftDelFlx, &
692 !$omp drho, DPRES_hyd, GradPhyd_x, GradPhyd_y, &
693 !$omp G11, G12, G22, Rgam2, GsqrtV, RGsqrtV, &
694 !$omp X, Y, twoOVdel2, &
695 !$omp OM1, OM2, OM3, DEL, R, ke, ke2D )
696
697 !$omp do
698 do ke2d = lmesh2d%NeS, lmesh2d%NeE
699 x2d(:,ke2d) = tan(lmesh2d%pos_en(:,ke2d,1))
700 y2d(:,ke2d) = tan(lmesh2d%pos_en(:,ke2d,2))
701 end do
702
703 !$omp do
704 do ke = lmesh%NeS, lmesh%NeE
705 !--
706 ke2d = lmesh%EMap3Dto2D(ke)
707 rgam2(:) = 1.0_rp / lmesh%gam(:,ke)**2
708 g11(:) = lmesh%GIJ(elem%IndexH2Dto3D,ke2d,1,1) * rgam2(:)
709 g12(:) = lmesh%GIJ(elem%IndexH2Dto3D,ke2d,1,2) * rgam2(:)
710 g22(:) = lmesh%GIJ(elem%IndexH2Dto3D,ke2d,2,2) * rgam2(:)
711 gsqrtv(:) = lmesh%Gsqrt(:,ke) * rgam2(:) / lmesh%GsqrtH(elem%IndexH2Dto3D,ke2d)
712 rgsqrtv(:) = 1.0_rp / gsqrtv(:)
713
714 !--
715 rhot_(:) = p0ovr * ( pres_hyd(:,ke) * rp0 )**rgamm + drhot_(:,ke)
716 ! DPRES_(:) = PRES00 * ( Rtot(:,ke) * rP0 * RHOT_(:) )**( CPtot(:,ke) / CVtot(:,ke) ) &
717 ! - PRES_hyd(:,ke)
718
719 rdens_(:) = 1.0_rp / ( ddens_(:,ke) + dens_hyd(:,ke) )
720 u_(:) = momx_(:,ke) * rdens_(:)
721 v_(:) = momy_(:,ke) * rdens_(:)
722 w_(:) = momz_(:,ke) * rdens_(:)
723 wt_(:) = w_(:) * rgsqrtv(:) + lmesh%GI3(:,ke,1) * u_(:) + lmesh%GI3(:,ke,2) * v_(:)
724
725 x(:) = x2d(elem%IndexH2Dto3D,ke2d)
726 y(:) = y2d(elem%IndexH2Dto3D,ke2d)
727 del(:) = sqrt( 1.0_rp + x(:)**2 + y(:)**2 )
728 twoovdel2(:) = 2.0_rp / ( 1.0_rp + x(:)**2 + y(:)**2 )
729
730 r(:) = rplanet * lmesh%gam(:,ke)
731
732 !- pnl=1~4: OM1: 0, OM2: s del / (r (1+Y^2)) , s OM3 : Y /del
733 !- pnl=5,6: OM1: - s X del/(r (1+X^2)), OM2: - s Y del/(r(1+Y^2)), OM3 : s/del
734 if ( is_panel1to4 ) then
735 om1(:) = 0.0_rp
736 om2(:) = s * del(:) / ( r(:) * ( 1.0_rp + y(:)**2 ) )
737 om3(:) = s * y(:) / del(:)
738 else
739 om1(:) = - s * x(:) * del(:) / ( r(:) * ( 1.0_rp + x(:)**2 ) )
740 om2(:) = - s * y(:) * del(:) / ( r(:) * ( 1.0_rp + y(:)**2 ) )
741 om3(:) = s / del(:)
742 end if
743
744 drho(:) = matmul(intrpmat_vpordm1, ddens_(:,ke))
745
746 !-- Gradient hydrostatic pressure
747
748 dpres_hyd(:) = pres_hyd(:,ke) - pres_hyd_ref(:,ke)
749
750 call sparsemat_matmul(dx, gsqrtv(:) * dpres_hyd(:), fx)
751 call sparsemat_matmul(dz, gsqrtv(:) * lmesh%GI3(:,ke,1) * dpres_hyd(:), fz)
752 call sparsemat_matmul(lift, lmesh%Fscale(:,ke) * del_flux_hyd(:,ke,1), liftdelflx)
753 gradphyd_x(:) = lmesh%Escale(:,ke,1,1) * fx(:) &
754 + lmesh%Escale(:,ke,3,3) * fz(:) &
755 + liftdelflx(:)
756
757 call sparsemat_matmul(dy, gsqrtv(:) * dpres_hyd(:), fy)
758 call sparsemat_matmul(dz, gsqrtv(:) * lmesh%GI3(:,ke,2) * dpres_hyd(:), fz)
759 call sparsemat_matmul(lift, lmesh%Fscale(:,ke) * del_flux_hyd(:,ke,2), liftdelflx)
760 gradphyd_y(:) = lmesh%Escale(:,ke,2,2) * fy(:) &
761 + lmesh%Escale(:,ke,3,3) * fz(:) &
762 + liftdelflx(:)
763
764 !-- DENS
765 call sparsemat_matmul(dx, lmesh%Gsqrt(:,ke) * momx_(:,ke), fx)
766 call sparsemat_matmul(dy, lmesh%Gsqrt(:,ke) * momy_(:,ke), fy)
767 call sparsemat_matmul(dz, lmesh%Gsqrt(:,ke) * ( ddens_(:,ke) + dens_hyd(:,ke) ) * wt_(:), fz)
768 call sparsemat_matmul(lift, lmesh%Fscale(:,ke) * del_flux(:,ke,dens_vid), liftdelflx)
769
770 dens_dt(:,ke) = - ( &
771 lmesh%Escale(:,ke,1,1) * fx(:) &
772 + lmesh%Escale(:,ke,2,2) * fy(:) &
773 + lmesh%Escale(:,ke,3,3) * fz(:) &
774 + liftdelflx(:) ) / lmesh%Gsqrt(:,ke)
775
776 !-- MOMX
777 call sparsemat_matmul(dx, lmesh%Gsqrt(:,ke) * ( u_(:) * momx_(:,ke) + g11(:) * dpres_(:,ke) ), fx)
778 call sparsemat_matmul(dy, lmesh%Gsqrt(:,ke) * ( v_(:) * momx_(:,ke) + g12(:) * dpres_(:,ke) ), fy)
779 call sparsemat_matmul(dz, lmesh%Gsqrt(:,ke) * ( wt_(:) * momx_(:,ke) &
780 + ( lmesh%GI3(:,ke,1) * g11(:) + lmesh%GI3(:,ke,2) * g12(:) ) * dpres_(:,ke) ), fz)
781 call sparsemat_matmul(lift, lmesh%Fscale(:,ke) * del_flux(:,ke,momx_vid), liftdelflx)
782
783 momx_dt(:,ke) = &
784 - ( lmesh%Escale(:,ke,1,1) * fx(:) &
785 + lmesh%Escale(:,ke,2,2) * fy(:) &
786 + lmesh%Escale(:,ke,3,3) * fz(:) &
787 + liftdelflx(:) ) / lmesh%Gsqrt(:,ke) &
788 - twoovdel2(:) * y(:) * & !-> metric terms
789 ( x(:) * y(:) * u_(:) - ( 1.0_rp + y(:)**2 ) * v_(:) ) * momx_(:,ke) & !
790 - 2.0_rp * u_(:) * momz_(:,ke) / r(:) & !<-
791 - ( g11(:) * gradphyd_x(:) + g12(:) * gradphyd_y(:) ) * rgsqrtv(:) & !-> gradient hydrostatic pressure
792 - lmesh%Gsqrt(:,ke) * ( g11(:) * ( om2(:) * momz_(:,ke) - om3(:) * momy_(:,ke) ) & !-> Coriolis term
793 - g12(:) * ( om1(:) * momz_(:,ke) - om3(:) * momx_(:,ke) ) )
794
795 !-- MOMY
796 call sparsemat_matmul(dx, lmesh%Gsqrt(:,ke) * ( u_(:) * momy_(:,ke) + g12(:) * dpres_(:,ke) ), fx)
797 call sparsemat_matmul(dy, lmesh%Gsqrt(:,ke) * ( v_(:) * momy_(:,ke) + g22(:) * dpres_(:,ke) ), fy)
798 call sparsemat_matmul(dz, lmesh%Gsqrt(:,ke) * ( wt_(:) * momy_(:,ke) &
799 + ( lmesh%GI3(:,ke,1) * g12(:) + lmesh%GI3(:,ke,2) * g22(:) ) * dpres_(:,ke) ), fz)
800 call sparsemat_matmul(lift, lmesh%Fscale(:,ke) * del_flux(:,ke,momy_vid), liftdelflx)
801
802 momy_dt(:,ke) = &
803 - ( lmesh%Escale(:,ke,1,1) * fx(:) &
804 + lmesh%Escale(:,ke,2,2) * fy(:) &
805 + lmesh%Escale(:,ke,3,3) * fz(:) &
806 + liftdelflx(:) ) / lmesh%Gsqrt(:,ke) &
807 - twoovdel2(:) * x(:) * & !-> metric terms
808 ( - (1.0_rp + x(:)**2) * u_(:) + x(:) * y(:) * v_(:) ) * momy_(:,ke) & !
809 - 2.0_rp * v_(:) * momz_(:,ke) / r(:) & !<-
810 - ( g12(:) * gradphyd_x(:) + g22(:) * gradphyd_y(:) ) * rgsqrtv(:) & !-> gradient hydrostatic pressure
811 - lmesh%Gsqrt(:,ke) * ( g12(:) * ( om2(:) * momz_(:,ke) - om3(:) * momy_(:,ke) ) & !-> Coriolis term
812 - g22(:) * ( om1(:) * momz_(:,ke) - om3(:) * momx_(:,ke) ) )
813
814 !-- MOMZ
815 call sparsemat_matmul(dx, lmesh%Gsqrt(:,ke) * u_(:) * momz_(:,ke), fx)
816 call sparsemat_matmul(dy, lmesh%Gsqrt(:,ke) * v_(:) * momz_(:,ke), fy)
817 call sparsemat_matmul(dz, lmesh%Gsqrt(:,ke) * ( wt_(:) * momz_(:,ke) + rgsqrtv(:) * dpres_(:,ke) ), fz)
818 call sparsemat_matmul(lift, lmesh%Fscale(:,ke) * del_flux(:,ke,momz_vid), liftdelflx)
819
820 momz_dt(:,ke) = &
821 - ( lmesh%Escale(:,ke,1,1) * fx(:) &
822 + lmesh%Escale(:,ke,2,2) * fy(:) &
823 + lmesh%Escale(:,ke,3,3) * fz(:) &
824 + liftdelflx(:) ) / lmesh%Gsqrt(:,ke) &
825 - 0.25_rp * r(:) * twoovdel2(:)**2 * ( 1.0_rp * x(:)**2 ) * ( 1.0_rp * y(:)**2 ) & !-> metric terms
826 * ( - ( 1.0_rp + x(:)**2 ) * momx_(:,ke) * u_(:) & !
827 + 2.0_rp * x(:) * y(:) * momx_(:,ke) * v_(:) & !
828 - ( 1.0_rp + y(:)**2 ) * momy_(:,ke) * v_(:) ) & !<-
829 + 2.0_rp * dpres_(:,ke) / r(:) & !-> metric term with gradient of pressure deviaition
830 - lmesh%Gsqrt(:,ke) * ( om1(:) * momy_(:,ke) - om2(:) * momx_(:,ke) ) & !-> Coriolis term
831 - grav * rgam2(:) * drho(:) !-> buoyancy term
832
833 !-- RHOT
834 call sparsemat_matmul(dx, lmesh%Gsqrt(:,ke) * u_(:) * rhot_(:), fx)
835 call sparsemat_matmul(dy, lmesh%Gsqrt(:,ke) * v_(:) * rhot_(:), fy)
836 call sparsemat_matmul(dz, lmesh%Gsqrt(:,ke) * wt_(:) * rhot_(:), fz)
837 call sparsemat_matmul(lift, lmesh%Fscale(:,ke) * del_flux(:,ke,rhot_vid), liftdelflx)
838
839 rhot_dt(:,ke) = &
840 - ( lmesh%Escale(:,ke,1,1) * fx(:) &
841 + lmesh%Escale(:,ke,2,2) * fy(:) &
842 + lmesh%Escale(:,ke,3,3) * fz(:) &
843 + liftdelflx(:) ) / lmesh%Gsqrt(:,ke)
844 end do
845 !$omp end do
846 !$omp end parallel
847 call prof_rapend('cal_dyn_tend_interior', 3)
848
849 return

References scale_atm_dyn_dgm_nonhydro3d_rhot_heve_numflux::atm_dyn_dgm_nonhydro3d_rhot_heve_numflux_get_generalhvc_asis(), and scale_atm_dyn_dgm_nonhydro3d_common::intrpmat_vpordm1.