Skip to content

GitLab

Commit 7f6c28a0 authored by ph290's avatar ph290
Browse files

adding more tidal components

parent fe0cb198
Hide whitespace changes
Inline Side-by-side
Showing with 2715 additions and 0 deletions
main/s2p3_rv2.0_more_tide.f90 0 → 100755
View file @ 7f6c28a0
!***************************************************************************************
!
! S2P3 v7.0
!
!------------------Shelf Sea Physics and Primary Production v7.0------------------------
!
! 1-D MODEL OF THE EQUATION OF MOTION USING THE Canuto k-e TURBULENCE CLOSURE SCHEME
! with a simple single-species model of carbon fixation in response to light and nitrate
!
!***************************************************************************************
!
! Jonathan Sharples
! University of Liverpool and NERC National Oceanography Centre
!
! In: J. H. Simpson and J. Sharples
! Introduction to the Physical and Biological Oceanography of Shelf Seas
! Cambridge University Press 2012
!
!***************************************************************************************
!
! Modified for regional application, compiled and executed in unix environments
! (Marsh, Hickman and Sharples, Geoscientific Model Development, submitted)
!
!***************************************************************************************!
!
MODULE physics
implicit none
double precision :: lon,lat,kb,f0,f1,lambda,heat_shade,vismax,qsdt,stressx,stressy,windspeed,polarisn(6),rad_out, &
dewT,stresss,stressb,Px,Py,slopex(5),slopey(5),tsteps,orient(5),semi_major(5),semi_minor(5), &
first_temp,sx,sy,radsum,omega(5),major_phase(5),background_mol,u_bed,v_bed,speed_mag(200),mean_speed, &
rad0,speed,speed1,speed2,speed3,time,depth,newdepth,grav_accel,radiation_zerocloud(40000),radiation_fullcloud(40000), &
uac,uc,gac1,gac2,gac,gc1,gc2,gc,strat_jul,rad_sum,acount
DOUBLE PRECISION :: wind_speed(40000),wind_dir(40000),std_dir(40000),std_speed(40000),radiation(40000),humid(40000),airP(40000),airT(40000), &
cloud(40000),declination(40000)
real :: tid,tid2,uamp(5),vamp(5),vphase(5),uphase(5)
integer :: idmet,iday_front1,iday_front2,iday_front3,iday_front4
data omega /1.405278d-4,1.454440d-4,1.378616d-4,6.759596d-5,7.293472e-5/
END MODULE
MODULE variables_all
implicit none
DOUBLE PRECISION :: c1,c2,c3,rad,time_step,dz,newdz
double precision :: cnpar=0.5
INTEGER :: N,newN,i,dataplot,ikey,icont,iquit,message_to_quit
INTEGER :: imode
END MODULE
MODULE turbulence
implicit none
double precision :: NN(0:200),SS(0:200),lscale(0:200),h(200),NN_mean(0:200),SS_mean(0:200),DD_mean(0:200),RN_mean(0:200)
double precision :: Kz(0:200),Nz(0:200),sm(0:200),KK_mean(0:200)
double precision :: sh(0:200),Ri(0:200),P(0:200),B(0:200)
double precision :: cmue1(0:200),cmue2(0:200),as(0:200),an(0:200)
double precision, allocatable :: au(:),bu(:),cu(:),du(:),ru(:),qu(:),tke(:),tkeold(:),eps(:)
double precision :: kappa=0.41
double precision :: Kz_bg,Nz_bg ! Used in simple IW model
END MODULE
MODULE biology
implicit none
double precision :: x_old(200),x_new(200),ni_old(200),ni_new(200),quo(200), &
net_grow(200),gross_grow(200),prod_net(200),prod_gross(200),prod_net_daily(200),prod_gross_daily(200),uptake(200),g0(40000)
double precision :: sub_quota,alpha,alpha_persec,chl_carbon,max_quota, &
Nrecycle,Pmax10,Pmax10_persec,Q10,gT0,uptake_max,uptake_max_persec, &
half_saturation,graze_thresh,respiration,respiration_persec,rT0,rQ10, &
swim_speed,swim_speed_persec,sink_speed, &
sink_speed_persec,graze_min,graze_min_persec,graze_max,graze_max_persec,graze_daymax, &
seed_biomass,chl_abscross,x_new_max,temp_x_new_max
DOUBLE PRECISION :: s_old(200),s_new(200),din_source(200),radbar(200),par_percent, &
din_rate,bed_din,par_atten,par_shade,graz(200),Pmax,Pm
integer :: growth_model
integer, parameter :: eppley=1
integer, parameter :: Qten=2
END MODULE
MODULE file_names
implicit none
character(len=100) :: metfile_in,surface_out,phys_profile_day,phys_profile_hour, &
bio_profile_day,bio_profile_hour,initialdata,monthly_out,user_metfile
END MODULE
!
MODULE VARIABLES
!
! Shared variables for any routine with 'USE VARIABLES'
!
IMPLICIT NONE
!
INTEGER, PARAMETER :: MAX_CHILD = 20 ! Max child windows
INTEGER :: ICHILD = 0 ! No. of children
INTEGER, DIMENSION(MAX_CHILD) :: CHILDREN ! Child handles
LOGICAL :: PLOTT = .FALSE. ! Graphic plotted?
!
DOUBLE PRECISION :: velx_old(200),velx_new(200),vely_old(200), &
vely_new(200),ri_mean(200),Kz_mean(200), &
tke_mean(200), diss_mean(200),rad_mean(200), &
temp_old(200),temp_new(200),density(200), &
vel_mean(200),u_mean(200),v_mean(200),c1flux_mean(200), &
n1flux_mean(200),sflux_mean(200),hourly_net, &
hourly_gross,daily_net,daily_gross,prod_net2, &
prod_gross2,tpn1,tpn2,tpg1,tpg2,total_flux, &
surf_gross1,surf_net1,total_uptake,surf_uptake
real :: height(200),xtotal(200),wind_distribution(100)
INTEGER :: n_rad_choice,nhr_out,ndays,infile_error,iyear,met_type,metfile_error
integer :: output_start,output_end, output_hr_start, output_hr_end,ibmp(5),irun
END MODULE VARIABLES
!
MODULE GRAPHICS_VARIABLES
implicit none
REAL, allocatable :: contour_data(:,:),xdata(:),ydata(:,:),profile(:,:)
REAL :: r_height(200),r_dz,parmax,velmax,wmax
real :: hu3a(5000),h_hu3(5000),t_hu3(5000),x1_hu3(5000),n1_hu3(5000),u1_hu3(5000),g1_hu3(5000), &
x2_hu3(5000),n2_hu3(5000),u2_hu3(5000),g2_hu3(5000),s_hu3(5000),hkz_hu3(5000),Kz_hu3(5000), &
hu3b(5000),sf_hu3(5000),c1f_hu3(5000),n1f_hu3(5000)
!
real :: originx(15),originy(15),startx(15),starty(15),deltax(15),deltay(15),xincr(15),yincr(15), &
endx(15),endy(15),xlen(15),ylen(15),ticklen,temp_min,temp_max,zmin,zmax1,zmax2,z_scale
real :: charhx=0.006
real :: charhy=0.0225
INTEGER ncolour(15),ig(15),nticksx(15),nticksy(15),irunplot,screenx,screeny,irenew,black,white,JD_front,ifront
character (len=15) :: labelx(15),labely(15),labelz(15),last_x(12),last_y(12)
integer red(29),green(29),blue(29)
data red /193,167,138,106,72,36,0,0,0,0,0,0,0,0,0,0,0,36,72,106,138,167,193,214,231,244,252, &
255,252/
data green /0,0,0,0,0,0,0,36,72,106,138,167,193,214,231,244,252,255,252,244,231,214,193,167,138, &
106,72,36,0/
data blue /138,167,193,214,231,244,252,255,252,244,231,214,193,167,138,106,72,36,0,0,0, &
0,0,0,0,0,0,0,0/
END MODULE GRAPHICS_VARIABLES
!
!*****************************************************************************
!
MODULE INTERFACES
IMPLICIT NONE
INTERFACE
SUBROUTINE ProcessMenu(IDENT,QUIT,unique_job_id)
IMPLICIT NONE
INTEGER, INTENT(IN ) :: IDENT
LOGICAL, INTENT(IN OUT) :: QUIT
character(len=36) :: unique_job_id
END SUBROUTINE ProcessMenu
!
SUBROUTINE About()
IMPLICIT NONE
END SUBROUTINE About
!
SUBROUTINE Redraw()
IMPLICIT NONE
END SUBROUTINE Redraw
!
SUBROUTINE Plotit(IPLT)
IMPLICIT NONE
INTEGER, INTENT (IN ) :: IPLT
END SUBROUTINE Plotit
SUBROUTINE WExit(QUIT,IWIN)
IMPLICIT NONE
LOGICAL, INTENT (IN OUT) :: QUIT
INTEGER, OPTIONAL, INTENT (IN ) :: IWIN
END SUBROUTINE WExit
END INTERFACE
END MODULE INTERFACES
!
!*****************************************************************************
!
PROGRAM phytowin
!
! Program generated by WiDE Wizard at 09:17 on 06 Jun 2001.
!
USE INTERFACES
use turbulence
use variables_all
use variables
use physics
use biology
use graphics_variables
use file_names
!
IMPLICIT NONE
!
! Variable declarations
!
LOGICAL :: QUIT = .FALSE.
INTEGER :: ITYPE
INTEGER :: IDENT
integer :: run_year,start_year,iline
character(len=12) :: lat_in_domain
character(len=12) :: lon_in_domain
character(len=36) :: unique_job_id
character(len=300) :: met_data_location
real :: woa_nutrient
!
! Initialise Winteracter
!
! print*, 'getting physics defaults'
call get_physics_defaults(lat_in_domain,lon_in_domain,run_year,start_year,unique_job_id,met_data_location,iline,woa_nutrient)
depth=newdepth; N=newN; dz=newdz
! print*, 'getting phyto defaults'
call get_phyto_defaults()
! print*, 'getting grazing defaults'
call get_grazing_defaults()
allocate (contour_data(1,1),xdata(1),ydata(1,1),profile(1,1))
allocate (au(1),bu(1),cu(1),du(1),ru(1),qu(1),tke(1),tkeold(1),eps(1))
! ifront=1; met_type=1; metfile_error=0
ifront=1; met_type=0; metfile_error=0
! print*, 'met_type = (1 for default)'
! read(5,*) met_type
! print*, 'getting met'
call get_met(lat_in_domain,lon_in_domain,run_year,start_year,unique_job_id,met_data_location,iline)
surface_out='surface'//unique_job_id//'.dat'; phys_profile_day='physday.dat'; bio_profile_day='biolday.dat'; monthly_out='monthly'//unique_job_id//'.dat'
phys_profile_hour='physhour.dat'; bio_profile_hour='biolhour.dat'
metfile_in='Model default meteorological data'
user_metfile=''
initialdata='Initialisation data not saved'
output_start=0; output_end=0; output_hr_start=0; output_hr_end=0; JD_front=182
temp_min=5.0; temp_max=20.0; zmin=0.0; zmax1=10.0; zmax2=10.0
irunplot=-1; irun=0
!!! new calls to ProcessMenu (rma, 21/10/11)
! print*, 'IDENT = '
IDENT = 1
! read(5,*) IDENT
CALL ProcessMenu(IDENT,QUIT,unique_job_id)
!!! end new calls to ProcessMenu (rma, 21/10/11)
STOP
END PROGRAM phytowin
!
!*****************************************************************************
!
SUBROUTINE ProcessMenu(IDENT,QUIT,unique_job_id)
!
! This subroutine processes the menu selections
!
USE INTERFACES
use turbulence
use biology
use physics
use graphics_variables
use variables_all
use variables
use file_names
!
IMPLICIT NONE
!
integer nhu3a,nhu3b
INTEGER, INTENT (IN) :: IDENT
LOGICAL, INTENT (IN OUT) :: QUIT
character(len=4) :: run_year,start_year,iline
character(len=12) :: lat_in_domain
character(len=12) :: lon_in_domain
character(len=36) :: unique_job_id
character(len=300) :: met_data_location
real :: woa_nutrient
!
! Branch depending on chosen menu item
!
call save_work(lat_in_domain,lon_in_domain,run_year,start_year,unique_job_id,met_data_location,iline,woa_nutrient)
deallocate (contour_data,xdata,ydata,profile,au,bu,cu,du,ru,qu,tke,tkeold,eps)
!
! Open required files for results, and write headers
!
open(21,file=surface_out,status='replace')
if(output_end.gt.0.0)then
open(22,file=bio_profile_day,status='replace')
open(23,file=phys_profile_day,status='replace')
end if
if(output_hr_end.gt.0.0)then
open(24,file=bio_profile_hour,status='replace')
open(25,file=phys_profile_hour,status='replace')
end if
open(26,file=monthly_out,status='replace')
if(ifront.eq.1)then
open(27,file='front_data'//unique_job_id//'.dat',status='replace')
write(27,fmt="(' hu3 height temp X1 N1 U1 G1 DIN hKz log10Kz')")
ifront=0
else
i=1
open(27,file='front_data'//unique_job_id//'.dat')
read(27,*,end=322)
320 read(27,*,end=322) hu3a(i),h_hu3(i),t_hu3(i),x1_hu3(i),n1_hu3(i),u1_hu3(i),g1_hu3(i), &
s_hu3(i),hkz_hu3(i),Kz_hu3(i)
i=i+1
goto 320
322 close(27)
nhu3a=i-1
open(27,file='front_data'//unique_job_id//'.dat',status='replace')
write(27,fmt="(' hu3 height temp X1 N1 U1 G1 DIN hKz log10Kz')")
if(nhu3a.gt.1)then
do i=1,nhu3a
write(27,fmt="(10f8.2)") hu3a(i),h_hu3(i),t_hu3(i),x1_hu3(i),n1_hu3(i),u1_hu3(i),g1_hu3(i),s_hu3(i),hkz_hu3(i),Kz_hu3(i)
end do
end if
end if
!
! Write headers to output files
WRITE(21,421)
if(output_end.gt.0)then
WRITE(22,422)
write(23,423)
end if
if(output_hr_end.gt.0)then
WRITE(24,424)
WRITE(25,425)
end if
write(26,426)
421 FORMAT(1x,' JD time SH Ts Tb Ts-Tb PHI &
& spd stress Qs Qflux CHLs CHLt &
& DINs netp grossp')
422 FORMAT(1x,' time hu3 height PAR chl1 n1 grow1 &
& uptake1 chl2 n2 grow2 uptake2 DIN')
423 FORMAT(1x,' time hu3 height temp sigmat u v &
& h_turb Ri logKz Logdiss Logtke')
424 FORMAT(1x,' time hu3 height PAR chl1 n1 grow1 &
& uptake1 chl2 n2 grow2 uptake2 DIN')
425 FORMAT(1x,' time hu3 height temp sigmat u v &
& h_turb N2 S2 Ri logKz Logdiss Logtke')
426 FORMAT(1x,' month SH Ts Tb delta-T Wstress Chlt Ct Cgross Cnet accumC')
!
CALL run_model(run_year,start_year,unique_job_id,iline,woa_nutrient)
surface_out='surface'//unique_job_id//'.dat'; phys_profile_day='physday.dat'; bio_profile_day='biolday.dat'; monthly_out='monthly'//unique_job_id//'.dat'
close(21)
if(output_end.gt.0.0)close(22)
if(output_end.gt.0.0)close(23)
if(output_hr_end.gt.0.0)close(24)
if(output_hr_end.gt.0.0)close(25)
close(26)
close(27)
close(28)
close(50)
close(55)
close(56)
close(57)
close(54)
RETURN
END SUBROUTINE ProcessMenu
!******************************************************************************************
!!!! SUBROUTINE SpawnChild(IHANDLE, in_flags,in_x,in_y,in_width,in_height,in_title)
!
!*****************************************************************************
!
Subroutine save_work(lat_in_domain,lon_in_domain,run_year,start_year,unique_job_id,met_data_location,iline,woa_nutrient)
USE VARIABLES
use physics
use turbulence
use variables_all
use biology
use file_names
implicit none
integer :: iflags,testdate(8),idate(8),iclock(6),int_date,int_clock
integer :: run_year,start_year,iline
real :: real_clock
real :: woa_nutrient_tmp, woa_nutrient
real :: smaj1, smin1, smaj2, smin2, smaj3, smin3, smaj4, smin4, smaj5, smin5, smaj6, smin6
character(len=6) :: met_time
character(len=30) :: filter
character (len=50) :: filein
character(len=8) :: date
character(len=10) :: clock,zone
character(len=100) :: initialdata2
character(len=1) :: ans
character(len=3) :: type
character(len=300) :: domain_file
character(len=300) :: nutrient_file
character(len=12) :: lat_in_domain
character(len=12) :: lon_in_domain
character(len=36) :: unique_job_id
character(len=300) :: met_data_location
if(initialdata.eq.'Initialisation data not saved')then
! initialdata2='initial.txt'
initialdata2='initial'//unique_job_id//'.txt'
else
initialdata2=initialdata
end if
initialdata=initialdata2
! initialdata = 'initial'//unique_job_id//'.txt'
open(20,file=initialdata,status='replace')
call date_and_time(date,clock,zone,testdate)
write(20,fmt="('Initialisation and driving parameters for S2P3 (Shelf Sea Physics and Primary Production) model')")
write(20,fmt="('In: Simpson & Sharples, Introduction to the Physical and Biological Oceanography of Shelf Seas')")
write(20,fmt="(' Cambridge University Press, 2012.')")
write(20,*)
write(20,fmt="('File generated:')")
write(20,*)
write(20,*)
write(20,fmt="('Physics parameters:')")
write(20,fmt="('Total depth (m) = ',f9.2)") depth
write(20,fmt="('Number of depth cells = ',i9)") N
write(20,fmt="('Depth cell thickness (m) = ',f9.2)") dz
write(20,fmt="('Time step (s) = ',f9.2)") time_step
write(20,fmt="('Longitude (degrees, positive east) = ',f9.2)") lon
write(20,fmt="('Latitude (degrees, positive north) = ',f9.2)") lat
write(20,fmt="('Bottom quadratic drag coefficient = ',f9.5)") kb
write(20,fmt="('Maximum diffusivity and viscosity (m2 s-1) = ',f9.6)") vismax
write(20,fmt="('Background viscosity (m2 s-1) = ',f9.6)") Nz_bg
write(20,fmt="('Background diffusivity (m2 s-1) = ',f9.6)") Kz_bg
write(20,fmt="('Initial water temperature (deg C) = ',f9.2)") first_temp
write(20,fmt="('Heat vertical attenuation (m-1) = ',f9.3)") lambda
write(20,fmt="('Chl effect on heat attenuation (m2 (mg Chl)-1) = ',f9.5)") heat_shade
write(20,fmt="('PAR vertical attenuation (m-1) = ',f9.3)") par_atten
write(20,fmt="('Fraction of surface radiation that is PAR = ',f9.2)") par_percent
write(20,fmt="('Maximum seabed dissolved inorganic N (mmol m-3) = ',f9.2)") bed_din
write(20,fmt="('Maximum flux of inorganic N from seabed (mmol m-2 d-1) = ',f9.2)") din_rate
write(20,*)
write(20,fmt="('Tidal parameters:')")
write(20,fmt="(' M2 S2 N2 O1 K1')")
write(20,fmt="('u amplitude (m s-1) : ',5f9.3)") (uamp(i),i=1,5)
write(20,fmt="('u phase (radians) : ',5f9.3)") (uphase(i),i=1,5)
write(20,fmt="('v amplitude (m s-1) : ',5f9.3)") (vamp(i),i=1,5)
write(20,fmt="('v phase (radians) : ',5f9.3)") (vphase(i),i=1,5)
! Calculate ellipse parameters
f0=(4.0*3.142/(24.0*3600.0))*dsin(lat*3.142/180.0); f1=-f0 ! calculate Coriolis parameter.
do i=1,5
uac=0.5*sqrt(uamp(i)**2.0+vamp(i)**2.0+(2.0*uamp(i)*vamp(i)*sin(vphase(i)-uphase(i))))
uc=0.5*sqrt(uamp(i)**2.0+vamp(i)**2.0-(2.0*uamp(i)*vamp(i)*sin(vphase(i)-uphase(i))))
gac1=(vamp(i)*COS(vphase(i)))-(uamp(i)*SIN(uphase(i)))
gac2=(uamp(i)*COS(uphase(i)))+(vamp(i)*SIN(vphase(i)))
gc1=(vamp(i)*COS(vphase(i)))+(uamp(i)*SIN(uphase(i)))
gc2=(uamp(i)*COS(uphase(i)))-(vamp(i)*SIN(vphase(i)))
gac=dATAN2(gac1,gac2)
gc=dATAN2(gc1,gc2)
if (uamp(i).lt.0.0001.and.vamp(i).lt.0.0001) then
orient(i)=0.0
major_phase(i)=0.0
semi_major(i)=0.0
semi_minor(i)=0.0
polarisn(i)=0.0
else
orient(i)=(gac+gc)/2.0
major_phase(i)=1.57-(gac-gc)/2.0
semi_major(i)=uac+uc
semi_minor(i)=uac-uc
polarisn(i)=semi_minor(i)/semi_major(i)
end if
END do
! print*, 'over-ride M2, S2, N2 tides (y/n)?'
! read(*,'(a1)') ans
! if(ans.eq.'y') then
ans = 'y'
! open(1,file='s12_m2_s2_n2_h.dat',status='old')
! open(1,file='s12_m2_s2_n2_h_sec.dat',status='old')
! open(1,file='s12_m2_s2_n2_h_tim.dat',status='old')
! get tidal current amplitudes
! first read whether map, sec or tim
read(5,'(i4)') start_year
read(5,'(i4)') run_year
read(5,'(a12)') lat_in_domain
read(5,'(a12)') lon_in_domain
read(5,'(a300)') domain_file
read(5,'(a300)') nutrient_file
read(5,'(a36)') unique_job_id
read(5,'(a300)') met_data_location
read(5,'(a3)') type
read(5,*) iline
open(1,file=TRIM(ADJUSTL(nutrient_file)),status='old')
do i = 1,iline
read(1,'(16x,f6.1)') woa_nutrient_tmp
if(i.eq.iline) then
woa_nutrient = woa_nutrient_tmp
endif
enddo
close (1)
open(1,file=TRIM(ADJUSTL(domain_file)),status='old')
read(1,'(i1)') imode
do i = 1,iline
read(1,'(16x,12f6.1)') smaj1, smin1, smaj2, smin2, smaj3, smin3, smaj4, smin4, smaj5, smin5, smaj6, smin6
!Note, reading in 1st 6 constituents of teh tides, but only using the first 5
! write(6,'(16x,6f6.1)') smaj1, smin1, smaj2, smin2, smaj3, smin3
if(i.eq.iline) then
semi_major(1)=smaj1*1e-2
semi_minor(1)=smin1*1e-2
semi_major(2)=smaj2*1e-2
semi_minor(2)=smin2*1e-2
semi_major(3)=smaj3*1e-2
semi_minor(3)=smin3*1e-2
semi_major(4)=smaj4*1e-2
semi_minor(4)=smin4*1e-2
semi_major(5)=smaj5*1e-2
semi_minor(5)=smin5*1e-2
endif
enddo
close (1)
do i=1,5
polarisn(i)=semi_minor(i)/semi_major(i)
enddo
write(20,*)
write(20,fmt="('tidal ellipse orientation (radians) : ',5f9.3)") (orient(i),i=1,5)
write(20,fmt="('tidal ellipse polarisation : ',5f9.3)") (polarisn(i),i=1,5)
write(20,fmt="('tidal ellipse semi-major axis (m s-1): ',5f9.3)") (semi_major(i),i=1,5)
write(20,*)
write(20,fmt="('Biology parameters:')")
write(20,fmt="('Growth model (1=Eppley, 2=Q10) = ',i4)") growth_model
if(growth_model.eq.2)then
write(20,*)
write(20,fmt="(' Parameters for Q10 growth model:')")
write(20,fmt="(' Reference maximum growth rate (d-1) = ',f9.2)") Pmax10
write(20,fmt="(' Reference temperature (deg C) = ',f9.2)") gT0
write(20,fmt="(' Q10 exponent for growth = ',f9.2)") Q10
write(20,*)
end if
write(20,fmt="('Max light utilisation coefficient (mg C (mg Chl)-1 d-1 (W m-2)-1) = ',f9.4)") alpha
write(20,fmt="('Reference respiration rate (mg C (mg Chl)-1 d-1) = ',f9.4)") respiration
write(20,fmt="('Reference temperature for respiration rate (deg C) = ',f9.4)") rT0
write(20,fmt="('Q10 exponent for respiration = ',f9.4)") rQ10
write(20,fmt="('Chl:carbon (mg Chl (mg C)-1) = ',f9.4)") chl_carbon
write(20,fmt="('Near-bed seed stock of phytoplankton (mg C m-3) = ',f9.4)") seed_biomass
write(20,fmt="('Pigment absorption cross-section (m2 (mg Chl)-1) = ',f9.4)") chl_abscross
write(20,fmt="('Maximum nitrate uptake rate (mmol (mg Chl)-1 d-1) = ',f9.4)") uptake_max
write(20,fmt="('Maximum cell nutrient quota (mmol N (mg Chl)-1) = ',f9.4)") max_quota
write(20,fmt="('Subsistence cell nutrient quota (mmol N (mg Chl)-1) = ',f9.4)") sub_quota
write(20,fmt="('Nitrate uptake half-saturation concentration (mmol m-3) = ',f9.4)") half_saturation
write(20,fmt="('Swimming speed (m d-1) = ',f9.4)") swim_speed
write(20,fmt="('Sinking speed (m d-1) = ',f9.4)") sink_speed
write(20,fmt="('Minimum grazing impact (d-1) = ',f9.4)") graze_min
write(20,fmt="('Amplitude of seasonal grazing impact (d-1) = ',f9.4)") Graze_max
write(20,fmt="('Year day on which maximum grazing impact is reached = ',i9)") int(graze_daymax)
write(20,fmt="('Biomass threshold for grazing (mg Chl m-3) = ',f9.4)") graze_thresh
write(20,fmt="('Proportion of grazed organic N recycled to ambient DIN (0.0-1.0) = ',f9.4)") Nrecycle
write(20,*)
write(20,fmt="('----End of initialisation data file----')")
close(20)
return
end subroutine save_work
!*****************************************************************************
!
Subroutine load_work()
USE VARIABLES
use physics
use turbulence
use variables_all
use biology
use file_names
implicit none
integer :: iflags,i_algae,inerror,npoints,igraze_daymax
real :: dummy
character(len=25) :: filter
character(len=6) :: met_time
character (len=8) :: ntxt
character(len=100) :: initialdata2