Changeset 2eefa58 in flexpart.git


Ignore:
Timestamp:
May 27, 2019, 3:28:44 PM (5 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
Children:
6741557
Parents:
f963113
Message:

Added Ronas changes for inversion output

Location:
src
Files:
3 added
16 edited

Legend:

Unmodified
Added
Removed
  • src/FLEXPART.f90

    r0a94e13 r2eefa58  
    137137    write(*,*) 'call readpaths'
    138138  endif
    139   call readpaths(pathfile)
     139  call readpaths
    140140 
    141141  if (verbosity.gt.1) then !show clock info
     
    452452  endif
    453453
     454  if (verbosity.gt.0) write (*,*) 'timemanager> call wetdepo'
    454455  call timemanager(metdata_format)
     456 
    455457
    456458  if (verbosity.gt.0) then
     
    467469      endif
    468470    end do
    469     write (*,*) 'timemanager> call wetdepo'
    470471  endif
    471472 
  • src/FLEXPART_MPI.f90

    r0c8c7f2 r2eefa58  
    146146    write(*,*) 'call readpaths'
    147147  endif
    148   call readpaths(pathfile)
     148  call readpaths
    149149 
    150150  if (verbosity.gt.1) then !show clock info
  • src/com_mod.f90

    r0a94e13 r2eefa58  
    124124  ! lnetcdfout   1 for netcdf grid output, 0 if not. Set in COMMAND (namelist input)
    125125
     126  integer :: linversionout
     127  ! linversionout 1 for one grid_time output file for each release containing all timesteps
     128
    126129  integer :: nageclass,lage(maxageclass)
    127130
     
    176179  real :: ri(5,numclass),rac(5,numclass),rcl(maxspec,5,numclass)
    177180  real :: rgs(maxspec,5,numclass),rlu(maxspec,5,numclass)
    178   real :: rm(maxspec),dryvel(maxspec)
     181  real :: rm(maxspec),dryvel(maxspec),kao(maxspec)
    179182  real :: ohcconst(maxspec),ohdconst(maxspec),ohnconst(maxspec)
    180183
     
    361364  real :: ciwc(0:nxmax-1,0:nymax-1,nzmax,numwfmem)=0.0 !ice      [kg/kg]
    362365  real :: clw(0:nxmax-1,0:nymax-1,nzmax,numwfmem)=0.0  !combined [m3/m3]
    363 
     366! RLT add pressure and dry air density
     367  real :: prs(0:nxmax-1,0:nymax-1,nzmax,numwfmem)
     368  real :: rho_dry(0:nxmax-1,0:nymax-1,nzmax,numwfmem)
    364369  real :: pv(0:nxmax-1,0:nymax-1,nzmax,numwfmem)
    365370  real :: rho(0:nxmax-1,0:nymax-1,nzmax,numwfmem)
     
    383388  ! uupol,vvpol [m/s]    wind components in polar stereographic projection
    384389  ! tt [K]               temperature data
     390  ! prs                  air pressure
    385391  ! qv                   specific humidity data
    386392  ! pv (pvu)             potential vorticity
  • src/concoutput.f90

    r20963b1 r2eefa58  
    7272  real :: sp_fact
    7373  real :: outnum,densityoutrecept(maxreceptor),xl,yl
     74! RLT
     75  real :: densitydryrecept(maxreceptor)
     76  real :: factor_dryrecept(maxreceptor)
    7477
    7578!real densityoutgrid(0:numxgrid-1,0:numygrid-1,numzgrid),
     
    106109  character(LEN=8),save :: file_stat='REPLACE'
    107110  logical :: ldates_file
     111  logical :: lexist
    108112  integer :: ierr
    109113  character(LEN=100) :: dates_char
     
    203207        densityoutgrid(ix,jy,kz)=(rho(iix,jjy,kzz,mind)*dz1+ &
    204208             rho(iix,jjy,kzz-1,mind)*dz2)/dz
     209! RLT
     210        densitydrygrid(ix,jy,kz)=(rho_dry(iix,jjy,kzz,mind)*dz1+ &
     211             rho_dry(iix,jjy,kzz-1,mind)*dz2)/dz 
    205212      end do
    206213    end do
     
    214221!densityoutrecept(i)=rho(iix,jjy,1,2)
    215222    densityoutrecept(i)=rho(iix,jjy,1,mind)
     223! RLT
     224    densitydryrecept(i)=rho_dry(iix,jjy,1,mind)
    216225  end do
    217226
     227! RLT
     228! conversion factor for output relative to dry air
     229  factor_drygrid=densityoutgrid/densitydrygrid
     230  factor_dryrecept=densityoutrecept/densitydryrecept
    218231
    219232! Output is different for forward and backward simulations
     
    353366! Concentration output
    354367!*********************
    355         if ((iout.eq.1).or.(iout.eq.3).or.(iout.eq.5).or.(iout.eq.6)) then
     368        if ((iout.eq.1).or.(iout.eq.3).or.(iout.eq.5)) then
    356369
    357370! Wet deposition
     
    614627  end do
    615628
     629! RLT Aug 2017
     630! Write out conversion factor for dry air
     631  inquire(file=path(2)(1:length(2))//'factor_drygrid',exist=lexist)
     632  if (lexist) then
     633    ! open and append
     634    open(unitoutfactor,file=path(2)(1:length(2))//'factor_drygrid',form='unformatted',&
     635            status='old',action='write',access='append')
     636  else
     637    ! create new
     638    open(unitoutfactor,file=path(2)(1:length(2))//'factor_drygrid',form='unformatted',&
     639            status='new',action='write')
     640  endif
     641  sp_count_i=0
     642  sp_count_r=0
     643  sp_fact=-1.
     644  sp_zer=.true.
     645  do kz=1,numzgrid
     646    do jy=0,numygrid-1
     647      do ix=0,numxgrid-1
     648        if (factor_drygrid(ix,jy,kz).gt.(1.+smallnum).or.factor_drygrid(ix,jy,kz).lt.(1.-smallnum)) then
     649          if (sp_zer.eqv..true.) then ! first value not equal to one
     650            sp_count_i=sp_count_i+1
     651            sparse_dump_i(sp_count_i)= &
     652                  ix+jy*numxgrid+kz*numxgrid*numygrid
     653            sp_zer=.false.
     654            sp_fact=sp_fact*(-1.)
     655          endif
     656          sp_count_r=sp_count_r+1
     657          sparse_dump_r(sp_count_r)= &
     658               sp_fact*factor_drygrid(ix,jy,kz)
     659        else ! factor is one
     660          sp_zer=.true.
     661        endif
     662      end do
     663    end do
     664  end do
     665  write(unitoutfactor) sp_count_i
     666  write(unitoutfactor) (sparse_dump_i(i),i=1,sp_count_i)
     667  write(unitoutfactor) sp_count_r
     668  write(unitoutfactor) (sparse_dump_r(i),i=1,sp_count_r)
     669  close(unitoutfactor)
     670
     671
    616672  if (gridtotal.gt.0.) gridtotalunc=gridsigmatotal/gridtotal
    617673  if (wetgridtotal.gt.0.) wetgridtotalunc=wetgridsigmatotal/ &
     
    640696  endif
    641697
    642 
     698! RLT Aug 2017
     699! Write out conversion factor for dry air
     700  if (numreceptor.gt.0) then
     701    inquire(file=path(2)(1:length(2))//'factor_dryreceptor',exist=lexist)
     702     if (lexist) then
     703     ! open and append
     704      open(unitoutfactor,file=path(2)(1:length(2))//'factor_dryreceptor',form='unformatted',&
     705              status='old',action='write',access='append')
     706    else
     707      ! create new
     708      open(unitoutfactor,file=path(2)(1:length(2))//'factor_dryreceptor',form='unformatted',&
     709              status='new',action='write')
     710    endif
     711    write(unitoutfactor) itime
     712    write(unitoutfactor) (factor_dryrecept(i),i=1,numreceptor)
     713    close(unitoutfactor)
     714  endif
    643715
    644716! Reinitialization of grid
  • src/concoutput_nest.f90

    r7d02c2f r2eefa58  
    7070  real :: sp_fact
    7171  real :: outnum,densityoutrecept(maxreceptor),xl,yl
     72! RLT
     73  real :: densitydryrecept(maxreceptor)
     74  real :: factor_dryrecept(maxreceptor)
    7275
    7376  !real densityoutgrid(0:numxgrid-1,0:numygrid-1,numzgrid),
     
    9699  character :: adate*8,atime*6
    97100  character(len=3) :: anspec
     101  logical :: lexist
    98102  integer :: mind
    99103! mind        eso:added to ensure identical results between 2&3-fields versions
     
    164168        densityoutgrid(ix,jy,kz)=(rho(iix,jjy,kzz,mind)*dz1+ &
    165169             rho(iix,jjy,kzz-1,mind)*dz2)/dz
     170! RLT
     171        densitydrygrid(ix,jy,kz)=(rho_dry(iix,jjy,kzz,mind)*dz1+ &
     172             rho_dry(iix,jjy,kzz-1,mind)*dz2)/dz
    166173      end do
    167174    end do
     
    175182    !densityoutrecept(i)=rho(iix,jjy,1,2)
    176183    densityoutrecept(i)=rho(iix,jjy,1,mind)
     184! RLT
     185    densitydryrecept(i)=rho_dry(iix,jjy,1,mind)
    177186  end do
    178187
     188! RLT
     189! conversion factor for output relative to dry air
     190  factor_drygrid=densityoutgrid/densitydrygrid
     191  factor_dryrecept=densityoutrecept/densitydryrecept
    179192
    180193  ! Output is different for forward and backward simulations
     
    551564  end do
    552565
     566! RLT Aug 2017
     567! Write out conversion factor for dry air
     568  inquire(file=path(2)(1:length(2))//'factor_drygrid_nest',exist=lexist)
     569  if (lexist) then
     570    ! open and append
     571    open(unitoutfactor,file=path(2)(1:length(2))//'factor_drygrid_nest',form='unformatted',&
     572            status='old',action='write',access='append')
     573  else
     574    ! create new
     575    open(unitoutfactor,file=path(2)(1:length(2))//'factor_drygrid_nest',form='unformatted',&
     576            status='new',action='write')
     577  endif
     578  sp_count_i=0
     579  sp_count_r=0
     580  sp_fact=-1.
     581  sp_zer=.true.
     582  do kz=1,numzgrid
     583    do jy=0,numygridn-1
     584      do ix=0,numxgridn-1
     585        if (factor_drygrid(ix,jy,kz).gt.(1.+smallnum).or.factor_drygrid(ix,jy,kz).lt.(1.-smallnum)) then
     586          if (sp_zer.eqv..true.) then ! first value not equal to one
     587            sp_count_i=sp_count_i+1
     588            sparse_dump_i(sp_count_i)= &
     589                  ix+jy*numxgridn+kz*numxgridn*numygridn
     590            sp_zer=.false.
     591            sp_fact=sp_fact*(-1.)
     592          endif
     593          sp_count_r=sp_count_r+1
     594          sparse_dump_r(sp_count_r)= &
     595               sp_fact*factor_drygrid(ix,jy,kz)
     596        else ! factor is one
     597          sp_zer=.true.
     598        endif
     599      end do
     600    end do
     601  end do
     602  write(unitoutfactor) sp_count_i
     603  write(unitoutfactor) (sparse_dump_i(i),i=1,sp_count_i)
     604  write(unitoutfactor) sp_count_r
     605  write(unitoutfactor) (sparse_dump_r(i),i=1,sp_count_r)
     606  close(unitoutfactor)
    553607
    554608
  • src/concoutput_surf.f90

    r16b61a5 r2eefa58  
    7272  real :: sp_fact
    7373  real :: outnum,densityoutrecept(maxreceptor),xl,yl
     74! RLT
     75  real :: densitydryrecept(maxreceptor)
     76  real :: factor_dryrecept(maxreceptor)
    7477
    7578!real densityoutgrid(0:numxgrid-1,0:numygrid-1,numzgrid),
     
    101104  character :: adate*8,atime*6
    102105  character(len=3) :: anspec
     106  logical :: lexist
    103107
    104108
     
    180184        densityoutgrid(ix,jy,kz)=(rho(iix,jjy,kzz,2)*dz1+ &
    181185             rho(iix,jjy,kzz-1,2)*dz2)/dz
     186! RLT
     187        densitydrygrid(ix,jy,kz)=(rho_dry(iix,jjy,kzz,2)*dz1+ &
     188             rho_dry(iix,jjy,kzz-1,2)*dz2)/dz
    182189      end do
    183190    end do
     
    190197    jjy=max(min(nint(yl),nymin1),0)
    191198    densityoutrecept(i)=rho(iix,jjy,1,2)
     199! RLT
     200    densitydryrecept(i)=rho_dry(iix,jjy,1,2)
    192201  end do
    193202
     203! RLT
     204! conversion factor for output relative to dry air
     205  factor_drygrid=densityoutgrid/densitydrygrid
     206  factor_dryrecept=densityoutrecept/densitydryrecept
    194207
    195208! Output is different for forward and backward simulations
     
    596609  end do
    597610
     611! RLT Aug 2017
     612! Write out conversion factor for dry air
     613  inquire(file=path(2)(1:length(2))//'factor_drygrid',exist=lexist)
     614  if (lexist) then
     615    ! open and append
     616    open(unitoutfactor,file=path(2)(1:length(2))//'factor_drygrid',form='unformatted',&
     617            status='old',action='write',access='append')
     618  else
     619    ! create new
     620    open(unitoutfactor,file=path(2)(1:length(2))//'factor_drygrid',form='unformatted',&
     621            status='new',action='write')
     622  endif
     623  sp_count_i=0
     624  sp_count_r=0
     625  sp_fact=-1.
     626  sp_zer=.true.
     627  do kz=1,1
     628    do jy=0,numygrid-1
     629      do ix=0,numxgrid-1
     630        if (factor_drygrid(ix,jy,kz).gt.(1.+smallnum).or.factor_drygrid(ix,jy,kz).lt.(1.-smallnum)) then
     631          if (sp_zer.eqv..true.) then ! first value not equal to one
     632            sp_count_i=sp_count_i+1
     633            sparse_dump_i(sp_count_i)= &
     634                  ix+jy*numxgrid+kz*numxgrid*numygrid
     635            sp_zer=.false.
     636            sp_fact=sp_fact*(-1.)
     637          endif
     638          sp_count_r=sp_count_r+1
     639          sparse_dump_r(sp_count_r)= &
     640               sp_fact*factor_drygrid(ix,jy,kz)
     641        else ! factor is one
     642          sp_zer=.true.
     643        endif
     644      end do
     645    end do
     646  end do
     647  write(unitoutfactor) sp_count_i
     648  write(unitoutfactor) (sparse_dump_i(i),i=1,sp_count_i)
     649  write(unitoutfactor) sp_count_r
     650  write(unitoutfactor) (sparse_dump_r(i),i=1,sp_count_r)
     651  close(unitoutfactor)
     652
     653
    598654  if (gridtotal.gt.0.) gridtotalunc=gridsigmatotal/gridtotal
    599655  if (wetgridtotal.gt.0.) wetgridtotalunc=wetgridsigmatotal/ &
     
    622678  endif
    623679
    624 
     680! RLT Aug 2017
     681! Write out conversion factor for dry air
     682  if (numreceptor.gt.0) then
     683    inquire(file=path(2)(1:length(2))//'factor_dryreceptor',exist=lexist)
     684     if (lexist) then
     685     ! open and append
     686      open(unitoutfactor,file=path(2)(1:length(2))//'factor_dryreceptor',form='unformatted',&
     687              status='old',action='write',access='append')
     688    else
     689      ! create new
     690      open(unitoutfactor,file=path(2)(1:length(2))//'factor_dryreceptor',form='unformatted',&
     691              status='new',action='write')
     692    endif
     693    write(unitoutfactor) itime
     694    write(unitoutfactor) (factor_dryrecept(i),i=1,numreceptor)
     695    close(unitoutfactor)
     696  endif
    625697
    626698! Reinitialization of grid
  • src/concoutput_surf_nest.f90

    r6a678e3 r2eefa58  
    7070  real :: sp_fact
    7171  real :: outnum,densityoutrecept(maxreceptor),xl,yl
     72! RLT
     73  real :: densitydryrecept(maxreceptor)
     74  real :: factor_dryrecept(maxreceptor)
    7275
    7376  !real densityoutgrid(0:numxgrid-1,0:numygrid-1,numzgrid),
     
    9699  character :: adate*8,atime*6
    97100  character(len=3) :: anspec
    98 
     101  logical :: lexist
    99102
    100103  ! Determine current calendar date, needed for the file name
     
    159162        densityoutgrid(ix,jy,kz)=(rho(iix,jjy,kzz,2)*dz1+ &
    160163             rho(iix,jjy,kzz-1,2)*dz2)/dz
     164! RLT
     165        densitydrygrid(ix,jy,kz)=(rho_dry(iix,jjy,kzz,2)*dz1+ &
     166             rho_dry(iix,jjy,kzz-1,2)*dz2)/dz
    161167      end do
    162168    end do
     
    169175      jjy=max(min(nint(yl),nymin1),0)
    170176      densityoutrecept(i)=rho(iix,jjy,1,2)
     177! RLT
     178    densitydryrecept(i)=rho_dry(iix,jjy,1,2)
    171179    end do
    172180
     181! RLT
     182! conversion factor for output relative to dry air
     183  factor_drygrid=densityoutgrid/densitydrygrid
     184  factor_dryrecept=densityoutrecept/densitydryrecept
    173185
    174186  ! Output is different for forward and backward simulations
     
    317329         write(unitoutgrid) sp_count_r
    318330         write(unitoutgrid) (sparse_dump_r(i),i=1,sp_count_r)
    319          write(unitoutgrid) sp_count_r
    320          write(unitoutgrid) (sparse_dump_u(i),i=1,sp_count_r)
     331!         write(unitoutgrid) sp_count_r
     332!         write(unitoutgrid) (sparse_dump_u(i),i=1,sp_count_r)
    321333
    322334  ! Dry deposition
     
    354366         write(unitoutgrid) sp_count_r
    355367         write(unitoutgrid) (sparse_dump_r(i),i=1,sp_count_r)
    356          write(unitoutgrid) sp_count_r
    357          write(unitoutgrid) (sparse_dump_u(i),i=1,sp_count_r)
     368!         write(unitoutgrid) sp_count_r
     369!         write(unitoutgrid) (sparse_dump_u(i),i=1,sp_count_r)
    358370
    359371
     
    399411         write(unitoutgrid) sp_count_r
    400412         write(unitoutgrid) (sparse_dump_r(i),i=1,sp_count_r)
    401          write(unitoutgrid) sp_count_r
    402          write(unitoutgrid) (sparse_dump_u(i),i=1,sp_count_r)
     413!         write(unitoutgrid) sp_count_r
     414!         write(unitoutgrid) (sparse_dump_u(i),i=1,sp_count_r)
    403415         else
    404416
     
    440452         write(unitoutgrid) sp_count_r
    441453         write(unitoutgrid) (sparse_dump_r(i),i=1,sp_count_r)
    442          write(unitoutgrid) sp_count_r
    443          write(unitoutgrid) (sparse_dump_u(i),i=1,sp_count_r)
     454!         write(unitoutgrid) sp_count_r
     455!         write(unitoutgrid) (sparse_dump_u(i),i=1,sp_count_r)
    444456         endif ! surf_only
    445457
     
    487499         write(unitoutgridppt) sp_count_r
    488500         write(unitoutgridppt) (sparse_dump_r(i),i=1,sp_count_r)
    489          write(unitoutgridppt) sp_count_r
    490          write(unitoutgridppt) (sparse_dump_u(i),i=1,sp_count_r)
     501!         write(unitoutgridppt) sp_count_r
     502!         write(unitoutgridppt) (sparse_dump_u(i),i=1,sp_count_r)
    491503
    492504
     
    526538         write(unitoutgridppt) sp_count_r
    527539         write(unitoutgridppt) (sparse_dump_r(i),i=1,sp_count_r)
    528          write(unitoutgridppt) sp_count_r
    529          write(unitoutgridppt) (sparse_dump_u(i),i=1,sp_count_r)
     540!         write(unitoutgridppt) sp_count_r
     541!         write(unitoutgridppt) (sparse_dump_u(i),i=1,sp_count_r)
    530542
    531543
     
    570582         write(unitoutgridppt) sp_count_r
    571583         write(unitoutgridppt) (sparse_dump_r(i),i=1,sp_count_r)
    572          write(unitoutgridppt) sp_count_r
    573          write(unitoutgridppt) (sparse_dump_u(i),i=1,sp_count_r)
     584!         write(unitoutgridppt) sp_count_r
     585!         write(unitoutgridppt) (sparse_dump_u(i),i=1,sp_count_r)
    574586         else
    575587
     
    611623         write(unitoutgridppt) sp_count_r
    612624         write(unitoutgridppt) (sparse_dump_r(i),i=1,sp_count_r)
    613          write(unitoutgridppt) sp_count_r
    614          write(unitoutgridppt) (sparse_dump_u(i),i=1,sp_count_r)
     625!         write(unitoutgridppt) sp_count_r
     626!         write(unitoutgridppt) (sparse_dump_u(i),i=1,sp_count_r)
    615627         endif ! surf_only
    616628
     
    624636
    625637  end do
     638
     639! RLT Aug 2017
     640! Write out conversion factor for dry air
     641  inquire(file=path(2)(1:length(2))//'factor_drygrid_nest',exist=lexist)
     642  if (lexist) then
     643    ! open and append
     644    open(unitoutfactor,file=path(2)(1:length(2))//'factor_drygrid_nest',form='unformatted',&
     645            status='old',action='write',access='append')
     646  else
     647    ! create new
     648    open(unitoutfactor,file=path(2)(1:length(2))//'factor_drygrid_nest',form='unformatted',&
     649            status='new',action='write')
     650  endif
     651  sp_count_i=0
     652  sp_count_r=0
     653  sp_fact=-1.
     654  sp_zer=.true.
     655  do kz=1,1
     656    do jy=0,numygridn-1
     657      do ix=0,numxgridn-1
     658        if (factor_drygrid(ix,jy,kz).gt.(1.+smallnum).or.factor_drygrid(ix,jy,kz).lt.(1.-smallnum)) then
     659          if (sp_zer.eqv..true.) then ! first value not equal to one
     660            sp_count_i=sp_count_i+1
     661            sparse_dump_i(sp_count_i)= &
     662                  ix+jy*numxgridn+kz*numxgridn*numygridn
     663            sp_zer=.false.
     664            sp_fact=sp_fact*(-1.)
     665          endif
     666          sp_count_r=sp_count_r+1
     667          sparse_dump_r(sp_count_r)= &
     668               sp_fact*factor_drygrid(ix,jy,kz)
     669        else ! factor is one
     670          sp_zer=.true.
     671        endif
     672      end do
     673    end do
     674  end do
     675  write(unitoutfactor) sp_count_i
     676  write(unitoutfactor) (sparse_dump_i(i),i=1,sp_count_i)
     677  write(unitoutfactor) sp_count_r
     678  write(unitoutfactor) (sparse_dump_r(i),i=1,sp_count_r)
     679  close(unitoutfactor)
    626680
    627681
  • src/getfields.f90

    r6ecb30a r2eefa58  
    8787  real :: pvhn(0:nxmaxn-1,0:nymaxn-1,nuvzmax,maxnests)
    8888  real :: wwhn(0:nxmaxn-1,0:nymaxn-1,nwzmax,maxnests)
     89! RLT added partial pressure water vapor
     90  real :: pwater(0:nxmax-1,0:nymax-1,nzmax,numwfmem)
     91  integer :: kz, ix
     92  character(len=100) :: rowfmt
    8993
    9094  integer :: indmin = 1
     
    15315740  indmin=indj
    154158
    155    if (WETBKDEP) then
    156         call writeprecip(itime,memind(1))
    157    endif
     159    if (WETBKDEP) then
     160      call writeprecip(itime,memind(1))
     161    endif
    158162
    159163  else
     
    20420860  indmin=indj
    205209
    206    if (WETBKDEP) then
    207         call writeprecip(itime,memind(1))
    208    endif
    209 
    210   endif
     210    if (WETBKDEP) then
     211      call writeprecip(itime,memind(1))
     212    endif
     213
     214  end if
     215
     216  ! RLT calculate dry air density
     217  pwater=qv*prs/((r_air/r_water)*(1.-qv)+qv)
     218  rho_dry=(prs-pwater)/(r_air*tt)
     219
     220  ! test density
     221!  write(rowfmt,'(A,I6,A)') '(',nymax,'(E11.4,1X))'
     222!  if(itime.eq.0) then
     223!    open(500,file=path(2)(1:length(2))//'rho_dry.txt',status='replace',action='write')
     224!    do kz=1,nzmax
     225!      do ix=1,nxmax
     226!        write(500,fmt=rowfmt) rho_dry(ix,:,kz,1)
     227!      end do
     228!    end do
     229!    close(500)
     230!    open(500,file=path(2)(1:length(2))//'rho.txt',status='replace',action='write')
     231!    do kz=1,nzmax
     232!      do ix=1,nxmax
     233!        write(500,fmt=rowfmt) rho(ix,:,kz,1)
     234!      end do
     235!    end do
     236!    close(500)
     237!  endif
    211238
    212239  lwindinterv=abs(memtime(2)-memtime(1))
  • src/makefile

    r0a94e13 r2eefa58  
    128128        redist.o                \
    129129        concoutput_surf.o       concoutput_surf_nest.o  \
     130        concoutput_inversion_nest.o     \
     131        concoutput_inversion.o \
    130132        getfields.o \
    131133        readwind_ecmwf.o
     
    200202drydepokernel_nest.o    zenithangle.o \
    201203ohreaction.o            getvdep_nests.o \
    202 initial_cond_calc.o     initial_cond_output.o \
     204initial_cond_calc.o     initial_cond_output.o initial_cond_output_inversion.o \
    203205dynamic_viscosity.o     get_settling.o  \
    204206initialize_cbl_vel.o    re_initialize_particle.o \
     
    272274conccalc_mpi.o: com_mod.o mpi_mod.o outg_mod.o par_mod.o unc_mod.o
    273275concoutput.o: com_mod.o outg_mod.o par_mod.o point_mod.o unc_mod.o mean_mod.o
     276concoutput_inversion.o: com_mod.o outg_mod.o par_mod.o point_mod.o unc_mod.o mean_mod.o
    274277concoutput_mpi.o: com_mod.o mpi_mod.o outg_mod.o par_mod.o point_mod.o \
    275278        unc_mod.o mean_mod.o
    276279concoutput_nest.o: com_mod.o outg_mod.o par_mod.o point_mod.o unc_mod.o mean_mod.o
     280concoutput_inversion_nest.o: com_mod.o outg_mod.o par_mod.o point_mod.o unc_mod.o mean_mod.o
    277281concoutput_nest_mpi.o: com_mod.o mpi_mod.o outg_mod.o par_mod.o point_mod.o \
    278282        unc_mod.o mean_mod.o
     
    321325initial_cond_calc.o: com_mod.o outg_mod.o par_mod.o unc_mod.o
    322326initial_cond_output.o: com_mod.o outg_mod.o par_mod.o point_mod.o unc_mod.o
     327initial_cond_output_inversion.o: com_mod.o outg_mod.o par_mod.o point_mod.o unc_mod.o
    323328initialize.o: com_mod.o hanna_mod.o interpol_mod.o par_mod.o random_mod.o
    324329initialize_cbl_vel.o: com_mod.o par_mod.o random_mod.o
  • src/outg_mod.f90

    r4c64400 r2eefa58  
    3737  real,allocatable, dimension (:,:,:) :: areanorth
    3838  real,allocatable, dimension (:,:,:) :: densityoutgrid
     39  real,allocatable, dimension (:,:,:) :: densitydrygrid ! added RLT
     40  real,allocatable, dimension (:,:,:) :: factor_drygrid ! added RLT
    3941  real,allocatable, dimension (:,:,:) :: factor3d
    4042  real,allocatable, dimension (:,:,:) :: grid
  • src/outgrid_init.f90

    rd2a5a83 r2eefa58  
    268268    if (stat.ne.0) write(*,*)'ERROR: could not allocate gridunc'
    269269  allocate(densityoutgrid(0:max(numxgrid,numxgridn)-1, &
     270       0:max(numygrid,numygridn)-1,numzgrid),stat=stat)
     271    if (stat.ne.0) write(*,*)'ERROR: could not allocate gridunc'
     272! RLT
     273  allocate(densitydrygrid(0:max(numxgrid,numxgridn)-1, &
     274       0:max(numygrid,numygridn)-1,numzgrid),stat=stat)
     275    if (stat.ne.0) write(*,*)'ERROR: could not allocate gridunc'
     276  allocate(factor_drygrid(0:max(numxgrid,numxgridn)-1, &
    270277       0:max(numygrid,numygridn)-1,numzgrid),stat=stat)
    271278    if (stat.ne.0) write(*,*)'ERROR: could not allocate gridunc'
  • src/par_mod.f90

    r0a94e13 r2eefa58  
    3030!        Update 15 August 2013 IP                                              *
    3131!                                                                              *
    32 !        ESO 2016:                                                             *
    33 !          GFS specific parameters moved to gfs_mod.f90                        *
    34 !          ECMWF specific parameters moved to ecmwf_mod.f90                    *
    3532!                                                                              *
    3633!*******************************************************************************
     
    8077  real,parameter :: pi=3.14159265, r_earth=6.371e6, r_air=287.05, ga=9.81
    8178  real,parameter :: cpa=1004.6, kappa=0.286, pi180=pi/180., vonkarman=0.4
     79  ! additional constants RLT Aug-2017
     80  real,parameter :: rgas=8.31447
     81  real,parameter :: r_water=461.495
    8282
    8383  ! pi                      number "pi"
     
    8989  ! kappa                   exponent of formula for potential temperature
    9090  ! vonkarman               von Karman constant
     91  ! rgas                    universal gas constant [J/mol/K]
     92  ! r_water                 specific gas constant for water vapor [J/kg/K]
    9193
    9294  real,parameter :: karman=0.40, href=15., convke=2.0
     
    290292  integer,parameter :: unitboundcond=89
    291293  integer,parameter :: unittmp=101
     294! RLT
     295  integer,parameter :: unitoutfactor=102
    292296
    293297!******************************************************
  • src/readcommand.f90

    r0a94e13 r2eefa58  
    114114  surf_only, &
    115115  cblflag, &
     116  linversionout, &
    116117  ohfields_path
    117118
     
    147148  surf_only=0
    148149  cblflag=0 ! if using old-style COMMAND file, set to 1 here to use mc cbl routine
     150  linversionout=0
    149151  ohfields_path="../../flexin/"
    150152
     
    414416  !**********************************************************************
    415417
    416   if ((iout.lt.1).or.(iout.gt.6)) then
     418  if ((iout.lt.1).or.(iout.gt.5)) then
    417419    write(*,*) ' #### FLEXPART MODEL ERROR! FILE COMMAND:     #### '
    418420    write(*,*) ' #### IOUT MUST BE 1, 2, 3, 4 OR 5 FOR        #### '
     
    474476      write(*,*) '#### FOR DOMAIN FILLING RUNS OUTPUT FOR      ####'
    475477      write(*,*) '#### EACH RELEASE IS FORBIDDEN !             ####'
     478      stop
     479  endif
     480
     481  ! Inversion output format only for backward runs
     482  !*****************************************************************************
     483 
     484  if ((linversionout.eq.1).and.(ldirect.eq.1)) then
     485      write(*,*) '#### FLEXPART MODEL ERROR! FILE COMMAND:     ####'
     486      write(*,*) '#### INVERSION OUTPUT FORMAT ONLY FOR        ####'
     487      write(*,*) '#### BACKWARD RUNS                           ####'
    476488      stop
    477489  endif
  • src/readpaths.f90

    r62e65c7 r2eefa58  
    2020!**********************************************************************
    2121
    22 subroutine readpaths !(pathfile)
     22subroutine readpaths
    2323
    2424  !*****************************************************************************
  • src/timemanager.f90

    r0a94e13 r2eefa58  
    408408#endif
    409409            else
     410              if (linversionout.eq.1) then
     411                call concoutput_inversion(itime,outnum,gridtotalunc,wetgridtotalunc,drygridtotalunc)
     412                if (verbosity.eq.1) then
     413                  print*,'called concoutput_inversion'
     414                  call system_clock(count_clock)
     415                  write(*,*) 'system clock',count_clock - count_clock0
     416                endif
     417              else
    410418              call concoutput_surf(itime,outnum,gridtotalunc,wetgridtotalunc,drygridtotalunc)
     419              endif
    411420              if (verbosity.eq.1) then
    412421                print*,'called concoutput_surf '
     
    422431                call concoutput_nest(itime,outnum)
    423432              else
     433                if(linversionout.eq.1) then
     434                  call concoutput_inversion_nest(itime,outnum)
     435                else
    424436                call concoutput_surf_nest(itime,outnum)
     437              endif
    425438              endif
    426439            else
     
    739752  if (ipout.eq.2) call partoutput(itime)     ! dump particle positions
    740753
    741   if (linit_cond.ge.1) call initial_cond_output(itime)   ! dump initial cond. field
     754  if (linit_cond.ge.1) then
     755    if(linversionout.eq.1) then
     756      call initial_cond_output_inversion(itime)   ! dump initial cond. field
     757    else
     758      call initial_cond_output(itime)   ! dump initial cond. fielf
     759    endif
     760  endif
    742761
    743762  !close(104)
  • src/verttransform_ecmwf.f90

    rd005a67 r2eefa58  
    8181  real,dimension(0:nxmax-1,0:nymax-1,nuvzmax) :: rhoh,uvzlev,wzlev
    8282  real,dimension(0:nxmax-1,0:nymax-1,nzmax) :: pinmconv
     83  ! RLT added pressure
     84  real,dimension(0:nxmax-1,0:nymax-1,nuvzmax) :: prsh
    8385  real,dimension(0:nxmax-1,0:nymax-1) ::  tvold,pold,pint,tv
    8486  real,dimension(0:nymax-1) :: cosf
     
    218220!*************************
    219221
    220 
    221222  do jy=0,nymin1
    222223    do ix=0,nxmin1
     
    229230  wzlev(:,:,1)=0.
    230231  rhoh(:,:,1)=pold/(r_air*tvold)
     232  ! RLT add pressure
     233  prsh(:,:,1)=ps(:,:,1,n)
    231234
    232235
     
    236239  do kz=2,nuvz
    237240    pint=akz(kz)+bkz(kz)*ps(:,:,1,n)
     241    ! RLT add pressure
     242    prsh(:,:,kz)=pint
    238243    tv=tth(:,:,kz,n)*(1.+0.608*qvh(:,:,kz,n))
    239244    rhoh(:,:,kz)=pint(:,:)/(r_air*tv)
     
    287292  pv(:,:,1,n)=pvh(:,:,1)
    288293  rho(:,:,1,n)=rhoh(:,:,1)
     294! RLT add pressure
     295  prs(:,:,1,n)=prsh(:,:,1)
    289296
    290297  uu(:,:,nz,n)=uuh(:,:,nuvz)
     
    300307  pv(:,:,nz,n)=pvh(:,:,nuvz)
    301308  rho(:,:,nz,n)=rhoh(:,:,nuvz)
    302 
     309! RLT
     310  prs(:,:,nz,n)=prsh(:,:,nuvz)
    303311
    304312  kmin=2
     
    320328          pv(ix,jy,iz,n)=pv(ix,jy,nz,n)
    321329          rho(ix,jy,iz,n)=rho(ix,jy,nz,n)
     330! RLT
     331          prs(ix,jy,iz,n)=prs(ix,jy,nz,n)
    322332        else
    323333          innuvz: do kz=idx(ix,jy),nuvz
     
    353363          pv(ix,jy,iz,n)=(pvh(ix,jy,kz-1)*dz2+pvh(ix,jy,kz)*dz1)/dz
    354364          rho(ix,jy,iz,n)=(rhoh(ix,jy,kz-1)*dz2+rhoh(ix,jy,kz)*dz1)/dz
     365! RLT add pressure
     366          prs(ix,jy,iz,n)=(prsh(ix,jy,kz-1)*dz2+prsh(ix,jy,kz)*dz1)/dz
    355367        endif
    356368      enddo
Note: See TracChangeset for help on using the changeset viewer.
hosted by ZAMG