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

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

Functions/Subroutines

subroutine, public atm_dyn_dgm_nonhydro3d_hevi_init (mesh)
 
subroutine, public atm_dyn_dgm_nonhydro3d_hevi_final ()
 
subroutine, public atm_dyn_dgm_nonhydro3d_hevi_cal_tend (dens_dt, momx_dt, momy_dt, momz_dt, rhot_dt, ddens_, momx_, momy_, momz_, drhot_, dens_hyd, pres_hyd, coriolis, gxu_, gyu_, gzu_, gxv_, gyv_, gzv_, gxw_, gyw_, gzw_, gxpt_, gypt_, gzpt_, visccoef_h, visccoef_v, diffcoef_h, diffcoef_v, dx, dy, dz, sx, sy, sz, lift, lmesh, elem, lmesh2d, elem2d)
 
subroutine, public atm_dyn_dgm_nonhydro3d_hevi_cal_grad_diffvars (gxu_, gyu_, gzu_, gxv_, gyv_, gzv_, gxw_, gyw_, gzw_, gxpt_, gypt_, gzpt_, ddens_, momx_, momy_, momz_, drhot_, dens_hyd, pres_hyd, dx, dy, dz, lift, lmesh, elem)
 
subroutine, public atm_dyn_dgm_nonhydro3d_hevi_cal_vi (dens_dt, momx_dt, momy_dt, momz_dt, rhot_dt, ddens_, momx_, momy_, momz_, drhot_, dens_hyd, pres_hyd, dz, lift, impl_fac, lmesh, elem, lmesh2d, elem2d)
 

Detailed Description

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

Description
HEVI DGM scheme for Atmospheric dynamical process in which GMRES is used for implicit solver. The governing equations is a fully compressibile nonhydrostic equations, which consist of mass, momentum, and thermodynamics (density * potential temperature conservation) equations.
Author
Yuta Kawai, Team SCALE

Function/Subroutine Documentation

◆ atm_dyn_dgm_nonhydro3d_hevi_init()

subroutine, public scale_atm_dyn_dgm_nonhydro3d_hevi_gmres::atm_dyn_dgm_nonhydro3d_hevi_init ( class(meshbase3d), intent(in) mesh)

Definition at line 93 of file scale_atm_dyn_dgm_nonhydro3d_hevi_gmres.F90.

94
95 implicit none
96 class(MeshBase3D), intent(in) :: mesh
97
98 integer :: p1, p2, p3, p_
99 type(ElementBase3D), pointer :: elem
100 real(RP) :: invV_POrdM1(mesh%refElem3D%Np,mesh%refElem3D%Np)
101 !--------------------------------------------
102
103 elem => mesh%refElem3D
104 allocate( intrpmat_vpordm1(elem%Np,elem%Np) )
105
106 invv_pordm1(:,:) = elem%invV
107 do p2=1, elem%Nnode_h1D
108 do p1=1, elem%Nnode_h1D
109 p_ = p1 + (p2-1)*elem%Nnode_h1D + (elem%Nnode_v-1)*elem%Nnode_h1D**2
110 invv_pordm1(p_,:) = 0.0_rp
111 end do
112 end do
113 intrpmat_vpordm1(:,:) = matmul(elem%V, invv_pordm1)
114
115 return

◆ atm_dyn_dgm_nonhydro3d_hevi_final()

subroutine, public scale_atm_dyn_dgm_nonhydro3d_hevi_gmres::atm_dyn_dgm_nonhydro3d_hevi_final

Definition at line 118 of file scale_atm_dyn_dgm_nonhydro3d_hevi_gmres.F90.

119 implicit none
120 !--------------------------------------------
121
122 deallocate( intrpmat_vpordm1 )
123
124 return

◆ atm_dyn_dgm_nonhydro3d_hevi_cal_tend()

subroutine, public scale_atm_dyn_dgm_nonhydro3d_hevi_gmres::atm_dyn_dgm_nonhydro3d_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) rhot_dt,
real(rp), dimension(elem%np,lmesh%nea), intent(in) ddens_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) momx_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) momy_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) momz_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) drhot_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) dens_hyd,
real(rp), dimension(elem%np,lmesh%nea), intent(in) pres_hyd,
real(rp), dimension(elem2d%np,lmesh2d%nea), intent(in) coriolis,
real(rp), dimension(elem%np,lmesh%nea), intent(in) gxu_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) gyu_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) gzu_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) gxv_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) gyv_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) gzv_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) gxw_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) gyw_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) gzw_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) gxpt_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) gypt_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) gzpt_,
real(rp), intent(in) visccoef_h,
real(rp), intent(in) visccoef_v,
real(rp), intent(in) diffcoef_h,
real(rp), intent(in) diffcoef_v,
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 129 of file scale_atm_dyn_dgm_nonhydro3d_hevi_gmres.F90.

