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

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

Functions/Subroutines

subroutine, public atm_dyn_dgm_globalnonhydro3d_etot_hevi_init (mesh)
subroutine, public atm_dyn_dgm_globalnonhydro3d_etot_hevi_final ()
subroutine, public atm_dyn_dgm_globalnonhydro3d_etot_hevi_cal_tend (dens_dt, momx_dt, momy_dt, momz_dt, entot_dt, ddens_, momx_, momy_, momz_, etot_, 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_etot_hevi_cal_vi (dens_dt, momx_dt, momy_dt, momz_dt, etot_dt, ddens_, momx_, momy_, momz_, etot_, dens_hyd, pres_hyd, ddens0_, momx0_, momy0_, momz0_, etot0_, rtot, cvtot, cptot, element3d_operation, dz, lift, impl_fac, dt, lmesh, elem, lmesh2d, elem2d)

Detailed Description

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

Description
HEVI DGM scheme for Global Atmospheric Dynamical process. The governing equations is a fully compressible nonhydrostatic equations, which consist of mass, momentum, and thermodynamics (total energy conservation) equations.
Author
Yuta Kawai, Team SCALE

Function/Subroutine Documentation

◆ atm_dyn_dgm_globalnonhydro3d_etot_hevi_init()

subroutine, public scale_atm_dyn_dgm_globalnonhydro3d_etot_hevi::atm_dyn_dgm_globalnonhydro3d_etot_hevi_init ( class(meshbase3d), intent(in) mesh)

Definition at line 75 of file scale_atm_dyn_dgm_globalnonhydro3d_etot_hevi.F90.

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

References scale_atm_dyn_dgm_nonhydro3d_common::atm_dyn_dgm_nonhydro3d_common_init().

◆ atm_dyn_dgm_globalnonhydro3d_etot_hevi_final()

subroutine, public scale_atm_dyn_dgm_globalnonhydro3d_etot_hevi::atm_dyn_dgm_globalnonhydro3d_etot_hevi_final

Definition at line 85 of file scale_atm_dyn_dgm_globalnonhydro3d_etot_hevi.F90.

86 implicit none
87 !--------------------------------------------
88
89 call atm_dyn_dgm_nonhydro3d_common_final()
90 return

References scale_atm_dyn_dgm_nonhydro3d_common::atm_dyn_dgm_nonhydro3d_common_final().

◆ atm_dyn_dgm_globalnonhydro3d_etot_hevi_cal_tend()

subroutine, public scale_atm_dyn_dgm_globalnonhydro3d_etot_hevi::atm_dyn_dgm_globalnonhydro3d_etot_hevi_cal_tend ( 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) entot_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) etot_,
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 96 of file scale_atm_dyn_dgm_globalnonhydro3d_etot_hevi.F90.

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

References scale_atm_dyn_dgm_nonhydro3d_etot_hevi_numflux::atm_dyn_dgm_nonhydro3d_etot_hevi_numflux_get_generalhvc().

◆ atm_dyn_dgm_globalnonhydro3d_etot_hevi_cal_vi()

subroutine, public scale_atm_dyn_dgm_globalnonhydro3d_etot_hevi::atm_dyn_dgm_globalnonhydro3d_etot_hevi_cal_vi ( 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) etot_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) etot_,
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) ddens0_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) momx0_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) momy0_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) momz0_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) etot0_,
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,
class(elementoperationbase3d), intent(in) element3d_operation,
class(sparsemat), intent(in) dz,
class(sparsemat), intent(in) lift,
real(rp), intent(in) impl_fac,
real(rp), intent(in) dt,
class(localmesh3d), intent(in) lmesh,
class(elementbase3d), intent(in) elem,
class(localmesh2d), intent(in) lmesh2d,
class(elementbase2d), intent(in) elem2d )

Definition at line 339 of file scale_atm_dyn_dgm_globalnonhydro3d_etot_hevi.F90.

