FE-Project
Loading...
Searching...
No Matches
Functions/Subroutines
scale_atm_dyn_dgm_trcadvect3d_heve Module Reference

module FElib / Fluid dyn solver / Atmosphere / Tracer advection More...

Functions/Subroutines

subroutine, public atm_dyn_dgm_trcadvect3d_heve_init (mesh, faceintmat)
subroutine, public atm_dyn_dgm_trcadvect3d_heve_final ()
subroutine, public atm_dyn_dgm_trcadvect3d_heve_cal_tend (qtrc_dt, qtrc_, momx_, momy_, momz_, alphdens_m, alphdens_p, fct_coef, rhoq_tp, element3d_operation, faceintmat, lmesh, elem, lmesh2d, elem2d)
subroutine, public atm_dyn_dgm_trcadvect3d_heve_calc_fct_coef (fct_coef, qtrc_, momx_, momy_, momz_, rhoq_tp_, alphdens_m, alphdens_p, dens_hyd, ddens_, ddens0_, rk_c_ssm1, dt, faceintmat, lmesh, elem, lmesh2d, elem2d, disable_limiter)
subroutine, public atm_dyn_dgm_trcadvect3d_tmar (qtrc_, dens_hyd, ddens_, lmesh, elem, lmesh2d, elem2d)
 Second Step of limiter in which nonlinear truncation and mass aware rescaling (TMAR)
subroutine, public atm_dyn_dgm_trcadvect3d_save_massflux (mflx_x_tavg, mflx_y_tavg, mflx_z_tavg, alph_dens_m, alph_dens_p, ddens, momx, momy, momz, dpres, dens_hyd, pres_hyd, rtot, cvtot, cptot, lmesh, elem, rkstage, tavg_weight_h, tavg_weight_v, is_hevi)
subroutine, public atm_dyn_dgm_trcadvect3d_heve_cal_alphdens_advtest (alph_dens_m, alph_dens_p, ddens_, momx_, momy_, momz_, dens_hyd, gsqrt, nx, ny, nz, vmapm, vmapp, lmesh, elem)

Detailed Description

module FElib / Fluid dyn solver / Atmosphere / Tracer advection

Description
HEVE DGM scheme for tracer advection.

To preserve nonnegativity, a limiter proposed by Light and Durran (2016, MWR) is used: we apply FCT for the lowest mode and the truncation and mass aware rescaling (TMAR).

Reference
  • Light and Durran 2016: Preserving Nonnegativity in Discontinuous Galerkin Approximations to Scalar Transport via Truncation and Mass Aware Rescaling (TMAR). Monthly Weather Review, 144(12), 4771–4786.
Author
Yuta Kawai, Team SCALE

Function/Subroutine Documentation

◆ atm_dyn_dgm_trcadvect3d_heve_init()

subroutine, public scale_atm_dyn_dgm_trcadvect3d_heve::atm_dyn_dgm_trcadvect3d_heve_init ( class(meshbase3d), intent(in), target mesh,
type(sparsemat), intent(inout) faceintmat )

Definition at line 75 of file scale_atm_dyn_dgm_trcadvect3d_heve.F90.