135
136 implicit none
137
138 class(LocalMesh3D), intent(in) :: lmesh
139 class(ElementBase3D), intent(in) :: elem
140 class(LocalMesh2D), intent(in) :: lmesh2D
141 class(ElementBase2D), intent(in) :: elem2D
142 type(SparseMat), intent(in) :: Dx, Dy, Dz, Sx, Sy, Sz, Lift
143 real(RP), intent(out) :: DENS_dt(elem%Np,lmesh%NeA)
144 real(RP), intent(out) :: MOMX_dt(elem%Np,lmesh%NeA)
145 real(RP), intent(out) :: MOMY_dt(elem%Np,lmesh%NeA)
146 real(RP), intent(out) :: MOMZ_dt(elem%Np,lmesh%NeA)
147 real(RP), intent(out) :: RHOT_dt(elem%Np,lmesh%NeA)
148 real(RP), intent(in) :: DDENS_(elem%Np,lmesh%NeA)
149 real(RP), intent(in) :: MOMX_(elem%Np,lmesh%NeA)
150 real(RP), intent(in) :: MOMY_(elem%Np,lmesh%NeA)
151 real(RP), intent(in) :: MOMZ_(elem%Np,lmesh%NeA)
152 real(RP), intent(in) :: DRHOT_(elem%Np,lmesh%NeA)
153 real(RP), intent(in) :: DENS_hyd(elem%Np,lmesh%NeA)
154 real(RP), intent(in) :: PRES_hyd(elem%Np,lmesh%NeA)
155 real(RP), intent(in) :: CORIOLIS(elem2D%Np,lmesh2D%NeA)
156 real(RP), intent(in) :: GxU_(elem%Np,lmesh%NeA)
157 real(RP), intent(in) :: GyU_(elem%Np,lmesh%NeA)
158 real(RP), intent(in) :: GzU_(elem%Np,lmesh%NeA)
159 real(RP), intent(in) :: GxV_(elem%Np,lmesh%NeA)
160 real(RP), intent(in) :: GyV_(elem%Np,lmesh%NeA)
161 real(RP), intent(in) :: GzV_(elem%Np,lmesh%NeA)
162 real(RP), intent(in) :: GxW_(elem%Np,lmesh%NeA)
163 real(RP), intent(in) :: GyW_(elem%Np,lmesh%NeA)
164 real(RP), intent(in) :: GzW_(elem%Np,lmesh%NeA)
165 real(RP), intent(in) :: GxPT_(elem%Np,lmesh%NeA)
166 real(RP), intent(in) :: GyPT_(elem%Np,lmesh%NeA)
167 real(RP), intent(in) :: GzPT_(elem%Np,lmesh%NeA)
168 real(RP), intent(in) :: viscCoef_h
169 real(RP), intent(in) :: viscCoef_v
170 real(RP), intent(in) :: diffCoef_h
171 real(RP), intent(in) :: diffCoef_v
172
173 real(RP) :: Fx(elem%Np), Fy(elem%Np), Fz(elem%Np), LiftDelFlx(elem%Np)
174 real(RP) :: del_flux(elem%NfpTot,lmesh%Ne,PROG_VARS_NUM)
175 real(RP) :: dens_(elem%Np), RHOT_hyd(elem%Np), RHOT_(elem%Np), dpres_(elem%Np)
176 real(RP) :: pres_(elem%Np), u_(elem%Np), v_(elem%Np), w_(elem%Np), pot_(elem%Np)
177 real(RP) :: Cori(elem%Np)
178
179 real(RP) :: tmp(elem%Np)
180 integer :: ke, ke2d
181 real(RP) :: gamm, rgamm
182 !------------------------------------------------------------------------
183
184 call prof_rapstart( 'cal_dyn_tend_bndflux', 3)
185 call cal_del_flux_dyn( del_flux, & ! (out)
186 ddens_, momx_, momy_, momz_, drhot_, dens_hyd, pres_hyd, & ! (in)
187 gxu_, gyu_, gzu_, gxv_, gyv_, gzv_, gxw_, gyw_, gzw_, & ! (in)
188 gxpt_, gypt_, gzpt_, & ! (in)
189 visccoef_h, visccoef_v, & ! (in)
190 diffcoef_h, diffcoef_v, & ! (in)
191 lmesh%normal_fn(:,:,1), lmesh%normal_fn(:,:,2), lmesh%normal_fn(:,:,3), & ! (in)
192 lmesh%vmapM, lmesh%vmapP, & ! (in)
193 lmesh, elem ) ! (in)
194 call prof_rapend( 'cal_dyn_tend_bndflux', 3)
195
196 !-----
197 call prof_rapstart( 'cal_dyn_tend_interior', 3)
198 gamm = cpdry / cvdry
199 rgamm = cvdry / cpdry
200
201 !$omp parallel do private(RHOT_hyd,RHOT_,pres_,dpres_,dens_,u_,v_,w_,pot_,ke2d,Cori,Fx,Fy,Fz,LiftDelFlx,tmp)
202 do ke = lmesh%NeS, lmesh%NeE
203 !--
204 rhot_hyd(:) = pres00/rdry * (pres_hyd(:,ke)/pres00)**rgamm
205 rhot_(:) = rhot_hyd(:) + drhot_(:,ke)
206 pres_(:) = pres_hyd(:,ke) * (1.0_rp + drhot_(:,ke)/rhot_hyd(:))**gamm
207 dpres_(:) = pres_(:) - pres_hyd(:,ke)
208 dens_(:) = ddens_(:,ke) + dens_hyd(:,ke)
209
210 u_(:) = momx_(:,ke)/dens_(:)
211 v_(:) = momy_(:,ke)/dens_(:)
212 w_(:) = momz_(:,ke)/dens_(:)
213 pot_(:) = rhot_(:)/dens_(:)
214
215 ke2d = lmesh%EMap3Dto2D(ke)
216 cori(:) = coriolis(elem%IndexH2Dto3D(:),ke2d)
217
218 !-- DENS
219 call sparsemat_matmul(dx, momx_(:,ke), fx)
220 call sparsemat_matmul(dy, momy_(:,ke), fy)
221 call sparsemat_matmul(lift, lmesh%Fscale(:,ke)*del_flux(:,ke,ddens_vid), liftdelflx)
222
223 dens_dt(:,ke) = - ( &
224 lmesh%Escale(:,ke,1,1) * fx(:) &
225 + lmesh%Escale(:,ke,2,2) * fy(:) &
226 + liftdelflx(:) )
227
228 !-- MOMX
229 call sparsemat_matmul(dx, u_(:)*momx_(:,ke) + pres_(:) - visccoef_h*dens_(:)*gxu_(:,ke), fx)
230 call sparsemat_matmul(dy, v_(:)*momx_(:,ke) - visccoef_h*dens_(:)*gyu_(:,ke), fy)
231 call sparsemat_matmul(dz, w_(:)*momx_(:,ke) - visccoef_v*dens_(:)*gzu_(:,ke) ,fz)
232 call sparsemat_matmul(lift, lmesh%Fscale(:,ke)*del_flux(:,ke,momx_vid), liftdelflx)
233
234 momx_dt(:,ke) = - ( &
235 lmesh%Escale(:,ke,1,1) * fx(:) &
236 + lmesh%Escale(:,ke,2,2) * fy(:) &
237 + lmesh%Escale(:,ke,3,3) * fz(:) &
238 + liftdelflx(:) &
239 - cori(:)*momy_(:,ke) &
240 )
241
242 !-- MOMY
243 call sparsemat_matmul(dx, u_(:)*momy_(:,ke) - visccoef_h*dens_(:)*gxv_(:,ke), fx)
244 call sparsemat_matmul(dy, v_(:)*momy_(:,ke) + pres_(:) - visccoef_h*dens_(:)*gyv_(:,ke), fy)
245 call sparsemat_matmul(dz, w_(:)*momy_(:,ke) - visccoef_v*dens_(:)*gzv_(:,ke) ,fz)
246 call sparsemat_matmul(lift, lmesh%Fscale(:,ke)*del_flux(:,ke,momy_vid), liftdelflx)
247
248 momy_dt(:,ke) = - ( &
249 lmesh%Escale(:,ke,1,1) * fx(:) &
250 + lmesh%Escale(:,ke,2,2) * fy(:) &
251 + lmesh%Escale(:,ke,3,3) * fz(:) &
252 + liftdelflx(:) &
253 + cori(:)*momx_(:,ke) &
254 )
255
256 !-- MOMZ
257 call sparsemat_matmul(dx, u_(:)*momz_(:,ke) - visccoef_h*dens_(:)*gxw_(:,ke), fx)
258 call sparsemat_matmul(dy, v_(:)*momz_(:,ke) - visccoef_h*dens_(:)*gyw_(:,ke), fy)
259 call sparsemat_matmul(dz, w_(:)*momz_(:,ke) - visccoef_v*dens_(:)*gzw_(:,ke), fz)
260 call sparsemat_matmul(lift, lmesh%Fscale(:,ke)*del_flux(:,ke,momz_vid), liftdelflx)
261
262 momz_dt(:,ke) = - ( &
263 lmesh%Escale(:,ke,1,1) * fx(:) &
264 + lmesh%Escale(:,ke,2,2) * fy(:) &
265 + lmesh%Escale(:,ke,3,3) * fz(:) &
266 + liftdelflx(:) )
267
268 !-- RHOT
269 call sparsemat_matmul(dx, pot_(:)*momx_(:,ke) - diffcoef_h*dens_(:)*gxpt_(:,ke), fx)
270 call sparsemat_matmul(dy, pot_(:)*momy_(:,ke) - diffcoef_h*dens_(:)*gypt_(:,ke), fy)
271 call sparsemat_matmul(dz, - diffcoef_v*dens_(:)*gzpt_(:,ke), fz)
272 call sparsemat_matmul(lift, lmesh%Fscale(:,ke)*del_flux(:,ke,drhot_vid), liftdelflx)
273
274 rhot_dt(:,ke) = - ( &
275 lmesh%Escale(:,ke,1,1) * fx(:) &
276 + lmesh%Escale(:,ke,2,2) * fy(:) &
277 + lmesh%Escale(:,ke,3,3) * fz(:) &
278 + liftdelflx(:) )
279 end do
280 call prof_rapend( 'cal_dyn_tend_interior', 3)
281
282 return

