Commit fd541384 authored by Espen Sollum's avatar Espen Sollum
Browse files

CTBTO vtables code merged with latest dev. Some other minor edits

parent 0eda0087
...@@ -138,7 +138,7 @@ program flexpart ...@@ -138,7 +138,7 @@ program flexpart
if (verbosity.gt.0) then if (verbosity.gt.0) then
write(*,*) 'call readpaths' write(*,*) 'call readpaths'
endif endif
call readpaths(pathfile) call readpaths
if (verbosity.gt.1) then !show clock info if (verbosity.gt.1) then !show clock info
!print*,'length(4)',length(4) !print*,'length(4)',length(4)
......
...@@ -150,7 +150,7 @@ program flexpart ...@@ -150,7 +150,7 @@ program flexpart
if (verbosity.gt.0) then if (verbosity.gt.0) then
write(*,*) 'call readpaths' write(*,*) 'call readpaths'
endif endif
call readpaths(pathfile) call readpaths
if (verbosity.gt.1) then !show clock info if (verbosity.gt.1) then !show clock info
!print*,'length(4)',length(4) !print*,'length(4)',length(4)
......
...@@ -132,7 +132,7 @@ subroutine conccalc(itime,weight) ...@@ -132,7 +132,7 @@ subroutine conccalc(itime,weight)
! Take density from 2nd wind field in memory (accurate enough, no time interpolation needed) ! Take density from 2nd wind field in memory (accurate enough, no time interpolation needed)
!***************************************************************************** !*****************************************************************************
do ind=indz,indzp do ind=indz,indzp
rhoprof(ind-indz+1)=p1*rho(ix ,jy ,ind,memind(2)) & rhoprof(ind-indz+1)=p1*rho(ix ,jy ,ind,2) &
+p2*rho(ixp,jy ,ind,2) & +p2*rho(ixp,jy ,ind,2) &
+p3*rho(ix ,jyp,ind,2) & +p3*rho(ix ,jyp,ind,2) &
+p4*rho(ixp,jyp,ind,2) +p4*rho(ixp,jyp,ind,2)
......
...@@ -337,14 +337,16 @@ k = current_grib_level ...@@ -337,14 +337,16 @@ k = current_grib_level
! Output of grid info ! Output of grid info
!******************** !********************
write(*,'(a,i2,a)') ' Nested domain ',l,':' if (lroot) then
write(*,'(a,f10.5,a,f10.5,a,f10.5)') ' Longitude range: ', & write(*,'(a,i2,a)') ' Nested domain ',l,':'
xlon0n(l),' to ',xlon0n(l)+(nxn(l)-1)*dxn(l), & write(*,'(a,f10.5,a,f10.5,a,f10.5)') ' Longitude range: ', &
' Grid distance: ',dxn(l) xlon0n(l),' to ',xlon0n(l)+(nxn(l)-1)*dxn(l), &
write(*,'(a,f10.5,a,f10.5,a,f10.5)') ' Latitude range : ', & ' Grid distance: ',dxn(l)
ylat0n(l),' to ',ylat0n(l)+(nyn(l)-1)*dyn(l), & write(*,'(a,f10.5,a,f10.5,a,f10.5)') ' Latitude range : ', &
' Grid distance: ',dyn(l) ylat0n(l),' to ',ylat0n(l)+(nyn(l)-1)*dyn(l), &
write(*,*) ' Grid distance: ',dyn(l)
write(*,*)
end if
! Determine, how much the resolutions in the nests are enhanced as ! Determine, how much the resolutions in the nests are enhanced as
! compared to the mother grid ! compared to the mother grid
......
...@@ -64,7 +64,7 @@ ifeq ($(gcc), 4.9) ...@@ -64,7 +64,7 @@ ifeq ($(gcc), 4.9)
INCPATH2 = ${ROOT_DIR}/include INCPATH2 = ${ROOT_DIR}/include
LIBPATH1 = ${ROOT_DIR}/lib LIBPATH1 = ${ROOT_DIR}/lib
else ifeq ($(gcc), 5.4) else
ROOT_DIR = /homevip/flexpart/ ROOT_DIR = /homevip/flexpart/
F90 = /usr/bin/gfortran F90 = /usr/bin/gfortran
MPIF90 = /usr/bin/mpif90 MPIF90 = /usr/bin/mpif90
...@@ -72,13 +72,6 @@ else ifeq ($(gcc), 5.4) ...@@ -72,13 +72,6 @@ else ifeq ($(gcc), 5.4)
INCPATH1 = ${ROOT_DIR}/gcc-5.4.0/include INCPATH1 = ${ROOT_DIR}/gcc-5.4.0/include
INCPATH2 = /usr/include INCPATH2 = /usr/include
LIBPATH1 = ${ROOT_DIR}/gcc-5.4.0/lib LIBPATH1 = ${ROOT_DIR}/gcc-5.4.0/lib
else
F90 = gfortran
MPIF90 = mpif90
INCPATH1 = /xnilu_wrk/projects/FLEXPART/flex_wrk/bin64/grib_api/include
INCPATH2 = /usr/include
LIBPATH1 = /xnilu_wrk/projects/FLEXPART/flex_wrk/bin64/grib_api/lib
LIBPATH2 = /xnilu_wrk/projects/FLEXPART/flex_wrk/bin64/grib_api/lib
endif endif
...@@ -112,8 +105,8 @@ FFLAGS = -I$(INCPATH1) -I$(INCPATH2) -O$(O_LEV) -g -cpp -m64 -mcmodel=medium - ...@@ -112,8 +105,8 @@ FFLAGS = -I$(INCPATH1) -I$(INCPATH2) -O$(O_LEV) -g -cpp -m64 -mcmodel=medium -
DBGFLAGS = -I$(INCPATH1) -I$(INCPATH2) -O$(O_LEV_DBG) -g3 -ggdb3 -cpp -m64 -mcmodel=medium -fconvert=little-endian -frecord-marker=4 -fmessage-length=0 -flto=jobserver -O$(O_LEV_DBG) $(NCOPT) -fbacktrace -Wall -fdump-core $(FUSER) # -ffpe-trap=invalid,overflow,denormal,underflow,zero -Warray-bounds -fcheck=all DBGFLAGS = -I$(INCPATH1) -I$(INCPATH2) -O$(O_LEV_DBG) -g3 -ggdb3 -cpp -m64 -mcmodel=medium -fconvert=little-endian -frecord-marker=4 -fmessage-length=0 -flto=jobserver -O$(O_LEV_DBG) $(NCOPT) -fbacktrace -Wall -fdump-core $(FUSER) # -ffpe-trap=invalid,overflow,denormal,underflow,zero -Warray-bounds -fcheck=all
LDFLAGS = $(FFLAGS) -L$(LIBPATH1) -Wl,-rpath,$(LIBPATH1) $(LIBS) -L$(LIBPATH2) LDFLAGS = $(FFLAGS) -L$(LIBPATH1) -Wl,-rpath,$(LIBPATH1) $(LIBS) # -L$(LIBPATH2)
LDDEBUG = $(DBGFLAGS) -L$(LIBPATH1) $(LIBS) -L$(LIBPATH2) LDDEBUG = $(DBGFLAGS) -L$(LIBPATH1) $(LIBS) # -L$(LIBPATH2)
# Default behaviour is to NOT use MPI-IN-PLACE. This may be safer. # Default behaviour is to NOT use MPI-IN-PLACE. This may be safer.
# DEFS = -DUSE_MPIINPLACE # DEFS = -DUSE_MPIINPLACE
......
...@@ -1960,7 +1960,7 @@ contains ...@@ -1960,7 +1960,7 @@ contains
! For now assume that data at all steps either have or do not have water ! For now assume that data at all steps either have or do not have water
if (readclouds) then if (readclouds) then
j=j+1 j=j+1
call MPI_Irecv(ctwc(:,:,mind),d2s1,mp_sp,id_read,MPI_ANY_TAG,& call MPI_Irecv(ctwc(:,:,mind),d2s1*5,mp_sp,id_read,MPI_ANY_TAG,&
&MPI_COMM_WORLD,reqs(j),mp_ierr) &MPI_COMM_WORLD,reqs(j),mp_ierr)
if (mp_ierr /= 0) goto 600 if (mp_ierr /= 0) goto 600
end if end if
...@@ -2325,7 +2325,7 @@ contains ...@@ -2325,7 +2325,7 @@ contains
! For now assume that data at all steps either have or do not have water ! For now assume that data at all steps either have or do not have water
if (readclouds) then if (readclouds) then
j=j+1 j=j+1
call MPI_Irecv(ctwcn(:,:,mind,k),d2s1,mp_sp,id_read,MPI_ANY_TAG,& call MPI_Irecv(ctwcn(:,:,mind,k),d2s1*5,mp_sp,id_read,MPI_ANY_TAG,&
&MPI_COMM_WORLD,reqs(j),mp_ierr) &MPI_COMM_WORLD,reqs(j),mp_ierr)
if (mp_ierr /= 0) goto 600 if (mp_ierr /= 0) goto 600
end if end if
...@@ -2461,6 +2461,19 @@ contains ...@@ -2461,6 +2461,19 @@ contains
& mp_comm_used, mp_ierr) & mp_comm_used, mp_ierr)
end if end if
!CGZ MOVED THIS PART HERE FROM OUTSIDE MPI_IN_PLACE (see below)
!**********************************************************
! Receptor concentrations
if (lroot) then
call MPI_Reduce(MPI_IN_PLACE,creceptor,rcpt_size,mp_sp,MPI_SUM,id_root, &
& mp_comm_used,mp_ierr)
if (mp_ierr /= 0) goto 600
else
call MPI_Reduce(creceptor,0,rcpt_size,mp_sp,MPI_SUM,id_root, &
& mp_comm_used,mp_ierr)
end if
!**********************************************************
#else #else
call MPI_Reduce(gridunc, gridunc0, grid_size3d, mp_sp, MPI_SUM, id_root, & call MPI_Reduce(gridunc, gridunc0, grid_size3d, mp_sp, MPI_SUM, id_root, &
...@@ -2481,15 +2494,18 @@ contains ...@@ -2481,15 +2494,18 @@ contains
if (mp_ierr /= 0) goto 600 if (mp_ierr /= 0) goto 600
end if end if
!CGZ MOVED THIS PART TO MPI_IN_PLACE (line 2467)
!**********************************************************
! Receptor concentrations ! Receptor concentrations
if (lroot) then ! if (lroot) then
call MPI_Reduce(MPI_IN_PLACE,creceptor,rcpt_size,mp_sp,MPI_SUM,id_root, & ! call MPI_Reduce(MPI_IN_PLACE,creceptor,rcpt_size,mp_sp,MPI_SUM,id_root, &
& mp_comm_used,mp_ierr) ! & mp_comm_used,mp_ierr)
if (mp_ierr /= 0) goto 600 ! if (mp_ierr /= 0) goto 600
else ! else
call MPI_Reduce(creceptor,0,rcpt_size,mp_sp,MPI_SUM,id_root, & ! call MPI_Reduce(creceptor,0,rcpt_size,mp_sp,MPI_SUM,id_root, &
& mp_comm_used,mp_ierr) ! & mp_comm_used,mp_ierr)
end if ! end if
!**********************************************************
if (mp_measure_time) call mpif_mtime('commtime',1) if (mp_measure_time) call mpif_mtime('commtime',1)
......
...@@ -159,7 +159,7 @@ module par_mod ...@@ -159,7 +159,7 @@ module par_mod
! Maximum dimensions of the nested input grids ! Maximum dimensions of the nested input grids
!********************************************* !*********************************************
integer,parameter :: maxnests=0,nxmaxn=0,nymaxn=0 integer,parameter :: maxnests=1,nxmaxn=451,nymaxn=226
! nxmax,nymax maximum dimension of wind fields in x and y ! nxmax,nymax maximum dimension of wind fields in x and y
! direction, respectively ! direction, respectively
...@@ -217,7 +217,7 @@ module par_mod ...@@ -217,7 +217,7 @@ module par_mod
! Maximum number of particles, species, and similar ! Maximum number of particles, species, and similar
!************************************************** !**************************************************
integer,parameter :: maxpart=100000 integer,parameter :: maxpart=1000000
integer,parameter :: maxspec=1 integer,parameter :: maxspec=1
real,parameter :: minmass=0.0001 real,parameter :: minmass=0.0001
...@@ -261,7 +261,7 @@ module par_mod ...@@ -261,7 +261,7 @@ module par_mod
! Dimension of random number field ! Dimension of random number field
!********************************* !*********************************
integer,parameter :: maxrand=1000000 integer,parameter :: maxrand=20000000
! maxrand number of random numbers used ! maxrand number of random numbers used
......
...@@ -339,7 +339,8 @@ subroutine readwind_ecmwf(indj,n,uuh,vvh,wwh) ...@@ -339,7 +339,8 @@ subroutine readwind_ecmwf(indj,n,uuh,vvh,wwh)
!! in snow depth, but I don't feel 100% good about this just yet. It may !! in snow depth, but I don't feel 100% good about this just yet. It may
!! need to be scrutinized more closely in the future. !! need to be scrutinized more closely in the future.
conversion_factor = 1.0 ! eso: reverted conversion factor to 1000.
conversion_factor = 1000.0
DO j=0,nymin1 DO j=0,nymin1
DO i=0,nxfield-1 DO i=0,nxfield-1
sd(i,j,1,n) = zsec4(nxfield*(ny-j-1)+i+1)/conversion_factor sd(i,j,1,n) = zsec4(nxfield*(ny-j-1)+i+1)/conversion_factor
...@@ -528,7 +529,7 @@ subroutine readwind_ecmwf(indj,n,uuh,vvh,wwh) ...@@ -528,7 +529,7 @@ subroutine readwind_ecmwf(indj,n,uuh,vvh,wwh)
END DO END DO
END DO END DO
ELSE IF(TRIM(fpname) .EQ. 'CICE') then !! CIWC Cloud ice water content ELSE IF(TRIM(fpname) .EQ. 'CIWC') then !! CIWC Cloud ice water content
DO j=0,nymin1 DO j=0,nymin1
DO i=0,nxfield-1 DO i=0,nxfield-1
ciwch(i,j,nlev_ec-k+2,n)=zsec4(nxfield*(ny-j-1)+i+1) ciwch(i,j,nlev_ec-k+2,n)=zsec4(nxfield*(ny-j-1)+i+1)
......
...@@ -82,7 +82,7 @@ subroutine readwind_nests(indj,n,uuhn,vvhn,wwhn) ...@@ -82,7 +82,7 @@ subroutine readwind_nests(indj,n,uuhn,vvhn,wwhn)
real,parameter :: eps=1.e-4 real,parameter :: eps=1.e-4
real :: ewss(0:nxmaxn-1,0:nymaxn-1),nsss(0:nxmaxn-1,0:nymaxn-1) real :: ewss(0:nxmaxn-1,0:nymaxn-1),nsss(0:nxmaxn-1,0:nymaxn-1)
real :: plev1,pmean,tv,fu,hlev1,ff10m,fflev1 real :: plev1,pmean,tv,fu,hlev1,ff10m,fflev1
real :: conversion_factor !added by mc to make it consistent with new gridchek.f90 real :: conversion_factor=1000. !added by mc to make it consistent with new gridchek.f90
logical :: hflswitch,strswitch logical :: hflswitch,strswitch
...@@ -124,8 +124,8 @@ subroutine readwind_nests(indj,n,uuhn,vvhn,wwhn) ...@@ -124,8 +124,8 @@ subroutine readwind_nests(indj,n,uuhn,vvhn,wwhn)
PRINT *, 'Loading Vtable: ', VTABLE_PATH PRINT *, 'Loading Vtable: ', VTABLE_PATH
call vtable_load_by_name(VTABLE_PATH, my_vtable) call vtable_load_by_name(VTABLE_PATH, my_vtable)
!! Debugging tool !! Debugging tool
PRINT *, 'Dump of Vtable...' !PRINT *, 'Dump of Vtable...'
call vtable_dump_records(my_vtable) !call vtable_dump_records(my_vtable)
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
...@@ -136,8 +136,8 @@ subroutine readwind_nests(indj,n,uuhn,vvhn,wwhn) ...@@ -136,8 +136,8 @@ subroutine readwind_nests(indj,n,uuhn,vvhn,wwhn)
! inventory of the FP-related messages, relative to the Vtable that's ! inventory of the FP-related messages, relative to the Vtable that's
! already been open. ! already been open.
CALL vtable_gribfile_inventory( path(numpath+2*(l-1)+1)(1:length(numpath+2*(l-1)+1))// trim(wfnamen(l, indj)), & ! CALL vtable_gribfile_inventory( path(numpath+2*(l-1)+1)(1:length(numpath+2*(l-1)+1))// trim(wfnamen(l, indj)), &
& my_vtable) ! & my_vtable)
!!!!!!!!!!!!!!!!!!! VTABLE code !!!!!!!!!!!!!!!!!!! VTABLE code
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
...@@ -177,7 +177,7 @@ subroutine readwind_nests(indj,n,uuhn,vvhn,wwhn) ...@@ -177,7 +177,7 @@ subroutine readwind_nests(indj,n,uuhn,vvhn,wwhn)
!!!!!!!!!!!!!!!!!!! VTABLE code !!!!!!!!!!!!!!!!!!! VTABLE code
! Get the fpname ! Get the fpname
fpname = vtable_get_fpname(igrib, my_vtable) fpname = vtable_get_fpname(igrib, my_vtable)
print *, 'fpname: ', trim(fpname) !print *, 'fpname: ', trim(fpname)
!!!!!!!!!!!!!!!!!!! VTABLE code !!!!!!!!!!!!!!!!!!! VTABLE code
...@@ -430,7 +430,8 @@ subroutine readwind_nests(indj,n,uuhn,vvhn,wwhn) ...@@ -430,7 +430,8 @@ subroutine readwind_nests(indj,n,uuhn,vvhn,wwhn)
!! in snow depth, but I don't feel 100% good about this just yet. It may !! in snow depth, but I don't feel 100% good about this just yet. It may
!! need to be scrutinized more closely in the future. !! need to be scrutinized more closely in the future.
conversion_factor = 1.0 ! eso: reverted conversion factor to 1000.
conversion_factor = 1000.0
DO j=0,nyn(l)-1 DO j=0,nyn(l)-1
DO i=0,nxn(l)-1 DO i=0,nxn(l)-1
sdn(i,j,1,n,l) = zsec4(nxn(l)*(nyn(l)-j-1)+i+1)/conversion_factor sdn(i,j,1,n,l) = zsec4(nxn(l)*(nyn(l)-j-1)+i+1)/conversion_factor
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment