Skip to content

Commit

Permalink
module_mp_nssl_2mom.F90: fix bug when nz > 128 where sedimentation d…
Browse files Browse the repository at this point in the history
…id not work for k > 128
  • Loading branch information
MicroTed committed Oct 29, 2023
1 parent 03772d6 commit 1abaff0
Showing 1 changed file with 115 additions and 37 deletions.
152 changes: 115 additions & 37 deletions physics/module_mp_nssl_2mom.F90
Original file line number Diff line number Diff line change
Expand Up @@ -4274,55 +4274,89 @@ subroutine sediment1d(dtp,nx,ny,nz,an,na,nor,norz,xfall,dn,dz3d,dz3dinv, &
! real gz(-nor+ng1:nz+nor),z1d(-nor+ng1:nz+nor,4)
real dtp
real xfall(nx,ny,na) ! array for stuff landing on the ground
real xfall0(nx,ny) ! dummy array
! real xfall0(nx,ny) ! dummy array
integer infdo
integer jslab ! which line of xfall to use

integer ix,jy,kz,ndfall,n,k,il,in
real tmp, vtmax, dtptmp, dtfrac
real, parameter :: dz = 200.

real :: xvt(nz+1,nx,3,lc:lhab) ! (nx,nz,2,lc:lhab) ! 1=mass-weighted, 2=number-weighted
real :: tmpn(-nor+ng1:nx+nor,-nor+ng1:ny+nor,-norz+ng1:nz+norz)
real :: tmpn2(-nor+ng1:nx+nor,-nor+ng1:ny+nor,-norz+ng1:nz+norz)
real :: z(-nor+ng1:nx+nor,-norz+ng1:nz+norz,lr:lhab)
real :: db1(nx,nz+1),dtz1(nz+1,nx,0:1),dz2dinv(nz+1,nx),db1inv(nx,nz+1)
! real :: xvt(nz+1,nx,3,lc:lhab) ! (nx,nz,2,lc:lhab) ! 1=mass-weighted, 2=number-weighted
! real :: tmpn(-nor+ng1:nx+nor,-nor+ng1:ny+nor,-norz+ng1:nz+norz)
! real :: tmpn2(-nor+ng1:nx+nor,-nor+ng1:ny+nor,-norz+ng1:nz+norz)
! real :: z(-nor+ng1:nx+nor,-norz+ng1:nz+norz,lr:lhab)
! real :: db1(nx,nz+1),dtz1(nz+1,nx,0:1),dz2dinv(nz+1,nx),db1inv(nx,nz+1)

real :: rhovtzx(nz,nx)
! real :: rhovtzx(nz,nx)

real, allocatable :: db1(:,:), dtz1(:,:,:),dz2dinv(:,:),db1inv(:,:) ! db1(nx,nz+1),dtz1(nz+1,nx,0:1),dz2dinv(nz+1,nx),db1inv(nx,nz+1)
real, allocatable :: rhovtzx(:,:)
real, allocatable :: xfall0(:,:), xvt(:,:,:,:),tmpn(:,:,:),tmpn2(:,:,:),z(:,:,:)

double precision :: timesed1,timesed2,timesed3, zmaxsed,timesetvt,dummy
double precision :: dt1,dt2,dt3,dt4

integer,parameter :: ngs = 128
integer :: ngs ! = 512
integer :: ngscnt,mgs,ipconc0

real :: qx(ngs,lv:lhab)
real :: qxw(ngs,ls:lhab)
real :: cx(ngs,lc:lhab)
real :: xv(ngs,lc:lhab)
real :: vtxbar(ngs,lc:lhab,3)
real :: xmas(ngs,lc:lhab)
real :: xdn(ngs,lc:lhab)
real :: xdia(ngs,lc:lhab,3)
real :: vx(ngs,li:lhab)
real :: alpha(ngs,lc:lhab)
real :: zx(ngs,lr:lhab)
logical :: hasmass(nx,lc+1:lhab)

integer igs(ngs),kgs(ngs)

real rho0(ngs),temcg(ngs)

real temg(ngs)

real rhovt(ngs)

real cwnc(ngs),cinc(ngs)
real fadvisc(ngs),cwdia(ngs),cipmas(ngs)

real cimasn,cimasx,cnina(ngs),cimas(ngs)

real cnostmp(ngs)
! real :: qx(ngs,lv:lhab)
! real :: qxw(ngs,ls:lhab)
! real :: cx(ngs,lc:lhab)
! real :: xv(ngs,lc:lhab)
! real :: vtxbar(ngs,lc:lhab,3)
! real :: xmas(ngs,lc:lhab)
! real :: xdn(ngs,lc:lhab)
! real :: xdia(ngs,lc:lhab,3)
! real :: vx(ngs,li:lhab)
! real :: alpha(ngs,lc:lhab)
! real :: zx(ngs,lr:lhab)
! logical :: hasmass(nx,lc+1:lhab)
!
! integer igs(ngs),kgs(ngs)
!
! real rho0(ngs),temcg(ngs)
!
! real temg(ngs)
!
! real rhovt(ngs)
!
! real cwnc(ngs),cinc(ngs)
! real fadvisc(ngs),cwdia(ngs),cipmas(ngs)
!
! real cimasn,cimasx,cnina(ngs),cimas(ngs)
!
! real cnostmp(ngs)

real, allocatable :: qx(:,:)
real, allocatable :: qxw(:,:)
real, allocatable :: cx(:,:)
real, allocatable :: xv(:,:)
real, allocatable :: vtxbar(:,:,:)
real, allocatable :: xmas(:,:)
real, allocatable :: xdn(:,:)
real, allocatable :: xdia(:,:,:)
real, allocatable :: vx(:,:)
real, allocatable :: alpha(:,:)
real, allocatable :: zx(:,:)
logical, allocatable :: hasmass(:,:)

integer, allocatable :: igs(:),kgs(:)

real, allocatable :: rho0(:),temcg(:)

real, allocatable :: temg(:)

real, allocatable :: rhovt(:)

real, allocatable :: cwnc(:),cinc(:)
real, allocatable :: fadvisc(:),cwdia(:),cipmas(:)

real, allocatable :: cnina(:),cimas(:)

real, allocatable :: cnostmp(:)

real :: cimasn,cimasx


!-----------------------------------------------------------------------------
Expand All @@ -4336,7 +4370,30 @@ subroutine sediment1d(dtp,nx,ny,nz,an,na,nor,norz,xfall,dn,dz3d,dz3dinv, &
! ###################################################################



allocate( db1(nx,nz+1),dtz1(nz+1,nx,0:1),dz2dinv(nz+1,nx),db1inv(nx,nz+1),rhovtzx(nz,nx) )
allocate( xfall0(nx,ny), xvt(nz+1,nx,3,lc:lhab), tmpn(-nor+ng1:nx+nor,-nor+ng1:ny+nor,-norz+ng1:nz+norz) )
allocate( tmpn2(-nor+ng1:nx+nor,-nor+ng1:ny+nor,-norz+ng1:nz+norz), z(-nor+ng1:nx+nor,-norz+ng1:nz+norz,lr:lhab))

ngs = nz+3

allocate( qx(ngs,lv:lhab), &
qxw(ngs,ls:lhab), &
cx(ngs,lc:lhab), &
xv(ngs,lc:lhab), &
vtxbar(ngs,lc:lhab,3), &
xmas(ngs,lc:lhab), &
xdn(ngs,lc:lhab), &
xdia(ngs,lc:lhab,3), &
vx(ngs,li:lhab), &
alpha(ngs,lc:lhab), &
zx(ngs,lr:lhab), &
hasmass(nx,lc+1:lhab), &
igs(ngs),kgs(ngs), &
rho0(ngs),temcg(ngs),temg(ngs), rhovt(ngs), &
cwnc(ngs),cinc(ngs), &
fadvisc(ngs),cwdia(ngs),cipmas(ngs), &
cnina(ngs),cimas(ngs), &
cnostmp(ngs) )

kzb = 1
kze = nz
Expand Down Expand Up @@ -4656,8 +4713,29 @@ subroutine sediment1d(dtp,nx,ny,nz,an,na,nor,norz,xfall,dn,dz3d,dz3dinv, &
ENDDO ! ix


deallocate( db1,dtz1,dz2dinv,db1inv,rhovtzx )
deallocate( xfall0, xvt, tmpn )
deallocate( tmpn2, z)

deallocate( qx, &
qxw, &
cx, &
xv, &
vtxbar, &
xmas, &
xdn, &
xdia, &
vx, &
alpha, &
zx, &
hasmass, &
igs,kgs, &
rho0,temcg,temg, rhovt, &
cwnc,cinc, &
fadvisc,cwdia,cipmas, &
cnina,cimas, &
cnostmp )


RETURN
END SUBROUTINE SEDIMENT1D

Expand Down

0 comments on commit 1abaff0

Please sign in to comment.