◆ atm_dyn_dgm_nonhydro3d_hevi_cal_grad_diffvars()

subroutine, public scale_atm_dyn_dgm_nonhydro3d_hevi_gmres::atm_dyn_dgm_nonhydro3d_hevi_cal_grad_diffvars ( real(rp), dimension(elem%np,lmesh%nea), intent(out) gxu_,
real(rp), dimension(elem%np,lmesh%nea), intent(out) gyu_,
real(rp), dimension(elem%np,lmesh%nea), intent(out) gzu_,
real(rp), dimension(elem%np,lmesh%nea), intent(out) gxv_,
real(rp), dimension(elem%np,lmesh%nea), intent(out) gyv_,
real(rp), dimension(elem%np,lmesh%nea), intent(out) gzv_,
real(rp), dimension(elem%np,lmesh%nea), intent(out) gxw_,
real(rp), dimension(elem%np,lmesh%nea), intent(out) gyw_,
real(rp), dimension(elem%np,lmesh%nea), intent(out) gzw_,
real(rp), dimension(elem%np,lmesh%nea), intent(out) gxpt_,
real(rp), dimension(elem%np,lmesh%nea), intent(out) gypt_,
real(rp), dimension(elem%np,lmesh%nea), intent(out) gzpt_,
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,
type(sparsemat), intent(in) dx,
type(sparsemat), intent(in) dy,
type(sparsemat), intent(in) dz,
type(sparsemat), intent(in) lift,
class(localmesh3d), intent(in) lmesh,
class(elementbase3d), intent(in) elem )

Definition at line 632 of file scale_atm_dyn_dgm_nonhydro3d_hevi_gmres.F90.

636
637 implicit none
638
639 class(LocalMesh3D), intent(in) :: lmesh
640 class(ElementBase3D), intent(in) :: elem
641 type(SparseMat), intent(in) :: Dx, Dy, Dz, Lift
642 real(RP), intent(out) :: GxU_(elem%Np,lmesh%NeA)
643 real(RP), intent(out) :: GyU_(elem%Np,lmesh%NeA)
644 real(RP), intent(out) :: GzU_(elem%Np,lmesh%NeA)
645 real(RP), intent(out) :: GxV_(elem%Np,lmesh%NeA)
646 real(RP), intent(out) :: GyV_(elem%Np,lmesh%NeA)
647 real(RP), intent(out) :: GzV_(elem%Np,lmesh%NeA)
648 real(RP), intent(out) :: GxW_(elem%Np,lmesh%NeA)
649 real(RP), intent(out) :: GyW_(elem%Np,lmesh%NeA)
650 real(RP), intent(out) :: GzW_(elem%Np,lmesh%NeA)
651 real(RP), intent(out) :: GxPT_(elem%Np,lmesh%NeA)
652 real(RP), intent(out) :: GyPT_(elem%Np,lmesh%NeA)
653 real(RP), intent(out) :: GzPT_(elem%Np,lmesh%NeA)
654 real(RP), intent(in) :: DDENS_(elem%Np,lmesh%NeA)
655 real(RP), intent(in) :: MOMX_(elem%Np,lmesh%NeA)
656 real(RP), intent(in) :: MOMY_(elem%Np,lmesh%NeA)
657 real(RP), intent(in) :: MOMZ_(elem%Np,lmesh%NeA)
658 real(RP), intent(in) :: DRHOT_(elem%Np,lmesh%NeA)
659 real(RP), intent(in) :: DENS_hyd(elem%Np,lmesh%NeA)
660 real(RP), intent(in) :: PRES_hyd(elem%Np,lmesh%NeA)
661
662 real(RP) :: DENS_(elem%Np), U_(elem%Np), V_(elem%Np), W_(elem%Np)
663 real(RP) :: DTHETA_(elem%Np), RHOT_(elem%Np), RHOT_hyd(elem%Np)
664 real(RP) :: Fx(elem%Np), Fy(elem%Np), Fz(elem%Np), LiftDelFlx(elem%Np)
665 real(RP) :: del_flux(elem%NfpTot,lmesh%Ne,AUX_DIFFVARS_NUM)
666
667 integer :: ke
668 !------------------------------------------------------------------------------
669
670 call cal_del_graddiffvar( del_flux, & ! (out)
671 ddens_, momx_, momy_, momz_, drhot_, dens_hyd, pres_hyd, & ! (in)
672 lmesh%normal_fn(:,:,1), lmesh%normal_fn(:,:,2), lmesh%normal_fn(:,:,3), & ! (in)
673 lmesh%vmapM, lmesh%vmapP, & ! (in)
674 lmesh, elem ) ! (in)
675
676 do ke=lmesh%NeS, lmesh%NeE
677 dens_(:) = ddens_(:,ke) + dens_hyd(:,ke)
678 rhot_hyd(:) = pres00/rdry * (pres_hyd(:,ke)/pres00)**(cvdry/cpdry)
679 rhot_(:) = rhot_hyd(:) + drhot_(:,ke)
680
681 u_(:) = momx_(:,ke)/dens_(:)
682 v_(:) = momy_(:,ke)/dens_(:)
683 w_(:) = momz_(:,ke)/dens_(:)
684 dtheta_(:) = rhot_(:)/dens_(:) - rhot_hyd(:)/dens_hyd(:,ke)
685
686 !- U
687
688 call sparsemat_matmul(dx, u_, fx)
689 call sparsemat_matmul(lift, lmesh%Fscale(:,ke)*del_flux(:,ke,vars_gxu_id), liftdelflx)
690 gxu_(:,ke) = lmesh%Escale(:,ke,1,1)*fx(:) + liftdelflx(:)
691
692 call sparsemat_matmul(dy, u_, fy)
693 call sparsemat_matmul(lift, lmesh%Fscale(:,ke)*del_flux(:,ke,vars_gyu_id), liftdelflx)
694 gyu_(:,ke) = lmesh%Escale(:,ke,2,2)*fx(:) + liftdelflx(:)
695
696 call sparsemat_matmul(dz, u_, fz)
697 call sparsemat_matmul(lift, lmesh%Fscale(:,ke)*del_flux(:,ke,vars_gzu_id), liftdelflx)
698 gzu_(:,ke) = lmesh%Escale(:,ke,3,3)*fz(:) + liftdelflx(:)
699
700 !- V
701
702 call sparsemat_matmul(dx, v_, fx)
703 call sparsemat_matmul(lift, lmesh%Fscale(:,ke)*del_flux(:,ke,vars_gxv_id), liftdelflx)
704 gxv_(:,ke) = lmesh%Escale(:,ke,1,1)*fx(:) + liftdelflx(:)
705
706 call sparsemat_matmul(dy, v_, fy)
707 call sparsemat_matmul(lift, lmesh%Fscale(:,ke)*del_flux(:,ke,vars_gyv_id), liftdelflx)
708 gyv_(:,ke) = lmesh%Escale(:,ke,2,2)*fx(:) + liftdelflx(:)
709
710 call sparsemat_matmul(dz, v_, fz)
711 call sparsemat_matmul(lift, lmesh%Fscale(:,ke)*del_flux(:,ke,vars_gzv_id), liftdelflx)
712 gzv_(:,ke) = lmesh%Escale(:,ke,3,3)*fz(:) + liftdelflx(:)
713
714 !- W
715
716 call sparsemat_matmul(dx, w_, fx)
717 call sparsemat_matmul(lift, lmesh%Fscale(:,ke)*del_flux(:,ke,vars_gxw_id), liftdelflx)
718 gxw_(:,ke) = lmesh%Escale(:,ke,1,1)*fx(:) + liftdelflx(:)
719
720 call sparsemat_matmul(dy, w_, fx)
721 call sparsemat_matmul(lift, lmesh%Fscale(:,ke)*del_flux(:,ke,vars_gyw_id), liftdelflx)
722 gyw_(:,ke) = lmesh%Escale(:,ke,2,2)*fy(:) + liftdelflx(:)
723
724 call sparsemat_matmul(dz, w_, fz)
725 call sparsemat_matmul(lift, lmesh%Fscale(:,ke)*del_flux(:,ke,vars_gzw_id), liftdelflx)
726 gzw_(:,ke) = lmesh%Escale(:,ke,3,3)*fz(:) + liftdelflx(:)
727
728 !- PT
729
730 call sparsemat_matmul(dx, dtheta_, fx)
731 call sparsemat_matmul(lift, lmesh%Fscale(:,ke)*del_flux(:,ke,vars_gxpt_id), liftdelflx)
732 gxpt_(:,ke) = lmesh%Escale(:,ke,1,1)*fx(:) + liftdelflx(:)
733
734 call sparsemat_matmul(dy, dtheta_, fy)
735 call sparsemat_matmul(lift, lmesh%Fscale(:,ke)*del_flux(:,ke,vars_gypt_id), liftdelflx)
736 gypt_(:,ke) = lmesh%Escale(:,ke,2,2)*fz(:) + liftdelflx(:)
737
738 call sparsemat_matmul(dz, dtheta_, fz)
739 call sparsemat_matmul(lift, lmesh%Fscale(:,ke)*del_flux(:,ke,vars_gzpt_id), liftdelflx)
740 gzpt_(:,ke) = lmesh%Escale(:,ke,3,3)*fz(:) + liftdelflx(:)
741
742 end do
743
744 return

