! source file: /nfs/tape_cache/smomw258/latin_hypercube_agg/template/updates/UVic_ESCM.F
program UVic_ESCM
!=======================================================================
! UNIVERSITY OF VICTORIA EARTH SYSTEM CLIMATE MODEL (UVic ESCM)
! The UVic ESCM is a climate model developed by researchers in the
! Climate Modelling Group, within the School of Earth and Ocean
! Sciences, located at the University of Victoria, Victoria,
! British Columbia, Canada.
! Many people have contributed to the development of this code.
! This is a collective effort and although individual contributions
! are appreciated, individual authorship is not indicated.
! Some of the people who have contributed to the development of the
! model are:
! D. Archer, C. Avis, A. Berger, R. Betts, A. Biastoch, C. Bitz,
! J. Brauch, C. Brennan, B. Bryan, J. Burton, S. Carto,
! M. Cottet-Puinel, M. Cox, P. Cox, G. Danabasoglu, K. Dixon,
! J. Dukowicz, M. Eby, R. Essery, T. Ewen, A. Fanning, J. Fyke,
! R. Gerdes, A. Gnanadesikan, C. Goldberg, D. Goldberg,
! D. Gregory, J. Gregory, S. Griffies, R. Hanson, R. Hetherington,
! H. Hickey, M. Holland, G. Holloway, T. Huck, T. Hughes,
! E. Hunke, W. Hurlin, W. Ingram, R. Key, E. Kluzek, C. Koeberle,
! K. Kvale, C. Lawson, J. Lewis, M. Loutre, R. Malone,
! D. Matthews, K. Meissner, A. Montenegro, A. Mouchet, T. Murdock,
! P. Myers, A. Oschlies, R. Pacanowski, M. Pahlow, P. Poussart
! S. Rahmstorf, R. Redler,, D. Robitaille, A. Rosati, M. Roth
! C. Sabine, O. Saenko, A. Schmittner, B. Semtner, W. Sijp
! T. Silva, H. Simmons, A. Skvortsov, R. Smith, P. Spence
! D. Stone, R. Tonkonojenkov, C. Tricot, S. Valcke, A. Weaver
! E. Wiebe, J. Willebrandt, M. Yoshimori, K. Zickfeld
! Please direct problems or questions to the code contact
! person at: http://climate.uvic.ca/model
! Requirements:
! Standard fortran 90 is used
! Disclaimer:
! The UVic Earth System Climate Model is a climate modelling
! research tool developed at the University of Victoria. Others may
! use it freely but we assume no responsibility for problems or
! incorrect use. It is left to the user to ensure that a particular
! configuration is working correctly.
!=======================================================================
implicit none
integer i, j, numots, numats, numseg, n, loop
include "size.h"
include "param.h"
include "pconst.h"
include "stdunits.h"
include "coord.h"
include "csbc.h"
include "iounit.h"
include "emode.h"
include "levind.h"
include "scalar.h"
include "switch.h"
include "tmngr.h"
include "cembm.h"
include "atm.h"
include "mw.h"
print*, '== UNIVERSITY OF VICTORIA EARTH SYSTEM CLIMATE MODEL =='
print*, ' '
!-----------------------------------------------------------------------
! initialize i/o units
!-----------------------------------------------------------------------
call ioinit
!-----------------------------------------------------------------------
! setup file renaming
!-----------------------------------------------------------------------
call file_names
!-----------------------------------------------------------------------
! Initialize S.B.C. indices
!-----------------------------------------------------------------------
call sbc_init
!-----------------------------------------------------------------------
! Initialize tracers
!-----------------------------------------------------------------------
call tracer_init
!-----------------------------------------------------------------------
! read namelist variables
!-----------------------------------------------------------------------
call read_namelist
!-----------------------------------------------------------------------
! read grid
!-----------------------------------------------------------------------
call grids
!-----------------------------------------------------------------------
! read topography
!-----------------------------------------------------------------------
call topog (kmt, kmu, map, xt, yt, zt, xu, yu, zw, imt, jmt, km
&, sg_bathy
& )
call isleperim (kmt, map, iperm, jperm, iofs, nippts, nisle, imt
&, jmt, km, mnisle, maxipp, xu, yu, zw)
!-----------------------------------------------------------------------
! common setup
!-----------------------------------------------------------------------
call setcom (1, imt, 1, jmt)
!-----------------------------------------------------------------------
! ocean setup
!-----------------------------------------------------------------------
call setmom (1, imt, 1, jmt)
!-----------------------------------------------------------------------
! atmosphere setup
!-----------------------------------------------------------------------
call setembm (1, imt, 1, jmt)
!-----------------------------------------------------------------------
! land setup (requires embm)
!-----------------------------------------------------------------------
call setmtlm (1, imt, 1, jmt)
!-----------------------------------------------------------------------
! initialize data reading routine
!-----------------------------------------------------------------------
call setdata (1, imt, 1, jmt)
!-----------------------------------------------------------------------
! compute the number of ocean time steps "numots" for this run and
! the number of ocean time steps per ocean segment "ntspos".
! compute the number of atmos time steps "numats" for this run and
! the number of atmos time steps per atmos segment "ntspas".
! divide the integration time "days" into "numseg" segments.
! each will be length "segtim" days. Surface boundary conditions
! are supplied every "segtim" days.
!-----------------------------------------------------------------------
numots = nint(rundays/(dtocn*secday))
ntspos = nint(segtim/(dtocn*secday))
numats = nint(rundays/(dtatm*secday))
ntspas = nint(segtim/(dtatm*secday))
numseg = numots/ntspos
if (segtim .gt. 1.) then
ntspls = nint(c1/(dtlnd*secday))
else
ntspls = nint(segtim/(dtlnd*secday))
endif
!-----------------------------------------------------------------------
! check for consistency in the S.B.C. setup
!-----------------------------------------------------------------------
call chkcpl
!-----------------------------------------------------------------------
! get global sums at the start of the run
!-----------------------------------------------------------------------
call globalsum (1)
!-----------------------------------------------------------------------
! S T A R T S E G M E N T L O O P
!-----------------------------------------------------------------------
tsicpl = .false.
tavgcpl = .false.
do n=1,numseg
cplts = .false.
!-----------------------------------------------------------------------
! get the atmospheric S.B.C.
!-----------------------------------------------------------------------
call gasbc (1, imt, 1, jmt)
!-----------------------------------------------------------------------
! call the atmospheric model once for each time step until one
! segment of "segtim" days is complete. hold atmos S.B.C. fixed
! during each segment and predict average S.B.C. for ocean
!-----------------------------------------------------------------------
do loop=1,ntspas
call embm (1, imt, 1, jmt)
enddo
!-----------------------------------------------------------------------
! get land S.B.C.s
!-----------------------------------------------------------------------
call glsbc (1, imt, 1, jmt)
!----------------------------------------------------------------------
! call the land-surface and vegetation model once for each time
! step until one segment of "segtim" days is complete.
!-----------------------------------------------------------------------
do loop=1,ntspls
call mtlm (1, imt, 1, jmt)
enddo
!-----------------------------------------------------------------------
! get ocean S.B.C.s
!-----------------------------------------------------------------------
call gosbc (1, imt, 1, jmt)
!-----------------------------------------------------------------------
! call the ocean model once for each time step until one
! segment of "segtim" days is complete. hold ocean S.B.C. fixed
! during each segment and predict average S.B.C. for atmos
!-----------------------------------------------------------------------
do loop=1,ntspos
if (loop .eq. ntspos) cplts = .true.
call mom
call embmout (1, imt, 1, jmt)
call mtlmout (1, imt, 1, jmt)
if (tsits .and. iotsi .ne. stdout .and. iotsi .gt. 0) then
write (*,'(1x, a3, i7, 1x, a32)') 'ts=',itt, stamp
endif
enddo
!-----------------------------------------------------------------------
! write change in global sums for heat and fresh water
!-----------------------------------------------------------------------
if (tsits) call globalsum (2)
!-----------------------------------------------------------------------
! close any open netcdf files to flush buffers
!-----------------------------------------------------------------------
call closeall
enddo
!-----------------------------------------------------------------------
! E N D S E G M E N T L O O P
!-----------------------------------------------------------------------
!-----------------------------------------------------------------------
! get global sums at the end of the run
!-----------------------------------------------------------------------
call globalsum (3)
print*, ' ==> UVIC_ESCM integration is complete.'
call closeall
call release_all
stop
end
subroutine chkcpl
implicit none
integer jrow, j, i, k
logical errorc
real critv, tmp, t1, r1, r2, r3, r4, r5, r6
include "size.h"
include "param.h"
include "pconst.h"
include "stdunits.h"
include "csbc.h"
include "cembm.h"
include "calendar.h"
include "switch.h"
include "tmngr.h"
!-----------------------------------------------------------------------
! adjust output intervals for forcing acceleration
!-----------------------------------------------------------------------
if (accel .le. 1) then
accel = 1.
else
print*, ' '
print*, '==> Warning: any forcing that changes on time scales'
& , ' longer than 1 year has been accelerated'
print*, ' acceleration factor: ', accel
print*, ' starts at relyr: ', accel_yr0
print*, ' current year: ', year0 + relyr
tmp = year0 + accel_yr0 + (relyr - accel_yr0)*accel
print*, ' current accelerated year: ', tmp
print*, ' '
print*, '==> Warning: adjusting output intervals for'
&, ' long term forcing acceleration'
print*, ' '
if (tsiint .ge. yrlen*accel) tsiint = tsiint/accel
if (tavgint .ge. yrlen*accel) tavgint = tavgint/accel
if (tmbint .ge. yrlen*accel) tmbint = tmbint/accel
if (stabint .ge. yrlen*accel) stabint = stabint/accel
if (zmbcint .ge. yrlen*accel) zmbcint = zmbcint/accel
if (glenint .ge. yrlen*accel) glenint = glenint/accel
if (trmbint .ge. yrlen*accel) trmbint = trmbint/accel
if (vmsfint .ge. yrlen*accel) vmsfint = vmsfint/accel
if (gyreint .ge. yrlen*accel) gyreint = gyreint/accel
if (extint .ge. yrlen*accel) extint = extint/accel
if (prxzint .ge. yrlen*accel) prxzint = prxzint/accel
if (exconvint .ge. yrlen*accel) exconvint = exconvint/accel
if (dspint .ge. yrlen*accel) dspint = dspint/accel
if (timavgint .ge. yrlen*accel) timavgint = timavgint/accel
if (cmixint .ge. yrlen*accel) cmixint = cmixint/accel
if (xbtint .ge. yrlen*accel) xbtint = xbtint/accel
if (crossint .ge. yrlen*accel) crossint = crossint/accel
if (densityint .ge. yrlen*accel) densityint = densityint/accel
if (fctint .ge. yrlen*accel) fctint = fctint/accel
if (tyzint .ge. yrlen*accel) tyzint = tyzint/accel
if (restint .ge. yrlen*accel) restint = restint/accel
if (tbtint .ge. yrlen*accel) tbtint = tbtint/accel
endif
if (timavgint .lt. timavgper) then
write (stdout,'(/,(1x,a))')
& '==> Warning: averaging period "timavgper" exceeds interval '
&, ' "timavgint". Setting timavgper = timavgint '
timavgper = timavgint
endif
write (stdout,'(/,1x,a,f10.2,a,/,1x,a,f10.2,a)')
& '==> Time averages written every ', timavgint, ' days, '
&, ' with an averaging period of ', timavgper, ' days. '
if (tsiint .lt. tsiper) then
write (stdout,'(/,(1x,a))')
& '==> Warning: averaging period "tsiper" exceeds interval '
&, ' "tsiint". Setting tsiper = tsiint '
tsiper = tsiint
endif
write (stdout,'(/,1x,a,f10.2,a,/,1x,a,f10.2,a)')
& '==> Time step integrals written every ', tsiint, ' days, '
&, ' with an averaging period of ', tsiper, ' days. '
if (mod(namix,2) .ne. mod(nint(segtim*daylen/dtatm),2)) then
write(*,*) '==> Error: time steps between mixing and coupling'
write(*,*) ' in the atmosphere must both be even to'
write(*,*) ' use even_fluxes.'
stop '=>UVic_ESCM'
endif
if (mod(nmix,2) .ne. mod(nint(segtim*daylen/dtocn),2)) then
write(*,*) '==> Error: time steps between mixing and coupling'
write(*,*) ' in the ocean must both be even to use'
write(*,*) ' even_fluxes.'
stop '=>UVic_ESCM'
endif
if (mod(namix,2) .ne. mod(nats,2)) then
write(*,*) '==> Warning: restart was not saved with even flux'
write(*,*) ' averaging. Starting with a mixing time'
write(*,*) ' step in atmosphere.'
nats = namix
endif
if (mod(nmix,2) .ne. mod(nots,2)) then
write(*,*) '==> Warning: restart was not saved with even flux'
write(*,*) ' averaging. Starting with a mixing time'
write(*,*) ' step in ocean.'
nots = nmix
endif
!-----------------------------------------------------------------------
! do consistency checks before allowing model to continue
!-----------------------------------------------------------------------
errorc = .false.
write (stdout,*) ' '
write (stdout,*) ' (checking setup)'
if (dtatm .eq. c0) then
write (stdout,9000)
& '==> Error: the atmospheric time step must be set in "setatm"'
errorc = .true.
dtatm = 1.e-6
endif
! critv = 1.e-6
critv = 1.e-4
if (segtim .ne. c0) then
r1 = rundays/segtim
else
r1 = 0.5
endif
if (segtim .eq. c0) then
write (stdout,9000)
& '==> Error: coupling period "segtim" must be specified'
errorc = .true.
endif
if (abs(r1-nint(r1)) .gt. critv) then
write (stdout,9000)
& '==> Error: there must be an integral number of segments '
&,'"segtim" within "rundays" (the length of the run)'
errorc = .true.
endif
r2 = segtim/(dtocn*secday)
if (abs(r2-nint(r2)) .gt. critv) then
write (stdout,9000)
& '==> Error: there must be an integral number of density time '
&,'steps "dtocn" within "segtim" (the segment time)'
errorc = .true.
endif
r3 = segtim/(dtatm*secday)
if (abs(r3-nint(r3)) .gt. critv) then
write (stdout,9000)
& '==> Error: there must be an integral number of atmos time '
&,'steps "dtatm" within "segtim" (the segment time)'
errorc = .true.
endif
r5 = segtim/(dtlnd*secday)
if (abs(r5-nint(r5)) .gt. critv) then
write (stdout,9000)
& '==> Error: there must be an integral number of land time '
&,'steps "dtlnd" within "segtim" (the segment time)'
errorc = .true.
endif
tmp = 0
& + 1
if (tmp .gt. 1) then
write (stdout,9000)
& '==> Error: use only one of O_co2ccn_data, O_co2emit_data, '
&,'or O_co2ccn_user'
errorc = .true.
endif
write (stdout,*) ' '
write (stdout,*) ' (End of checks) '
write (stdout,*) ' '
if (errorc) stop '=> ERRORS found in chkcpl'
9000 format (/,(1x,a))
return
end
subroutine set (index, num, name, text, inc)
!-----------------------------------------------------------------------
! increment counter, set index and text
!-----------------------------------------------------------------------
character(*) :: name, text
name = text
index = num
inc = index + 1
return
end
subroutine getst (jrow, ocnout, ntabc)
!-----------------------------------------------------------------------
! read surface tracers from disk row "jrow"
!-----------------------------------------------------------------------
implicit none
integer i, jrow, ntabc
include "size.h"
include "param.h"
include "pconst.h"
include "stdunits.h"
include "iounit.h"
include "mw.h"
include "tmngr.h"
real ocnout(imt,jmt)
call getrow (latdisk(taup1disk), nslab, jrow
&, u(1,1,jmw,1,taup1), t(1,1,jmw,1,taup1))
do i=1,imt
ocnout(i,jrow) = t(i,1,jmw,ntabc,taup1)
enddo
return
end
subroutine sbc_init
!-----------------------------------------------------------------------
! Initialize S.B.C. indices
!-----------------------------------------------------------------------
implicit none
integer m
include "size.h"
include "param.h"
include "pconst.h"
include "stdunits.h"
include "csbc.h"
sbc(:,:,:) = 0.0
mapsbc(:) = " "
itaux = 0
itauy = 0
iws = 0
iaca = 0
isca = 0
ihflx = 0
isflx = 0
isst = 0
isss = 0
iro = 0
iwa = 0
iwxq = 0
iwyq = 0
iwxt = 0
iwyt = 0
iwxc = 0
iwyc = 0
ipsw = 0
isu = 0
isv = 0
igu = 0
igv = 0
issdic = 0
idicflx = 0
issalk = 0
ialkflx = 0
isso2 = 0
io2flx = 0
isspo4 = 0
ipo4flx = 0
issphyt = 0
iphytflx = 0
isszoop = 0
izoopflx = 0
issdetr = 0
idetrflx = 0
issmp = 0
impflx = 0
issmpa = 0
impaflx = 0
issmpp = 0
imppflx = 0
issdetrz = 0
idetrzflx = 0
issno3 = 0
ino3flx = 0
issdiaz = 0
idiazflx = 0
isscocc = 0
icoccflx = 0
icaco3flx = 0
isscaco3 = 0
issdetr_B = 0
idetrflx_B = 0
issc14 = 0
ic14flx = 0
isscfc11 = 0
icfc11flx = 0
isscfc12 = 0
icfc12flx = 0
iat = 0
irh = 0
ipr = 0
ips = 0
iaws = 0
iswr = 0
ilwr = 0
isens = 0
ievap = 0
idtr = 0
isr = 0
inpp = 0
iburn = 0
ibtemp = 0
ibsalt = 0
ibdic = 0
ibdicfx = 0
ibalk = 0
ibalkfx = 0
ibo2 = 0
ircal = 0
irorg = 0
m = 1
call set (itaux, m, mapsbc(m), 'taux', m)
call set (itauy, m, mapsbc(m), 'tauy', m)
call set (iws, m, mapsbc(m), 'ws', m)
call set (iaca, m, mapsbc(m), 'a_calb', m)
call set (isca, m, mapsbc(m), 's_calb', m)
call set (ihflx, m, mapsbc(m), 'hflx', m)
call set (isflx, m, mapsbc(m), 'sflx', m)
call set (isst, m, mapsbc(m), 'sst', m)
call set (isss, m, mapsbc(m), 'sss', m)
call set (iro, m, mapsbc(m), 'ro', m)
call set (iwxq, m, mapsbc(m), 'wx_q', m)
call set (iwyq, m, mapsbc(m), 'wy_q', m)
call set (iwxt, m, mapsbc(m), 'wx_t', m)
call set (iwyt, m, mapsbc(m), 'wy_t', m)
call set (isu, m, mapsbc(m), 'su', m)
call set (isv, m, mapsbc(m), 'sv', m)
call set (igu, m, mapsbc(m), 'gu', m)
call set (igv, m, mapsbc(m), 'gv', m)
call set (issdic, m, mapsbc(m), 'ssdic', m)
call set (idicflx, m, mapsbc(m), 'dicflx', m)
call set (issc14, m, mapsbc(m), 'ssc14', m)
call set (ic14flx, m, mapsbc(m), 'c14flx', m)
call set (issmp, m, mapsbc(m), 'ssmp', m)
call set (impflx, m, mapsbc(m), 'mpflx', m)
call set (issmpa, m, mapsbc(m), 'ssmpa', m)
call set (impaflx, m, mapsbc(m), 'mpaflx', m)
call set (issmpp, m, mapsbc(m), 'ssmpp', m)
call set (imppflx, m, mapsbc(m), 'mppflx', m)
call set (issalk, m, mapsbc(m), 'ssalk', m)
call set (ialkflx, m, mapsbc(m), 'alkflx', m)
call set (isso2, m, mapsbc(m), 'sso2', m)
call set (io2flx, m, mapsbc(m), 'o2flx', m)
call set (isspo4, m, mapsbc(m), 'sspo4', m)
call set (ipo4flx, m, mapsbc(m), 'po4flx', m)
call set (issno3, m, mapsbc(m), 'ssno3', m)
call set (ino3flx, m, mapsbc(m), 'no3flx', m)
call set (iat, m, mapsbc(m), 'at', m)
call set (irh, m, mapsbc(m), 'rh', m)
call set (ipr, m, mapsbc(m), 'pr', m)
call set (ips, m, mapsbc(m), 'ps', m)
call set (iaws, m, mapsbc(m), 'aws', m)
call set (iswr, m, mapsbc(m), 'swr', m)
call set (ilwr, m, mapsbc(m), 'lwr', m)
call set (isens, m, mapsbc(m), 'sens', m)
call set (ievap, m, mapsbc(m), 'evap', m)
call set (idtr, m, mapsbc(m), 'dtr', m)
call set (isr, m, mapsbc(m), 'sr', m)
call set (inpp, m, mapsbc(m), 'npp', m)
call set (iburn, m, mapsbc(m), 'burn', m)
if ( m-1 .gt. numsbc) then
print*, '==> Error: increase numsbc in csbc.h to ', m-1
stop '=>UVic_ESCM'
endif
return
end
subroutine tracer_init
!-----------------------------------------------------------------------
! Initialize ocean tracer names
!-----------------------------------------------------------------------
implicit none
integer m
include "size.h"
include "param.h"
include "pconst.h"
include "stdunits.h"
include "atm.h"
include "mw.h"
!-----------------------------------------------------------------------
! Initialize ocean tracer names
!-----------------------------------------------------------------------
mapt(:) = " "
itemp = 0
isalt = 0
iidealage = 0
idic = 0
ic14 = 0
icfc11 = 0
icfc12 = 0
ialk = 0
io2 = 0
ipo4 = 0
ino3 = 0
iphyt = 0
izoop = 0
idetr = 0
imp = 0
impa = 0
impp = 0
idetrz = 0
idiaz = 0
icocc = 0
icaco3 = 0
idetr_B = 0
m = 1
call set (itemp, m, mapt(m), 'temp', m)
call set (isalt, m, mapt(m), 'salt', m)
call set (idic, m, mapt(m), 'dic', m)
call set (ic14, m, mapt(m), 'c14', m)
call set (imp, m, mapt(m), 'mp', m)
call set (impa, m, mapt(m), 'mpa', m)
call set (impp, m, mapt(m), 'mpp', m)
call set (ialk, m, mapt(m), 'alk', m)
call set (io2, m, mapt(m), 'o2', m)
call set (ipo4, m, mapt(m), 'po4', m)
call set (iphyt, m, mapt(m), 'phyt', m)
call set (izoop, m, mapt(m), 'zoop', m)
call set (idetr, m, mapt(m), 'detr', m)
call set (idetrz, m, mapt(m), 'detrz', m)
call set (ino3, m, mapt(m), 'no3', m)
call set (idiaz, m, mapt(m), 'diaz', m)
call set (icocc, m, mapt(m), 'cocc', m)
call set (icaco3, m, mapt(m), 'caco3', m)
call set (idetr_B, m, mapt(m), 'detr_B', m)
if ( m-1 .gt. nt) then
print*, '==> Error: increase nt for tracers in size.h'
stop '=>UVic_ESCM'
endif
!-----------------------------------------------------------------------
! Initialize ocean biology tracer names
! Added by KK
!-----------------------------------------------------------------------
mapb(:) = " "
ibion = 0
ibiop = 0
ibioz = 0
ibiod = 0
ibiomp = 0.
ibiompa = 0.
ibiompp = 0.
ibiodz = 0.
ibiono3 = 0
ibiodiaz = 0
ibioc = 0
ibiocaco3 = 0
ibiod_B = 0
m = 1
call set (ibion, m, mapb(m), 'po4', m)
call set (ibioc, m, mapb(m), 'cocc', m)
call set (ibiop, m, mapb(m), 'phyt', m)
call set (ibioz, m, mapb(m), 'zoop', m)
call set (ibiod, m, mapb(m), 'detr', m)
call set (ibiomp, m, mapb(m), 'mp', m)
call set (ibiompa, m, mapb(m), 'mpa', m)
call set (ibiompp, m, mapb(m), 'mpp', m)
call set (ibiodz, m, mapb(m), 'detrz', m)
call set (ibiono3, m, mapb(m), 'no3', m)
call set (ibiodiaz, m, mapb(m), 'diaz', m)
call set (ibiocaco3, m, mapb(m), 'caco3', m)
call set (ibiod_B, m, mapb(m), 'detr_B', m)
if ( m-1 .gt. ntnpzd) then
print*, '==> Error: increase ntnpzd for tracers in size.h'
stop '=>UVic_ESCM'
endif
!-----------------------------------------------------------------------
! Initialize ocean tracer source names, must have equivalent tracer
!-----------------------------------------------------------------------
mapst(:) = " "
itrc(:) = 0
m = 1
call set (isdic, m, mapst(m), 'sdic', m)
itrc(idic) = m-1
call set (isc14, m, mapst(m), 'sc14', m)
itrc(ic14) = m-1
call set (ismp, m, mapst(m), 'smp', m)
itrc(imp) = m-1
call set (ismpa, m, mapst(m), 'smpa', m)
itrc(impa) = m-1
call set (ismpp, m, mapst(m), 'smpp', m)
itrc(impp) = m-1
call set (isalk, m, mapst(m), 'salk', m)
itrc(ialk) = m-1
call set (iso2, m, mapst(m), 'so2', m)
itrc(io2) = m-1
call set (ispo4, m, mapst(m), 'spo4', m)
itrc(ipo4) = m-1
call set (isphyt, m, mapst(m), 'sphyt', m)
itrc(iphyt) = m-1
call set (iszoop, m, mapst(m), 'szoop', m)
itrc(izoop) = m-1
call set (isdetr, m, mapst(m), 'sdetr', m)
itrc(idetr) = m-1
call set (isdetrz, m, mapst(m), 'sdetrz', m)
itrc(idetrz) = m-1
call set (isno3, m, mapst(m), 'sno3', m)
itrc(ino3) = m-1
call set (isdiaz, m, mapst(m), 'sdiaz', m)
itrc(idiaz) = m-1
call set (iscocc, m, mapst(m), 'scocc', m)
itrc(icocc) = m-1
call set (iscaco3, m, mapst(m), 'scaco3', m)
itrc(icaco3) = m-1
call set (isdetr_B, m, mapst(m), 'sdetr_B', m)
itrc(idetr_B) = m-1
if ( m-1 .gt. nt) then
print*, '==> Error: increase nsrc for tracer sources in size.h'
stop '=>UVic_ESCM'
endif
!-----------------------------------------------------------------------
! Initialize atmosphere tracer names
!-----------------------------------------------------------------------
mapat(:) = " "
isat = 0
ishum = 0
ico2 = 0
m = 1
call set (isat, m, mapat(m), 'sat', m)
call set (ishum, m, mapat(m), 'shum', m)
if ( m-1 .gt. nat) then
print*, '==> Error: increase nat in size.h'
stop '=>UVic_ESCM'
endif
return
end
subroutine read_namelist
!-----------------------------------------------------------------------
! read all model namelist variables
!-----------------------------------------------------------------------
implicit none
character (120) :: fname, new_file_name, logfile
integer iotraj, i, j, k, ip, kr, jq, n, ioun, num_processors
logical exists, annlevobc, annlev, initpt
real dtatms, snapint, snapls, snaple, snapde, ahbkg, runstep
real slmx, s_dm, afkph, dfkph, sfkph, zfkph, ahs, ahb, trajint
real crops_yr
include "size.h"
include "param.h"
include "pconst.h"
include "stdunits.h"
include "accel.h"
include "calendar.h"
include "cembm.h"
include "coord.h"
include "cnep.h"
include "csbc.h"
include "cprnts.h"
include "diag.h"
include "emode.h"
include "fwa.h"
include "hmixc.h"
include "insolation.h"
include "iounit.h"
include "isopyc.h"
include "mtlm.h"
include "npzd.h"
include "scalar.h"
include "sed.h"
include "stab.h"
include "switch.h"
include "tmngr.h"
include "veg.h"
include "vmixc.h"
!-----------------------------------------------------------------------
! read all model namelist variables
!-----------------------------------------------------------------------
namelist /contrl/ init, runlen, rununits, restrt, initpt
&, num_processors, runstep, initbgc
namelist /tsteps/ dtts, dtuv, dtsf, dtatm, dtatms, namix, segtim
&, daylen
namelist /riglid/ mxscan, sor, tolrsf, tolrsp, tolrfs
namelist /mixing/ am, ah, ahbkg, ambi, ahbi, kappa_m, kappa_h
&, cdbot, spnep, senep, aidif, ncon, nmix, eb
&, acor, dampts, dampdz, annlev, annlevobc
namelist /diagn/ tsiint, tsiper, tavgint, itavg, tmbint, tmbper
&, itmb, stabint, zmbcint, glenint, trmbint, itrmb
&, vmsfint, gyreint,igyre, extint, prxzint, trajint
&, exconvint, dspint, dspper, snapint, snapls
&, snaple, snapde, timavgint, timavgper, cmixint
&, prlat, prslon, prelon, prsdpt, predpt
&, slatxy, elatxy, slonxy, elonxy
&, cflons, cflone, cflats, cflate, cfldps, cfldpe
&, maxcfl, xbtint, xbtper, crossint, densityint
&, fctint, tyzint, restint, tbtint, tbtper
&, accel, accel_yr0
namelist /io/ expnam, logfile, runstamp, iotavg, iotmb, iotrmb
&, ioglen, iovmsf, iogyre, ioprxz, ioext, iodsp
&, iotsi, iozmbc, iotraj, ioxbt, isot1, ieot1
&, isot2, ieot2, jsot, jeot, ksot, keot, mrot
namelist /ictime/ year0, month0, day0, hour0, min0, sec0, ryear
&, rmonth, rday, rhour, rmin, rsec, refrun, refinit
&, refuser, eqyear, eqmon, monlen, init_time
&, init_time_in, init_time_out, omega
namelist /npzd/ alpha, kw, kc, abio, bbio, cbio, k1n, nup
&, gamma1, gbio0, nuz, nud0, wdd, alpha_D
&, wmp0, pinp, pbc, pagg, prise, kp, ping
&, wdz0, rpoo
&, wd0, par, dtnpzd, redctn, ki, redptn, capr
&, dcaco3, nupt0, redotn, jdiar, kzoo, geZ
&, zprefP, zprefZ, zprefDet, zprefDiaz, kfe, kfe_D
&, zprefMP, zprefDiat, zprefC, ro2n, o2t
&, abioc, k1n_C, nuct0,nuc, alpha_C
&, kfe_C
&, bapr
&, wc0, dissk0, wdc
&, kc_c
namelist /fwa/ isfwa1, iefwa1, isfwa2, iefwa2, jsfwa, jefwa
&, mrfwa, fwaflxi, fwayri, fwayrf, fwarate
&, compensate
namelist /blmix/ Ahv, Ahh
namelist /hlmix/ hl_depth, hl_back, hl_max
namelist /isopyc/ slmx, ahisop, athkdf, del_dm, s_dm
namelist /ppmix/ wndmix, fricmx, diff_cbt_back, visc_cbu_back
&, visc_cbu_limit, diff_cbt_limit
namelist /smagnl/ diff_c_back
namelist /sed/ dtsedyr, nsedacc, weath, kmin
namelist /embm/ rlapse, rf1, rf2, scatter, rhmax, vcsref, vcsfac
&, vcsyri, aggfor_os, adiff, cropf, pastf
namelist /carbon/ co2ccn, co2emit, co2for, co2_yr, c14_yr, dc14ccn
&, c14prod
namelist /paleo/ orbit_yr, pyear, eccen, obliq, mvelp, sealev
&, sealev_yr, volcano_yr, sulph_yr, aggfor_yr
&, cfcs_yr
namelist /ice/ niats, nivts, dampice, tsno, hsno_max, ice_yr
&, landice_yr, ice_calb, sno_calb
namelist /veg/ veg_alb, veg_rl, veg_rs, veg_smd, agric_yr
&, crops_yr, iagric, icrops, idesert, iice
namelist /solar/ solarconst, solar_yr
namelist /mtlm/ TIMESTEP, INT_VEG, VEG_EQUIL, DAY_TRIF, DAY_PHEN
&, BF
! physical constants
rho0 = 1.035
rho0r = c1/rho0
grav = 980.6
radius = 6370.0e5
pi = atan(1.0) * 4.0
! set defaults for namelist contrl
init = .true.
runlen = 365.0
rununits = 'days'
restrt = .true.
initpt = .false.
num_processors = 1
runstep = -1.0
initbgc = .true.
! set defaults for namelist tsteps
dtts = 43200.0
dtuv = 600.0
dtsf = 600.0
dtatm = 43200.0
dtatms = 43200.0
namix = 16
segtim = 1.0
daylen = 86400.0
! set defaults for namelist riglid
mxscan = 300
sor = 1.60
tolrsf = 5.0e8
tolrsp = 1.0e4
tolrfs = 1.0e4
! set defaults for namelist mixing
am = 2.0e9
ah = 2.0e7
ahbkg = 0.
ambi = 1.0e23
ahbi = 5.0e22
kappa_m = 10.0
kappa_h = 1.0
cdbot = 1.3e-3
spnep = 3.0e5
senep = 12.0e5
aidif = 1.0
ncon = 1
nmix = 16
eb = .false.
acor = 0.0
do n=1,nt
dampts(n) = 50.0
dampdz(n) = 26.575e2
enddo
annlev = .false.
annlevobc = .false.
! set defaults for namelist diagn
tsiint = 1.0
tsiper = 1.0
tavgint = -36500.0
itavg = .true.
tmbint = -36500.0
tmbper = 365.0
itmb = .true.
stabint = -36500.0
zmbcint = -36500.0
glenint = -36500.0
trmbint = -36500.0
itrmb = .true.
vmsfint = -36500.0
gyreint = -36500.0
igyre = .true.
extint = -36500.0
prxzint = -36500.0
trajint = 0.0
exconvint = -36500.0
dspint = -36500.0
dspper = -365.0
snapint = 0.0
snapls = 0.0
snaple = 0.0
snapde = 0.0
timavgint = 36500.0
timavgper = 365.0
cmixint = -36500.0
do n=1, nlatpr
prlat(n) = 100.0
prslon(n) = 0.0
prelon(n) = 0.0
prsdpt(n) = 0.0
predpt(n) = 6000.0e2
if (n. le. 4) then
prslon(n) = 180.0
prelon(n) = 250.0
endif
enddo
prlat(1) = -60.0
prlat(2) = 0.0
prlat(3) = 27.0
prlat(4) = 55.0
slatxy = -90.0
elatxy = 90.0
slonxy = 3.0
elonxy = 357.0
cflons = 0.0
cflone = 360.0
cflats = -90.0
cflate = 90.0
cfldps = 0.0
cfldpe = 6000.0e2
maxcfl = 3
xbtint = -36500.0
xbtper = -365.0
crossint = 365000.0
densityint = -36500.0
fctint = -36500.0
tyzint = -36500.0
restint = 36500.0
tbtint = -36500.0
tbtper = -365.0
accel = 1.
accel_yr0 = 0.
! set defaults for namelist io
expnam = ' '
logfile = 'machine.log'
runstamp = ' '
restrt = .false.
iotavg = -1
iotmb = -1
iotrmb = -1
ioglen = -1
iovmsf = -1
iogyre = -1
ioprxz = -1
ioext = -1
iodsp = -1
iotsi = -1
iozmbc = -1
iotraj = -1
ioxbt = -1
isot1 = 2
ieot1 = imtm1
isot2 = 2
ieot2 = 1
jsot = 2
jeot = jmtm1
ksot = 1
keot = km
mrot = 0
! set defaults for namelist ictime
year0 = 1
month0 = 1
day0 = 1
hour0 = 0
min0 = 0
sec0 = 0
ryear = 1
rmonth = 1
rday = 1
rhour = 0
rmin = 0
rsec = 0
refrun = .false.
refinit = .true.
refuser = .false.
eqyear = .true.
eqmon = .false.
monlen = 30
init_time = .false.
init_time_in = .false.
init_time_out = .false.
omega = pi/43082.0
! set defaults for namelist npzd
alpha = 0.1 ! Initial slope P-I curve [(W/m^2)^(-1)/day]
alpha_D = 0.38
kw = 0.04 ! Light attenuation due to water [1/m]
kc = 0.07 ! Light atten. by phytoplankton [1/(m*mmol/m^3)]
ki = 5.0 ! Light attenuation through sea ice & snow
alpha_C = 0.1 ! Initial slope P-I curve [(W/m^2)^(-1)/day]
abioc = 1.03 ! a; Maximum growth rate parameter coccolithophore [1/day]
k1n_C = 0.00074 ! Half sat const for N uptake (coccolithophores) [mmol/m^3]
nuc = 0.03 ! Specific mortality rate (coccolithophores) [1/day]
nuct0 = 0.015 ! Specific mortality rate (coccolithophores) [1/day]
bapr = 0.05 !+-0.011 detritus to carbonate ratio [mg POC/mg PIC]
kc_c = 0.2 ! Light atten. by calcite [1/(m*mmol/m^3)]
wc0 = 35.0 ! constant calcite sinking speed [m/day]
wdc = 5.0e-2 ! increase in calcite sinking speed with depth
! Gehlen et al 2007 use a surface sink speed of 50 m/day inc w/ depth
! dissk0 = 0.013 !dissolution rate parameter [1/day] from Gehlen et al 2007
dissk0 = 0.
kcal = 100. !mol m^-3 half saturation for PIC dissl from PlankTOM10
abio = 1.0 ! a; Maximum growth rate parameter [1/day]
bbio = 1.066 ! b
cbio = 1.0 ! c [1/deg_C]
k1n = 0.001 ! Half saturation constant for N uptake [mmol/m^3]
nup = 0.03 ! Specific mortality rate (Phytoplankton) [1/day]
nupt0 = 0.015 ! Specific mortality rate (Phytoplankton) [1/day]
gamma1 = 0.70 ! gama1; Assimilation efficiency (zpk)
gbio0 = 0.4 ! Maximum grazing rate at 0 deg C [1/day]
nuz = 0.06 ! Quadratic mortality (zpk) [(mmol/m^3)^(-2)day^(-1)]
nud0 = 0.07 ! remineralization rate [1/day]
wd0 = 12. ! Sinking speed of detritus at surface [m/day]
wmp0 = 0.019 !rising plastic rate m/s
pinp = 0.2 !portion of total plastic production into ocean
pbc = 0.2 !portion of bottom amount returned to water colum
pagg = 0.03 !aggregation rate
prise = 0.02 ! portion of free plastic that rises
ping = 1.0 ! replacement of food by plastic
kp= 10. !half saturation constant for free plastic uptake by mar snow
wdz0 = 24. !sinking poo rate m/d
rpoo = 0.5 ! ratio of zoop detr going in to poo
wdd = 6.0e-2 ! increase in sinking speed with depth
par = 0.43 ! fraction of photosythetically active radiation
dtnpzd = dtts ! time step of biology [s]
redctn = 6.625 ! C/N Redfield ratio
redptn = 1./16. ! P/N Redfield ratio
capr = 0.04 ! carbonate to carbon production ratio
dcaco3 = 650000.0 ! remineralisation depth of calcite [cm]
redotn = 10. ! O2/N Redfield ratio
jdiar = 0.4 ! factor reducing the growth rate of diazotrophs
kzoo = 0.15 ! half sat. constant for Z grazing [mmol N m^3]
geZ = 0.40 ! Zooplankton growth efficiency
! graz prefs needed to be restructured- too confusing w/so many choices! KK
zprefMP = 0. !zooplankton preference for microplastic
zprefP = 0. ! Zooplankton preference for P
zprefC = 0. ! Zooplankton preference for C
zprefZ = 0. ! Zooplankton preference for other Z
zprefDet = 0. ! Zooplankton preference for Detritus
zprefDiaz = 0. ! Zooplankton preference for diazotrophs
zprefDiat = 0. ! Zooplankton preference for diatoms
ro2n = 800. !ratio of o to n for denitrification, default is 800
o2t = 5. !threshold for denitrification, default is 5 mmol O2
kmfe = 3 ! number of depth levels used for Fe limitation
kfe = 0.12 ! Half saturation constant for Fe limitation
kfe_D = 0.1 ! Half saturation constant for Diaz Fe limitation
kfe_C = 0.06
! set defaults for namelist fwa
isfwa1 = 1
iefwa1 = imt
isfwa2 = 0
iefwa2 = 0
jsfwa = 1
jefwa = jmt
mrfwa = 1
fwaflxi = 0.
fwayri = -1.e20
fwayrf = 1.e20
fwarate = 0.
compensate = .false.
! set defaults for namelist blmix
! Reference:
! A Water Mass Model of the World Ocean K. Bryan, L.J. Lewis
! JGR, vol 84, No. C5, May 20, 1979
afkph = 0.8
dfkph = 1.05
sfkph = 4.5e-5
zfkph = 2500.0e2
! Use Bryan & Lewis values for vertical tracer diffusion
! Ahv range of 0.3 to 1.3, crossover at 2500m.
! compute depth dependent vertical diffusion coefficients for
! tracers using the relationship of Bryan and Lewis
pi = 4.0 * atan(1.0)
! compute depth dependent horizontal diffusion coefficients for
! tracers using the relationship of Bryan and Lewis
ahs = 5.0e+3
ahb = 1.0e+3
! set defaults for namelist hlmix
hl_depth = 500.0e2
hl_back = 1.e4
hl_max = 1.e9
! set defaults for namelist isopyc
slmx = 1.0/100.0 ! maximum isopycnal slope
ahisop = 1.e7 ! isopycnal diffusion coefficient
athkdf = 1.0e7 ! isopycnal thickness diffusion coefficient
del_dm = 4.0/1000.0 ! transition for scaling diffusion coefficient
s_dm = 1.0/1000.0 ! half width scaling for diffusion coefficient
! set defaults for namelist ppmix
wndmix = 10.0
fricmx = 50.0
diff_cbt_back = 0.1
visc_cbu_back = 1.0
! in regions of gravitational instability set mixing limits to the
! maximum consistant with the "cfl" criterion. convective adjustment
! will also act on the instability.
visc_cbu_limit = fricmx
diff_cbt_limit = fricmx
! set defaults for namelist smagnl
diff_c_back = c0
! set defaults for namelist sed
dtsedyr = 1. ! time step of sediment model in years
nsedacc = 1 ! number of steps for accelerating sediments
weath = 2.e20 ! weathering flux in Pg year-1
! set high to check for constant namelist input
kmin = 8 ! minimum ocean level for sediments
! set defaults namelist embm
rlapse = 5.e-5
rf1 = 0.3
rf2 = 3.e5
scatter = 0.23
rhmax = 0.85
vcsref = 0.
vcsfac = 0.
vcsyri = 0.
aggfor_os = 0.
adiff = 0.03
cropf = 1.
pastf = 0.5
! set defaults namelist carbon
co2ccn = 280.
co2emit = 0.
co2for = 5.35e03
co2_yr = 2.e20
c14_yr = 2.e20
dc14ccn = 0.
c14prod = 0.
! set defaults namelist paleo
pyear = 2.e20
orbit_yr = 2.e20
! default user orbit is for 1950
eccen = 1.672393E-02
obliq = 23.44627
mvelp = 102.0391
sealev = 0.
sealev_yr = 2.e20
volcano_yr = 2.e20
sulph_yr = 2.e20
aggfor_yr = 2.e20
cfcs_yr = 2.e20
! set defaults namelist ice
niats = 1
nivts = 1
dampice = 5.
tsno = 0.
hsno_max = 1000.
ice_yr = 2.e20
landice_yr = 2.e20
ice_calb = 0.25
sno_calb = 0.2
! set defaults namelist veg
veg_alb(1:7) = (/0.17, 0.17, 0.22, 0.22, 0.22, 0.30, 0.60/)
veg_rl(1:7) = (/1.5, 1.0, 0.1, 0.3, 0.01, 0.01, 0.005/)
veg_rs(1:7) = (/0.75, 1.0, 1.0, 1.0, 1.0, 1.0, 1.0/)
veg_smd(1:7) = (/2.0, 3.0, 0.1, 0.1, 0.5, 0.01, 0.01/)
agric_yr = 2.e20
crops_yr = 2.e20
iagric = 3
icrops = 3
idesert = 6
iice = 7
! set defaults namelist solar
solarconst = 1.368e6 ! mw m-2
solar_yr = 2.e20
! set defaults namelist mtlm
TIMESTEP = 3600.
INT_VEG = .true.
VEG_EQUIL = .false.
DAY_TRIF = 30
DAY_PHEN = 1
BF = 1.
!-----------------------------------------------------------------------
! provide for change in above presets using "namelist"
!-----------------------------------------------------------------------
call getunit (ioun, 'control.in', 'f s r')
read (ioun, contrl, end=101)
101 continue
runstep = float(int(runstep))
write (stdout, contrl)
call relunit (ioun)
call getunit (ioun, 'control.in', 'f s r')
read (ioun, tsteps, end=102)
102 continue
write (stdout, tsteps)
call relunit (ioun)
call getunit (ioun, 'control.in', 'f s r')
read (ioun, riglid, end=103)
103 continue
write (stdout, riglid)
call relunit (ioun)
call getunit (ioun, 'control.in', 'f s r')
read (ioun, mixing, end=104)
104 continue
write (stdout, mixing)
call relunit (ioun)
call getunit (ioun, 'control.in', 'f s r')
read (ioun, diagn, end=105)
105 continue
write (stdout, diagn)
call relunit (ioun)
call getunit (ioun, 'control.in', 'f s r')
read (ioun, io, end=106)
106 continue
write (stdout, io)
call relunit (ioun)
call getunit (ioun, 'control.in', 'f s r')
read (ioun, ictime, end=107)
107 continue
write (stdout, ictime)
call relunit (ioun)
call getunit (ioun, 'control.in', 'f s r')
read (ioun, npzd, end=108)
108 continue
write (stdout, npzd)
call relunit (ioun)
call getunit (ioun, 'control.in', 'f s r')
read (ioun, fwa, end=109)
109 continue
write (stdout, fwa)
call relunit (ioun)
call getunit (ioun, 'control.in', 'f s r')
read (ioun, blmix, end=110)
110 continue
write (stdout, blmix)
call relunit (ioun)
call getunit (ioun, 'control.in', 'f s r')
read (ioun, hlmix, end=111)
write (stdout, hlmix)
call relunit (ioun)
111 continue
call getunit (ioun, 'control.in', 'f s r')
read (ioun, isopyc, end=112)
112 continue
write (stdout, isopyc)
call relunit (ioun)
call getunit (ioun, 'control.in', 'f s r')
read (ioun, ppmix, end=113)
113 continue
write (stdout, ppmix)
call relunit (ioun)
call getunit (ioun, 'control.in', 'f s r')
read (ioun, smagnl, end=114)
114 continue
write (stdout, smagnl)
call relunit (ioun)
call getunit (ioun, 'control.in', 'f s r')
read (ioun, sed, end=115)
115 continue
write (stdout, sed)
call relunit (ioun)
call getunit (ioun, 'control.in', 'f s r')
read (ioun, embm, end=116)
116 continue
write (stdout, embm)
call relunit (ioun)
call getunit (ioun, 'control.in', 'f s r')
read (ioun, carbon, end=117)
117 continue
write (stdout, carbon)
call relunit (ioun)
call getunit (ioun, 'control.in', 'f s r')
read (ioun, paleo, end=118)
118 continue
write (stdout, paleo)
call relunit (ioun)
call getunit (ioun, 'control.in', 'f s r')
read (ioun, ice, end=119)
119 continue
write (stdout, ice)
call relunit (ioun)
call getunit (ioun, 'control.in', 'f s r')
read (ioun, veg, end=120)
120 continue
write (stdout, veg)
call relunit (ioun)
call getunit (ioun, 'control.in', 'f s r')
read (ioun, solar, end=121)
121 continue
write (stdout, solar)
call relunit (ioun)
call getunit (ioun, 'control.in', 'f s r')
read (ioun, mtlm, end=122)
122 continue
write (stdout, mtlm)
call relunit (ioun)
! limit or convert some variables
! ensure pass is between zero and one.
pass = min(max((1. - scatter), 0.), 1.)
dtocn = dtts
nots = nmix
isot1 = max(isot1, 1)
ieot1 = min(ieot1, imt)
isot2 = max(isot1, 1)
ieot2 = min(ieot1, imt)
jsot = max(jsot, 1)
jeot = min(jeot, jmt)
ksot = max(ksot, 1)
keot = min(keot, km)
do k=1,km
dtxcel(k) = 1.0
enddo
ah = ahbkg
slmxr = 0.
if (slmx .gt. 0.) slmxr = c1/slmx
s_dm = 0.
if (s_dm .gt. 0.) s_dmr = c1/s_dm
! sort out forcing years
if (pyear .gt. 1.e20) pyear = 1800.
if (c14_yr .gt. 1.e20) c14_yr = pyear
if (orbit_yr .gt. 1.e20) orbit_yr = pyear
if (volcano_yr .gt. 1.e20) volcano_yr = pyear
if (sulph_yr .gt. 1.e20) sulph_yr = pyear
if (aggfor_yr .gt. 1.e20) aggfor_yr = pyear
if (cfcs_yr .gt. 1.e20) cfcs_yr = pyear
if (sealev_yr .gt. 1.e20) sealev_yr = pyear
if (agric_yr .gt. 1.e20 .and. crops_yr .lt. 1.e20)
& agric_yr = crops_yr
if (landice_yr .gt. 1.e20 .and. ice_yr .lt. 1.e20)
& landice_yr = ice_yr
ice_yr = landice_yr
tsiper = min(tsiper, tsiint)
tbtper = min(tbtper, tbtint)
if (accel .le. 1) accel = 1.
runlen = runlen/accel
runstep = runstep/accel
if (runstep .gt. 0.0) runlen = min(runstep, runlen)
if (init_time) then
init_time_in = .true.
init_time_out = .true.
endif
monlen = 30
eqyear = .true.
eqmon = .false.
calendar = 'noleap'
!-----------------------------------------------------------------------
! set runstamp
!-----------------------------------------------------------------------
if (runstamp .eq. ' ') then
! read runstamp from last line of logfile if it exists
fname = new_file_name (logfile)
inquire (file=trim(fname), exist=exists)
if (exists) then
call getunit (ioun, fname, 'f s r')
do while (exists)
read (ioun, '(a130)', end=130) runstamp
enddo
130 continue
call relunit (ioun)
print*, " "
print*, "runstamp: ", trim(runstamp)
endif
endif
return
end