Changeset 4c64400 in flexpart.git


Ignore:
Timestamp:
Nov 8, 2016, 4:42:27 PM (7 years ago)
Author:
Espen Sollum ATMOS <eso@…>
Branches:
master, 10.4.1_pesei, GFS_025, bugfixes+enhancements, dev, release-10, release-10.4.1, scaling-bug, univie
Children:
c8fc724
Parents:
16b61a5
Message:

Bugfix for double precision dry deposition calculation. Added (hardcoded) option to not use output kernel. Version/date string updated.

Location:
src
Files:
20 edited

Legend:

Unmodified
Added
Removed

  • src/FLEXPART.f90

    rdb712a8 r4c64400  
    6767  call gasdev1(idummy,rannumb(maxrand),rannumb(maxrand-1))
    6868
     69
    6970  ! FLEXPART version string
    7071  flexversion_major = '10' ! Major version number, also used for species file names
    71   flexversion='Version '//trim(flexversion_major)//'.0beta (2015-05-01)'
     72  flexversion='Version '//trim(flexversion_major)//'.1beta (2016-11-02)'
    7273  verbosity=0
    7374
    … …  
    384385  end do
    385386
     387  ! Inform whether output kernel is used or not
     388  !*********************************************
     389  if (lroot) then
     390    if (lnokernel) then
     391      write(*,*) "Concentrations are calculated without using kernel"
     392    else
     393      write(*,*) "Concentrations are calculated using kernel"
     394    end if
     395  end if
     396
     397
    386398  ! Calculate particle trajectories
    387399  !********************************

  • src/FLEXPART_MPI.f90

    r16b61a5 r4c64400  
    6363  ! Initialize arrays in com_mod
    6464  !*****************************
    65   if (.not.lmpreader) call com_mod_allocate_part(maxpart_mpi)
     65
     66  if(.not.(lmpreader.and.lmp_use_reader)) call com_mod_allocate_part(maxpart_mpi)
     67
    6668
    6769  ! Generate a large number of random numbers
    … …  
    7981  flexversion_major = '10' ! Major version number, also used for species file names
    8082!  flexversion='Ver. 10 Beta MPI (2015-05-01)'
    81   flexversion='Ver. '//trim(flexversion_major)//' Beta MPI (2015-05-01)'
     83  flexversion='Ver. '//trim(flexversion_major)//'.1beta MPI (2016-11-02)'
    8284  verbosity=0
    8385
    … …  
    306308  endif
    307309
    308   if (.not.lmpreader) then
     310  if (.not.(lmpreader.and.lmp_use_reader)) then
    309311    do j=1, size(itra1) ! maxpart_mpi
    310312      itra1(j)=-999999999
    … …  
    320322    endif
    321323    ! readwind process skips this step
    322     if (lmp_use_reader.and..not.lmpreader) call readpartpositions
     324    if (.not.(lmpreader.and.lmp_use_reader)) call readpartpositions
    323325  else
    324326    if (verbosity.gt.0 .and. lroot) then
    … …  
    427429  end do
    428430
     431  ! Inform whether output kernel is used or not
     432  !*********************************************
     433
     434  if (lroot) then
     435    if (lnokernel) then
     436      write(*,*) "Concentrations are calculated without using kernel"
     437    else
     438      write(*,*) "Concentrations are calculated using kernel"
     439    end if
     440  end if
    429441
    430442! Calculate particle trajectories

  • src/com_mod.f90

    r62e65c7 r4c64400  
    132132
    133133!ZHG SEP 2015 wheather or not to read clouds from GRIB
    134   logical :: readclouds
     134  logical :: readclouds=.false.
    135135!ESO DEC 2015 whether or not both clwc and ciwc are present (if so they are summed)
    136   logical :: sumclouds
     136  logical :: sumclouds=.false.
    137137
    138138  logical,dimension(maxnests) :: readclouds_nest, sumclouds_nest
    … …  
    743743
    744744  ! These variables are used to avoid having separate versions of
    745   ! files in cases where differences with MPI version is minor (eso)
     745  ! files in cases where differences with MPI version are minor (eso)
    746746  !*****************************************************************
    747747  integer :: mpi_mode=0 ! .gt. 0 if running MPI version

  • src/conccalc.f90

    r16b61a5 r4c64400  
    181181  !*****************************************************************************
    182182
    183       if ((itage.lt.10800).or.(xl.lt.0.5).or.(yl.lt.0.5).or. &
     183      if (lnokernel.or.(itage.lt.10800).or.(xl.lt.0.5).or.(yl.lt.0.5).or. &
    184184           (xl.gt.real(numxgrid-1)-0.5).or. &
    185185           (yl.gt.real(numygrid-1)-0.5)) then             ! no kernel, direct attribution to grid cell

  • src/conccalc_mpi.f90

    r5f9d14a r4c64400  
    194194  !*****************************************************************************
    195195
    196       if ((itage.lt.10800).or.(xl.lt.0.5).or.(yl.lt.0.5).or. &
     196      if (lnokernel.or.(itage.lt.10800).or.(xl.lt.0.5).or.(yl.lt.0.5).or. &
    197197           (xl.gt.real(numxgrid-1)-0.5).or. &
    198198           (yl.gt.real(numygrid-1)-0.5)) then             ! no kernel, direct attribution to grid cell

  • src/drydepokernel.f90

    re200b7a r4c64400  
    4040  !                                                                            *
    4141  !*****************************************************************************
     42  ! Changes:
     43  ! eso 10/2016: Added option to disregard kernel
     44  !
     45  !*****************************************************************************
     46
    4247
    4348  use unc_mod
    … …  
    4752  implicit none
    4853
    49   real :: x,y,deposit(maxspec),ddx,ddy,xl,yl,wx,wy,w
     54  real(dep_prec), dimension(maxspec) :: deposit
     55  real :: x,y,ddx,ddy,xl,yl,wx,wy,w
    5056  integer :: ix,jy,ixp,jyp,ks,nunc,nage,kp
    5157
    … …  
    7480  endif
    7581
     82  ! If no kernel is used, direct attribution to grid cell
     83  !******************************************************
     84
     85  if (lnokernel) then
     86    do ks=1,nspec
     87      if ((abs(deposit(ks)).gt.0).and.DRYDEPSPEC(ks)) then
     88        if ((ix.ge.0).and.(jy.ge.0).and.(ix.le.numxgrid-1).and. &
     89             (jy.le.numygrid-1)) then
     90          drygridunc(ix,jy,ks,kp,nunc,nage)= &
     91               drygridunc(ix,jy,ks,kp,nunc,nage)+deposit(ks)
     92        end if
     93      end if
     94    end do
     95  else ! use kernel
     96
    7697
    7798  ! Determine mass fractions for four grid points
    7899  !**********************************************
    79     do ks=1,nspec
     100  do ks=1,nspec
    80101
    81     if ((abs(deposit(ks)).gt.0).and.DRYDEPSPEC(ks)) then
     102   if ((abs(deposit(ks)).gt.0).and.DRYDEPSPEC(ks)) then
    82103
    83104   if ((ix.ge.0).and.(jy.ge.0).and.(ix.le.numxgrid-1).and. &
    84105        (jy.le.numygrid-1)) then
    85     w=wx*wy
    86       drygridunc(ix,jy,ks,kp,nunc,nage)= &
    87            drygridunc(ix,jy,ks,kp,nunc,nage)+deposit(ks)*w
    88       continue
    89   endif
     106     w=wx*wy
     107     drygridunc(ix,jy,ks,kp,nunc,nage)= &
     108          drygridunc(ix,jy,ks,kp,nunc,nage)+deposit(ks)*w
     109     continue
     110   endif
    90111
    91112  if ((ixp.ge.0).and.(jyp.ge.0).and.(ixp.le.numxgrid-1).and. &
    … …  
    112133  endif
    113134
    114     end do
     135  end do
     136end if
    115137
    116138end subroutine drydepokernel

  • src/drydepokernel_nest.f90

    re200b7a r4c64400  
    5050  implicit none
    5151
    52   real :: x,y,deposit(maxspec),ddx,ddy,xl,yl,wx,wy,w
     52  real(dep_prec), dimension(maxspec) :: deposit
     53  real :: x,y,ddx,ddy,xl,yl,wx,wy,w
    5354  integer :: ix,jy,ixp,jyp,ks,kp,nunc,nage
    5455

  • src/gfs_mod.f90

    r62e65c7 r4c64400  
    4242  !*********************************************
    4343
    44   integer,parameter :: nxmax=721,nymax=361,nuvzmax=64,nwzmax=64,nzmax=64
     44!  integer,parameter :: nxmax=721,nymax=361,nuvzmax=64,nwzmax=64,nzmax=64
     45  integer,parameter :: nxmax=361,nymax=181,nuvzmax=64,nwzmax=64,nzmax=64
    4546  integer,parameter :: nxshift=0     ! for GFS or FNL
    4647

  • src/gridcheck_gfs.f90

    r5f9d14a r4c64400  
    421421  !********************
    422422
    423   write(*,*)
    424   write(*,*)
    425   write(*,'(a,2i7)') 'Vertical levels in NCEP data: ', &
    426        nuvz,nwz
    427   write(*,*)
    428   write(*,'(a)') 'Mother domain:'
    429   write(*,'(a,f10.2,a1,f10.2,a,f10.2)') '  Longitude range: ', &
    430        xlon0,' to ',xlon0+(nx-1)*dx,'   Grid distance: ',dx
    431   write(*,'(a,f10.2,a1,f10.2,a,f10.2)') '  Latitude range : ', &
    432        ylat0,' to ',ylat0+(ny-1)*dy,'   Grid distance: ',dy
    433   write(*,*)
    434 
     423  if (lroot) then
     424    write(*,*)
     425    write(*,*)
     426    write(*,'(a,2i7)') 'Vertical levels in NCEP data: ', &
     427         nuvz,nwz
     428    write(*,*)
     429    write(*,'(a)') 'Mother domain:'
     430    write(*,'(a,f10.2,a1,f10.2,a,f10.2)') '  Longitude range: ', &
     431         xlon0,' to ',xlon0+(nx-1)*dx,'   Grid distance: ',dx
     432    write(*,'(a,f10.2,a1,f10.2,a,f10.2)') '  Latitude range : ', &
     433         ylat0,' to ',ylat0+(ny-1)*dy,'   Grid distance: ',dy
     434    write(*,*)
     435  end if
    435436
    436437  ! CALCULATE VERTICAL DISCRETIZATION OF ECMWF MODEL

  • src/makefile

    r341f4b7 r4c64400  
    419419        par_mod.o point_mod.o unc_mod.o xmass_mod.o
    420420timemanager_mpi.o: com_mod.o flux_mod.o mpi_mod.o oh_mod.o outg_mod.o \
    421         par_mod.o point_mod.o unc_mod.o xmass_mod.o
     421        par_mod.o point_mod.o unc_mod.o xmass_mod.o netcdf_output_mod.o
    422422unc_mod.o: par_mod.o
    423423verttransform.o: cmapf_mod.o com_mod.o par_mod.o

  • src/mean_mod.f90

    r16b61a5 r4c64400  
    2929    module procedure mean_sp
    3030    module procedure mean_dp
    31     module procedure mean_mixed_prec
     31    module procedure mean_mixed_dss
     32    module procedure mean_mixed_dsd
     33
    3234  end interface mean
    3335
    … …  
    142144  end subroutine mean_dp
    143145
    144   subroutine mean_mixed_prec(x_dp,xm,xs,number)
    145 
    146 !*****************************************************************************
    147 !                                                                            *
    148 !  This subroutine calculates mean and standard deviation of a given element.*
    149 !                                                                            *
    150 !      AUTHOR: Andreas Stohl, 25 January 1994                                *
    151 !                                                                            *
    152 !      Mixed precision version ESO 2016 (dp input, sp output)                *
     146  subroutine mean_mixed_dss(x_dp,xm,xs,number)
     147
     148!*****************************************************************************
     149!                                                                            *
     150!  This subroutine calculates mean and standard deviation of a given element.*
     151!                                                                            *
     152!      AUTHOR: Andreas Stohl, 25 January 1994                                *
     153!                                                                            *
     154!      Mixed precision version ESO 2016 (dp in, sp out, sp out)              *
    153155!*****************************************************************************
    154156!                                                                            *
    … …  
    192194    endif
    193195
    194   end subroutine mean_mixed_prec
     196  end subroutine mean_mixed_dss
     197
     198  subroutine mean_mixed_dsd(x_dp,xm,xs_dp,number)
     199
     200!*****************************************************************************
     201!                                                                            *
     202!  This subroutine calculates mean and standard deviation of a given element.*
     203!                                                                            *
     204!      AUTHOR: Andreas Stohl, 25 January 1994                                *
     205!                                                                            *
     206!      Mixed precision version ESO 2016 (dp in, sp out, dp out)              *
     207!*****************************************************************************
     208!                                                                            *
     209! Variables:                                                                 *
     210! x_dp(number)        field of input data                                    *
     211! xm                  mean                                                   *
     212! xs_dp               standard deviation                                     *
     213! number              number of elements of field x_dp                       *
     214!                                                                            *
     215! Constants:                                                                 *
     216! eps                 tiny number                                            *
     217!                                                                            *
     218!*****************************************************************************
     219
     220    use par_mod, only: sp,dp
     221
     222    implicit none
     223
     224    integer,intent(in) :: number
     225    real(dp), intent(in) :: x_dp(number)
     226    real(sp), intent(out) ::xm
     227    real(dp), intent(out) ::xs_dp
     228    real(dp) :: xl,xq,xaux
     229    real(dp),parameter :: eps=1.0e-30_dp
     230    integer :: i
     231
     232    xl=0._dp
     233    xq=0._dp
     234    do i=1,number
     235      xl=xl+x_dp(i)
     236      xq=xq+x_dp(i)*x_dp(i)
     237    end do
     238
     239    xm=xl/real(number,kind=sp)
     240
     241    xaux=xq-xl*xl/real(number,kind=dp)
     242
     243    if (xaux.lt.eps) then
     244      xs_dp=0._dp
     245    else
     246      xs_dp=sqrt(xaux/real(number-1,kind=dp))
     247    endif
     248
     249  end subroutine mean_mixed_dsd
     250
    195251end module mean_mod

  • src/mpi_mod.f90

    r16b61a5 r4c64400  
    331331! Set maxpart per process
    332332! eso 08/2016: Increase maxpart per process, in case of unbalanced distribution
    333     maxpart_mpi=int(mp_maxpart_factor*maxpart/mp_partgroup_np)
     333    maxpart_mpi=int(mp_maxpart_factor*real(maxpart)/real(mp_partgroup_np))
    334334
    335335! Do not allocate particle data arrays for readwind process
    … …  
    832832      numpart=max(numpart,ipart)
    833833
    834       deallocate(itra1_tmp,npoint_tmp,nclass_tmp,idt_tmp,itramem_tmp,itrasplit_tmp,&
    835            & xtra1_tmp,ytra1_tmp,ztra1_tmp,xmass1_tmp)
     834      deallocate(itra1_tmp,npoint_tmp,nclass_tmp,idt_tmp,itramem_tmp,&
     835           &itrasplit_tmp,xtra1_tmp,ytra1_tmp,ztra1_tmp,xmass1_tmp)
    836836
    837837      if (mp_dev_mode) then
    … …  
    23402340!
    23412341! TODO
    2342 !   take into account nested wind fields
     2342!   NB: take into account nested wind fields by using a separate
     2343!   subroutine mpif_gf_request_nest (and separate request arrays for the
     2344!   nested variables)
    23432345!
    23442346!*******************************************************************************

  • src/netcdf_output_mod.f90

    rf28aa0a r4c64400  
    1919!**********************************************************************
    2020
     21
    2122  !*****************************************************************************
    2223  !                                                                            *
    … …  
    2425  !  residence time and wet and dry deposition output.                         *
    2526  !                                                                            *
    26   !  - writeheader_netcdf generates file including all information previously    *
     27  !  - writeheader_netcdf generates file including all information previously  *
    2728  !    stored in separate header files                                         *
    28   !  - concoutput_netcdf write concentration output and wet and dry deposition   *
     29  !  - concoutput_netcdf write concentration output and wet and dry deposition *
    2930  !                                                                            *
    3031  !     Author: D. Brunner                                                     *
    … …  
    893894  gridsigmatotal=0.
    894895  gridtotalunc=0.
    895   wetgridtotal=0.
    896   wetgridsigmatotal=0.
    897   wetgridtotalunc=0.
    898   drygridtotal=0.
    899   drygridsigmatotal=0.
    900   drygridtotalunc=0.
     896  wetgridtotal=0._dep_prec
     897  wetgridsigmatotal=0._dep_prec
     898  wetgridtotalunc=0._dep_prec
     899  drygridtotal=0._dep_prec
     900  drygridsigmatotal=0._dep_prec
     901  drygridtotalunc=0._dep_prec
    901902
    902903  do ks=1,nspec
    … …  
    922923                   wetgridsigma(ix,jy),nclassunc)
    923924              ! Multiply by number of classes to get total concentration
    924               wetgrid(ix,jy)=wetgrid(ix,jy)*real(nclassunc,kind=dep_prec)
     925              wetgrid(ix,jy)=wetgrid(ix,jy)*real(nclassunc,kind=sp)
    925926              wetgridtotal=wetgridtotal+wetgrid(ix,jy)
    926927              ! Calculate standard deviation of the mean
    … …  
    946947                   drygridsigma(ix,jy),nclassunc)
    947948              ! Multiply by number of classes to get total concentration
    948               drygrid(ix,jy)=drygrid(ix,jy)*real(nclassunc)
     949              drygrid(ix,jy)=drygrid(ix,jy)*real(nclassunc,kind=sp)
    949950              drygridtotal=drygridtotal+drygrid(ix,jy)
    950951              ! Calculate standard deviation of the mean
    951952              drygridsigma(ix,jy)= &
    952953                   drygridsigma(ix,jy)* &
    953                    sqrt(real(nclassunc))
     954                   sqrt(real(nclassunc, kind=dep_prec))
    954955              drygridsigmatotal=drygridsigmatotal+ &
    955956                   drygridsigma(ix,jy)
    … …  
    10541055  if (wetgridtotal.gt.0.) wetgridtotalunc=wetgridsigmatotal/ &
    10551056       wetgridtotal
    1056   if (drygridtotal.gt.0.) drygridtotalunc=drygridsigmatotal/ &
    1057        drygridtotal
     1057  if (drygridtotal.gt.0.) drygridtotalunc=real(drygridsigmatotal/ &
     1058       drygridtotal, kind=dep_prec)
    10581059
    10591060  ! Dump of receptor concentrations
    … …  
    12981299              wetgridsigma(ix,jy)= &
    12991300                   wetgridsigma(ix,jy)* &
    1300                    sqrt(real(nclassunc))
     1301                   sqrt(real(nclassunc,kind=dep_prec))
    13011302            endif
    13021303
    … …  
    13191320              drygridsigma(ix,jy)= &
    13201321                   drygridsigma(ix,jy)* &
    1321                    sqrt(real(nclassunc))
     1322                   sqrt(real(nclassunc,kind=dep_prec))
    13221323            endif
    13231324

  • src/outg_mod.f90

    r6a678e3 r4c64400  
    2222module outg_mod
    2323
    24   use par_mod, only: dep_prec
     24  use par_mod, only: dep_prec, sp
    2525
    2626  implicit none
    … …  
    3939  real,allocatable, dimension (:,:,:) :: factor3d
    4040  real,allocatable, dimension (:,:,:) :: grid
    41   real(dep_prec),allocatable, dimension (:,:) :: wetgrid
    42   real(dep_prec),allocatable, dimension (:,:) :: drygrid
     41  real(sp),allocatable, dimension (:,:) :: wetgrid
     42  real(sp),allocatable, dimension (:,:) :: drygrid
    4343  real,allocatable, dimension (:,:,:) :: gridsigma
    4444  real(dep_prec),allocatable, dimension (:,:) :: drygridsigma

  • src/par_mod.f90

    r3b80e98 r4c64400  
    5858
    5959  integer,parameter :: dep_prec=sp
     60
     61  !****************************************************************
     62  ! Set to T to disable use of kernel for concentrations/deposition
     63  !****************************************************************
     64
     65  logical, parameter :: lnokernel=.false.
    6066
    6167  !***********************************************************

  • src/readavailable.f90

    rdb712a8 r4c64400  
    236236  do i=2,numbwf
    237237    idiff=abs(wftime(i)-wftime(i-1))
    238     if (idiff.gt.idiffmax) then
     238    if (idiff.gt.idiffmax.and.lroot) then
    239239      write(*,*) 'FLEXPART WARNING: TIME DIFFERENCE BETWEEN TWO'
    240240      write(*,*) 'WIND FIELDS IS TOO BIG FOR TRANSPORT CALCULATION.&
    241241           &'
    242242      write(*,*) 'THEREFORE, TRAJECTORIES HAVE TO BE SKIPPED.'
    243     else if (idiff.gt.idiffnorm) then
     243    else if (idiff.gt.idiffnorm.and.lroot) then
    244244      write(*,*) 'FLEXPART WARNING: TIME DIFFERENCE BETWEEN TWO'
    245245      write(*,*) 'WIND FIELDS IS BIG. THIS MAY CAUSE A DEGRADATION'

  • src/readcommand.f90

    rc04b739 r4c64400  
    170170    if (index(line,'LDIRECT') .eq. 0) then
    171171      old = .false.
    172       write(*,*) 'COMMAND in old short format, please update to namelist format'
     172      if (lroot) write(*,*) 'COMMAND in old short format, &
     173           &please update to namelist format'
    173174    else
    174175      old = .true.
    175       write(*,*) 'COMMAND in old long format, please update to namelist format'
     176      if (lroot) write(*,*) 'COMMAND in old long format, &
     177           &please update to namelist format'
    176178    endif
    177179    rewind(unitcommand)

  • src/readreleases.f90

    r05cf28d r4c64400  
    218218  rewind(unitreleases)
    219219
     220  if (nspec.gt.maxspec) goto 994
     221
    220222  ! allocate arrays of matching size for number of species (namelist output)
    221223  deallocate(mass)
    … …  
    574576  endif
    575577
     578
    576579  ! Determine the release rate (particles per second) and total number
    577580  ! of particles released during the simulation
    … …  
    633636  endif
    634637
     638  if (lroot) then
     639    write(*,FMT='(A,ES14.7)') ' Total mass released:', sum(xmass(1:numpoint,1:nspec))
     640  end if
     641
    635642  return
    636643

  • src/unc_mod.f90

    r6a678e3 r4c64400  
    2626  implicit none
    2727
    28   real,allocatable ,dimension (:,:,:,:,:,:,:) :: gridunc
     28  real,allocatable, dimension (:,:,:,:,:,:,:) :: gridunc
    2929  real,allocatable, dimension (:,:,:,:,:,:,:) :: griduncn
    3030  real(dep_prec),allocatable, dimension (:,:,:,:,:,:) :: drygridunc

  • src/wetdepokernel.f90

    re200b7a r4c64400  
    4040  !                                                                            *
    4141  !*****************************************************************************
     42  ! Changes:
     43  ! eso 10/2016: Added option to disregard kernel
     44  !
     45  !*****************************************************************************
    4246
    4347  use unc_mod
    … …  
    7377  endif
    7478
     79  ! If no kernel is used, direct attribution to grid cell
     80  !******************************************************
    7581
     82  if (lnokernel) then
     83    do ks=1,nspec
     84      if ((ix.ge.0).and.(jy.ge.0).and.(ix.le.numxgrid-1).and. &
     85           (jy.le.numygrid-1)) then
     86        wetgridunc(ix,jy,ks,kp,nunc,nage)= &
     87             wetgridunc(ix,jy,ks,kp,nunc,nage)+deposit(ks)
     88      end if
     89    end do
     90  else ! use kernel
     91   
    7692  ! Determine mass fractions for four grid points
    7793  !**********************************************
    … …  
    107123  endif
    108124  end do
     125  end if
    109126
    110127end subroutine wetdepokernel
Note: See TracChangeset for help on using the changeset viewer.
hosted by ZAMG