◆ atm_dyn_dgm_nonhydro3d_hevi_cal_vi()

subroutine, public scale_atm_dyn_dgm_nonhydro3d_hevi_gmres::atm_dyn_dgm_nonhydro3d_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) rhot_dt,
real(rp), dimension(elem%np,lmesh%nea), intent(in) ddens_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) momx_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) momy_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) momz_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) drhot_,
real(rp), dimension(elem%np,lmesh%nea), intent(in) dens_hyd,
real(rp), dimension(elem%np,lmesh%nea), intent(in) pres_hyd,
class(sparsemat), intent(in) dz,
class(sparsemat), intent(in) lift,
real(rp), intent(in) impl_fac,
class(localmesh3d), intent(in) lmesh,
class(elementbase3d), intent(in) elem,
class(localmesh2d), intent(in) lmesh2d,
class(elementbase2d), intent(in) elem2d )

Definition at line 818 of file scale_atm_dyn_dgm_nonhydro3d_hevi_gmres.F90.

822
823 implicit none
824
825 class(LocalMesh3D), intent(in) :: lmesh
826 class(ElementBase3D), intent(in) :: elem
827 class(LocalMesh2D), intent(in) :: lmesh2D
828 class(ElementBase2D), intent(in) :: elem2D
829 real(RP), intent(out) :: DENS_dt(elem%Np,lmesh%NeA)
830 real(RP), intent(out) :: MOMX_dt(elem%Np,lmesh%NeA)
831 real(RP), intent(out) :: MOMY_dt(elem%Np,lmesh%NeA)
832 real(RP), intent(out) :: MOMZ_dt(elem%Np,lmesh%NeA)
833 real(RP), intent(out) :: RHOT_dt(elem%Np,lmesh%NeA)
834 real(RP), intent(in) :: DDENS_(elem%Np,lmesh%NeA)
835 real(RP), intent(in) :: MOMX_(elem%Np,lmesh%NeA)
836 real(RP), intent(in) :: MOMY_(elem%Np,lmesh%NeA)
837 real(RP), intent(in) :: MOMZ_(elem%Np,lmesh%NeA)
838 real(RP), intent(in) :: DRHOT_(elem%Np,lmesh%NeA)
839 real(RP), intent(in) :: DENS_hyd(elem%Np,lmesh%NeA)
840 real(RP), intent(in) :: PRES_hyd(elem%Np,lmesh%NeA)
841 class(SparseMat), intent(in) :: Dz, Lift
842 real(RP), intent(in) :: impl_fac
843
844 real(RP) :: PROG_VARS(elem%Np,PROG_VARS_NUM,lmesh%NeZ)
845 real(RP) :: PROG_VARS0(elem%Np,PROG_VARS_NUM,lmesh%NeZ)
846 real(RP) :: PROG_VARS00(elem%Np,PROG_VARS_NUM,lmesh%NeZ)
847 real(RP) :: PROG_DEL(elem%Np,PROG_VARS_NUM,lmesh%NeZ)
848 real(RP) :: b(elem%Np,PROG_VARS_NUM,lmesh%NeZ)
849 real(RP) :: Ax(elem%Np,PROG_VARS_NUM,lmesh%NeZ)
850 real(RP) :: tend(elem%Np,PROG_VARS_NUM,lmesh%NeZ)
851 real(RP) :: DENS_hyd_z(elem%Np,lmesh%NeZ)
852 real(RP) :: PRES_hyd_z(elem%Np,lmesh%NeZ)
853 real(RP) :: nz(elem%NfpTot,lmesh%NeZ)
854 integer :: vmapM(elem%NfpTot,lmesh%NeZ)
855 integer :: vmapP(elem%NfpTot,lmesh%NeZ)
856 integer :: ke_x, ke_y, ke_z, ke, p, v
857 integer :: itr_lin, itr_nlin
858 integer :: m, N
859 integer :: f, vs, ve
860 logical :: is_converged
861
862
863 type(GMRES) :: gmres_hevi
864 real(RP), allocatable :: wj(:), pinv_v(:)
865 real(RP) :: PmatDlu(elem%Np*PROG_VARS_NUM,elem%Np*PROG_VARS_NUM,lmesh%NeZ)
866 integer :: PmatDlu_ipiv(elem%Np*PROG_VARS_NUM,lmesh%NeZ)
867 real(RP) :: PmatL(elem%Np,elem%Np,PROG_VARS_NUM,PROG_VARS_NUM,lmesh%NeZ)
868 real(RP) :: PmatU(elem%Np,elem%Np,PROG_VARS_NUM,PROG_VARS_NUM,lmesh%NeZ)
869 real(RP), parameter :: EPS0 = 1.0e-16_rp
870 real(RP), parameter :: EPS = 1.0e-16_rp
871
872 !------------------------------------------------------------------------
873
874 call prof_rapstart( 'hevi_cal_vi', 3)
875 call prof_rapstart( 'hevi_cal_vi_prep', 3)
876
877 n = elem%Np * prog_vars_num * lmesh%NeZ
878 m = n / prog_vars_num / elem%Nnode_h1D**2
879
880 call gmres_hevi%Init( n, m, eps, eps )
881 allocate( wj(n), pinv_v(n) )
882
883 do ke_z=1, lmesh%NeZ
884 do f=1, elem%Nfaces_h
885 vs = 1 + (f-1)*elem%Nfp_h
886 ve = vs + elem%Nfp_h - 1
887 vmapm(vs:ve,ke_z) = elem%Fmask_h(:,f) + (ke_z-1)*elem%Np
888 end do
889 do f=1, elem%Nfaces_v
890 vs = elem%Nfp_h*elem%Nfaces_h + 1 + (f-1)*elem%Nfp_v
891 ve = vs + elem%Nfp_v - 1
892 vmapm(vs:ve,ke_z) = elem%Fmask_v(:,f) + (ke_z-1)*elem%Np
893 end do
894 vmapp(:,ke_z) = vmapm(:,ke_z)
895 end do
896
897 do ke_z=1, lmesh%NeZ
898 vs = elem%Nfp_h*elem%Nfaces_h + 1
899 ve = vs + elem%Nfp_v - 1
900 if (ke_z > 1) &
901 vmapp(vs:ve,ke_z) = elem%Fmask_v(:,2) + (ke_z-2)*elem%Np
902 ! if (ke_z > 1) then
903 ! vmapP(vs:ve,ke_z) = elem%Fmask_v(:,2) + (ke_z-2)*elem%Np
904 ! else
905 ! vmapP(vs:ve,ke_z) = elem%Fmask_v(:,2) + (lmesh%NeZ-1)*elem%Np
906 ! end if
907
908 vs = elem%Nfp_h*elem%Nfaces_h + elem%Nfp_v + 1
909 ve = vs + elem%Nfp_v - 1
910 if (ke_z < lmesh%NeZ) &
911 vmapp(vs:ve,ke_z) = elem%Fmask_v(:,1) + ke_z*elem%Np
912 ! if (ke_z < lmesh%NeZ) then
913 ! vmapP(vs:ve,ke_z) = elem%Fmask_v(:,1) + ke_z*elem%Np
914 ! else
915 ! vmapP(vs:ve,ke_z) = elem%Fmask_v(:,1)
916 ! end if
917 end do
918 call prof_rapend( 'hevi_cal_vi_prep', 3)
919
920 do ke_y=1, lmesh%NeY
921 do ke_x=1, lmesh%NeX
922
923 call prof_rapstart( 'hevi_cal_vi_get_var', 3)
924
925 do ke_z=1, lmesh%NeZ
926 ke = ke_x + (ke_y-1)*lmesh%NeX + (ke_z-1)*lmesh%NeX*lmesh%NeY
927
928 prog_vars(:,ddens_vid,ke_z) = ddens_(:,ke)
929 prog_vars(:,momx_vid,ke_z) = momx_(:,ke)
930 prog_vars(:,momy_vid,ke_z) = momy_(:,ke)
931 prog_vars(:,momz_vid,ke_z) = momz_(:,ke)
932 prog_vars(:,drhot_vid,ke_z) = drhot_(:,ke)
933 dens_hyd_z(:,ke_z) = dens_hyd(:,ke)
934 pres_hyd_z(:,ke_z) = pres_hyd(:,ke)
935
936 prog_vars0(:,:,ke_z) = prog_vars(:,:,ke_z)
937 prog_vars00(:,:,ke_z) = prog_vars(:,:,ke_z)
938 nz(:,ke_z) = lmesh%normal_fn(:,ke,3)
939 end do
940 call prof_rapend( 'hevi_cal_vi_get_var', 3)
941
942 if ( impl_fac > 0.0_rp ) then
943 call prof_rapstart( 'hevi_cal_vi_itr', 3)
944
945 do itr_nlin = 1, 1
946
947 do ke_z=1, lmesh%NeZ
948 prog_del(:,:,ke_z) = 0.0_rp
949 end do
950
951 call vi_eval_ax( ax(:,:,:), & ! (out)
952 prog_vars, prog_vars0, dens_hyd_z, pres_hyd_z, & ! (in)
953 dz, lift, impl_fac, lmesh, elem, & ! (in)
954 nz, vmapm, vmapp, ke_x, ke_y, .false. )
955
956 do ke_z=1, lmesh%NeZ
957 b(:,:,ke_z) = - ax(:,:,ke_z) + prog_vars00(:,:,ke_z)
958 end do
959
960 call prof_rapstart( 'hevi_cal_vi_pmatinv', 3)
961
962 call vi_construct_pmatinv( pmatdlu, pmatdlu_ipiv, pmatl, pmatu, & ! (out)
963 prog_vars0, dens_hyd_z, pres_hyd_z, & ! (in)
964 dz, lift, impl_fac, lmesh, elem, & ! (in)
965 nz, vmapm, vmapp, ke_x, ke_y )
966
967 call prof_rapend( 'hevi_cal_vi_pmatinv', 3)
968
969 call prof_rapstart( 'hevi_cal_vi_itr_lin', 3)
970 do itr_lin=1, 2*int(n/m)
971
972 call vi_gmres_core( gmres_hevi, prog_del(:,:,:), wj, & ! (inout)
973 is_converged, & ! (out)
974 prog_vars(:,:,:), b(:,:,:), n, m, pmatdlu, pmatdlu_ipiv, pmatl, pmatu, pinv_v, & ! (in)
975 dens_hyd_z, pres_hyd_z, & ! (in)
976 dz, lift, impl_fac, lmesh, elem, & ! (in)
977 nz, vmapm, vmapp, ke_x, ke_y )
978
979 if (is_converged) exit
980 end do ! itr lin
981! write(*,*) lmesh%tileID, ke_x, ke_y, itr_lin, gmres_hevi%m_out
982 call prof_rapend( 'hevi_cal_vi_itr_lin', 3)
983
984 do ke_z=1, lmesh%NeZ
985 prog_vars(:,:,ke_z) = prog_vars(:,:,ke_z) + prog_del(:,:,ke_z)
986 prog_vars0(:,:,ke_z) = prog_vars(:,:,ke_z)
987 end do
988 end do ! itr nlin
989
990 call prof_rapend( 'hevi_cal_vi_itr', 3)
991 end if
992
993 call prof_rapstart( 'hevi_cal_vi_retrun_var', 3)
994 call vi_eval_ax( tend(:,:,:), & ! (out)
995 prog_vars, prog_vars, dens_hyd_z, pres_hyd_z, & ! (in)
996 dz, lift, impl_fac, lmesh, elem, & ! (in)
997 nz, vmapm, vmapp, ke_x, ke_y, .true. )
998
999 !tend(:,:,:) = 0.0_RP
1000 do ke_z=1, lmesh%NeZ
1001 ke = ke_x + (ke_y-1)*lmesh%NeX + (ke_z-1)*lmesh%NeX*lmesh%NeY
1002 dens_dt(:,ke) = - tend(:,ddens_vid,ke_z)
1003 momx_dt(:,ke) = - tend(:,momx_vid,ke_z)
1004 momy_dt(:,ke) = - tend(:,momy_vid,ke_z)
1005 momz_dt(:,ke) = - tend(:,momz_vid,ke_z)
1006 rhot_dt(:,ke) = - tend(:,drhot_vid,ke_z)
1007 end do
1008 !tend(:,:,:) = 0.0_RP
1009 call prof_rapend( 'hevi_cal_vi_retrun_var', 3)
1010 end do
1011 end do
1012
1013 call gmres_hevi%Final()
1014
1015 call prof_rapend( 'hevi_cal_vi', 3)
1016
1017 return

References get_pmatd_lu(), and scale_linalgebra::linalgebra_lu().

Generated by doxygen 1.10.0