scale_atm_dyn_dgm_nonhydro3d_rhot_heve_splitform Module Reference

module Atmosphere / Dynamics HEVE More...

Functions/Subroutines
subroutine, public	atm_dyn_dgm_nonhydro3d_rhot_heve_splitform_init (mesh)
subroutine, public	atm_dyn_dgm_nonhydro3d_rhot_heve_splitform_final ()
subroutine, public	atm_dyn_dgm_nonhydro3d_rhot_heve_splitform_cal_tend (dens_dt, momx_dt, momy_dt, momz_dt, rhot_dt, ddens_, momx_, momy_, momz_, drhot_, 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)

Detailed Description

module Atmosphere / Dynamics HEVE

Description

HEVE DGM scheme for Atmospheric dynamical process. To improve the numerical instability due to the aliasing errors, the split form based on Gassner et al. (2016, JCP) is used for advection terms.

Author

Yuta Kawai, Team SCALE

Function/Subroutine Documentation

atm_dyn_dgm_nonhydro3d_rhot_heve_splitform_init()

subroutine, public scale_atm_dyn_dgm_nonhydro3d_rhot_heve_splitform::atm_dyn_dgm_nonhydro3d_rhot_heve_splitform_init

(

class(meshbase3d), intent(in)

mesh

)

Definition at line 79 of file scale_atm_dyn_dgm_nonhydro3d_rhot_heve_splitform.F90.

    implicit none
    class(MeshBase3D), intent(in) :: mesh
 
    integer :: p1, p2, p3, p_
    type(ElementBase3D), pointer :: elem
    !--------------------------------------------
 
    call atm_dyn_dgm_nonhydro3d_common_init( mesh )
    elem => mesh%refElem3D
 
    allocate( dxt1d_(elem%Nnode_h1D,elem%Nnode_h1D) )
    allocate( dyt1d_(elem%Nnode_h1D,elem%Nnode_h1D) )
    allocate( dzt1d_(elem%Nnode_v,elem%Nnode_v) )
 
    do p1=1, elem%Nnode_h1D
      dxt1d_(:,p1) = elem%Dx1(p1,1:elem%Nnode_h1D)
    end do
 
    do p2=1, elem%Nnode_h1D
    do p_=1, elem%Nnode_h1D
      dyt1d_(p_,p2) = elem%Dx2(1+(p2-1)*elem%Nnode_h1D,1+(p_-1)*elem%Nnode_h1D)
    end do
    end do
 
    do p3=1, elem%Nnode_v
      dzt1d_(:,p3) = elem%Dx3(elem%Colmask(p3,1),elem%Colmask(:,1))
    end do
 
    return

References scale_atm_dyn_dgm_nonhydro3d_common::atm_dyn_dgm_nonhydro3d_common_init().

atm_dyn_dgm_nonhydro3d_rhot_heve_splitform_final()

subroutine, public scale_atm_dyn_dgm_nonhydro3d_rhot_heve_splitform::atm_dyn_dgm_nonhydro3d_rhot_heve_splitform_final

Definition at line 112 of file scale_atm_dyn_dgm_nonhydro3d_rhot_heve_splitform.F90.

    implicit none
    !--------------------------------------------
 
    deallocate( dxt1d_, dyt1d_, dzt1d_ )
    call atm_dyn_dgm_nonhydro3d_common_final()
    
    return

References scale_atm_dyn_dgm_nonhydro3d_common::atm_dyn_dgm_nonhydro3d_common_final().

atm_dyn_dgm_nonhydro3d_rhot_heve_splitform_cal_tend()

subroutine, public scale_atm_dyn_dgm_nonhydro3d_rhot_heve_splitform::atm_dyn_dgm_nonhydro3d_rhot_heve_splitform_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)	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 124 of file scale_atm_dyn_dgm_nonhydro3d_rhot_heve_splitform.F90.

 
    use scale_atm_dyn_dgm_nonhydro3d_rhot_heve_numflux, only: &
      get_ebnd_flux => atm_dyn_dgm_nonhydro3d_rhot_heve_numflux_get_generalvc_asis
    
    implicit none
 
    class(LocalMesh3D), intent(in) :: lmesh
    class(ElementBase3D), intent(in) :: elem
    class(LocalMesh2D), intent(in) :: lmesh2D
    class(ElementBase2D), intent(in) :: elem2D
    class(ElementOperationBase3D), intent(in) :: element3D_operation
    type(SparseMat), intent(in) :: Dx, Dy, Dz, Sx, Sy, Sz, Lift
    real(RP), intent(out) :: DENS_dt(elem%Np,lmesh%NeA)
    real(RP), intent(out) :: MOMX_dt(elem%Np,lmesh%NeA)
    real(RP), intent(out) :: MOMY_dt(elem%Np,lmesh%NeA)
    real(RP), intent(out) :: MOMZ_dt(elem%Np,lmesh%NeA)
    real(RP), intent(out) :: RHOT_dt(elem%Np,lmesh%NeA)
    real(RP), intent(in)  :: DDENS_(elem%Np,lmesh%NeA)
    real(RP), intent(in)  :: MOMX_(elem%Np,lmesh%NeA)
    real(RP), intent(in)  :: MOMY_(elem%Np,lmesh%NeA)
    real(RP), intent(in)  :: MOMZ_(elem%Np,lmesh%NeA)
    real(RP), intent(in)  :: DRHOT_(elem%Np,lmesh%NeA)
    real(RP), intent(in)  :: DPRES_(elem%Np,lmesh%NeA)
    real(RP), intent(in)  :: DENS_hyd(elem%Np,lmesh%NeA)
    real(RP), intent(in)  :: PRES_hyd(elem%Np,lmesh%NeA)
    real(RP), intent(in)  :: PRES_hyd_ref(elem%Np,lmesh%NeA)
    real(RP), intent(in)  :: THERM_hyd(elem%Np,lmesh%NeA)
    real(RP), intent(in)  :: CORIOLIS(elem2D%Np,lmesh2D%NeA)
    real(RP), intent(in)  :: Rtot(elem%Np,lmesh%NeA)
    real(RP), intent(in)  :: CVtot(elem%Np,lmesh%NeA)
    real(RP), intent(in)  :: CPtot(elem%Np,lmesh%NeA)
    real(RP), intent(in) :: DPhydDx(elem%Np,lmesh%NeA)
    real(RP), intent(in) :: DPhydDy(elem%Np,lmesh%NeA)
 
    real(RP) :: Fx(elem%Np), Fy(elem%Np), Fz(elem%Np), LiftDelFlx(elem%Np)
    real(RP) :: Fx_sp(elem%Np), Fy_sp(elem%Np), Fz_sp(elem%Np)  
    real(RP) :: DPRES_hyd(elem%Np), GradPhyd_x(elem%Np), GradPhyd_y(elem%Np)
    real(RP) :: del_flux(elem%NfpTot,lmesh%Ne,PRGVAR_NUM)
    real(RP) :: del_flux_hyd(elem%NfpTot,lmesh%Ne,2)
    real(RP) :: GsqrtDens_(elem%Np), rdens_(elem%Np), RHOT_(elem%Np)
    real(RP) :: u_(elem%Np), v_(elem%Np), w_(elem%Np), wt_(elem%Np), pot_(elem%Np)
    real(RP) :: drho(elem%Np), Cori(elem%Np)
    real(RP) :: GsqrtV(elem%Np), RGsqrtV(elem%Np)
 
    integer :: ke, ke2d
 
    real(RP) :: gamm, rgamm    
    real(RP) :: rP0
    real(RP) :: RovP0, P0ovR       
    !------------------------------------------------------------------------
 
    call prof_rapstart( 'cal_dyn_tend_bndflux', 3)
    call get_ebnd_flux( &
      del_flux, del_flux_hyd,                                                 & ! (out)
      ddens_, momx_, momy_, momz_, drhot_, dpres_, dens_hyd, pres_hyd,        & ! (in)
      rtot, cvtot, cptot,                                                     & ! (in)
      lmesh%Gsqrt, lmesh%GI3(:,:,1), lmesh%GI3(:,:,2),                        & ! (in)    
      lmesh%normal_fn(:,:,1), lmesh%normal_fn(:,:,2), lmesh%normal_fn(:,:,3), & ! (in)
      lmesh%vmapM, lmesh%vmapP,                                               & ! (in)
      lmesh, elem, lmesh2d, elem2d )                                            ! (in)
    call prof_rapend( 'cal_dyn_tend_bndflux', 3)
 
    !-----
    call prof_rapstart( 'cal_dyn_tend_interior', 3)
    gamm  = cpdry / cvdry
    rgamm = cvdry / cpdry
    rp0   = 1.0_rp / pres00
    rovp0 = rdry * rp0
    p0ovr = pres00 / rdry
 
    !$omp parallel do private( ke2d, Cori,                          &
    !$omp RHOT_, GsqrtDens_, rdens_, u_, v_, w_, wt_, pot_,         &
    !$omp DPRES_hyd, GradPhyd_x, GradPhyd_y, drho,                  &
    !$omp GsqrtV, RGsqrtV,                                          &
    !$omp Fx, Fy, Fz, Fx_sp, Fy_sp, Fz_sp, LiftDelFlx               )
    do ke = lmesh%NeS, lmesh%NeE
      !--
      ke2d = lmesh%EMap3Dto2D(ke)
      cori(:) = coriolis(elem%IndexH2Dto3D(:),ke2d)
 
      gsqrtv(:)  = lmesh%Gsqrt(:,ke) / lmesh%GsqrtH(elem%IndexH2Dto3D,ke2d)
      rgsqrtv(:) = 1.0_rp / gsqrtv(:)
 
      !--
      rhot_(:) = p0ovr * (pres_hyd(:,ke) * rp0)**rgamm + drhot_(:,ke)
      ! DPRES_(:) = PRES00 * ( Rtot(:,ke) * rP0 * RHOT_(:) )**( CPtot(:,ke) / CVtot(:,ke) ) &
      !           - PRES_hyd(:,ke)
 
      gsqrtdens_(:) = lmesh%Gsqrt(:,ke) * ( ddens_(:,ke) + dens_hyd(:,ke) )
      rdens_(:) = 1.0_rp / gsqrtdens_(:)
      u_(:) = momx_(:,ke) * rdens_(:)
      v_(:) = momy_(:,ke) * rdens_(:)
      w_(:) = momz_(:,ke) * rdens_(:)
      wt_(:) = w_(:) * rgsqrtv(:) + lmesh%GI3(:,ke,1) * u_(:) + lmesh%GI3(:,ke,2) * v_(:) 
      pot_(:) = rhot_(:) * rdens_(:)
 
      ke2d = lmesh%EMap3Dto2D(ke)
      cori(:) = coriolis(elem%IndexH2Dto3D(:),ke2d)
 
      drho(:) = matmul(intrpmat_vpordm1, ddens_(:,ke))
 
      !-- Gradient hydrostatic pressure
      
      dpres_hyd(:) = pres_hyd(:,ke) - pres_hyd_ref(:,ke)
 
      call sparsemat_matmul(dx, gsqrtv(:) * dpres_hyd(:), fx)
      call sparsemat_matmul(dz, gsqrtv(:) * lmesh%GI3(:,ke,1) * dpres_hyd(:), fz)
      call sparsemat_matmul(lift, lmesh%Fscale(:,ke) * del_flux_hyd(:,ke,1), liftdelflx)
      gradphyd_x(:) = lmesh%Escale(:,ke,1,1) * fx(:) &
                    + lmesh%Escale(:,ke,3,3) * fz(:) &
                    + liftdelflx(:)
 
      call sparsemat_matmul(dy, gsqrtv(:) * dpres_hyd(:), fy)
      call sparsemat_matmul(dz, gsqrtv(:) * lmesh%GI3(:,ke,2) * dpres_hyd(:), fz)
      call sparsemat_matmul(lift, lmesh%Fscale(:,ke) * del_flux_hyd(:,ke,2), liftdelflx)
      gradphyd_y(:) = lmesh%Escale(:,ke,2,2) * fy(:) &
                    + lmesh%Escale(:,ke,3,3) * fz(:) &
                    + liftdelflx(:)
 
      !-- DENS
      call dx_ab( dxt1d_, gsqrtdens_(:),     u_, elem%Nnode_h1D, elem%Nnode_v, fx_sp )
      call dy_ab( dyt1d_, gsqrtdens_(:),     v_, elem%Nnode_h1D, elem%Nnode_v, fy_sp )
      call dz_ab( dzt1d_, gsqrtdens_(:), wt_(:), elem%Nnode_h1D, elem%Nnode_v, fz_sp )
      call sparsemat_matmul(lift, lmesh%Fscale(:,ke) * del_flux(:,ke,dens_vid), liftdelflx)
 
      dens_dt(:,ke) = - ( &
            lmesh%Escale(:,ke,1,1) * fx_sp(:) &
          + lmesh%Escale(:,ke,2,2) * fy_sp(:) &
          + lmesh%Escale(:,ke,3,3) * fz_sp(:) &
          + liftdelflx(:) )
      
      !-- MOMX
      call dx_abc( dxt1d_, gsqrtdens_, u_,  u_, elem%Nnode_h1D, elem%Nnode_v, fx_sp )
      call dy_abc( dyt1d_, gsqrtdens_, u_,  v_, elem%Nnode_h1D, elem%Nnode_v, fy_sp )
      call dz_abc( dzt1d_, gsqrtdens_, u_, wt_, elem%Nnode_h1D, elem%Nnode_v, fz_sp )
      call sparsemat_matmul(dx, lmesh%Gsqrt(:,ke) * dpres_(:,ke)              , fx)
      call sparsemat_matmul(lift, lmesh%Fscale(:,ke) * del_flux(:,ke,momx_vid), liftdelflx)
 
      momx_dt(:,ke) = &
          - ( lmesh%Escale(:,ke,1,1) * ( fx_sp(:) + fx(:) )  &
            + lmesh%Escale(:,ke,2,2) * fy_sp(:)              &
            + lmesh%Escale(:,ke,3,3) * fz_sp(:)              &
            + liftdelflx(:) ) / lmesh%Gsqrt(:,ke)            &
          - gradphyd_x(:) * rgsqrtv(:)                       &
          + cori(:) * momy_(:,ke)
 
      !-- MOMY
      call dx_abc( dxt1d_, gsqrtdens_, v_,  u_, elem%Nnode_h1D, elem%Nnode_v, fx_sp )
      call dy_abc( dyt1d_, gsqrtdens_, v_,  v_, elem%Nnode_h1D, elem%Nnode_v, fy_sp )  
      call dz_abc( dzt1d_, gsqrtdens_, v_, wt_, elem%Nnode_h1D, elem%Nnode_v, fz_sp )
      call sparsemat_matmul(dy, lmesh%Gsqrt(:,ke) * dpres_(:,ke)              , fy)
      call sparsemat_matmul(lift, lmesh%Fscale(:,ke) * del_flux(:,ke,momy_vid), liftdelflx)
 
      momy_dt(:,ke) = &
          - ( lmesh%Escale(:,ke,1,1) * fx_sp(:)             &
            + lmesh%Escale(:,ke,2,2) * ( fy_sp(:) + fy(:) ) &
            + lmesh%Escale(:,ke,3,3) * fz_sp(:)             &
            + liftdelflx(:) ) / lmesh%Gsqrt(:,ke)           &
          - gradphyd_y(:) * rgsqrtv(:)                      &
          - cori(:) * momx_(:,ke)
 
      !-- MOMZ
      call dx_abc( dxt1d_, gsqrtdens_, w_,  u_, elem%Nnode_h1D, elem%Nnode_v, fx_sp )
      call dy_abc( dyt1d_, gsqrtdens_, w_,  v_, elem%Nnode_h1D, elem%Nnode_v, fy_sp )
      call dz_abc( dzt1d_, gsqrtdens_, w_, wt_, elem%Nnode_h1D, elem%Nnode_v, fz_sp )
      call sparsemat_matmul(dz, lmesh%Gsqrt(:,ke) * rgsqrtv(:) * dpres_(:,ke) , fz)
      call sparsemat_matmul(lift, lmesh%Fscale(:,ke) * del_flux(:,ke,momz_vid), liftdelflx)
      
      momz_dt(:,ke) = - ( &
            lmesh%Escale(:,ke,1,1) * fx_sp(:)             &
          + lmesh%Escale(:,ke,2,2) * fy_sp(:)             &
          + lmesh%Escale(:,ke,3,3) * ( fz_sp(:) + fz(:) ) &
          + liftdelflx(:)                               ) &
          - grav * drho(:)  
 
      !-- RHOT
      call dx_abc( dxt1d_, gsqrtdens_, pot_,      u_, elem%Nnode_h1D, elem%Nnode_v, fx_sp )
      call dy_abc( dyt1d_, gsqrtdens_, pot_,      v_, elem%Nnode_h1D, elem%Nnode_v, fy_sp )        
      call dz_abc( dzt1d_, gsqrtdens_, pot_,  wt_(:), elem%Nnode_h1D, elem%Nnode_v, fz_sp )        
      call sparsemat_matmul(lift, lmesh%Fscale(:,ke) * del_flux(:,ke,rhot_vid), liftdelflx)
      
      rhot_dt(:,ke) =  - (  &
            lmesh%Escale(:,ke,1,1) * fx_sp(:) &
          + lmesh%Escale(:,ke,2,2) * fy_sp(:) &
          + lmesh%Escale(:,ke,3,3) * fz_sp(:) &
          + liftdelflx(:)                     ) 
    end do
    call prof_rapend( 'cal_dyn_tend_interior', 3)
 
    return

References scale_atm_dyn_dgm_nonhydro3d_rhot_heve_numflux::atm_dyn_dgm_nonhydro3d_rhot_heve_numflux_get_generalvc_asis(), and scale_atm_dyn_dgm_nonhydro3d_common::intrpmat_vpordm1.