76 use scale_polynominal, only: polynominal_gengausslobattoptintweight
77 implicit none
78 class(MeshBase3D), intent(in), target :: mesh
79 type(SparseMat), intent(inout) :: FaceIntMat
80
81 class(LocalMesh3D), pointer :: lcmesh
82 class(ElementBase3D), pointer :: elem
83 real(RP), allocatable :: intWeight_lgl1DPts_h(:)
84 real(RP), allocatable :: intWeight_lgl1DPts_v(:)
85 real(RP), allocatable :: intWeight_h(:)
86 real(RP), allocatable :: intWeight_v(:)
87
88 integer :: f
89 integer :: i, j, k, l
90 integer :: is, ie
91
92 real(RP), allocatable :: IntWeight(:,:)
93 !-------------------------------------------------------------
94
95 lcmesh => mesh%lcmesh_list(1)
96 elem => lcmesh%refElem3D
97 allocate( intweight(elem%Nfaces,elem%NfpTot) )
98 intweight(:,:) = 0.0_rp
99
100 allocate( intweight_lgl1dpts_h(elem%Nnode_h1D) )
101 allocate( intweight_lgl1dpts_v(elem%Nnode_v) )
102 allocate( intweight_h(elem%Nnode_h1D*elem%Nnode_v) )
103 allocate( intweight_v(elem%Nnode_h1D**2) )
104
105 intweight_lgl1dpts_h(:) = polynominal_gengausslobattoptintweight(elem%PolyOrder_h)
106 intweight_lgl1dpts_v(:) = polynominal_gengausslobattoptintweight(elem%PolyOrder_v)
107
108 do f=1, elem%Nfaces_h
109 do k=1, elem%Nnode_v
110 do i=1, elem%Nnode_h1D
111 l = i + (k-1)*elem%Nnode_h1D
112 intweight_h(l) = intweight_lgl1dpts_h(i) * intweight_lgl1dpts_v(k)
113 end do
114 end do
115
116 is = (f-1)*elem%Nfp_h + 1
117 ie = is + elem%Nfp_h - 1
118 intweight(f,is:ie) = intweight_h(:)
119 end do
120
121 do f=1, elem%Nfaces_v
122 do j=1, elem%Nnode_h1D
123 do i=1, elem%Nnode_h1D
124 l = i + (j-1)*elem%Nnode_h1D
125 intweight_v(l) = intweight_lgl1dpts_h(i) * intweight_lgl1dpts_h(j)
126 end do
127 end do
128
129 is = elem%Nfaces_h*elem%Nfp_h + (f-1)*elem%Nfp_v + 1
130 ie = is + elem%Nfp_v - 1
131 intweight(elem%Nfaces_h+f,is:ie) = intweight_v(:)
132 end do
133
134 call faceintmat%Init( intweight )
135
136 return
scale_polynominal
Module common / Polynominal.
Definition scale_polynominal.F90:11
scale_polynominal::polynominal_gengausslobattoptintweight
real(rp) function, dimension(nord+1), public polynominal_gengausslobattoptintweight(nord)
A function to calculate the Gauss-Lobbato weights.
Definition scale_polynominal.F90:272

References scale_polynominal::polynominal_gengausslobattoptintweight().

◆ atm_dyn_dgm_trcadvect3d_heve_final()

subroutine, public scale_atm_dyn_dgm_trcadvect3d_heve::atm_dyn_dgm_trcadvect3d_heve_final

Definition at line 140 of file scale_atm_dyn_dgm_trcadvect3d_heve.F90.

141 implicit none
142 !--------------------------------------------
143
144 return

◆ atm_dyn_dgm_trcadvect3d_heve_cal_tend()

subroutine, public scale_atm_dyn_dgm_trcadvect3d_heve::atm_dyn_dgm_trcadvect3d_heve_cal_tend ( real(rp), dimension(elem%np,lmesh%nea), intent(out) qtrc_dt,
real(rp), dimension(elem%np,lmesh%nea), intent(in) qtrc_,
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%nfptot,lmesh%ne), intent(in) alphdens_m,
real(rp), dimension(elem%nfptot,lmesh%ne), intent(in) alphdens_p,
real(rp), dimension(elem%np,lmesh%nea), intent(in) fct_coef,
real(rp), dimension(elem%np,lmesh%nea), intent(in) rhoq_tp,
class(elementoperationbase3d), intent(in) element3d_operation,
type(sparsemat), intent(in) faceintmat,
class(localmesh3d), intent(in) lmesh,
class(elementbase3d), intent(in) elem,
class(localmesh2d), intent(in) lmesh2d,
class(elementbase2d), intent(in) elem2d )

Definition at line 149 of file scale_atm_dyn_dgm_trcadvect3d_heve.F90.

156
157 class(LocalMesh3D), intent(in) :: lmesh
158 class(ElementBase3D), intent(in) :: elem
159 class(LocalMesh2D), intent(in) :: lmesh2D
160 class(ElementBase2D), intent(in) :: elem2D
161
162 real(RP), intent(out) :: QTRC_dt(elem%Np,lmesh%NeA)
163 real(RP), intent(in) :: QTRC_(elem%Np,lmesh%NeA)
164 real(RP), intent(in) :: MOMX_(elem%Np,lmesh%NeA)
165 real(RP), intent(in) :: MOMY_(elem%Np,lmesh%NeA)
166 real(RP), intent(in) :: MOMZ_(elem%Np,lmesh%NeA)
167 real(RP), intent(in) :: alphDENS_M(elem%NfpTot,lmesh%Ne)
168 real(RP), intent(in) :: alphDENS_P(elem%NfpTot,lmesh%Ne)
169 real(RP), intent(in) :: fct_coef(elem%Np,lmesh%NeA)
170 real(RP), intent(in) :: RHOQ_tp(elem%Np,lmesh%NeA)
171 class(ElementOperationBase3D), intent(in) :: element3D_operation
172 type(SparseMat), intent(in) :: FaceIntMat
173
174 real(RP) :: Flux(elem%Np,3), DFlux(elem%Np,4)
175 real(RP) :: del_flux(elem%NfpTot,lmesh%Ne)
176 real(RP) :: momwt_
177 real(RP) :: Gsqrt_
178 real(RP) :: RGsqrt(elem%Np), RGsqrtV(elem%Np)
179
180 integer :: ke, ke2d
181
182 real(RP) :: Q0, Q1, vol
183 integer :: p
184 !---------------------------------------------------------------------------
185
186 call prof_rapstart('cal_trcadv_tend_bndflux', 3)
187 call atm_dyn_dgm_trcadvect3d_heve_get_delflux_generalhvc( &
188 del_flux, & ! (out)
189 qtrc_, momx_, momy_, momz_, alphdens_m, alphdens_p, fct_coef, & ! (in)
190 lmesh%Gsqrt, lmesh%GsqrtH, lmesh%GI3(:,:,1), lmesh%GI3(:,:,2), & ! (in)
191 lmesh%normal_fn(:,:,1), lmesh%normal_fn(:,:,2), lmesh%normal_fn(:,:,3), & ! (in)
192 lmesh%J(:,:), lmesh%Fscale(:,:), & ! (in)
193 lmesh%vmapM, lmesh%vmapP, elem%IndexH2Dto3D_bnd, faceintmat, & ! (in)
194 lmesh, elem, lmesh2d, elem2d ) ! (in)
195 call prof_rapend('cal_trcadv_tend_bndflux', 3)
196
197 call prof_rapstart('cal_trcadv_tend_interior', 3)
198 !$omp parallel do private( ke, ke2D, &
199 !$omp momwt_, Flux, DFlux, Gsqrt_, RGsqrt, RGsqrtV )
200 do ke=lmesh%NeS, lmesh%NeE
201 ke2d = lmesh%EMap3Dto2D(ke)
202
203 do p=1, elem%Np
204 rgsqrt(p) = 1.0_rp / lmesh%Gsqrt(p,ke)
205 rgsqrtv(p) = lmesh%GsqrtH(elem%IndexH2Dto3D(p),ke2d) * rgsqrt(p)
206 end do
207
208 do p=1, elem%Np
209 gsqrt_ = lmesh%Gsqrt(p,ke)
210 flux(p,1) = gsqrt_ * momx_(p,ke) * qtrc_(p,ke)
211 flux(p,2) = gsqrt_ * momy_(p,ke) * qtrc_(p,ke)
212 flux(p,3) = gsqrt_ * ( &
213 momz_(p,ke) * rgsqrtv(p) &
214 + lmesh%GI3(p,ke,1) * momx_(p,ke) &
215 + lmesh%GI3(p,ke,2) * momy_(p,ke) ) * qtrc_(p,ke)
216 end do
217 call element3d_operation%Div( &
218 flux, del_flux(:,ke), & ! (in)
219 dflux ) ! (out)
220
221 qtrc_dt(:,ke) = - ( &
222 lmesh%Escale(:,ke,1,1) * dflux(:,1) &
223 + lmesh%Escale(:,ke,2,2) * dflux(:,2) &
224 + lmesh%Escale(:,ke,3,3) * dflux(:,3) &
225 + dflux(:,4) ) / lmesh%Gsqrt(:,ke) &
226 + rhoq_tp(:,ke)
227 end do
228 call prof_rapend('cal_trcadv_tend_interior', 3)
229
230 return

◆ atm_dyn_dgm_trcadvect3d_heve_calc_fct_coef()

subroutine, public scale_atm_dyn_dgm_trcadvect3d_heve::atm_dyn_dgm_trcadvect3d_heve_calc_fct_coef ( real(rp), dimension(elem%np,lmesh%nea), intent(out) fct_coef,
real(rp), dimension(elem%np,lmesh%nea), intent(in) qtrc_,
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) rhoq_tp_,
real(rp), dimension(elem%nfptot,lmesh%ne), intent(in) alphdens_m,
real(rp), dimension(elem%nfptot,lmesh%ne), intent(in) alphdens_p,
real(rp), dimension(elem%np,lmesh%nea), intent(in) dens_hyd,
real(rp), dimension (elem%np,lmesh%nea), intent(in) ddens_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) ddens0_,
real(rp), intent(in) rk_c_ssm1,
real(rp), intent(in) dt,
type(sparsemat), intent(in) faceintmat,
class(localmesh3d), intent(in) lmesh,
class(elementbase3d), intent(in) elem,
class(localmesh2d), intent(in) lmesh2d,
class(elementbase2d), intent(in) elem2d,
logical, intent(in), optional disable_limiter )

Definition at line 234 of file scale_atm_dyn_dgm_trcadvect3d_heve.F90.

240
241 implicit none
242 class(LocalMesh3D), intent(in) :: lmesh
243 class(ElementBase3D), intent(in) :: elem
244 class(LocalMesh2D), intent(in) :: lmesh2D
245 class(ElementBase2D), intent(in) :: elem2D
246 real(RP), intent(out) :: fct_coef(elem%Np,lmesh%NeA)
247 real(RP), intent(in) :: QTRC_(elem%Np,lmesh%NeA)
248 real(RP), intent(in) :: MOMX_(elem%Np,lmesh%NeA)
249 real(RP), intent(in) :: MOMY_(elem%Np,lmesh%NeA)
250 real(RP), intent(in) :: MOMZ_(elem%Np,lmesh%NeA)
251 real(RP), intent(in) :: RHOQ_tp_(elem%Np,lmesh%NeA)
252 real(RP), intent(in) :: alphDENS_M(elem%NfpTot,lmesh%Ne)
253 real(RP), intent(in) :: alphDENS_P(elem%NfpTot,lmesh%Ne)
254 real(RP), intent(in) :: DENS_hyd(elem%Np,lmesh%NeA)
255 real(RP), intent(in) :: DDENS_ (elem%Np,lmesh%NeA)
256 real(RP), intent(in) :: DDENS0_(elem%Np,lmesh%NeA)
257 real(RP), intent(in) :: rk_c_ssm1
258 real(RP), intent(in) :: dt
259 type(SparseMat), intent(in) :: FaceIntMat
260 logical, intent(in), optional :: disable_limiter
261
262 real(RP) :: netOutwardFlux(lmesh%Ne)
263 real(RP) :: momwt_(elem%Np)
264
265 integer :: ke
266 real(RP) :: Q
267 real(RP) :: dens_ssm1(elem%Np)
268 !---------------------------------------------------------------------------
269
270 if ( present(disable_limiter) ) then
271 if ( disable_limiter ) then
272 !$omp parallel do
273 do ke=lmesh%NeS, lmesh%NeE
274 fct_coef(:,ke) = 1.0_rp
275 end do
276 return
277 end if
278 end if
279
280 call prof_rapstart('cal_trcadv_fct_coef_bndflux', 3)
281 call atm_dyn_dgm_trcadvect3d_heve_get_netoutwardflux_generalhvc( &
282 netoutwardflux, & ! (out)
283 qtrc_, momx_, momy_, momz_, alphdens_m, alphdens_p, & ! (in)
284 lmesh%Gsqrt, lmesh%GsqrtH, lmesh%GI3(:,:,1), lmesh%GI3(:,:,2), & ! (in)
285 lmesh%normal_fn(:,:,1), lmesh%normal_fn(:,:,2), lmesh%normal_fn(:,:,3), & ! (in)
286 lmesh%J(:,:), lmesh%Fscale(:,:), lmesh%vmapM, lmesh%vmapP, elem%IndexH2Dto3D_bnd, & ! (in)
287 faceintmat, & ! (in)
288 lmesh, elem, lmesh2d, elem2d ) ! (in)
289 call prof_rapend('cal_trcadv_fct_coef_bndflux', 3)
290
291 call prof_rapstart('cal_trcadv_fct_coef', 3)
292
293 !$omp parallel do private( ke, Q, dens_ssm1 )
294 do ke=lmesh%NeS, lmesh%NeE
295
296 dens_ssm1(:) = dens_hyd(:,ke) &
297 + ( 1.0_rp - rk_c_ssm1 ) * ddens0_(:,ke) + rk_c_ssm1 * ddens_(:,ke)
298 q = sum( lmesh%Gsqrt(:,ke) * lmesh%J(:,ke) * elem%IntWeight_lgl(:) * ( dens_ssm1(:) * qtrc_(:,ke) / dt + rhoq_tp_(:,ke) ) )
299
300 fct_coef(:,ke) = max( 0.0_rp, min( 1.0_rp, q / ( netoutwardflux(ke) + 1.0e-10_rp ) ) )
301 end do
302
303 call prof_rapend('cal_trcadv_fct_coef', 3)
304
305 return

◆ atm_dyn_dgm_trcadvect3d_tmar()

subroutine, public scale_atm_dyn_dgm_trcadvect3d_heve::atm_dyn_dgm_trcadvect3d_tmar ( real(rp), dimension(elem%np,lmesh%nea), intent(inout) qtrc_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) dens_hyd,
real(rp), dimension (elem%np,lmesh%nea), intent(in) ddens_,
class(localmesh3d), intent(in) lmesh,
class(elementbase3d), intent(in) elem,
class(localmesh2d), intent(in) lmesh2d,
class(elementbase2d), intent(in) elem2d )

Second Step of limiter in which nonlinear truncation and mass aware rescaling (TMAR)

Definition at line 311 of file scale_atm_dyn_dgm_trcadvect3d_heve.F90.

313
314 implicit none
315 class(LocalMesh3D), intent(in) :: lmesh
316 class(ElementBase3D), intent(in) :: elem
317 class(LocalMesh2D), intent(in) :: lmesh2D
318 class(ElementBase2D), intent(in) :: elem2D
319 real(RP), intent(inout) :: QTRC_(elem%Np,lmesh%NeA)
320 real(RP), intent(in) :: DENS_hyd(elem%Np,lmesh%NeA)
321 real(RP), intent(in) :: DDENS_ (elem%Np,lmesh%NeA)
322
323 integer :: ke
324 real(RP) :: Q0, Q1
325 real(RP) :: Q(elem%Np)
326 real(RP) :: dens(elem%Np)
327 !--------------------------------------------------------
328
329 !$omp parallel do private( dens, Q, Q0, Q1 )
330 do ke=lmesh%NeS, lmesh%NeE
331 dens(:) = dens_hyd(:,ke) + ddens_(:,ke)
332 q(:) = max( 0.0_rp, qtrc_(:,ke) )
333
334 q0 = sum( lmesh%Gsqrt(:,ke) * lmesh%J(:,ke) * elem%IntWeight_lgl(:) * dens(:) * qtrc_(:,ke) )
335 q1 = sum( lmesh%Gsqrt(:,ke) * lmesh%J(:,ke) * elem%IntWeight_lgl(:) * dens(:) * q(:) )
336 qtrc_(:,ke) = q0 / ( q1 + 1.0e-32_rp ) * q(:)
337 end do
338
339 return

◆ atm_dyn_dgm_trcadvect3d_save_massflux()

subroutine, public scale_atm_dyn_dgm_trcadvect3d_heve::atm_dyn_dgm_trcadvect3d_save_massflux ( real(rp), dimension(elem%np,lmesh%nea), intent(inout) mflx_x_tavg,
real(rp), dimension(elem%np,lmesh%nea), intent(inout) mflx_y_tavg,
real(rp), dimension(elem%np,lmesh%nea), intent(inout) mflx_z_tavg,
real(rp), dimension(elem%nfptot,lmesh%ne), intent(inout) alph_dens_m,
real(rp), dimension(elem%nfptot,lmesh%ne), intent(inout) alph_dens_p,
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) 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) rtot,
real(rp), dimension(elem%np,lmesh%nea), intent(in) cvtot,
real(rp), dimension(elem%np,lmesh%nea), intent(in) cptot,
class(localmesh3d), intent(in) lmesh,
class(elementbase3d), intent(in) elem,
integer, intent(in) rkstage,
real(rp), intent(in) tavg_weight_h,
real(rp), intent(in) tavg_weight_v,
logical, intent(in) is_hevi )

Definition at line 343 of file scale_atm_dyn_dgm_trcadvect3d_heve.F90.

348
349 implicit none
350 class(LocalMesh3D), intent(in) :: lmesh
351 class(ElementBase3D), intent(in) :: elem
352 real(RP), intent(inout) :: MFLX_x_tavg(elem%Np,lmesh%NeA)
353 real(RP), intent(inout) :: MFLX_y_tavg(elem%Np,lmesh%NeA)
354 real(RP), intent(inout) :: MFLX_z_tavg(elem%Np,lmesh%NeA)
355 real(RP), intent(inout) :: alph_dens_M(elem%NfpTot,lmesh%Ne)
356 real(RP), intent(inout) :: alph_dens_P(elem%NfpTot,lmesh%Ne)
357 real(RP), intent(in) :: DDENS(elem%Np,lmesh%NeA)
358 real(RP), intent(in) :: MOMX(elem%Np,lmesh%NeA)
359 real(RP), intent(in) :: MOMY(elem%Np,lmesh%NeA)
360 real(RP), intent(in) :: MOMZ(elem%Np,lmesh%NeA)
361 real(RP), intent(in) :: DPRES(elem%Np,lmesh%NeA)
362 real(RP), intent(in) :: DENS_hyd(elem%Np,lmesh%NeA)
363 real(RP), intent(in) :: PRES_hyd(elem%Np,lmesh%NeA)
364 real(RP), intent(in) :: Rtot(elem%Np,lmesh%NeA)
365 real(RP), intent(in) :: CVtot(elem%Np,lmesh%NeA)
366 real(RP), intent(in) :: CPtot(elem%Np,lmesh%NeA)
367 integer, intent(in) :: rkstage
368 real(RP), intent(in) :: tavg_weight_h
369 real(RP), intent(in) :: tavg_weight_v
370 logical, intent(in) :: is_hevi
371
372 integer :: ke
373 !--------------------------------------------------------------
374
375 !$omp parallel do private(ke)
376 do ke=lmesh%NeS, lmesh%NeE
377 if (rkstage == 1) then
378 mflx_x_tavg(:,ke) = tavg_weight_h * momx(:,ke)
379 mflx_y_tavg(:,ke) = tavg_weight_h * momy(:,ke)
380 mflx_z_tavg(:,ke) = tavg_weight_v * momz(:,ke)
381 alph_dens_m(:,ke) = 0.0_rp
382 alph_dens_p(:,ke) = 0.0_rp
383 else
384 mflx_x_tavg(:,ke) = mflx_x_tavg(:,ke) + tavg_weight_h * momx(:,ke)
385 mflx_y_tavg(:,ke) = mflx_y_tavg(:,ke) + tavg_weight_h * momy(:,ke)
386 mflx_z_tavg(:,ke) = mflx_z_tavg(:,ke) + tavg_weight_v * momz(:,ke)
387 end if
388 end do
389
390 call atm_dyn_dgm_trcadvect3d_heve_cal_alphdens_dyn( &
391 alph_dens_m, alph_dens_p, & ! (out)
392 ddens, momx, momy, momz, dpres, dens_hyd, pres_hyd, & ! (in)
393 rtot, cvtot, cptot, & ! (in)
394 lmesh%Gsqrt, lmesh%GIJ(:,:,1,1), lmesh%GIJ(:,:,1,2), lmesh%GIJ(:,:,2,2), & ! (in)
395 lmesh%normal_fn(:,:,1), lmesh%normal_fn(:,:,2), lmesh%normal_fn(:,:,3), & ! (in)
396 lmesh%vmapM, lmesh%vmapP, elem%IndexH2Dto3D_bnd, & ! (in)
397 lmesh, lmesh%lcmesh2D, lmesh%refElem3D, lmesh%lcmesh2D%refElem2D, & ! (in)
398 tavg_weight_h, tavg_weight_v, is_hevi ) ! (in)
399
400 return

◆ atm_dyn_dgm_trcadvect3d_heve_cal_alphdens_advtest()

subroutine, public scale_atm_dyn_dgm_trcadvect3d_heve::atm_dyn_dgm_trcadvect3d_heve_cal_alphdens_advtest ( real(rp), dimension(elem%nfptot*lmesh%ne), intent(inout) alph_dens_m,
real(rp), dimension(elem%nfptot*lmesh%ne), intent(inout) alph_dens_p,
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) dens_hyd,
real(rp), dimension(elem%np*lmesh%nea), intent(in) gsqrt,
real(rp), dimension(elem%nfptot*lmesh%ne), intent(in) nx,
real(rp), dimension(elem%nfptot*lmesh%ne), intent(in) ny,
real(rp), dimension(elem%nfptot*lmesh%ne), intent(in) nz,
integer, dimension(elem%nfptot*lmesh%ne), intent(in) vmapm,
integer, dimension(elem%nfptot*lmesh%ne), intent(in) vmapp,
class(localmesh3d), intent(in) lmesh,
class(elementbase3d), intent(in) elem )

Definition at line 404 of file scale_atm_dyn_dgm_trcadvect3d_heve.F90.

407
408 use scale_const, only: &
409 grav => const_grav, &
410 rdry => const_rdry, &
411 cpdry => const_cpdry, &
412 cvdry => const_cvdry, &
413 pres00 => const_pre00
414
415 implicit none
416
417 class(LocalMesh3D), intent(in) :: lmesh
418 class(ElementBase3D), intent(in) :: elem
419 real(RP), intent(inout) :: alph_dens_M(elem%NfpTot*lmesh%Ne)
420 real(RP), intent(inout) :: alph_dens_P(elem%NfpTot*lmesh%Ne)
421 real(RP), intent(in) :: DDENS_(elem%Np*lmesh%NeA)
422 real(RP), intent(in) :: MOMX_(elem%Np*lmesh%NeA)
423 real(RP), intent(in) :: MOMY_(elem%Np*lmesh%NeA)
424 real(RP), intent(in) :: MOMZ_(elem%Np*lmesh%NeA)
425 real(RP), intent(in) :: DENS_hyd(elem%Np*lmesh%NeA)
426 real(RP), intent(in) :: Gsqrt(elem%Np*lmesh%NeA)
427 real(RP), intent(in) :: nx(elem%NfpTot*lmesh%Ne)
428 real(RP), intent(in) :: ny(elem%NfpTot*lmesh%Ne)
429 real(RP), intent(in) :: nz(elem%NfpTot*lmesh%Ne)
430 integer, intent(in) :: vmapM(elem%NfpTot*lmesh%Ne)
431 integer, intent(in) :: vmapP(elem%NfpTot*lmesh%Ne)
432
433 integer :: i, iP, iM
434 real(RP) :: VelP, VelM, alpha, densM, densP
435 !------------------------------------------------------------------------
436
437
438 !$omp parallel do private( &
439 !$omp iM, iP, VelP, VelM, alpha, densM, densP )
440 do i=1, elem%NfpTot*lmesh%Ne
441 im = vmapm(i); ip = vmapp(i)
442
443 densm = ddens_(im) + dens_hyd(im)
444 densp = ddens_(ip) + dens_hyd(ip)
445
446 velm = ( momx_(im) * nx(i) + momy_(im) * ny(i) + momz_(im) * nz(i) ) / densm
447 velp = ( momx_(ip) * nx(i) + momy_(ip) * ny(i) + momz_(ip) * nz(i) ) / densp
448
449 alpha = max( abs(velm), abs(velp) )
450 alph_dens_m(i) = alpha * densm * gsqrt(im)
451 alph_dens_p(i) = alpha * densp * gsqrt(ip)
452 end do
453
454 return
Generated by doxygen 1.14.0