347
353
354 implicit none
355
356 class(LocalMesh3D), intent(in) :: lmesh
357 class(ElementBase3D), intent(in) :: elem
358 class(LocalMesh2D), intent(in) :: lmesh2D
359 class(ElementBase2D), intent(in) :: elem2D
360 real(RP), intent(out) :: DENS_dt(elem%Np,lmesh%NeA)
361 real(RP), intent(out) :: MOMX_dt(elem%Np,lmesh%NeA)
362 real(RP), intent(out) :: MOMY_dt(elem%Np,lmesh%NeA)
363 real(RP), intent(out) :: MOMZ_dt(elem%Np,lmesh%NeA)
364 real(RP), intent(out) :: ETOT_dt(elem%Np,lmesh%NeA)
365 real(RP), intent(in) :: DDENS_(elem%Np,lmesh%NeA)
366 real(RP), intent(in) :: MOMX_(elem%Np,lmesh%NeA)
367 real(RP), intent(in) :: MOMY_(elem%Np,lmesh%NeA)
368 real(RP), intent(in) :: MOMZ_(elem%Np,lmesh%NeA)
369 real(RP), intent(in) :: ETOT_(elem%Np,lmesh%NeA)
370 real(RP), intent(in) :: DENS_hyd(elem%Np,lmesh%NeA)
371 real(RP), intent(in) :: PRES_hyd(elem%Np,lmesh%NeA)
372 real(RP), intent(in) :: DDENS0_(elem%Np,lmesh%NeA)
373 real(RP), intent(in) :: MOMX0_(elem%Np,lmesh%NeA)
374 real(RP), intent(in) :: MOMY0_(elem%Np,lmesh%NeA)
375 real(RP), intent(in) :: MOMZ0_(elem%Np,lmesh%NeA)
376 real(RP), intent(in) :: ETOT0_(elem%Np,lmesh%NeA)
377 real(RP), intent(in) :: Rtot(elem%Np,lmesh%NeA)
378 real(RP), intent(in) :: CVtot(elem%Np,lmesh%NeA)
379 real(RP), intent(in) :: CPtot(elem%Np,lmesh%NeA)
380 class(ElementOperationBase3D), intent(in) :: element3D_operation
381 class(SparseMat), intent(in) :: Dz, Lift
382 real(RP), intent(in) :: impl_fac
383 real(RP), intent(in) :: dt
384
385 real(RP) :: PROG_VARS (elem%Np,lmesh%NeZ,PRGVAR_NUM,lmesh%NeX*lmesh%NeY)
386 real(RP) :: DPRES (elem%Np,lmesh%NeZ,lmesh%NeX*lmesh%NeY)
387 real(RP) :: PROG_VARS0(elem%Np,lmesh%NeZ,PRGVAR_NUM,lmesh%NeX*lmesh%NeY)
388 real(RP) :: DPRES0 (elem%Np,lmesh%NeZ,lmesh%NeX*lmesh%NeY)
389 real(RP) :: b1D(3,elem%Nnode_v,lmesh%NeZ,elem%Nnode_h1D**2,lmesh%NeX*lmesh%NeY)
390 real(RP) :: GeoPot (elem%Np,lmesh%NeZ,lmesh%NeX*lmesh%NeY)
391 real(RP) :: KinHovDENS(elem%Np,lmesh%NeZ,lmesh%NeX*lmesh%NeY)
392 integer :: ipiv(elem%Nnode_v*3*lmesh%NeZ,elem%Nnode_h1D**2)
393 real(RP) :: b1D_uv(elem%Nnode_v,lmesh%NeZ,2,elem%Nnode_h1D**2,lmesh%NeX*lmesh%NeY)
394 integer :: ipiv_uv(elem%Nnode_v*1*lmesh%NeZ,elem%Nnode_h1D**2)
395 real(RP) :: alph(elem%NfpTot,lmesh%NeZ,lmesh%NeX*lmesh%NeY)
396 real(RP) :: Rtot_z(elem%Np,lmesh%NeZ,lmesh%NeX*lmesh%NeY)
397 real(RP) :: CPtot_ov_CVtot(elem%Np,lmesh%NeZ,lmesh%NeX*lmesh%NeY)
398 real(RP) :: DENS_hyd_z(elem%Np,lmesh%NeZ,lmesh%NeX*lmesh%NeY)
399 real(RP) :: PRES_hyd_z(elem%Np,lmesh%NeZ,lmesh%NeX*lmesh%NeY)
400 real(RP) :: GnnM_z(elem%Np,lmesh%NeZ,lmesh%NeX*lmesh%NeY)
401 real(RP) :: G13_z(elem%Np,lmesh%NeZ,lmesh%NeX*lmesh%NeY)
402 real(RP) :: G23_z(elem%Np,lmesh%NeZ,lmesh%NeX*lmesh%NeY)
403 real(RP) :: GsqrtV_z(elem%Np,lmesh%NeZ,lmesh%NeX*lmesh%NeY)
404 real(RP) :: nz(elem%NfpTot,lmesh%NeZ,lmesh%NeX*lmesh%NeY)
405 integer :: vmapM(elem%NfpTot,lmesh%NeZ)
406 integer :: vmapP(elem%NfpTot,lmesh%NeZ)
407 integer :: ColMask(elem%Nnode_v)
408 integer :: ke_xy, ke_z, ke, ke2D, v
409 integer :: itr_nlin
410 integer :: kl, ku, nz_1D
411 integer :: kl_uv, ku_uv, nz_1D_uv
412 integer :: ij
413 logical :: is_converged
414
415 real(RP), allocatable :: PmatBnd(:,:,:)
416 real(RP), allocatable :: PmatBnd_uv(:,:,:)
417 integer :: info, info_uv
418
419 real(RP) :: DENS_(elem%Np)
420 !------------------------------------------------------------------------
421
422 call prof_rapstart( 'hevi_cal_vi_prep', 3)
423
424 nz_1d = elem%Nnode_v * 3 * lmesh%NeZ
425 kl = ( elem%Nnode_v + 1 ) * 3 - 1
426 ku = kl
427 nz_1d_uv = elem%Nnode_v * 1 * lmesh%NeZ
428 kl_uv = elem%Nnode_v
429 ku_uv = kl_uv
430 allocate( pmatbnd(2*kl+ku+1,nz_1d,elem%Nnode_h1D**2) )
431 allocate( pmatbnd_uv(2*kl_uv+ku_uv+1,nz_1d_uv,elem%Nnode_h1D**2) )
432
433 call lmesh%GetVmapZ1D( vmapm, vmapp ) ! (out)
434
435 !-
436
437 !$omp parallel private( ke_xy, ke_z, ke, ke2D, DENS_ )
438 !$omp do collapse(2)
439 do ke_xy=1, lmesh%NeX*lmesh%NeY
440 do ke_z=1, lmesh%NeZ
441 ke = ke_xy + (ke_z-1)*lmesh%NeX*lmesh%NeY
442 ke2d = lmesh%EMap3Dto2D(ke)
443
444 prog_vars(:,ke_z,dens_vid,ke_xy) = ddens0_(:,ke)
445 prog_vars(:,ke_z,momx_vid,ke_xy) = momx0_(:,ke)
446 prog_vars(:,ke_z,momy_vid,ke_xy) = momy0_(:,ke)
447 prog_vars(:,ke_z,momz_vid,ke_xy) = momz0_(:,ke)
448 prog_vars(:,ke_z,etot_vid,ke_xy) = etot0_(:,ke)
449
450 dens_hyd_z(:,ke_z,ke_xy) = dens_hyd(:,ke)
451 pres_hyd_z(:,ke_z,ke_xy) = pres_hyd(:,ke)
452 geopot(:,ke_z,ke_xy) = grav * lmesh%zlev(:,ke)
453
454 dens_(:) = dens_hyd(:,ke) + ddens0_(:,ke)
455 kinhovdens(:,ke_z,ke_xy) = 0.5_rp * ( &
456 momx0_(:,ke) * ( lmesh%G_ij(elem%IndexH2Dto3D,ke2d,1,1) * momx0_(:,ke) + lmesh%G_ij(elem%IndexH2Dto3D,ke2d,2,1) * momy0_(:,ke) ) &
457 + momy0_(:,ke) * ( lmesh%G_ij(elem%IndexH2Dto3D,ke2d,2,1) * momx0_(:,ke) + lmesh%G_ij(elem%IndexH2Dto3D,ke2d,2,2) * momy0_(:,ke) ) &
458 ) / dens_(:)**2
459
460 rtot_z(:,ke_z,ke_xy) = rtot(:,ke)
461 cptot_ov_cvtot(:,ke_z,ke_xy) = cptot(:,ke) / cvtot(:,ke)
462
463 dpres(:,ke_z,ke_xy) = &
464 ( cptot_ov_cvtot(:,ke_z,ke_xy) - 1.0_rp ) &
465 * ( etot0_(:,ke) - ( dens_(:) * ( kinhovdens(:,ke_z,ke_xy) + geopot(:,ke_z,ke_xy) ) + 0.5_rp * momz0_(:,ke)**2 / dens_(:) ) ) &
466 - pres_hyd(:,ke)
467
468 nz(:,ke_z,ke_xy) = lmesh%normal_fn(:,ke,3)
469 g13_z(:,ke_z,ke_xy) = lmesh%GI3(:,ke,1)
470 g23_z(:,ke_z,ke_xy) = lmesh%GI3(:,ke,2)
471 gsqrtv_z(:,ke_z,ke_xy) = lmesh%Gsqrt(:,ke) / lmesh%GsqrtH(elem%IndexH2Dto3D,ke2d)
472
473 gnnm_z(:,ke_z,ke_xy) = ( &
474 1.0_rp / gsqrtv_z(:,ke_z,ke_xy)**2 &
475 + g13_z(:,ke_z,ke_xy) * ( lmesh%GIJ(elem%IndexH2Dto3D,ke2d,1,1) * g13_z(:,ke_z,ke_xy) &
476 + lmesh%GIJ(elem%IndexH2Dto3D,ke2d,1,2) * g23_z(:,ke_z,ke_xy) ) &
477 + g23_z(:,ke_z,ke_xy) * ( lmesh%GIJ(elem%IndexH2Dto3D,ke2d,1,2) * g13_z(:,ke_z,ke_xy) &
478 + lmesh%GIJ(elem%IndexH2Dto3D,ke2d,2,2) * g23_z(:,ke_z,ke_xy) ) )
479 end do
480 end do
481 !$omp end do
482 !$omp workshare
483 prog_vars0(:,:,:,:) = prog_vars(:,:,:,:)
484 dpres0(:,:,:) = dpres(:,:,:)
485 !$omp end workshare
486 !$omp end parallel
487
488 call prof_rapend( 'hevi_cal_vi_prep', 3)
489
490 !--
491
492 if ( abs(impl_fac) > 0.0_rp ) then
493 call prof_rapstart( 'hevi_cal_vi_itr', 3)
494
495 ! G = (q^n+1 - q^n*) + impl_fac * A(q^n+1) = 0
496 ! dG/dq^n+1 del[q] = - G(q^n*)
497 do itr_nlin = 1, 1
498 call prof_rapstart( 'hevi_cal_vi_ax', 3)
499
500 call vi_eval_ax_uv( &
501 momx_dt(:,:), momy_dt(:,:), alph(:,:,:), & ! (out, dummy)
502 prog_vars, dpres, prog_vars0, dpres0, & ! (in)
503 ddens_, momx_, momy_, momz_, etot_, & ! (in)
504 dens_hyd_z, pres_hyd_z, & ! (in)
505 rtot_z, cptot_ov_cvtot, & ! (in)
506 dz, lift, intrpmat_vpordm1, & ! (in)
507 gnnm_z, g13_z, g23_z, gsqrtv_z, & ! (in)
508 impl_fac, dt, & ! (in)
509 lmesh, elem, nz, vmapm, vmapp, & ! (in)
510 b1d_uv(:,:,:,:,:) ) ! (out)
511
512 call prof_rapend( 'hevi_cal_vi_ax', 3)
513
514 do ke_xy=1, lmesh%NeX * lmesh%NeY
515 call prof_rapstart( 'hevi_cal_vi_matbnd', 3)
516 call vi_construct_matbnd_uv( pmatbnd_uv(:,:,:), & ! (out)
517 kl_uv, ku_uv, nz_1d_uv, & ! (in)
518 prog_vars(:,:,:,ke_xy), kinhovdens(:,:,ke_xy), & ! (in)
519 dens_hyd_z(:,:,ke_xy), pres_hyd_z(:,:,ke_xy), & ! (in)
520 g13_z(:,:,ke_xy), g23_z(:,:,ke_xy), gsqrtv_z(:,:,ke_xy), & ! (in)
521 alph(:,:,ke_xy), & ! (in)
522 rtot_z(:,:,ke_xy), cptot_ov_cvtot(:,:,ke_xy), & ! (in)
523 geopot(:,:,ke_xy), & ! (in)
524 dz, lift, intrpmat_vpordm1, & ! (in)
525 impl_fac, dt, & ! (in)
526 lmesh, elem, nz(:,:,ke_xy), vmapm, vmapp, ke_xy, 1 ) ! (in)
527 call prof_rapend( 'hevi_cal_vi_matbnd', 3)
528
529 call prof_rapstart( 'hevi_cal_vi_lin', 3)
530 !$omp parallel private(ij, v, ke_z, info, ColMask)
531 !$omp do
532 do ij=1, elem%Nnode_h1D**2
533 call dgbsv( nz_1d_uv, kl_uv, ku_uv, 2, pmatbnd_uv(:,:,ij), 2*kl_uv+ku_uv+1, ipiv_uv(:,ij), b1d_uv(:,:,:,ij,ke_xy), nz_1d_uv, info)
534
535 colmask(:) = elem%Colmask(:,ij)
536 do ke_z=1, lmesh%NeZ
537 prog_vars(colmask(:),ke_z,momx_vid,ke_xy) = prog_vars(colmask(:),ke_z,momx_vid,ke_xy) + b1d_uv(:,ke_z,1,ij,ke_xy)
538 prog_vars(colmask(:),ke_z,momy_vid,ke_xy) = prog_vars(colmask(:),ke_z,momy_vid,ke_xy) + b1d_uv(:,ke_z,2,ij,ke_xy)
539 end do
540 end do ! for ij
541 !$omp end do
542 !$omp end parallel
543 call prof_rapend( 'hevi_cal_vi_lin', 3)
544
545 end do ! for ke_xy
546
547 call prof_rapstart( 'hevi_cal_vi_ax', 3)
548 call vi_eval_ax( &
549 dens_dt(:,:), momz_dt(:,:), etot_dt(:,:), & ! (out, dummy)
550 alph(:,:,:), & ! (in)
551 prog_vars, dpres, prog_vars0, dpres0, & ! (in)
552 ddens_, momx_, momy_, momz_, etot_, & ! (in)
553 dens_hyd_z, pres_hyd_z, & ! (in)
554 rtot_z, cptot_ov_cvtot, & ! (in)
555 dz, lift, intrpmat_vpordm1, & ! (in)
556 gnnm_z, g13_z, g23_z, gsqrtv_z, & ! (in)
557 impl_fac, dt, & ! (in)
558 lmesh, elem, nz, vmapm, vmapp, & ! (in)
559 b1d(:,:,:,:,:) ) ! (out)
560 call prof_rapend( 'hevi_cal_vi_ax', 3)
561
562 do ke_xy=1, lmesh%NeX * lmesh%NeY
563 call prof_rapstart( 'hevi_cal_vi_matbnd', 3)
564 call vi_construct_matbnd( pmatbnd(:,:,:), & ! (out)
565 kl, ku, nz_1d, & ! (in)
566 prog_vars(:,:,:,ke_xy), kinhovdens(:,:,ke_xy), & ! (in)
567 dens_hyd_z(:,:,ke_xy), pres_hyd_z(:,:,ke_xy), & ! (in)
568 g13_z(:,:,ke_xy), g23_z(:,:,ke_xy), gsqrtv_z(:,:,ke_xy), & ! (in)
569 alph(:,:,ke_xy), & ! (in)
570 rtot_z(:,:,ke_xy), cptot_ov_cvtot(:,:,ke_xy), & ! (in)
571 geopot(:,:,ke_xy), & ! (in)
572 dz, lift, intrpmat_vpordm1, & ! (in)
573 impl_fac, dt, & ! (in)
574 lmesh, elem, nz(:,:,ke_xy), vmapm, vmapp, ke_xy, 1 ) ! (in)
575 call prof_rapend( 'hevi_cal_vi_matbnd', 3)
576
577 call prof_rapstart( 'hevi_cal_vi_lin', 3)
578 !$omp parallel private(ij, v, ke_z, info, ColMask)
579 !$omp do
580 do ij=1, elem%Nnode_h1D**2
581 call dgbsv( nz_1d, kl, ku, 1, pmatbnd(:,:,ij), 2*kl+ku+1, ipiv(:,ij), b1d(:,:,:,ij,ke_xy), nz_1d, info)
582
583 colmask(:) = elem%Colmask(:,ij)
584 do ke_z=1, lmesh%NeZ
585 prog_vars(colmask(:),ke_z,dens_vid,ke_xy) = prog_vars(colmask(:),ke_z,dens_vid,ke_xy) + b1d(1,:,ke_z,ij,ke_xy)
586 prog_vars(colmask(:),ke_z,momz_vid,ke_xy) = prog_vars(colmask(:),ke_z,momz_vid,ke_xy) + b1d(2,:,ke_z,ij,ke_xy)
587 prog_vars(colmask(:),ke_z,etot_vid,ke_xy) = prog_vars(colmask(:),ke_z,etot_vid,ke_xy) + b1d(3,:,ke_z,ij,ke_xy)
588 end do
589 end do ! for ij
590 !$omp end do
591 !$omp do
592 do ke_z=1, lmesh%NeZ
593 dens_(:) = dens_hyd_z(:,ke_z,ke_xy) + prog_vars(:,ke_z,dens_vid,ke_xy)
594 dpres(:,ke_z,ke_xy) = &
595 ( cptot_ov_cvtot(:,ke_z,ke_xy) - 1.0_rp ) &
596 * ( prog_vars(:,ke_z,etot_vid,ke_xy) &
597 - ( dens_(:) * ( kinhovdens(:,ke_z,ke_xy) + geopot(:,ke_z,ke_xy) ) + 0.5_rp * prog_vars(:,ke_z,momz_vid,ke_xy)**2 / dens_(:) ) ) &
598 - pres_hyd_z(:,ke_z,ke_xy)
599 end do
600 !$omp end parallel
601 call prof_rapend( 'hevi_cal_vi_lin', 3)
602
603 end do ! for ke_xy
604 end do ! itr nlin
605
606 call prof_rapend( 'hevi_cal_vi_itr', 3)
607 end if
608
609 call prof_rapstart( 'hevi_cal_vi_retrun_var', 3)
610 if ( abs(impl_fac) > 0.0_rp) then
611 !$omp parallel do collapse(2) private(ke_xy, ke_z, ke)
612 do ke_xy=1, lmesh%NeX * lmesh%NeY
613 do ke_z=1, lmesh%NeZ
614 ke = ke_xy + (ke_z-1)*lmesh%NeX*lmesh%NeY
615 dens_dt(:,ke) = ( prog_vars(:,ke_z,dens_vid,ke_xy) - ddens_(:,ke) ) / impl_fac
616 momx_dt(:,ke) = ( prog_vars(:,ke_z,momx_vid,ke_xy) - momx_(:,ke) ) / impl_fac
617 momy_dt(:,ke) = ( prog_vars(:,ke_z,momy_vid,ke_xy) - momy_(:,ke) ) / impl_fac
618 momz_dt(:,ke) = ( prog_vars(:,ke_z,momz_vid,ke_xy) - momz_(:,ke) ) / impl_fac
619 etot_dt(:,ke) = ( prog_vars(:,ke_z,etot_vid,ke_xy) - etot_(:,ke) ) / impl_fac
620 end do
621 end do
622 else
623 call vi_eval_ax_uv( &
624 momx_dt(:,:), momy_dt(:,:), alph(:,:,:), & ! (out)
625 prog_vars, dpres, prog_vars0, dpres0, & ! (in)
626 ddens_, momx_, momy_, momz_, etot_, & ! (in)
627 dens_hyd_z, pres_hyd_z, & ! (in)
628 rtot_z, cptot_ov_cvtot, & ! (in)
629 dz, lift, intrpmat_vpordm1, & ! (in)
630 gnnm_z, g13_z, g23_z, gsqrtv_z, & ! (in)
631 impl_fac, dt, & ! (in)
632 lmesh, elem, nz, vmapm, vmapp ) ! (in)
633
634 call vi_eval_ax( &
635 dens_dt(:,:), momz_dt(:,:), etot_dt(:,:), & ! (out)
636 alph(:,:,:), & ! (in)
637 prog_vars, dpres, prog_vars0, dpres0, & ! (in)
638 ddens_, momx_, momy_, momz_, etot_, & ! (in)
639 dens_hyd_z, pres_hyd_z, & ! (in)
640 rtot_z, cptot_ov_cvtot, & ! (in)
641 dz, lift, intrpmat_vpordm1, & ! (in)
642 gnnm_z, g13_z, g23_z, gsqrtv_z, & ! (in)
643 impl_fac, dt, & ! (in)
644 lmesh, elem, nz, vmapm, vmapp ) ! (in)
645 end if
646 call prof_rapend( 'hevi_cal_vi_retrun_var', 3)
647
648
649 return
module FElib / Fluid dyn solver / Atmosphere / Nonhydrostatic model / HEVI / Common
subroutine, public atm_dyn_dgm_nonhydro3d_etot_hevi_common_construct_matbnd(pmatbnd, kl, ku, nz_1d, prog_vars0, kinhovdens00, dens_hyd, pres_hyd, g13, g23, gsqrtv, alph, rtot, cptot_ov_cvtot, geopot, dz, lift, intrpmat_vpordm1, impl_fac, dt, lmesh, elem, nz, vmapm, vmapp, ke_x, ke_y)
subroutine, public atm_dyn_dgm_nonhydro3d_etot_hevi_common_eval_ax(dens_t, momz_t, etot_t, alph, prog_vars, dpres, prog_vars0, dpres0, ddens00, momx00, momy00, momz00, entot00, dens_hyd, pres_hyd, rtot, cptot_ov_cvtot, dz, lift, intrpmat_vpordm1, gnnm, g13, g23, gsqrtv, impl_fac, dt, lmesh, elem, nz, vmapm, vmapp, b1d_ij)
subroutine, public atm_dyn_dgm_nonhydro3d_etot_hevi_common_eval_ax_uv(momx_t, momy_t, alph, prog_vars, dpres, prog_vars0, dpres0, ddens00, momx00, momy00, momz00, entot00, dens_hyd, pres_hyd, rtot, cptot_ov_cvtot, dz, lift, intrpmat_vpordm1, gnnm, g13, g23, gsqrtv, impl_fac, dt, lmesh, elem, nz, vmapm, vmapp, b1d_ij_uv)
subroutine, public atm_dyn_dgm_nonhydro3d_etot_hevi_common_construct_matbnd_uv(pmatbnd_uv, kl_uv, ku_uv, nz_1d_uv, prog_vars0, kinhovdens00, dens_hyd, pres_hyd, g13, g23, gsqrtv, alph, rtot, cptot_ov_cvtot, geopot, dz, lift, intrpmat_vpordm1, impl_fac, dt, lmesh, elem, nz, vmapm, vmapp, ke_x, ke_y)

References scale_atm_dyn_dgm_nonhydro3d_etot_hevi_common::atm_dyn_dgm_nonhydro3d_etot_hevi_common_construct_matbnd(), scale_atm_dyn_dgm_nonhydro3d_etot_hevi_common::atm_dyn_dgm_nonhydro3d_etot_hevi_common_construct_matbnd_uv(), scale_atm_dyn_dgm_nonhydro3d_etot_hevi_common::atm_dyn_dgm_nonhydro3d_etot_hevi_common_eval_ax(), scale_atm_dyn_dgm_nonhydro3d_etot_hevi_common::atm_dyn_dgm_nonhydro3d_etot_hevi_common_eval_ax_uv(), and scale_atm_dyn_dgm_nonhydro3d_common::intrpmat_vpordm1.