FE-Project
Loading...
Searching...
No Matches
Functions/Subroutines
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, 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(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.

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

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