Commit db712a8f authored by Espen Sollum's avatar Espen Sollum

Completed handling of nested wind fields with cloud water (for wet deposition).

parent b0434e19
......@@ -314,11 +314,6 @@ program flexpart
call readOHfield
endif
!! testing !!
! open(999,file=trim(path(1))//'OH_FIELDS/jscalar_50N.txt',action='write',status='new')
! open(998,file=trim(path(1))//'OH_FIELDS/jscalar_50S.txt',action='write',status='new')
! Write basic information on the simulation to a file "header"
! and open files that are to be kept open throughout the simulation
!******************************************************************
......
......@@ -347,11 +347,6 @@ program flexpart
call readOHfield
endif
!! testing !!
! open(999,file=trim(path(1))//'OH_FIELDS/jscalar_50N.txt',action='write',status='new')
! open(998,file=trim(path(1))//'OH_FIELDS/jscalar_50S.txt',action='write',status='new')
! Write basic information on the simulation to a file "header"
! and open files that are to be kept open throughout the simulation
!******************************************************************
......
......@@ -102,6 +102,7 @@ subroutine calcpar_nests(n,uuhn,vvhn,pvhn)
ustarn(ix,jy,1,n,l)=scalev(psn(ix,jy,1,n,l),tt2n(ix,jy,1,n,l), &
td2n(ix,jy,1,n,l),surfstrn(ix,jy,1,n,l))
if (ustarn(ix,jy,1,n,l).le.1.e-8) ustarn(ix,jy,1,n,l)=1.e-8
! 2) Calculation of inverse Obukhov length scale
!***********************************************
......@@ -134,7 +135,7 @@ subroutine calcpar_nests(n,uuhn,vvhn,pvhn)
if(lsubgrid.eq.1) then
subsceff=min(excessoron(ix,jy,l),hmixplus)
else
subsceff=0
subsceff=0.0
endif
!
! CALCULATE HMIX EXCESS ACCORDING TO SUBGRIDSCALE VARIABILITY AND STABILITY
......@@ -148,8 +149,9 @@ subroutine calcpar_nests(n,uuhn,vvhn,pvhn)
!********************************************
if (DRYDEP) then
z0(4)=0.016*ustarn(ix,jy,1,n,l)*ustarn(ix,jy,1,n,l)/ga
z0(9)=0.016*ustarn(ix,jy,1,n,l)*ustarn(ix,jy,1,n,l)/ga
! z0(4)=0.016*ustarn(ix,jy,1,n,l)*ustarn(ix,jy,1,n,l)/ga
! z0(9)=0.016*ustarn(ix,jy,1,n,l)*ustarn(ix,jy,1,n,l)/ga
z0(7)=0.016*ustarn(ix,jy,1,n,l)*ustarn(ix,jy,1,n,l)/ga
! Calculate relative humidity at surface
!***************************************
......@@ -227,8 +229,10 @@ subroutine calcpar_nests(n,uuhn,vvhn,pvhn)
end do
end do
! Calculation of potential vorticity on 3-d grid
!***********************************************
call calcpv_nests(l,n,uuhn,vvhn,pvhn)
call calcpv_nests(l,n,uuhn,vvhn,pvhn)
end do
......
......@@ -69,10 +69,8 @@ subroutine calcpv(n,uuh,vvh,pvh)
enddo
enddo
! :DBG: halts with SIGFPE if compiling with -ffpe-trap
! where(ppml.le.0) ppml=10000.0
! ppmk(:,:,1:nuvz)=(100000./ppml(:,:,1:nuvz))**kappa
ppmk=(100000./ppml)**kappa
ppmk(0:nxmin1,0:nymin1,1:nuvz)=(100000./ppml(0:nxmin1,0:nymin1,1:nuvz))**kappa
do jy=0,nymin1
if (sglobal.and.jy.eq.0) goto 10
......
......@@ -67,7 +67,8 @@ subroutine calcpv_nests(l,n,uuhn,vvhn,pvhn)
enddo
enddo
enddo
ppmk=(100000./ppml)**kappa
! ppmk=(100000./ppml)**kappa
ppmk(0:nxn(l)-1,0:nyn(l)-1,1:nuvz)=(100000./ppml(0:nxn(l)-1,0:nyn(l)-1,1:nuvz))**kappa
do jy=0,nyn(l)-1
phi = (ylat0n(l) + jy * dyn(l)) * pi / 180.
......
......@@ -133,7 +133,9 @@ module com_mod
!ZHG SEP 2015 wheather or not to read clouds from GRIB
logical :: readclouds
!ESO DEC 2015 whether or not both clwc and ciwc are present (if so they are summed)
logical :: sumclouds
logical :: sumclouds
logical,dimension(maxnests) :: readclouds_nest, sumclouds_nest
!NIK 16.02.2015
......@@ -369,7 +371,6 @@ module com_mod
!ZHG Sep 2015
real :: icloud_stats(0:nxmax-1,0:nymax-1,5,numwfmem)
real :: clh(0:nxmax-1,0:nymax-1,nuvzmax,numwfmem)
real :: clw4(0:nxmax-1,0:nymax-1,numwfmem) ! eso: =icloud_stats(:,:,4,:)
......@@ -486,8 +487,9 @@ module com_mod
!*****************
real,allocatable,dimension(:,:,:,:,:) :: uun, vvn, wwn, ttn, qvn, pvn,&
& rhon, drhodzn, tthn, qvhn
integer,allocatable,dimension(:,:,:,:) :: cloudsnh
& rhon, drhodzn, tthn, qvhn, clwcn, ciwcn, clwn, clwchn, ciwchn
real,allocatable,dimension(:,:,:,:) :: clw4n
integer,allocatable,dimension(:,:,:,:) :: cloudshn
integer(kind=1),allocatable,dimension(:,:,:,:,:) :: cloudsn
! 2d nested fields
......@@ -787,13 +789,29 @@ contains
allocate(ttn(0:nxmaxn-1,0:nymaxn-1,nzmax,numwfmem,numbnests))
allocate(qvn(0:nxmaxn-1,0:nymaxn-1,nzmax,numwfmem,numbnests))
allocate(pvn(0:nxmaxn-1,0:nymaxn-1,nzmax,numwfmem,numbnests))
allocate(cloudsn(0:nxmaxn-1,0:nymaxn-1,0:nzmax,numwfmem,numbnests))
allocate(cloudsnh(0:nxmaxn-1,0:nymaxn-1,numwfmem,numbnests))
allocate(clwcn(0:nxmaxn-1,0:nymaxn-1,nzmax,numwfmem,numbnests))
allocate(ciwcn(0:nxmaxn-1,0:nymaxn-1,nzmax,numwfmem,numbnests))
allocate(clwn(0:nxmaxn-1,0:nymaxn-1,nzmax,numwfmem,numbnests))
allocate(cloudsn(0:nxmaxn-1,0:nymaxn-1,nzmax,numwfmem,numbnests))
allocate(cloudshn(0:nxmaxn-1,0:nymaxn-1,numwfmem,numbnests))
allocate(rhon(0:nxmaxn-1,0:nymaxn-1,nzmax,numwfmem,numbnests))
allocate(drhodzn(0:nxmaxn-1,0:nymaxn-1,nzmax,numwfmem,numbnests))
allocate(tthn(0:nxmaxn-1,0:nymaxn-1,nuvzmax,numwfmem,numbnests))
allocate(qvhn(0:nxmaxn-1,0:nymaxn-1,nuvzmax,numwfmem,numbnests))
allocate(clwchn(0:nxmaxn-1,0:nymaxn-1,nuvzmax,numwfmem,numbnests))
allocate(ciwchn(0:nxmaxn-1,0:nymaxn-1,nuvzmax,numwfmem,numbnests))
allocate(clw4n(0:nxmax-1,0:nymax-1,numwfmem,numbnests))
! clw4n(:,:,:,:)=0.
clwcn(:,:,:,:,:)=0.
ciwcn(:,:,:,:,:)=0.
clwchn(:,:,:,:,:)=0.
ciwchn(:,:,:,:,:)=0.
! clwn(:,:,:,:,:)=0.
! cloudsn(:,:,:,:,:)=0
! cloudshn(:,:,:,:)=0
end subroutine com_mod_allocate_nests
......
......@@ -625,24 +625,26 @@ subroutine concoutput(itime,outnum,gridtotalunc,wetgridtotalunc, &
! Reinitialization of grid
!*************************
do ks=1,nspec
do kp=1,maxpointspec_act
do i=1,numreceptor
creceptor(i,ks)=0.
end do
do jy=0,numygrid-1
do ix=0,numxgrid-1
do l=1,nclassunc
do nage=1,nageclass
do kz=1,numzgrid
gridunc(ix,jy,kz,ks,kp,l,nage)=0.
end do
end do
end do
end do
end do
end do
end do
! do ks=1,nspec
! do kp=1,maxpointspec_act
! do i=1,numreceptor
! creceptor(i,ks)=0.
! end do
! do jy=0,numygrid-1
! do ix=0,numxgrid-1
! do l=1,nclassunc
! do nage=1,nageclass
! do kz=1,numzgrid
! gridunc(ix,jy,kz,ks,kp,l,nage)=0.
! end do
! end do
! end do
! end do
! end do
! end do
! end do
creceptor(:,:)=0.
gridunc(:,:,:,:,:,:,:)=0.
end subroutine concoutput
......@@ -544,24 +544,26 @@ subroutine concoutput_nest(itime,outnum)
! Reinitialization of grid
!*************************
do ks=1,nspec
do kp=1,maxpointspec_act
do i=1,numreceptor
creceptor(i,ks)=0.
end do
do jy=0,numygridn-1
do ix=0,numxgridn-1
do l=1,nclassunc
do nage=1,nageclass
do kz=1,numzgrid
griduncn(ix,jy,kz,ks,kp,l,nage)=0.
end do
end do
end do
end do
end do
end do
end do
! do ks=1,nspec
! do kp=1,maxpointspec_act
! do i=1,numreceptor
! creceptor(i,ks)=0.
! end do
! do jy=0,numygridn-1
! do ix=0,numxgridn-1
! do l=1,nclassunc
! do nage=1,nageclass
! do kz=1,numzgrid
! griduncn(ix,jy,kz,ks,kp,l,nage)=0.
! end do
! end do
! end do
! end do
! end do
! end do
! end do
creceptor(:,:)=0.
griduncn(:,:,:,:,:,:,:)=0.
end subroutine concoutput_nest
......
......@@ -40,15 +40,19 @@ module wind_mod
! Maximum dimensions of the input mother grids
!*********************************************
integer,parameter :: nxmax=361,nymax=181,nuvzmax=152,nwzmax=152,nzmax=152 !ECMWF new
! integer,parameter :: nxmax=361,nymax=181,nuvzmax=92,nwzmax=92,nzmax=92 !ECMWF new
integer,parameter :: nxmax=361,nymax=181,nuvzmax=138,nwzmax=138,nzmax=138 !ECMWF new
integer,parameter :: nxshift=359
! integer,parameter :: nxmax=361,nymax=181,nuvzmax=138,nwzmax=138,nzmax=138 !test BUG
! integer,parameter :: nxshift=0
!
!*********************************************
! Maximum dimensions of the nested input grids
!*********************************************
integer,parameter :: maxnests=1,nxmaxn=361,nymaxn=181 !ECMWF
integer,parameter :: maxnests=1,nxmaxn=361,nymaxn=181
! integer,parameter :: maxnests=1,nxmaxn=86,nymaxn=31
end module wind_mod
......@@ -193,16 +193,12 @@ subroutine gridcheck
!ZHG FOR CLOUDS FROM GRIB
elseif ((parCat.eq.1).and.(parNum.eq.83).and.(typSurf.eq.105)) then ! clwc
isec1(6)=246 ! indicatorOfParameter
! readclouds=.true.
! sumclouds=.false.
elseif ((parCat.eq.1).and.(parNum.eq.84).and.(typSurf.eq.105)) then ! ciwc
isec1(6)=247 ! indicatorOfParameter
!ZHG end
! ESO qc(=clwc+ciwc)
elseif ((parCat.eq.201).and.(parNum.eq.31).and.(typSurf.eq.105)) then ! qc
isec1(6)=201031 ! indicatorOfParameter
! readclouds=.true.
! sumclouds=.true.
elseif ((parCat.eq.3).and.(parNum.eq.0).and.(typSurf.eq.1)) then !SP
isec1(6)=134 ! indicatorOfParameter
elseif ((parCat.eq.2).and.(parNum.eq.32)) then ! W, actually eta dot
......@@ -365,6 +361,22 @@ subroutine gridcheck
if (nxshift.ge.nxfield) stop 'nxshift (par_mod) too large'
endif ! gotGrid
if (nx.gt.nxmax) then
write(*,*) 'FLEXPART error: Too many grid points in x direction.'
write(*,*) 'Reduce resolution of wind fields.'
write(*,*) 'Or change parameter settings in file par_mod.'
write(*,*) nx,nxmax
stop
endif
if (ny.gt.nymax) then
write(*,*) 'FLEXPART error: Too many grid points in y direction.'
write(*,*) 'Reduce resolution of wind fields.'
write(*,*) 'Or change parameter settings in file par_mod.'
write(*,*) ny,nymax
stop
endif
k=isec1(8)
if(isec1(6).eq.131) iumax=max(iumax,nlev_ec-k+1)
if(isec1(6).eq.135) iwmax=max(iwmax,nlev_ec-k+1)
......@@ -410,22 +422,6 @@ subroutine gridcheck
nwz =iwmax
if(nuvz.eq.nlev_ec) nwz=nlev_ec+1
if (nx.gt.nxmax) then
write(*,*) 'FLEXPART error: Too many grid points in x direction.'
write(*,*) 'Reduce resolution of wind fields.'
write(*,*) 'Or change parameter settings in file par_mod.'
write(*,*) nx,nxmax
stop
endif
if (ny.gt.nymax) then
write(*,*) 'FLEXPART error: Too many grid points in y direction.'
write(*,*) 'Reduce resolution of wind fields.'
write(*,*) 'Or change parameter settings in file par_mod.'
write(*,*) ny,nymax
stop
endif
if (nuvz+1.gt.nuvzmax) then
write(*,*) 'FLEXPART error: Too many u,v grid points in z '// &
'direction.'
......
......@@ -171,22 +171,18 @@ subroutine gridcheck_nests
isec1(6)=133 ! indicatorOfParameter
elseif ((parCat.eq.1).and.(parNum.eq.83).and.(typSurf.eq.105)) then ! clwc
isec1(6)=246 ! indicatorOfParameter
! readclouds=.true.
! sumclouds=.false.
elseif ((parCat.eq.1).and.(parNum.eq.84).and.(typSurf.eq.105)) then ! ciwc
isec1(6)=247 ! indicatorOfParameter
!ZHG end
! ESO qc(=clwc+ciwc)
elseif ((parCat.eq.201).and.(parNum.eq.31).and.(typSurf.eq.105)) then ! qc
isec1(6)=201031 ! indicatorOfParameter
! readclouds=.true.
! sumclouds=.true.
elseif ((parCat.eq.3).and.(parNum.eq.0).and.(typSurf.eq.1)) then !SP
isec1(6)=134 ! indicatorOfParameter
elseif ((parCat.eq.2).and.(parNum.eq.32)) then ! W, actually eta dot
isec1(6)=135 ! indicatorOfParameter
elseif ((parCat.eq.128).and.(parNum.eq.77)) then ! W, actually eta dot !added bymc to make it consistent with new gridchek.f90
isec1(6)=135 ! indicatorOfParameter ! ! " " " " " " " " " " " " " " " " " " " " " " " " " " "
elseif ((parCat.eq.128).and.(parNum.eq.77)) then ! W, actually eta dot !added bymc to make it consistent with new gridcheck.f90
isec1(6)=135 ! indicatorOfParameter !
elseif ((parCat.eq.3).and.(parNum.eq.0).and.(typSurf.eq.101)) then !SLP
isec1(6)=151 ! indicatorOfParameter
elseif ((parCat.eq.2).and.(parNum.eq.2).and.(typSurf.eq.103)) then ! 10U
......@@ -257,6 +253,24 @@ subroutine gridcheck_nests
nlev_ecn=isec2(12)/2-1
endif ! ifield
if (nxn(l).gt.nxmaxn) then
write(*,*) 'FLEXPART error: Too many grid points in x direction.'
write(*,*) 'Reduce resolution of wind fields (file GRIDSPEC)'
write(*,*) 'for nesting level ',l
write(*,*) 'Or change parameter settings in file par_mod.'
write(*,*) nxn(l),nxmaxn
stop
endif
if (nyn(l).gt.nymaxn) then
write(*,*) 'FLEXPART error: Too many grid points in y direction.'
write(*,*) 'Reduce resolution of wind fields (file GRIDSPEC)'
write(*,*) 'for nesting level ',l
write(*,*) 'Or change parameter settings in file par_mod.'
write(*,*) nyn(l),nymaxn
stop
endif
!HSO get the second part of the grid dimensions only from GRiB1 messages
if (isec1(6) .eq. 167 .and. (gotGrib.eq.0)) then !added by mc to make it consistent with new gridchek.f90 note that gotGrid must be changed in gotGrib!!
call grib_get_real8(igrib,'longitudeOfFirstGridPointInDegrees', & !comment by mc: note that this was in the (if (ifield.eq.1) ..end above in gridchek.f90 see line 257
......@@ -332,24 +346,6 @@ subroutine gridcheck_nests
nwzn=iwmax
if(nuvzn.eq.nlev_ec) nwzn=nlev_ecn+1
if (nxn(l).gt.nxmaxn) then
write(*,*) 'FLEXPART error: Too many grid points in x direction.'
write(*,*) 'Reduce resolution of wind fields (file GRIDSPEC)'
write(*,*) 'for nesting level ',l
write(*,*) 'Or change parameter settings in file par_mod.'
write(*,*) nxn(l),nxmaxn
stop
endif
if (nyn(l).gt.nymaxn) then
write(*,*) 'FLEXPART error: Too many grid points in y direction.'
write(*,*) 'Reduce resolution of wind fields (file GRIDSPEC)'
write(*,*) 'for nesting level ',l
write(*,*) 'Or change parameter settings in file par_mod.'
write(*,*) nyn(l),nymaxn
stop
endif
if ((nuvzn.gt.nuvzmax).or.(nwzn.gt.nwzmax)) then
write(*,*) 'FLEXPART error: Nested wind fields have too many'// &
'vertical levels.'
......
......@@ -98,7 +98,7 @@ subroutine initialize(itime,ldt,up,vp,wp, &
ixp=ix+1
jyp=jy+1
h=max(hmix(ix ,jy ,1,memind(1)), &
h=max(hmix(ix ,jy,1,memind(1)), &
hmix(ixp,jy ,1,memind(1)), &
hmix(ix ,jyp,1,memind(1)), &
hmix(ixp,jyp,1,memind(1)), &
......@@ -156,15 +156,15 @@ subroutine initialize(itime,ldt,up,vp,wp, &
vp=rannumb(nrand+1)*sigv
wp=rannumb(nrand+2)
if (.not.turbswitch) then ! modified by mc
wp=wp*sigw
wp=wp*sigw
else if (cblflag.eq.1) then ! modified by mc
if(-h/ol.gt.5) then
!if (ol.lt.0.) then
!if (ol.gt.0.) then !by mc : only for test correct is lt.0
call initialize_cbl_vel(idummy,zt,ust,wst,h,sigw,wp,ol)
else
wp=wp*sigw
end if
if(-h/ol.gt.5) then
!if (ol.lt.0.) then
!if (ol.gt.0.) then !by mc : only for test correct is lt.0
call initialize_cbl_vel(idummy,zt,ust,wst,h,sigw,wp,ol)
else
wp=wp*sigw
end if
end if
......
......@@ -87,10 +87,14 @@ subroutine interpol_rain_nests(yy1,yy2,yy3,nxmaxn,nymaxn,nzmax, &
! If point at border of grid -> small displacement into grid
!***********************************************************
if (xt.ge.(real(nxn(ngrid)-1)-0.0001)) &
xt=real(nxn(ngrid)-1)-0.0001
if (yt.ge.(real(nyn(ngrid)-1)-0.0001)) &
yt=real(nyn(ngrid)-1)-0.0001
! if (xt.ge.(real(nxn(ngrid)-1)-0.0001)) &
! xt=real(nxn(ngrid)-1)-0.0001
! if (yt.ge.(real(nyn(ngrid)-1)-0.0001)) &
! yt=real(nyn(ngrid)-1)-0.0001
! ESO make it consistent with interpol_rain
if (xt.ge.(real(nxn(ngrid)-1))) xt=real(nxn(ngrid)-1)-0.00001
if (yt.ge.(real(nyn(ngrid)-1))) yt=real(nyn(ngrid)-1)-0.00001
......@@ -104,6 +108,7 @@ subroutine interpol_rain_nests(yy1,yy2,yy3,nxmaxn,nymaxn,nzmax, &
ix=int(xt)
jy=int(yt)
ixp=ix+1
jyp=jy+1
ddx=xt-real(ix)
......
......@@ -387,7 +387,7 @@ readavailable.o: com_mod.o par_mod.o
readcommand.o: com_mod.o par_mod.o
readdepo.o: com_mod.o par_mod.o
readlanduse.o: com_mod.o par_mod.o
readlanduse_int1.o: com_mod.o par_mod.o
#readlanduse_int1.o: com_mod.o par_mod.o
readOHfield.o: com_mod.o oh_mod.o par_mod.o
readoutgrid.o: com_mod.o outg_mod.o par_mod.o
readoutgrid_nest.o: com_mod.o outg_mod.o par_mod.o
......
......@@ -859,8 +859,6 @@ contains
! MPI_Bcast is used, so implicitly all processes are synchronized at this
! step
!
! TODO
! GFS version
!
!*******************************************************************************
use com_mod
......@@ -921,7 +919,6 @@ contains
!**********************************************************************
! The non-reader processes need to know if cloud water were read.
! TODO: only at first step or always?
call MPI_Bcast(readclouds,1,MPI_LOGICAL,id_read,MPI_COMM_WORLD,mp_ierr)
if (mp_ierr /= 0) goto 600
......@@ -1016,9 +1013,6 @@ contains
call MPI_Bcast(oli(:,:,:,li:ui),d2s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
if (mp_ierr /= 0) goto 600
call MPI_Bcast(z0,numclass,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
if (mp_ierr /= 0) goto 600
if (mp_measure_time) call mpif_mtime('commtime',1)
goto 601
......@@ -1044,8 +1038,6 @@ contains
! MPI_Bcast is used, so implicitly all processes are synchronized at this
! step
!
! TODO
! Transfer cloud water/ice if and when available for nested
!
!***********************************************************************
use com_mod
......@@ -1059,11 +1051,11 @@ contains
integer :: li,ui
! Common array sizes used for communications
integer :: d3_size1 = nxmaxn*nymaxn*nzmax*maxnests
integer :: d3_size2 = nxmaxn*nymaxn*nuvzmax*maxnests
integer :: d2_size1 = nxmaxn*nymaxn*maxnests
integer :: d2_size2 = nxmaxn*nymaxn*maxspec*maxnests
integer :: d2_size3 = nxmaxn*nymaxn*maxnests
integer :: d3_size1 = nxmaxn*nymaxn*nzmax
integer :: d3_size2 = nxmaxn*nymaxn*nuvzmax
integer :: d2_size1 = nxmaxn*nymaxn
integer :: d2_size2 = nxmaxn*nymaxn*maxspec
integer :: d2_size3 = nxmaxn*nymaxn
integer :: d3s1,d3s2,d2s1,d2s2
......@@ -1105,6 +1097,10 @@ contains
! processes
!**********************************************************************
! The non-reader processes need to know if cloud water were read.
call MPI_Bcast(readclouds_nest,maxnests,MPI_LOGICAL,id_read,MPI_COMM_WORLD,mp_ierr)
if (mp_ierr /= 0) goto 600
! Static fields/variables sent only at startup
if (first_call) then
call MPI_Bcast(oron(:,:,:),d2_size3,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
......@@ -1119,8 +1115,8 @@ contains
end if
! MPI prefers contiguous arrays for sending (else a buffer is created),
! hence the loop
! hence the loop over nests
!**********************************************************************
do i=1, numbnests
! 3D fields
call MPI_Bcast(uun(:,:,:,li:ui,i),d3s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
......@@ -1144,10 +1140,18 @@ contains
call MPI_Bcast(pvn(:,:,:,li:ui,i),d3s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
if (mp_ierr /= 0) goto 600
call MPI_Bcast(cloudsn(:,:,:,li:ui,i),d3s1,MPI_INTEGER1,id_read,MPI_COMM_WORLD,mp_ierr)
if (mp_ierr /= 0) goto 600
if (mp_ierr /= 0) goto 600
! cloud water/ice:
if (readclouds_nest(i)) then
! call MPI_Bcast(icloud_stats(:,:,:,li:ui),d2s1*5,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
! if (mp_ierr /= 0) goto 600
call MPI_Bcast(clw4n(:,:,li:ui,i),d2s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
if (mp_ierr /= 0) goto 600
end if
! 2D fields
call MPI_Bcast(cloudsnh(:,:,li:ui,i),d2s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
call MPI_Bcast(cloudshn(:,:,li:ui,i),d2s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
if (mp_ierr /= 0) goto 600
call MPI_Bcast(vdepn(:,:,:,li:ui,i),d2s2,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
if (mp_ierr /= 0) goto 600
......@@ -1205,7 +1209,6 @@ contains
! memstat -- input, for resolving pointer to windfield index being read
! mind -- index where to place new fields
!
! TODO
!
!*******************************************************************************
use com_mod
......@@ -1391,7 +1394,6 @@ contains
! VARIABLES
! memstat -- input, used to resolve windfield index being received
!
! TODO
!
!*******************************************************************************
use com_mod
......@@ -2018,7 +2020,6 @@ contains
! j=mp_pid*nvar_async
! In the implementation with 3 fields, the processes may have posted
! MPI_Irecv requests that should be cancelled here
!! TODO:
! if (.not.lmp_sync) then
! r=mp_pid*nvar_async
! do j=r,r+nvar_async-1
......@@ -2104,7 +2105,6 @@ contains
! VARIABLES
!
!
! TODO
!
!*******************************************************************************
use com_mod
......@@ -2131,7 +2131,7 @@ contains
! Time-varying fields:
uu(:,:,:,li:ui) = 10.0
vv(:,:,:,li:ui) = 0.0
uupol(:,:,:,li:ui) = 10.0 ! TODO check if ok
uupol(:,:,:,li:ui) = 10.0
vvpol(:,:,:,li:ui)=0.0
ww(:,:,:,li:ui)=0.
tt(:,:,:,li:ui)=300.
......@@ -2163,7 +2163,6 @@ contains
hmix(:,:,:,li:ui)=10000.
tropopause(:,:,:,li:ui)=10000.
oli(:,:,:,li:ui)=0.01
z0=1.0
end subroutine set_fields_synthetic
......
......@@ -78,6 +78,11 @@ module netcdf_output_mod
implicit none
private
public :: writeheader_netcdf,concoutput_surf_nest_netcdf,concoutput_netcdf,&
&concoutput_nest_netcdf,concoutput_surf_netcdf
! include 'netcdf.inc'
! parameter for data compression (1-9, 9 = most aggressive)
......@@ -96,7 +101,7 @@ module netcdf_output_mod
integer, dimension(5) :: depdimids, depdimidsn
real,parameter :: eps=nxmax/3.e5
private:: writemetadata, output_units, nf90_err
! private:: writemetadata, output_units, nf90_err
! switch output of release point information on/off
logical, parameter :: write_releases = .true.
......@@ -627,10 +632,18 @@ subroutine writeheader_netcdf(lnest)
call nf90_err(nf90_put_var(ncid, levID, outheight(1:numzgrid)))
! volume
call nf90_err(nf90_put_var(ncid, volID, volume(:,:,:)))
if (lnest) then
call nf90_err(nf90_put_var(ncid, volID, volumen(:,:,:)))
else
call nf90_err(nf90_put_var(ncid, volID, volume(:,:,:)))
end if
! area
call nf90_err(nf90_put_var(ncid, areaID, area(:,:)))
if (lnest) then
call nf90_err(nf90_put_var(ncid, areaID, arean(:,:)))
else
call nf90_err(nf90_put_var(ncid, areaID, area(:,:)))
end if
if (write_releases.eqv..true.) then
! release point information
......
......@@ -138,9 +138,7 @@ module par_mod
!integer,parameter :: nxmax=361,nymax=181,nuvzmax=92,nwzmax=92,nzmax=92 !FNL XF
!integer,parameter :: nxmax=361,nymax=181,nuvzmax=152,nwzmax=152,nzmax=152 !ECMWF new
!integer,parameter :: nxmax=361,nymax=181,nuvzmax=92,nwzmax=92,nzmax=92 !ECMWF
!integer,parameter :: nxmax=361,nymax=181,nuvzmax=26,nwzmax=26,nzmax=26
!integer,parameter :: nxmax=721,nymax=361,nuvzmax=64,nwzmax=64,nzmax=64
!integer,parameter :: nxmax=1201,nymax=235,nuvzmax=58,nwzmax=58,nzmax=58
! integer,parameter :: nxshift=359 ! for ECMWF
!integer,parameter :: nxshift=0 ! for GFS or FNL
......@@ -215,9 +213,9 @@ module par_mod
! Maximum number of particles, species, and similar
!**************************************************
integer,parameter :: maxpart=40000000
integer,parameter :: maxpart=400000
integer,parameter :: maxspec=6
integer,parameter :: minmass=0.0 !0.0001
real,parameter :: minmass=0.0 !0.0001
! maxpart Maximum number of particles
! maxspec Maximum number of chemical species per release
......
......@@ -282,7 +282,7 @@ subroutine readavailable
write(*,*) ' #### CANNOT BE OPENED #### '
stop
999 write(*,*) ' #### FLEXPART MODEL ERROR! AVAILABLE IILE #### '
999 write(*,*) ' #### FLEXPART MODEL ERROR! AVAILABLE FILE #### '
write(*,'(a)') ' '//path(4)(1:length(4))
write(*,*) ' #### CANNOT BE OPENED #### '
stop
......
......@@ -346,22 +346,20 @@ subroutine readwind_nests(indj,n,uuhn,vvhn,wwhn)
! ESO TODO:
! -add check for if one of clwc/ciwc missing (error),
! also if all 3 cw fields present, use qc and disregard the others
! -use same flags readclouds/sumclouds as in mother grid? this assumes
! that both the nested and mother grids contain CW in same format
if(isec1(6).eq.246) then !! CLWC Cloud liquid water content [kg/kg]
clwch(i,j,nlev_ec-k+2,n)=zsec4(nxfield*(ny-j-1)+i+1)
! readclouds=.true.
! sumclouds=.false.
clwchn(i,j,nlev_ec-k+2,n,l)=zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
readclouds_nest(l)=.true.
sumclouds_nest(l)=.false.
endif
if(isec1(6).eq.247) then !! CIWC Cloud ice water content
ciwch(i,j,nlev_ec-k+2,n)=zsec4(nxfield*(ny-j-1)+i+1)
ciwchn(i,j,nlev_ec-k+2,n,l)=zsec4(nxn(l)*(nyn(l)-j-1)+i+1)
endif