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

module FElib / Atmosphere / Physics turbulence More...

Functions/Subroutines

subroutine, public atm_phy_tb_dgm_dns_init (mesh)
subroutine, public atm_phy_tb_dgm_dns_final ()
subroutine, public atm_phy_tb_dgm_dns_cal_grad (t11, t12, t13, t21, t22, t23, t31, t32, t33, df1, df2, df3, tke, nu, kh, ddens_, momx_, momy_, momz_, drhot_, dens_hyd, pres_hyd, pres, pt, dx, dy, dz, sx, sy, sz, lift, lmesh, elem, lmesh2d, elem2d, is_bound)
 Calculate parameterized stress tensor and eddy heat flux with turbulent model.

Detailed Description

module FElib / Atmosphere / Physics turbulence

Description
Turbulence process for DNS
Author
Yuta Kawai, Team SCALE
NAMELIST
  • PARAM_ATMOS_PHY_TB_DGM_DNS
    nametypedefault valuecomment
    DNS_NUreal(RP)1.512E-5_RPkinematic viscosity coefficient [m2/s] for air at 20degC
    DNS_MUreal(RP)1.8E-5_RPmolecular diffusive coefficient [m2/s] for air at 20degC

History Output
No history output

Function/Subroutine Documentation

◆ atm_phy_tb_dgm_dns_init()

subroutine, public scale_atm_phy_tb_dgm_dns::atm_phy_tb_dgm_dns_init ( class(meshbase3d), intent(in) mesh)

Definition at line 69 of file scale_atm_phy_tb_dgm_dns.F90.

70 implicit none
71 class(MeshBase3D), intent(in) :: mesh
72
73 logical :: consistent_tke = .true.
74
75 namelist / param_atmos_phy_tb_dgm_dns / &
76 dns_nu, dns_mu
77
78 integer :: ierr
79 !--------------------------------------------------------------------
80
81 log_newline
82 log_info("ATMOS_PHY_TB_dgm_dns_setup",*) 'Setup'
83 log_info("ATMOS_PHY_TB_dgm_dns_setup",*) 'Eddy Viscocity Model for DNS'
84
85 !--- read namelist
86 rewind(io_fid_conf)
87 read(io_fid_conf,nml=param_atmos_phy_tb_dgm_dns,iostat=ierr)
88 if( ierr < 0 ) then !--- missing
89 log_info("ATMOS_PHY_TB_dgm_dns_setup",*) 'Not found namelist. Default used.'
90 elseif( ierr > 0 ) then !--- fatal error
91 log_error("ATMOS_PHY_TB_dgm_dns_setup",*) 'Not appropriate names in namelist PARAM_ATMOS_PHY_TB_DGM_DNS. Check!'
92 call prc_abort
93 endif
94 log_nml(param_atmos_phy_tb_dgm_dns)
95
96 return

◆ atm_phy_tb_dgm_dns_final()

subroutine, public scale_atm_phy_tb_dgm_dns::atm_phy_tb_dgm_dns_final

Definition at line 100 of file scale_atm_phy_tb_dgm_dns.F90.

101 implicit none
102 !--------------------------------------------------------------------
103
104 return

◆ atm_phy_tb_dgm_dns_cal_grad()

subroutine, public scale_atm_phy_tb_dgm_dns::atm_phy_tb_dgm_dns_cal_grad ( real(rp), dimension(elem%np,lmesh%nea), intent(out) t11,
real(rp), dimension(elem%np,lmesh%nea), intent(out) t12,
real(rp), dimension(elem%np,lmesh%nea), intent(out) t13,
real(rp), dimension(elem%np,lmesh%nea), intent(out) t21,
real(rp), dimension(elem%np,lmesh%nea), intent(out) t22,
real(rp), dimension(elem%np,lmesh%nea), intent(out) t23,
real(rp), dimension(elem%np,lmesh%nea), intent(out) t31,
real(rp), dimension(elem%np,lmesh%nea), intent(out) t32,
real(rp), dimension(elem%np,lmesh%nea), intent(out) t33,
real(rp), dimension(elem%np,lmesh%nea), intent(out) df1,
real(rp), dimension(elem%np,lmesh%nea), intent(out) df2,
real(rp), dimension(elem%np,lmesh%nea), intent(out) df3,
real(rp), dimension(elem%np,lmesh%nea), intent(out) tke,
real(rp), dimension(elem%np,lmesh%nea), intent(out) nu,
real(rp), dimension(elem%np,lmesh%nea), intent(out) kh,
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) dens_hyd,
real(rp), dimension(elem%np,lmesh%nea), intent(in) pres_hyd,
real(rp), dimension(elem%np,lmesh%nea), intent(in) pres,
real(rp), dimension(elem%np,lmesh%nea), intent(in) pt,
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,
logical, dimension(elem%nfptot,lmesh%ne), intent(in) is_bound )

Calculate parameterized stress tensor and eddy heat flux with turbulent model.

Parameters
[out]t11(1,1) component of stress tensor
[out]t12(1,2) component of stress tensor
[out]t13(1,3) component of stress tensor
[out]t21(2,1) component of stress tensor
[out]t22(2,2) component of stress tensor
[out]t23(2,3) component of stress tensor
[out]t31(3,1) component of stress tensor
[out]t32(3,2) component of stress tensor
[out]t33(3,3) component of stress tensor
[out]df1Diffusive heat flux in x1 direction / density
[out]df2Diffusive heat flux in x2 direction / density
[out]df3Diffusive heat flux in x3 direction / density
[out]tkeParameterized turbulent kinetic energy
[out]nuEddy viscosity
[out]khEddy diffusivity
[in]ddens_Density perturbation
[in]momx_Momentum in x1 direction
[in]momy_Momentum in x2 direction
[in]momz_Momentum in x3 direction
[in]drhot_Density x potential temperature perturbation
[in]dens_hydReference pressure in hydrostatic balance
[in]pres_hydReference density in hydrostatic balance
[in]presPressure
[in]ptPotential temperature
[in]dzDifferential matrix managed by sparse matrix type
[in]szStiffness matrix managed by sparse matrix type
[in]liftLifting matrix managed by sparse matrix type
[in]is_boundFlag whether nodes are located at domain boundaries

Definition at line 110 of file scale_atm_phy_tb_dgm_dns.F90.

