Changeset d6a0977 in flexpart.git for src


Ignore:
Timestamp:
Dec 14, 2015, 3:10:04 PM (8 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:
f75967d
Parents:
88d8c3d
Message:

Updates to Henrik's wet depo scheme

Location:
src
Files:
12 edited

Legend:

Unmodified
Added
Removed

  • src/FLEXPART_MPI.f90

    rdf967a99 rd6a0977  
    165165    write(*,FMT='(80("#"))')
    166166    write(*,*) '#### FLEXPART_MPI> ERROR: ', &
    167          & 'MPI version not (yet) working with backward runs. ',&
    168          & 'Use the serial version instead.'
     167         & 'MPI version not (yet) working with backward runs. '
     168    write(*,*) '#### Use the serial version instead.'
    169169    write(*,FMT='(80("#"))')
    170     stop
     170    ! call mpif_finalize
     171    ! stop
    171172  end if
    172173

  • src/com_mod.f90

    r6b22af9 rd6a0977  
    130130  ! gdomainfill             .T., if domain-filling is global, .F. if not
    131131
    132 !hg
     132!ZHG SEP 2015 wheather or not to read clouds from GRIB
    133133  logical :: readclouds
    134134
    … …  
    347347  real :: tt(0:nxmax-1,0:nymax-1,nzmax,numwfmem)
    348348  real :: qv(0:nxmax-1,0:nymax-1,nzmax,numwfmem)
    349 !hg adding cloud water
    350   real :: clwc(0:nxmax-1,0:nymax-1,nzmax,numwfmem)
    351   real :: ciwc(0:nxmax-1,0:nymax-1,nzmax,numwfmem)
     349!ZHG adding cloud water
     350  real :: clwc(0:nxmax-1,0:nymax-1,nzmax,numwfmem) !liquid   [kg/kg]
     351  real :: ciwc(0:nxmax-1,0:nymax-1,nzmax,numwfmem) !ice      [kg/kg]
     352  real :: clw(0:nxmax-1,0:nymax-1,nzmax,numwfmem)  !combined [m3/m3]
     353
    352354  real :: pv(0:nxmax-1,0:nymax-1,nzmax,numwfmem)
    353355  real :: rho(0:nxmax-1,0:nymax-1,nzmax,numwfmem)
    … …  
    362364  integer(kind=1) :: clouds(0:nxmax-1,0:nymax-1,nzmax,numwfmem)
    363365  integer :: cloudsh(0:nxmax-1,0:nymax-1,numwfmem)
    364   ! integer :: icloudbot(0:nxmax-1,0:nymax-1,numwfmem)
    365   ! integer :: icloudthck(0:nxmax-1,0:nymax-1,numwfmem)
     366
     367!ZHG Sep 2015 
     368   real :: icloud_stats(0:nxmax-1,0:nymax-1,5,numwfmem)
     369   real :: clh(0:nxmax-1,0:nymax-1,nuvzmax,numwfmem)
    366370
    367371
    … …  
    405409  real :: tropopause(0:nxmax-1,0:nymax-1,1,numwfmem)
    406410  real :: oli(0:nxmax-1,0:nymax-1,1,numwfmem)
    407 ! real :: diffk(0:nxmax-1,0:nymax-1,1,numwfmem) this is not in use?
    408 
     411! real :: diffk(0:nxmax-1,0:nymax-1,1,numwfmem) ESO: this is not in use?
     412! logical :: beneath_cloud=.true.
    409413  ! ps                   surface pressure
    410414  ! sd                   snow depth
    … …  
    667671
    668672
     673  !CGZ-lifetime
     674  real, allocatable, dimension(:,:) ::checklifetime, species_lifetime
     675  !CGZ-lifetime
     676
    669677  ! numpart                 actual number of particles in memory
    670678  ! itra1 (maxpart) [s]     temporal positions of the particles
    … …  
    759767         & idt(nmpart),itramem(nmpart),itrasplit(nmpart),&
    760768         & xtra1(nmpart),ytra1(nmpart),ztra1(nmpart),&
    761          & xmass1(nmpart, maxspec))
     769         & xmass1(nmpart, maxspec),&
     770         & checklifetime(nmpart,maxspec), species_lifetime(maxspec,2))!CGZ-lifetime
     771
    762772
    763773    allocate(uap(nmpart),ucp(nmpart),uzp(nmpart),us(nmpart),&

  • src/gridcheck.f90

    r5f9d14a rd6a0977  
    195195  elseif ((parCat.eq.1).and.(parNum.eq.0).and.(typSurf.eq.105)) then ! Q
    196196    isec1(6)=133         ! indicatorOfParameter
    197 !hg
     197!ZHG FOR CLOUDS FROM GRIB
    198198  elseif ((parCat.eq.1).and.(parNum.eq.83).and.(typSurf.eq.105)) then ! clwc
    199199    isec1(6)=246         ! indicatorOfParameter
    200200  elseif ((parCat.eq.1).and.(parNum.eq.84).and.(typSurf.eq.105)) then ! ciwc
    201201    isec1(6)=247         ! indicatorOfParameter
    202 !hg end
     202!ZHG end
    203203  elseif ((parCat.eq.3).and.(parNum.eq.0).and.(typSurf.eq.1)) then !SP
    204204    isec1(6)=134         ! indicatorOfParameter

  • src/mpi_mod.f90

    ra1f4dd6 rd6a0977  
    8989! MPI tags/requests for send/receive operation
    9090  integer :: tm1
    91   integer, parameter :: nvar_async=27 !29 :DBG:
     91  integer, parameter :: nvar_async=26 !27 !29 :DBG:
    9292!integer, dimension(:), allocatable :: tags
    9393  integer, dimension(:), allocatable :: reqs
    … …  
    119119  logical, parameter :: mp_dbg_out = .false.
    120120  logical, parameter :: mp_time_barrier=.true.
    121   logical, parameter :: mp_measure_time=.true.
     121  logical, parameter :: mp_measure_time=.false.
    122122  logical, parameter :: mp_exact_numpart=.true.
    123123
    … …  
    207207      write(*,FMT='(80("#"))')
    208208! Force "syncronized" version if all processes will call getfields
    209     else if (.not.lmp_sync.and.mp_np.lt.read_grp_min) then
     209    else if ((.not.lmp_sync.and.mp_np.lt.read_grp_min).or.(mp_np.eq.1)) then
    210210      if (lroot) then
    211211        write(*,FMT='(80("#"))')
    … …  
    963963! cloud water/ice:
    964964    if (readclouds) then
    965       call MPI_Bcast(clwc(:,:,:,li:ui),d3s1,mp_pp,id_read,MPI_COMM_WORLD,mp_ierr)
    966       if (mp_ierr /= 0) goto 600
    967       call MPI_Bcast(ciwc(:,:,:,li:ui),d3s1,mp_pp,id_read,MPI_COMM_WORLD,mp_ierr)
    968       if (mp_ierr /= 0) goto 600
     965      call MPI_Bcast(icloud_stats(:,:,:,li:ui),d2_size1*5,mp_pp,id_read,MPI_COMM_WORLD,mp_ierr)
     966      if (mp_ierr /= 0) goto 600
     967
     968      ! call MPI_Bcast(clwc(:,:,:,li:ui),d3s1,mp_pp,id_read,MPI_COMM_WORLD,mp_ierr)
     969      ! if (mp_ierr /= 0) goto 600
     970      ! call MPI_Bcast(ciwc(:,:,:,li:ui),d3s1,mp_pp,id_read,MPI_COMM_WORLD,mp_ierr)
     971      ! if (mp_ierr /= 0) goto 600
    969972    end if
    970973
    … …  
    11891192!
    11901193! TODO
    1191 !   Transfer cloud water/ice
    11921194!
    11931195!*******************************************************************************
    … …  
    13341336      if (readclouds) then
    13351337        i=i+1
    1336         call MPI_Isend(clwc(:,:,:,mind),d3s1,mp_pp,dest,tm1,&
     1338        call MPI_Isend(icloud_stats(:,:,:,mind),d2s1*5,mp_pp,dest,tm1,&
    13371339             &MPI_COMM_WORLD,reqs(i),mp_ierr)
    13381340        if (mp_ierr /= 0) goto 600
    1339         i=i+1
    1340 
    1341         call MPI_Isend(ciwc(:,:,:,mind),d3s1,mp_pp,dest,tm1,&
    1342              &MPI_COMM_WORLD,reqs(i),mp_ierr)
    1343         if (mp_ierr /= 0) goto 600
    1344 ! else
    1345 !   i=i+2
     1341
     1342        ! call MPI_Isend(clwc(:,:,:,mind),d3s1,mp_pp,dest,tm1,&
     1343        !      &MPI_COMM_WORLD,reqs(i),mp_ierr)
     1344        ! if (mp_ierr /= 0) goto 600
     1345        ! i=i+1
     1346
     1347        ! call MPI_Isend(ciwc(:,:,:,mind),d3s1,mp_pp,dest,tm1,&
     1348        !      &MPI_COMM_WORLD,reqs(i),mp_ierr)
     1349        ! if (mp_ierr /= 0) goto 600
     1350
    13461351      end if
    1347 
    13481352    end do
    13491353
    … …  
    13711375!
    13721376! TODO
    1373 !   Transfer cloud water/ice
    13741377!
    13751378!*******************************************************************************
    … …  
    13901393!    integer :: d3s1,d3s2,d2s1,d2s2
    13911394!*******************************************************************************
    1392 
    1393 ! :TODO: don't need these
    1394 ! d3s1=d3_size1
    1395 ! d3s2=d3_size2
    1396 ! d2s1=d2_size1
    1397 ! d2s2=d2_size2
    13981395
    13991396! At the time this immediate receive is posted, memstat is the state of
    … …  
    15451542    if (readclouds) then
    15461543      j=j+1
    1547       call MPI_Irecv(clwc(:,:,:,mind),d3s1,mp_pp,id_read,MPI_ANY_TAG,&
    1548            &MPI_COMM_WORLD,reqs(j),mp_ierr)   
    1549       if (mp_ierr /= 0) goto 600
    1550       j=j+1
    1551       call MPI_Irecv(ciwc(:,:,:,mind),d3s1,mp_pp,id_read,MPI_ANY_TAG,&
    1552            &MPI_COMM_WORLD,reqs(j),mp_ierr)   
    1553       if (mp_ierr /= 0) goto 600
     1544
     1545      call MPI_Irecv(icloud_stats(:,:,:,mind),d2s1*5,mp_pp,id_read,MPI_ANY_TAG,&
     1546           &MPI_COMM_WORLD,reqs(j),mp_ierr)
     1547      if (mp_ierr /= 0) goto 600
     1548
     1549      ! call MPI_Irecv(clwc(:,:,:,mind),d3s1,mp_pp,id_read,MPI_ANY_TAG,&
     1550      !      &MPI_COMM_WORLD,reqs(j),mp_ierr)   
     1551      ! if (mp_ierr /= 0) goto 600
     1552      ! j=j+1
     1553      ! call MPI_Irecv(ciwc(:,:,:,mind),d3s1,mp_pp,id_read,MPI_ANY_TAG,&
     1554      !      &MPI_COMM_WORLD,reqs(j),mp_ierr)   
     1555      ! if (mp_ierr /= 0) goto 600
     1556
    15541557    end if
    15551558
    … …  
    20892092    implicit none
    20902093
    2091     integer, parameter :: li=1, ui=3 ! wfmem indices (i.e, operate on entire field)
     2094    integer :: li=1, ui=2 ! wfmem indices (i.e, operate on entire field)
     2095
     2096    if (.not.lmp_sync) ui=3
    20922097
    20932098

  • src/par_mod.f90

    r4c4261c rd6a0977  
    7676  real,parameter :: hmixmin=100., hmixmax=4500., turbmesoscale=0.16
    7777  real,parameter :: d_trop=50., d_strat=0.1
    78 
     78  real,parameter :: rho_water=1000 !ZHG 2015 [kg/m3]
    7979  ! karman                  Karman's constant
    8080  ! href [m]                Reference height for dry deposition
    … …  
    226226  ! ---------
    227227  integer,parameter :: maxwf=50000, maxtable=1000, numclass=13, ni=11
    228   !integer,parameter :: numwfmem=2 ! Serial version/MPI with 2 fields
    229   integer,parameter :: numwfmem=3 ! MPI with 3 fields
     228  integer,parameter :: numwfmem=2 ! Serial version/MPI with 2 fields
     229  !integer,parameter :: numwfmem=3 ! MPI with 3 fields
    230230
    231231  ! maxwf                   maximum number of wind fields to be used for simulation

  • src/pathnames

    rde7afbd rd6a0977  
    22../output/
    33/
    4 /xnilu_wrk/flex_wrk/WIND_FIELDS/AVAILABLE_ECMWF_OPER_fields_global/
     4/xnilu_wrk/flex_wrk/WIND_FIELDS/AVAILABLE_ECMWF_OPER_fields_global
    55============================================
     6
     7/xnilu_wrk/flex_wrk/zhg/ECMWF_CLWC/Availables
    68
    79

  • src/readwind.f90

    r5f9d14a rd6a0977  
    106106  hflswitch=.false.
    107107  strswitch=.false.
    108 !hg
     108!ZHG test the grib fields that have lcwc without using them
    109109!  readcloud=.false.
    110 !hg end
     110
    111111  levdiff2=nlev_ec-nwz+1
    112112  iumax=0
    … …  
    190190    elseif ((parCat.eq.1).and.(parNum.eq.0).and.(typSurf.eq.105)) then ! Q
    191191      isec1(6)=133         ! indicatorOfParameter
    192 !hg
     192!ZHG ! READ CLOUD FIELD : assume these to be added together to one variable
    193193    elseif ((parCat.eq.1).and.(parNum.eq.83).and.(typSurf.eq.105)) then ! clwc
    194194      isec1(6)=246         ! indicatorOfParameter
    195     elseif ((parCat.eq.1).and.(parNum.eq.84).and.(typSurf.eq.105)) then ! ciwc
    196       isec1(6)=247         ! indicatorOfParameter
    197  !hg end
     195! ICE AND WATER IS ADDED TOGETHER IN NEW WINDFIELDS
     196!    elseif ((parCat.eq.1).and.(parNum.eq.84).and.(typSurf.eq.105)) then ! ciwc
     197!      isec1(6)=247         ! indicatorOfParameter
     198!ZHG end
    198199    elseif ((parCat.eq.3).and.(parNum.eq.0).and.(typSurf.eq.1)) then !SP
    199200      isec1(6)=134         ! indicatorOfParameter
    … …  
    226227    elseif ((parCat.eq.4).and.(parNum.eq.9).and.(typSurf.eq.1)) then ! SR
    227228      isec1(6)=176         ! indicatorOfParameter
    228     elseif ((parCat.eq.2).and.(parNum.eq.17) .or. parId .eq. 180) then ! EWSS
     229!    elseif ((parCat.eq.2).and.(parNum.eq.17) .or. parId .eq. 180) then ! EWSS --wrong
     230    elseif ((parCat.eq.2).and.(parNum.eq.38) .or. parId .eq. 180) then ! EWSS --correct
    229231      isec1(6)=180         ! indicatorOfParameter
    230     elseif ((parCat.eq.2).and.(parNum.eq.18) .or. parId .eq. 181) then ! NSSS
     232!    elseif ((parCat.eq.2).and.(parNum.eq.18) .or. parId .eq. 181) then ! NSSS --wrong
     233    elseif ((parCat.eq.2).and.(parNum.eq.37) .or. parId .eq. 181) then ! NSSS --correct
    231234      isec1(6)=181         ! indicatorOfParameter
    232235    elseif ((parCat.eq.3).and.(parNum.eq.4)) then ! ORO
    … …  
    356359      if(isec1(6).eq.131) iumax=max(iumax,nlev_ec-k+1)
    357360      if(isec1(6).eq.135) iwmax=max(iwmax,nlev_ec-k+1)
    358 !hg READING CLOUD FIELDS ASWELL
    359       if(isec1(6).eq.246) then  !! CLWC  Cloud liquid water content
     361!ZHG READING CLOUD FIELDS ASWELL
     362      if(isec1(6).eq.246) then  !! CLWC  Cloud liquid water content [kg/kg]
    360363        clwch(i,j,nlev_ec-k+2,n)=zsec4(nxfield*(ny-j-1)+i+1)
    361364        readclouds = .true.
    362         !write(*,*) 'found water!'
     365!if (clwch(i,j,nlev_ec-k+2,n) .gt. 0)        write(*,*) 'readwind: found water!', clwch(i,j,nlev_ec-k+2,n)
    363366      endif
    364       if(isec1(6).eq.247) then  !! CIWC  Cloud ice water content
    365         ciwch(i,j,nlev_ec-k+2,n)=zsec4(nxfield*(ny-j-1)+i+1)
     367!      if(isec1(6).eq.247) then  !! CIWC  Cloud ice water content
     368!        ciwch(i,j,nlev_ec-k+2,n)=zsec4(nxfield*(ny-j-1)+i+1)
    366369        !write(*,*) 'found ice!'
    367       endif
    368 !hg end
     370!      endif
     371!ZHG end
    369372
    370373    end do
    … …  
    422425    call shift_field(vvh,nxfield,ny,nuvzmax,nuvz,1,1)
    423426    call shift_field(wwh,nxfield,ny,nwzmax,nwz,1,1)
    424 !hg
     427!ZHG
    425428    call shift_field(clwch,nxfield,ny,nuvzmax,nuvz,2,n)
    426429    call shift_field(ciwch,nxfield,ny,nuvzmax,nuvz,2,n)
    427 !hg end
     430!ZHG end
    428431
    429432  endif

  • src/readwind_mpi.f90

    r5f9d14a rd6a0977  
    200200    elseif ((parCat.eq.1).and.(parNum.eq.0).and.(typSurf.eq.105)) then ! Q
    201201      isec1(6)=133         ! indicatorOfParameter
    202 !hg
     202!ZHG ! READ CLOUD FIELD : assume these to be added together to one variable
    203203    elseif ((parCat.eq.1).and.(parNum.eq.83).and.(typSurf.eq.105)) then ! clwc
    204204      isec1(6)=246         ! indicatorOfParameter
    205     elseif ((parCat.eq.1).and.(parNum.eq.84).and.(typSurf.eq.105)) then ! ciwc
    206       isec1(6)=247         ! indicatorOfParameter
    207  !hg end
     205! ICE AND WATER IS ADDED TOGETHER IN NEW WINDFIELDS
     206!    elseif ((parCat.eq.1).and.(parNum.eq.84).and.(typSurf.eq.105)) then ! ciwc
     207!      isec1(6)=247         ! indicatorOfParameter
     208!ZHG end
    208209    elseif ((parCat.eq.3).and.(parNum.eq.0).and.(typSurf.eq.1)) then !SP
    209210      isec1(6)=134         ! indicatorOfParameter
    … …  
    368369      if(isec1(6).eq.131) iumax=max(iumax,nlev_ec-k+1)
    369370      if(isec1(6).eq.135) iwmax=max(iwmax,nlev_ec-k+1)
    370 !hg READING CLOUD FIELDS ASWELL
    371       if(isec1(6).eq.246) then  !! CLWC  Cloud liquid water content
     371!ZHG READING CLOUD FIELDS ASWELL
     372      if(isec1(6).eq.246) then  !! CLWC  Cloud liquid water content [kg/kg]
    372373        clwch(i,j,nlev_ec-k+2,n)=zsec4(nxfield*(ny-j-1)+i+1)
    373374        readclouds = .true.
    374         !write(*,*) 'found water!'
     375!if (clwch(i,j,nlev_ec-k+2,n) .gt. 0)        write(*,*) 'readwind: found water!', clwch(i,j,nlev_ec-k+2,n)
    375376      endif
    376 
    377       if(isec1(6).eq.247) then  !! CIWC  Cloud ice water content
    378         ciwch(i,j,nlev_ec-k+2,n)=zsec4(nxfield*(ny-j-1)+i+1)
     377!      if(isec1(6).eq.247) then  !! CIWC  Cloud ice water content
     378!        ciwch(i,j,nlev_ec-k+2,n)=zsec4(nxfield*(ny-j-1)+i+1)
    379379        !write(*,*) 'found ice!'
    380       endif
    381 !hg end
     380!      endif
     381!ZHG end
    382382
    383383    end do

  • src/timemanager.f90

    r78e62dc rd6a0977  
    136136  !**********************************************************************
    137137
     138!ZHG 2015
     139!CGZ-lifetime: set lifetime to 0
     140  checklifetime(:,:)=0
     141  species_lifetime(:,:)=0
     142  print*, 'Initialized lifetime'
     143!CGZ-lifetime: set lifetime to 0
     144 
     145
    138146
    139147  if (verbosity.gt.0) then
    … …  
    411419        if (iflux.eq.1) call fluxoutput(itime)
    412420        write(*,45) itime,numpart,gridtotalunc,wetgridtotalunc,drygridtotalunc
     421 
     422        !CGZ-lifetime: output species lifetime
     423!ZHG
     424        write(*,*) 'Overview species lifetime in days', &
     425             real((species_lifetime(:,1)/species_lifetime(:,2))/real(3600.0*24.0))
     426        write(*,*) 'all info:',species_lifetime
     427!ZHG
     428        !CGZ-lifetime: output species lifetime
     429
    413430        !write(*,46) float(itime)/3600,itime,numpart
    41443145      format(i9,' SECONDS SIMULATED: ',i8, ' PARTICLES:    Uncertainty: ',3f7.3)
    … …  
    575592                   xmassfract=max(xmassfract,real(npart(npoint(j)))* &
    576593                   xmass1(j,ks)/xmass(npoint(j),ks))
     594!ZHG 2015
     595                  !CGZ-lifetime: Check mass fraction left/save lifetime
     596                   if(real(npart(npoint(j)))*xmass1(j,ks)/xmass(npoint(j),ks).lt.0.01.and.checklifetime(j,ks).eq.0.)then
     597                       !Mass below 1% of initial >register lifetime
     598                       checklifetime(j,ks)=abs(itra1(j)-itramem(j))
     599                       species_lifetime(ks,1)=species_lifetime(ks,1)+abs(itra1(j)-itramem(j))
     600                       species_lifetime(ks,2)= species_lifetime(ks,2)+1
     601                   endif
     602                   !CGZ-lifetime: Check mass fraction left/save lifetime
     603!ZHG 2015
    577604            else
    578605              xmassfract=1.

  • src/timemanager_mpi.f90

    ra1f4dd6 rd6a0977  
    139139!**********************************************************************
    140140
     141
     142!  itime=0
     143  if (lroot) then
     144  !  write(*,45) itime,numpart*mp_partgroup_np,gridtotalunc,wetgridtotalunc,drygridtotalunc
     145    write(*,46) float(itime)/3600,itime,numpart*mp_partgroup_np
     146   
     147    if (verbosity.gt.0) then
     148      write (*,*) 'timemanager> starting simulation'
     149    end if
     150  end if ! (lroot)
     151
     152!CGZ-lifetime: set lifetime to 0
     153  checklifetime(:,:)=0
     154  species_lifetime(:,:)=0
     155  print*, 'Initialized lifetime'
     156!CGZ-lifetime: set lifetime to 0
     157
     158
     159
    141160  do itime=0,ideltas,lsynctime
    142161
    … …  
    544563        endif
    545564       
    546         if (lroot) then
     565        !CGZ-lifetime: output species lifetime
     566        ! if (lroot) then
     567        !   write(*,*) 'Overview species lifetime in days', &
     568        !        real((species_lifetime(:,1)/species_lifetime(:,2))/real(3600.0*24.0))
     569        !   write(*,*) 'all info:',species_lifetime
    547570          write(*,45) itime,numpart_tot_mpi,gridtotalunc,&
    548571               &wetgridtotalunc,drygridtotalunc
    549           if (verbosity.gt.0) then
    550             write (*,*) 'timemanager> starting simulation'
    551           end if
    552         end if
     572        !   if (verbosity.gt.0) then
     573        !     write (*,*) 'timemanager> starting simulation'
     574        !   end if
     575        ! end if
    553576
    55457745      format(i13,' SECONDS SIMULATED: ',i13, ' PARTICLES:    Uncertainty: ',3f7.3)
    … …  
    729752
    730753            if (mdomainfill.eq.0) then
    731               if (xmass(npoint(j),ks).gt.0.) &
     754              if (xmass(npoint(j),ks).gt.0.)then
    732755                   xmassfract=max(xmassfract,real(npart(npoint(j)))* &
    733756                   xmass1(j,ks)/xmass(npoint(j),ks))
     757                   
     758                   !CGZ-lifetime: Check mass fraction left/save lifetime
     759                   if(lroot.and.real(npart(npoint(j)))*xmass1(j,ks)/xmass(npoint(j),ks).lt.0.01.and.checklifetime(j,ks).eq.0.)then
     760                       !Mass below 1% of initial >register lifetime
     761                       checklifetime(j,ks)=abs(itra1(j)-itramem(j))
     762
     763                       species_lifetime(ks,1)=species_lifetime(ks,1)+abs(itra1(j)-itramem(j))
     764                       species_lifetime(ks,2)= species_lifetime(ks,2)+1
     765                   endif
     766                   !CGZ-lifetime: Check mass fraction left/save lifetime
     767                   
     768              endif
    734769            else
    735770              xmassfract=1.
    … …  
    737772          end do
    738773
    739           if (xmassfract.lt.0.0001) then   ! terminate all particles carrying less mass
     774          if (xmassfract.lt.0.00005 .and. sum(real(npart(npoint(j)))*xmass1(j,:)).lt.1.0) then   ! terminate all particles carrying less mass
     775          !            print*,'terminated particle ',j,' for small mass (', sum(real(npart(npoint(j)))* &
     776          !         xmass1(j,:)), ' of ', sum(xmass(npoint(j),:)),')'
    740777            itra1(j)=-999999999
    741778          endif
    … …  
    758795                 call initial_cond_calc(itime+lsynctime,j)
    759796            itra1(j)=-999999999
     797            !print*, 'terminated particle ',j,'for age'
    760798          endif
    761799        endif
    … …  
    828866  endif
    829867  deallocate(gridunc)
    830   deallocate(xpoint1,xpoint2,ypoint1,ypoint2,zpoint1,zpoint2,xmass)
     868  deallocate(xpoint1,xpoint2,ypoint1,ypoint2,zpoint1,zpoint2,xmass, checklifetime)
    831869  deallocate(ireleasestart,ireleaseend,npart,kindz)
    832870  deallocate(xmasssave)

  • src/verttransform.f90

    r4d45639 rd6a0977  
    9797  logical :: init = .true.
    9898
    99 !hg
     99  !ZHG SEP 2014 tests
    100100  integer :: cloud_ver,cloud_min, cloud_max
     101  integer ::teller(5), convpteller=0, lspteller=0
    101102  real :: cloud_col_wat, cloud_water
    102 !hg temporary variables for testing
     103  !ZHG 2015 temporary variables for testing
    103104  real :: rcw(0:nxmax-1,0:nymax-1)
    104105  real :: rpc(0:nxmax-1,0:nymax-1)
    105 !hg
     106  character(len=60) :: zhgpath='/xnilu_wrk/flex_wrk/zhg/'
     107  character(len=60) :: fnameA,fnameB,fnameC,fnameD,fnameE,fnameF,fnameG,fnameH
     108  CHARACTER(LEN=3)  :: aspec
     109  integer :: virr=0
     110  real :: tot_cloud_h
     111!ZHG
    106112
    107113!*************************************************************************
    … …  
    561567
    562568
     569!*********************************************************************************** 
     570if (readclouds) then !HG METHOD
     571! The method is loops all grids vertically and constructs the 3D matrix for clouds
     572! Cloud top and cloud bottom gid cells are assigned as well as the total column
     573! cloud water. For precipitating grids, the type and whether it is in or below
     574! cloud scavenging are assigned with numbers 2-5 (following the old metod).
     575! Distinction is done for lsp and convp though they are treated the same in regards
     576! to scavenging. Also clouds that are not precipitating are defined which may be
     577! to include future cloud processing by non-precipitating-clouds.
     578!***********************************************************************************
    563579  if (readclouds) write(*,*) 'using cloud water from ECMWF'
    564 !  if (.not.readclouds) write(*,*) 'using cloud water from parameterization'
    565 
    566   rcw(:,:)=0
    567   rpc(:,:)=0
    568 
     580  clw(:,:,:,n)=0
     581  icloud_stats(:,:,:,n)=0
     582  clouds(:,:,:,n)=0
     583  do jy=0,nymin1
     584   do ix=0,nxmin1
     585    lsp=lsprec(ix,jy,1,n)
     586    convp=convprec(ix,jy,1,n)
     587    prec=lsp+convp
     588    tot_cloud_h=0
     589    ! Find clouds in the vertical
     590    do kz=1, nz-1 !go from top to bottom
     591     if (clwc(ix,jy,kz,n).gt.0) then     
     592      ! assuming rho is in kg/m3 and hz in m gives: kg/kg * kg/m3 *m3/kg /m = m2/m3
     593      clw(ix,jy,kz,n)=(clwc(ix,jy,kz,n)*rho(ix,jy,kz,n))*(height(kz+1)-height(kz))
     594      tot_cloud_h=tot_cloud_h+(height(kz+1)-height(kz))
     595      icloud_stats(ix,jy,4,n)= icloud_stats(ix,jy,4,n)+clw(ix,jy,kz,n)          ! Column cloud water [m3/m3]
     596      icloud_stats(ix,jy,3,n)= min(height(kz+1),height(kz))                     ! Cloud BOT height stats      [m]
     597!ZHG 2015 extra for testing
     598      clh(ix,jy,kz,n)=height(kz+1)-height(kz)
     599      icloud_stats(ix,jy,1,n)=icloud_stats(ix,jy,1,n)+(height(kz+1)-height(kz)) ! Cloud total vertical extent [m]
     600      icloud_stats(ix,jy,2,n)= max(icloud_stats(ix,jy,2,n),height(kz))          ! Cloud TOP height            [m]
     601!ZHG
     602    endif
     603    end do
     604
     605   ! If Precipitation. Define removal type in the vertical
     606  if ((lsp.gt.0.01).or.(convp.gt.0.01)) then ! cloud and precipitation
     607
     608    do kz=nz,1,-1 !go Bottom up!
     609     if (clw(ix,jy,kz,n).gt. 0) then ! is in cloud
     610        cloudsh(ix,jy,n)=cloudsh(ix,jy,n)+height(kz)-height(kz-1)
     611        clouds(ix,jy,kz,n)=1                               ! is a cloud
     612        if (lsp.ge.convp) then
     613           clouds(ix,jy,kz,n)=3                            ! lsp in-cloud
     614        else
     615          clouds(ix,jy,kz,n)=2                             ! convp in-cloud
     616        endif                                              ! convective or large scale
     617     elseif((clw(ix,jy,kz,n).le.0) .and. (icloud_stats(ix,jy,3,n).ge.height(kz)) ) then ! is below cloud
     618        if (lsp.ge.convp) then
     619          clouds(ix,jy,kz,n)=5                             ! lsp dominated washout
     620        else
     621          clouds(ix,jy,kz,n)=4                             ! convp dominated washout
     622        endif                                              ! convective or large scale
     623     endif
     624
     625     if (height(kz).ge. 19000) then                        ! set a max height for removal
     626        clouds(ix,jy,kz,n)=0
     627     endif !clw>0
     628    end do !nz
     629   endif ! precipitation
     630  end do
     631 end do
     632!**************************************************************************
     633 else       ! use old definitions
     634!**************************************************************************
     635!   create a cloud and rainout/washout field, clouds occur where rh>80%
     636!   total cloudheight is stored at level 0
     637 if (.not.readclouds) write(*,*) 'using cloud water from Parameterization'
    569638  do jy=0,nymin1
    570639    do ix=0,nxmin1
    … …  
    603672      end do
    604673!END OLD METHOD
     674      end do
     675     end do
     676endif !readclouds
     677
     678     !********* TEST ***************
     679     ! WRITE OUT SOME TEST VARIABLES
     680     !********* TEST ************'**
     681!teller(:)=0
     682!virr=virr+1
     683!WRITE(aspec, '(i3.3)'), virr
     684
     685!if (readclouds) then
     686!fnameH=trim(zhgpath)//trim(aspec)//'Vertical_placement.txt'
     687!else
     688!fnameH=trim(zhgpath)//trim(aspec)//'Vertical_placement_old.txt'
     689!endif
     690!
     691!OPEN(UNIT=118, FILE=fnameH,FORM='FORMATTED',STATUS = 'UNKNOWN')
     692!do kz_inv=1,nz-1
     693!  kz=nz-kz_inv+1
     694!  !kz=91
     695!  do jy=0,nymin1
     696!     do ix=0,nxmin1
     697!          if (clouds(ix,jy,kz,n).eq.1) teller(1)=teller(1)+1 ! no precipitation cloud
     698!          if (clouds(ix,jy,kz,n).eq.2) teller(2)=teller(2)+1 ! convp dominated rainout
     699!          if (clouds(ix,jy,kz,n).eq.3) teller(3)=teller(3)+1 ! lsp dominated rainout
     700!          if (clouds(ix,jy,kz,n).eq.4) teller(4)=teller(4)+1 ! convp dominated washout
     701!          if (clouds(ix,jy,kz,n).eq.5) teller(5)=teller(5)+1 ! lsp dominated washout
     702!         
     703!        !  write(*,*) height(kz),teller
     704!     end do
     705!  end do
     706!  write(118,*) height(kz),teller
     707!  teller(:)=0
     708!end do
     709!teller(:)=0
     710!write(*,*) teller
     711!write(*,*) aspec
     712!
     713!fnameA=trim(zhgpath)//trim(aspec)//'cloudV.txt'
     714!fnameB=trim(zhgpath)//trim(aspec)//'cloudT.txt'
     715!fnameC=trim(zhgpath)//trim(aspec)//'cloudB.txt'
     716!fnameD=trim(zhgpath)//trim(aspec)//'cloudW.txt'
     717!fnameE=trim(zhgpath)//trim(aspec)//'old_cloudV.txt'
     718!fnameF=trim(zhgpath)//trim(aspec)//'lsp.txt'
     719!fnameG=trim(zhgpath)//trim(aspec)//'convp.txt'
     720!if (readclouds) then
     721!OPEN(UNIT=111, FILE=fnameA,FORM='FORMATTED',STATUS = 'UNKNOWN')
     722!OPEN(UNIT=112, FILE=fnameB,FORM='FORMATTED',STATUS = 'UNKNOWN')
     723!OPEN(UNIT=113, FILE=fnameC,FORM='FORMATTED',STATUS = 'UNKNOWN')
     724!OPEN(UNIT=114, FILE=fnameD,FORM='FORMATTED',STATUS = 'UNKNOWN')
     725!else
     726!OPEN(UNIT=115, FILE=fnameE,FORM='FORMATTED',STATUS = 'UNKNOWN')
     727!OPEN(UNIT=116, FILE=fnameF,FORM='FORMATTED',STATUS = 'UNKNOWN')
     728!OPEN(UNIT=117, FILE=fnameG,FORM='FORMATTED',STATUS = 'UNKNOWN')
     729!endif
     730!
     731!do ix=0,nxmin1
     732!if (readclouds) then
     733!write(111,*) (icloud_stats(ix,jy,1,n),jy=0,nymin1)
     734!write(112,*) (icloud_stats(ix,jy,2,n),jy=0,nymin1)
     735!write(113,*) (icloud_stats(ix,jy,3,n),jy=0,nymin1)
     736!write(114,*) (icloud_stats(ix,jy,4,n),jy=0,nymin1)
     737!else
     738!write(115,*) (cloudsh(ix,jy,n),jy=0,nymin1)    !integer
     739!write(116,*) (lsprec(ix,jy,1,n),jy=0,nymin1)   !7.83691406E-02
     740!write(117,*) (convprec(ix,jy,1,n),jy=0,nymin1) !5.38330078E-02
     741!endif
     742!end do
     743!
     744!if (readclouds) then
     745!CLOSE(111)
     746!CLOSE(112)
     747!CLOSE(113)
     748!CLOSE(114)
     749!else
     750!CLOSE(115)
     751!CLOSE(116)
     752!CLOSE(117)
     753!endif
     754!
     755!END ********* TEST *************** END
     756! WRITE OUT SOME TEST VARIABLES
     757!END ********* TEST *************** END
     758
    605759
    606760! PS 2012
    … …  
    613767!            lconvectprec = .true.
    614768!           endif
    615 !HG METHOD
    616 !readclouds =.true.
    617 !      if (readclouds) then
    618 !hg added APR 2014  Cloud Water=clwc(ix,jy,kz,n)  Cloud Ice=ciwc(ix,jy,kz,n)
    619 !hg Use the cloud water variables to determine existence of clouds. This makes the PS code obsolete
    620 !        cloud_min=99999
    621 !        cloud_max=-1
    622 !        cloud_col_wat=0
    623 
    624 !        do kz=1, nz
    625 !         !clw & ciw are given in kg/kg 
    626           ! cloud_water=clwc(ix,jy,kz,n)+ciwc(ix,jy,kz,n)
    627           ! if (cloud_water .gt. 0) then
    628           !   cloud_min=min(nint(height(kz)),cloud_min) !hg needs reset each grid
    629           !   cloud_max=max(nint(height(kz)),cloud_max) !hg needs reset each grid
    630           !   cloud_col_wat=cloud_col_wat+cloud_water !hg needs reset each grid
    631           ! endif
    632           ! cloud_ver=max(0,cloud_max-cloud_min)
    633 
    634           ! icloudbot(ix,jy,n)=cloud_min
    635           ! icloudthck(ix,jy,n)=cloud_ver
    636           ! rcw(ix,jy)=cloud_col_wat
    637           ! rpc(ix,jy)=prec
    638 !write(*,*) 'Using clouds from ECMWF' !hg END Henrik Code
    639 !END HG METHOD
    640 
    641 
    642 
    643769!      else ! windfields does not contain cloud data
    644770!          rhmin = 0.90 ! standard condition for presence of clouds
    … …  
    698824!      endif !hg read clouds
    699825
    700     end do
    701   end do
    702 
    703 !do 102 kz=1,nuvz
    704 !write(an,'(i02)') kz+10
    705 !write(*,*) nuvz,nymin1,nxmin1,'--',an,'--'
    706 !open(4,file='/nilu_wrk2/sec/cloudtest/cloud'//an,form='formatted')
    707 !do 101 jy=0,nymin1
    708 !    write(4,*) (clouds(ix,jy,kz,n),ix=1,nxmin1)
    709 !101   continue
    710 ! close(4)
    711 !102   continue
    712 
    713 ! open(4,file='/nilu_wrk2/sec/cloudtest/height',form='formatted')
    714 ! do 103 jy=0,nymin1
    715 !     write (4,*)
    716 !+       (height(kz),kz=1,nuvz)
    717 !103   continue
    718 ! close(4)
    719 
    720 !open(4,file='/nilu_wrk2/sec/cloudtest/p',form='formatted')
    721 ! do 104 jy=0,nymin1
    722 !     write (4,*)
    723 !+       (r_air*tt(ix,jy,1,n)*rho(ix,jy,1,n),ix=1,nxmin1)
    724 !104   continue
    725 ! close(4)
     826
     827
    726828
    727829!eso measure CPU time

  • src/wetdepo.f90

    r0f20c31 rd6a0977  
    9090
    9191  integer :: blc_count, inc_count
     92  real    :: Si_dummy, wetscav_dummy
    9293
    9394
    … …  
    190191      clouds_h=cloudsh(ix,jy,n)
    191192    else
     193      ! new removal not implemented for nests yet
    192194      clouds_v=cloudsn(ix,jy,hz,n,ngrid)
    193195      clouds_h=cloudsnh(ix,jy,n,ngrid)
    … …  
    231233
    232234
    233   !ZHG calculated for 1) both 2) lsp 3) convp
     235  !ZHG oct 2014 : Calculated for 1) both 2) lsp 3) convp
     236  ! Tentatively differentiate the grfraction for lsp and convp for treating differently the two forms
     237  ! for now they are treated the same
    234238    grfraction(1)=max(0.05,cc*(lsp*lfr(i)+convp*cfr(j))/(lsp+convp))
    235239    grfraction(2)=max(0.05,cc*(lfr(i)))
    … …  
    252256      wetscav=0.   
    253257
     258      !ZHG test if it nested?
    254259      if (ngrid.gt.0) then
    255260        act_temp=ttn(ix,jy,hz,n,ngrid)
    … …  
    259264     
    260265
    261   !****i*******************
     266  !***********************
    262267  ! BELOW CLOUD SCAVENGING
    263268  !*********************** 
    … …  
    267272          blc_count=blc_count+1
    268273
    269 
    270274!GAS
    271275          if (dquer(ks) .le. 0.) then  !is gas
    272   ! Gas scavenging coefficient based on Hertel et al 1995 using the below-cloud scavenging parameters A (=weta) and B (=wetb) from SPECIES file
    273             wetscav=weta(ks)*prec(1)**wetb(ks)
    274 
    275   !AEROSOL
     276            wetscav=weta(ks)*prec(1)**wetb(ks)
     277           
     278!AEROSOL
    276279          else !is particle
    277   !NIK 17.02.2015
    278   ! For the calculation here particle size needs to be in meter and not um as dquer is changed to in readreleases
    279             dquer_m=dquer(ks)/1000000. !conversion from um to m
    280            
    281   !ZHG snow or rain removal is applied based on the temperature.
     280!NIK 17.02.2015
     281! For the calculation here particle size needs to be in meter and not um as dquer is changed to in readreleases
     282! for particles larger than 10 um use the largest size defined in the parameterizations (10um)
     283            dquer_m=min(10.,dquer(ks))/1000000. !conversion from um to m
    282284            if (act_temp .ge. 273 .and. weta(ks).gt.0.)  then !Rain
    283 
    284   !Particle RAIN scavenging coefficient based on Laakso et al 2003, the below-cloud scavenging (rain efficienty) parameter A (=weta) from SPECIES file
     285              ! ZHG 2014 : Particle RAIN scavenging coefficient based on Laakso et al 2003,
     286              ! the below-cloud scavenging (rain efficienty)
     287              ! parameter A (=weta) from SPECIES file
    285288              wetscav= weta(ks)*10**(bclr(1)+ (bclr(2)*(log10(dquer_m))**(-4))+(bclr(3)*(log10(dquer_m))**(-3))+ (bclr(4)* &
    286289                   (log10(dquer_m))**(-2))+ (bclr(5)*(log10(dquer_m))**(-1))+ bclr(6)* (prec(1))**(0.5))
    287290
    288             elseif (act_temp .lt. 273 .and. wetb(ks).gt.0.)  then !snow
    289 
    290   !Particle SNOW scavenging coefficient based on Kyro et al 2009, the below-cloud scavenging (Snow efficiency) parameter B (=wetb) from SPECIES file
     291            elseif (act_temp .lt. 273 .and. wetb(ks).gt.0.)  then ! Snow
     292              ! ZHG 2014 : Particle SNOW scavenging coefficient based on Kyro et al 2009,
     293              ! the below-cloud scavenging (Snow efficiency)
     294              ! parameter B (=wetb) from SPECIES file
    291295              wetscav= wetb(ks)*10**(bcls(1)+ (bcls(2)*(log10(dquer_m))**(-4))+(bcls(3)*(log10(dquer_m))**(-3))+ (bcls(4)* &
    292296                   (log10(dquer_m))**(-2))+ (bcls(5)*(log10(dquer_m))**(-1))+ bcls(6)* (prec(1))**(0.5))
    … …  
    298302          endif !gas or particle
    299303        endif ! positive below-cloud scavenging parameters given in Species file
    300       endif !end below-cloud
     304      endif !end BELOW
    301305
    302306
    303307  !********************
    304308  ! IN CLOUD SCAVENGING
    305   !********************
    306       if (clouds_v.lt.4) then !in-cloud
    307 
     309      !******************************************************
     310      if (clouds_v.lt.4) then ! In-cloud
    308311! NIK 13 may 2015: only do incloud if positive in-cloud scavenging parameters are given in species file
    309         if (weta_in(ks).gt.0. .or. wetb_in(ks).gt.0.) then !if positive in-cloud parameters given in SPECIES file (either Ai or Bi)
    310 
     312        if (weta_in(ks).gt.0. .or. wetb_in(ks).gt.0.) then
    311313! if negative coefficients (turned off) set to zero for use in equation
    312314          if (weta_in(ks).lt.0.) weta_in(ks)=0.
    313315          if (wetb_in(ks).lt.0.) wetb_in(ks)=0.
    314316
    315           inc_count=inc_count+1
    316 
    317   !ZHG liquid water parameterization (CLWC+CIWC)
    318           if (readclouds) then !get cloud water clwc & ciwc units Kg/Kg
    319             cl=clwc(ix,jy,hz,n)+ciwc(ix,jy,hz,n)
    320           else !parameterize cloudwater
    321             cl=2E-7*prec(1)**0.36
     317          !ZHG 2015 Cloud liquid & ice water (CLWC+CIWC) from ECMWF
     318          if (readclouds) then                  !icloud_stats(ix,jy,4,n) has units kg/m2
     319            cl =icloud_stats(ix,jy,4,n)*(grfraction(1)/cc)
     320          else                                  !parameterize cloudwater m2/m3
     321            !ZHG updated parameterization of cloud water to better reproduce the values coming from ECMWF
     322            cl=1.6E-6*prec(1)**0.36
    322323          endif
    323324
     325            !ZHG: Calculate the partition between liquid and water phase water.
    324326            if (act_temp .le. 253) then
    325327              liq_frac=0
    … …  
    329331              liq_frac =((act_temp-273)/(273-253))**2
    330332            endif
    331 
    332 ! ZHG  calculate the activated fraction based on the In-cloud scavenging parameters Ai (=weta_in) and Bi (=wetb_in) from SPECIES file
    333 ! frac_act is the combined IN and CCN efficiency
    334 ! The default values are 0.9 for CCN and 0.1 IN
    335 ! This parameterization is based on Verheggen et al. (2007) & Cozich et al. (2006)
     333           ! ZHG: Calculate the aerosol partition based on cloud phase and Ai and Bi
    336334            frac_act = liq_frac*weta_in(ks) +(1-liq_frac)*wetb_in(ks)
    337335 
    … …  
    339337
    340338  ! AEROSOL
     339          !**************************************************
    341340          if (dquer(ks).gt. 0.) then ! is particle
    342341
    343342            S_i= frac_act/cl
    344343
     344          !*********************
    345345  ! GAS
    346346          else ! is gas
    347347               
    348348            cle=(1-cl)/(henry(ks)*(r_air/3500.)*act_temp)+cl
    349             S_i=frac_act/cle
    350 
     349            !REPLACE to switch old/ new scheme
     350            ! S_i=frac_act/cle
     351            S_i=1/cle
    351352          endif ! gas or particle
    352353
    353354  ! scavenging coefficient based on Hertel et al 1995 - using the S_i for either gas or aerosol
     355           !OLD
     356          if (readclouds) then
     357           wetscav=S_i*(prec(1)/3.6E6)
     358          else
    354359          wetscav=S_i*(prec(1)/3.6E6)/clouds_h
    355 
    356 !          write(*,*) 'in-cloud, act_temp=',act_temp,',prec=',prec(1),',wetscav=',wetscav,',jpart=',jpart,',clouds_h=,', &
    357 !          clouds_h,',cl=',cl, 'diff to old scheme=', cl-2E-7*prec(1)**0.36
     360          endif
     361
     362!ZHG 2015 TEST
     363!          Si_dummy=frac_act/2E-7*prec(1)**0.36
     364!           wetscav_dummy=Si_dummy*(prec(1)/3.6E6)/clouds_h
     365!           if (clouds_v.lt.4) then
     366!           talltest=talltest+1
     367!if(talltest .eq. 1) OPEN(UNIT=199, FILE=utfil,FORM='FORMATTED',STATUS = 'UNKNOWN')
     368!if(talltest .lt. 100001)  write(199,*) prec(1)/3.6E6, cl, clouds_h*2E-7*prec(1)**0.36,clouds_v,ytra1(jpart)-90
     369!if(talltest .lt. 100001)  write(199,*) wetscav, wetscav_dummy,prec(1),ytra1(jpart)-90,clouds_v,cl
     370!if(talltest .eq. 100001) CLOSE(199)
     371!if(talltest .eq. 100001) STOP
     372!
     373!write(*,*)  'PREC kg/m2s CLOUD kg/m2', (prec(1)/3.6E6), cl !, '2E-7*prec(1)**0.36',  2E-7*prec(1)**0.36,'2E-7*prec(1)**0.36*clouds_h',2E-7*prec(1)**0.36*clouds_h
     374!write(*,*)  'PREC kg/m2s LSP+convp kg/m2', prec(1), convp+lsp
     375!write(*,*)  wetscav, wetscav_dummy
     376!write(*,*) cc, grfraction(1), cc/grfraction(1)
     377
     378!write(*,*)  'Lmbda_old', (prec(1)/3.6E6)/(clouds_h*2E-7*prec(1)**0.36)
     379
     380
     381!write(*,*) '**************************************************'
     382!write(*,*)  'clouds_h', clouds_h, 'clouds_v',clouds_v,'abs(ltsample)', abs(ltsample)
     383!write(*,*)  'readclouds', readclouds, 'wetscav',wetscav, 'wetscav_dummy', wetscav_dummy
     384!write(*,*)  'S_i', S_i , 'Si_dummy', Si_dummy, 'prec(1)', prec(1)
     385
     386
     387!           write(*,*) 'PRECIPITATION ,cl  ECMWF , cl PARAMETIZED, clouds_v, lat' &
     388!                      ,prec(1)/3.6E6, cl, clouds_h*2E-7*prec(1)**0.36,clouds_v,ytra1(jpart)-90
     389
     390!endif
     391
    358392        endif ! positive in-cloud scavenging parameters given in Species file
    359393      endif !incloud
     394!END ZHG TEST
    360395     
    361396  !**************************************************
    … …  
    368403        wetdeposit(ks)=xmass1(jpart,ks)* &
    369404             (1.-exp(-wetscav*abs(ltsample)))*grfraction(1)  ! wet deposition
     405!write(*,*) 'MASS DEPOSITED: PREC, WETSCAV, WETSCAVP', prec(1), wetdeposit(ks), xmass1(jpart,ks)* &
     406!             (1.-exp(-wetscav_dummy*abs(ltsample)))*grfraction(1), clouds_v
     407
     408
    370409      else ! if no scavenging
    371410        wetdeposit(ks)=0.
Note: See TracChangeset for help on using the changeset viewer.
hosted by ZAMG