118
119 implicit none
120
121 class(LocalMesh3D), intent(in) :: lmesh
122 class(ElementBase3D), intent(in) :: elem
123 class(LocalMesh2D), intent(in) :: lmesh2D
124 class(ElementBase2D), intent(in) :: elem2D
125 real(RP), intent(out) :: T11(elem%Np,lmesh%NeA) !< (1,1) component of stress tensor
126 real(RP), intent(out) :: T12(elem%Np,lmesh%NeA) !< (1,2) component of stress tensor
127 real(RP), intent(out) :: T13(elem%Np,lmesh%NeA) !< (1,3) component of stress tensor
128 real(RP), intent(out) :: T21(elem%Np,lmesh%NeA) !< (2,1) component of stress tensor
129 real(RP), intent(out) :: T22(elem%Np,lmesh%NeA) !< (2,2) component of stress tensor
130 real(RP), intent(out) :: T23(elem%Np,lmesh%NeA) !< (2,3) component of stress tensor
131 real(RP), intent(out) :: T31(elem%Np,lmesh%NeA) !< (3,1) component of stress tensor
132 real(RP), intent(out) :: T32(elem%Np,lmesh%NeA) !< (3,2) component of stress tensor
133 real(RP), intent(out) :: T33(elem%Np,lmesh%NeA) !< (3,3) component of stress tensor
134 real(RP), intent(out) :: DF1(elem%Np,lmesh%NeA) !< Diffusive heat flux in x1 direction / density
135 real(RP), intent(out) :: DF2(elem%Np,lmesh%NeA) !< Diffusive heat flux in x2 direction / density
136 real(RP), intent(out) :: DF3(elem%Np,lmesh%NeA) !< Diffusive heat flux in x3 direction / density
137 real(RP), intent(out) :: TKE(elem%Np,lmesh%NeA) !< Parameterized turbulent kinetic energy
138 real(RP), intent(out) :: Nu(elem%Np,lmesh%NeA) !< Eddy viscosity
139 real(RP), intent(out) :: Kh(elem%Np,lmesh%NeA) !< Eddy diffusivity
140 real(RP), intent(in) :: DDENS_(elem%Np,lmesh%NeA) !< Density perturbation
141 real(RP), intent(in) :: MOMX_ (elem%Np,lmesh%NeA) !< Momentum in x1 direction
142 real(RP), intent(in) :: MOMY_ (elem%Np,lmesh%NeA) !< Momentum in x2 direction
143 real(RP), intent(in) :: MOMZ_ (elem%Np,lmesh%NeA) !< Momentum in x3 direction
144 real(RP), intent(in) :: DRHOT_(elem%Np,lmesh%NeA) !< Density x potential temperature perturbation
145 real(RP), intent(in) :: DENS_hyd(elem%Np,lmesh%NeA) !< Reference pressure in hydrostatic balance
146 real(RP), intent(in) :: PRES_hyd(elem%Np,lmesh%NeA) !< Reference density in hydrostatic balance
147 real(RP), intent(in) :: PRES(elem%Np,lmesh%NeA) !< Pressure
148 real(RP), intent(in) :: PT(elem%Np,lmesh%NeA) !< Potential temperature
149 type(SparseMat), intent(in) :: Dx, Dy, Dz !< Differential matrix managed by sparse matrix type
150 type(SparseMat), intent(in) :: Sx, Sy, Sz !< Stiffness matrix managed by sparse matrix type
151 type(SparseMat), intent(in) :: Lift !< Lifting matrix managed by sparse matrix type
152 logical, intent(in) :: is_bound(elem%NfpTot,lmesh%Ne) !< Flag whether nodes are located at domain boundaries
153
154 real(RP) :: Fx(elem%Np), Fy(elem%Np), Fz(elem%Np), LiftDelFlx(elem%Np)
155 real(RP) :: DENS(elem%Np), RDENS(elem%Np), RHOT(elem%Np), Q(elem%Np)
156 real(RP) :: DdensDxi(elem%Np,3)
157 real(RP) :: DVelDxi(elem%Np,3,3)
158 real(RP) :: del_flux_rho (elem%NfpTot,lmesh%Ne,3)
159 real(RP) :: del_flux_mom (elem%NfpTot,lmesh%Ne,3,3)
160 real(RP) :: del_flux_rhot(elem%NfpTot,lmesh%Ne,3)
161
162 real(RP) :: S11(elem%Np), S12(elem%Np), S22(elem%Np), S23(elem%Np), S31(elem%Np), S33(elem%Np)
163 real(RP) :: SkkOvThree
164 real(RP) :: coef
165
166 integer :: ke
167 integer :: p
168 !--------------------------------------------------------------------
169
170 call cal_del_flux_grad( del_flux_rho, del_flux_mom, del_flux_rhot, & ! (out)
171 ddens_, momx_, momy_, momz_, drhot_, dens_hyd, pres_hyd, pt, & ! (in)
172 lmesh%normal_fn(:,:,1), lmesh%normal_fn(:,:,2), lmesh%normal_fn(:,:,3), & ! (in)
173 lmesh%vmapM, lmesh%vmapP, & ! (in)
174 lmesh, elem, is_bound ) ! (in)
175
176 !$omp parallel do private( &
177 !$omp Fx, Fy, Fz, LiftDelFlx, &
178 !$omp DENS, RHOT, RDENS, Q, DdensDxi, DVelDxi, &
179 !$omp S11, S12, S22, S23, S31, S33, &
180 !$omp coef, SkkOvThree, &
181 !$omp p )
182 do ke=lmesh%NeS, lmesh%NeE
183 !---
184 dens(:) = dens_hyd(:,ke) + ddens_(:,ke)
185 rdens(:) = 1.0_rp / dens(:)
186 rhot(:) = dens(:) * pt(:,ke)
187
188 ! gradient of density
189 call sparsemat_matmul( dx, dens, fx )
190 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_rho(:,ke,1), liftdelflx )
191 ddensdxi(:,1) = lmesh%Escale(:,ke,1,1) * fx(:) + liftdelflx(:)
192
193 call sparsemat_matmul( dy, dens, fy )
194 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_rho(:,ke,2), liftdelflx )
195 ddensdxi(:,2) = lmesh%Escale(:,ke,2,2) * fy(:) + liftdelflx(:)
196
197 call sparsemat_matmul( dz, dens, fz )
198 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_rho(:,ke,3), liftdelflx )
199 ddensdxi(:,3) = lmesh%Escale(:,ke,3,3) * fz(:) + liftdelflx(:)
200
201 ! gradient of u
202 q(:) = momx_(:,ke) * rdens(:)
203
204 call sparsemat_matmul( dx, momx_(:,ke), fx )
205 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_mom(:,ke,1,1), liftdelflx )
206 dveldxi(:,1,1) = ( lmesh%Escale(:,ke,1,1) * fx(:) + liftdelflx(:) - q(:) * ddensdxi(:,1) ) * rdens(:)
207
208 call sparsemat_matmul( dy, momx_(:,ke), fy )
209 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_mom(:,ke,2,1), liftdelflx )
210 dveldxi(:,2,1) = ( lmesh%Escale(:,ke,2,2) * fy(:) + liftdelflx(:) - q(:) * ddensdxi(:,2) ) * rdens(:)
211
212 call sparsemat_matmul( dz, momx_(:,ke), fz )
213 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_mom(:,ke,3,1), liftdelflx )
214 dveldxi(:,3,1) = ( lmesh%Escale(:,ke,3,3) * fz(:) + liftdelflx(:) - q(:) * ddensdxi(:,3) ) * rdens(:)
215
216 ! gradient of v
217 q(:) = momy_(:,ke) * rdens(:)
218
219 call sparsemat_matmul( dx, momy_(:,ke), fx )
220 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_mom(:,ke,1,2), liftdelflx )
221 dveldxi(:,1,2) = ( lmesh%Escale(:,ke,1,1) * fx(:) + liftdelflx(:) - q(:) * ddensdxi(:,1) ) * rdens(:)
222
223 call sparsemat_matmul( dy, momy_(:,ke), fy )
224 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_mom(:,ke,2,2), liftdelflx )
225 dveldxi(:,2,2) = ( lmesh%Escale(:,ke,2,2) * fy(:) + liftdelflx(:) - q(:) * ddensdxi(:,2) ) * rdens(:)
226
227 call sparsemat_matmul( dz, momy_(:,ke), fz )
228 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_mom(:,ke,3,2), liftdelflx )
229 dveldxi(:,3,2) = ( lmesh%Escale(:,ke,3,3) * fz(:) + liftdelflx(:) - q(:) * ddensdxi(:,3) ) * rdens(:)
230
231 ! gradient of w
232 q(:) = momz_(:,ke) * rdens(:)
233
234 call sparsemat_matmul( dx, momz_(:,ke), fx )
235 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_mom(:,ke,1,3), liftdelflx )
236 dveldxi(:,1,3) = ( lmesh%Escale(:,ke,1,1) * fx(:) + liftdelflx(:) - q(:) * ddensdxi(:,1) ) * rdens(:)
237
238 call sparsemat_matmul( dy, momz_(:,ke), fy )
239 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_mom(:,ke,2,3), liftdelflx )
240 dveldxi(:,2,3) = ( lmesh%Escale(:,ke,2,2) * fy(:) + liftdelflx(:) - q(:) * ddensdxi(:,2) ) * rdens(:)
241
242 call sparsemat_matmul( dz, momz_(:,ke), fz )
243 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_mom(:,ke,3,3), liftdelflx )
244 dveldxi(:,3,3) = ( lmesh%Escale(:,ke,3,3) * fz(:) + liftdelflx(:) - q(:) * ddensdxi(:,3) ) * rdens(:)
245
246 ! gradient of pt
247 q(:) = rhot(:) * rdens(:)
248
249 call sparsemat_matmul( dx, rhot, fx )
250 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_rhot(:,ke,1), liftdelflx )
251 df1(:,ke) = ( lmesh%Escale(:,ke,1,1) * fx(:) + liftdelflx(:) - q(:) * ddensdxi(:,1) ) * rdens(:)
252
253 call sparsemat_matmul( dy, rhot, fy )
254 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_rhot(:,ke,2), liftdelflx )
255 df2(:,ke) = ( lmesh%Escale(:,ke,2,2) * fy(:) + liftdelflx(:) - q(:) * ddensdxi(:,2) ) * rdens(:)
256
257 call sparsemat_matmul( dz, rhot, fz )
258 call sparsemat_matmul( lift, lmesh%Fscale(:,ke) * del_flux_rhot(:,ke,3), liftdelflx )
259 df3(:,ke) = ( lmesh%Escale(:,ke,3,3) * fz(:) + liftdelflx(:) - q(:) * ddensdxi(:,3) ) * rdens(:)
260
261
262 ! Calculate the component of strain velocity tensor
263 s11(:) = dveldxi(:,1,1)
264 s12(:) = 0.5_rp * ( dveldxi(:,1,2) + dveldxi(:,2,1) )
265 s22(:) = dveldxi(:,2,2)
266 s23(:) = 0.5_rp * ( dveldxi(:,2,3) + dveldxi(:,3,2) )
267 s31(:) = 0.5_rp * ( dveldxi(:,1,3) + dveldxi(:,3,1) )
268 s33(:) = dveldxi(:,3,3)
269
270 ! Caclulate eddy viscosity & eddy diffusivity
271
272 nu(:,ke) = dns_nu
273 kh(:,ke) = dns_mu
274
275 !---
276
277 do p=1, elem%Np
278 skkovthree = ( s11(p) + s22(p) + s33(p) ) * oneoverthree
279 coef = 2.0_rp * nu(p,ke)
280
281 t11(p,ke) = dens(p) * ( coef * ( s11(p) - skkovthree ) )
282 t12(p,ke) = dens(p) * coef * s12(p)
283 t13(p,ke) = dens(p) * coef * s31(p)
284
285 t21(p,ke) = dens(p) * coef * s12(p)
286 t22(p,ke) = dens(p) * ( coef * ( s22(p) - skkovthree ) )
287 t23(p,ke) = dens(p) * coef * s23(p)
288
289 t31(p,ke) = dens(p) * coef * s31(p)
290 t32(p,ke) = dens(p) * coef * s23(p)
291 t33(p,ke) = dens(p) * ( coef * ( s33(p) - skkovthree ) )
292 end do
293
294 df1(:,ke) = kh(:,ke) * df1(:,ke)
295 df2(:,ke) = kh(:,ke) * df2(:,ke)
296 df3(:,ke) = kh(:,ke) * df3(:,ke)
297 end do
298
299 return