Merge pull request #247 from derpycode/master.daylength

master.daylength -> master

genie-embm/src/fortran/embm.cmn

@@ -160,9 +160,17 @@ c adjustable freshwater forcing parameters
       real extra0,range0,extra1a,extra1b,extra1c
       common /embm_fwf_real/ extra0,range0,extra1a,extra1b,extra1c
 
+c CL (01/15/24) : seconds per solar day (necessary for GENIE)
+      real sodaylen
+      common /embm_sodaylen/sodaylen
+
+c CL (01/15/24) : seconds per sidereal day (necessary for GENIE)
+      real sidaylen
+      common /embm_sidaylen/sidaylen
+
 c AY (08/04/04) : days per year (necessary for GENIE)
       real yearlen
-      common /ocn_yearlen/yearlen
+      common /embm_yearlen/yearlen
 
 c v2 seasonal diagnostics
 c AP (03/08/06) : Addition of avg precipitated atm. humidity

genie-embm/src/fortran/initialise_embm.F

@@ -214,7 +214,7 @@ subroutine initialise_embm(alon1,alat1,alon2,alat2,
       NAMELIST /ini_embm_nml/xu_wstress,yu_wstress,xv_wstress
       NAMELIST /ini_embm_nml/yv_wstress,u_wspeed,v_wspeed
       NAMELIST /ini_embm_nml/npstp,iwstp,itstp,ianav,ans
-      NAMELIST /ini_embm_nml/yearlen,nyear,ndta,scf
+      NAMELIST /ini_embm_nml/sodaylen,sidaylen,yearlen,nyear,ndta,scf
       NAMELIST /ini_embm_nml/diffamp,diffwid,difflin
       NAMELIST /ini_embm_nml/betaz,betam
       NAMELIST /ini_embm_nml/radfor_scl_co2,radfor_pc_co2_rise
@@ -430,6 +430,16 @@ subroutine initialise_embm(alon1,alat1,alon2,alat2,
       if (debug_init) print*,'new or continuing run ?'
       if (debug_init) print*,ans
 
+c CL (01/15/24) : number of seconds per solar day
+c     sodaylen = 86400.0
+      if (debug_init) print*,'number of seconds per solar day'
+      if (debug_init) print*,sodaylen
+
+c CL (01/15/24) : number of seconds per sidereal day
+c     sidaylen = 86164.0
+      if (debug_init) print*,'number of seconds per sidereal day'
+      if (debug_init) print*,sidaylen
+
 c AY (05/05/04) : number of days per EMBM year (usually 365.25)
 c     yearlen = 365.25
       if (debug_init) print*,'number of days per EMBM year'
@@ -562,8 +572,8 @@ subroutine initialise_embm(alon1,alat1,alon2,alat2,
       if (debug_init) print*,'timesteps per year and A/O dt ratio'
       if (debug_init) print*,nyear,ndta
       if(nyear.gt.maxnyr)stop 'embm : nyear > maxnyr' 
-      tv = 86400.0*yearlen/(nyear*tsc)
-      ryear = 1.0/(yearlen*86400)
+      tv = sodaylen*yearlen/(nyear*tsc)
+      ryear = 1.0/(yearlen*sodaylen)
 
       dtatm = tv/ndta
       if (debug_init) print*, 'embm timestep (s) =',dtatm*tsc
@@ -577,7 +587,7 @@ subroutine initialise_embm(alon1,alat1,alon2,alat2,
 c        dzu(1,k) = 0
 c        dzu(2,k) = 0
       enddo
-      if (debug_init) print*,'dimensional ocean timestep',tv*tsc/86400
+      if (debug_init) print*,'dimensional ocean timestep',tv*tsc/sodaylen
       if (debug_init) print*,'dimensionless O/A timesteps',tv,dtatm
 
       rdtdim = 1.0/(tsc*dt(kmax))
@@ -1067,8 +1077,8 @@ subroutine initialise_embm(alon1,alat1,alon2,alat2,
 
 c     diffmod0 = 60e6
       diffmod0 = 0.
-      ppmin = 2./(yearlen*86400.)
-      ppmax = 4./(yearlen*86400.)
+      ppmin = 2./(yearlen*sodaylen)
+      ppmax = 4./(yearlen*sodaylen)
 
 c nre simpler diffusivity
 
@@ -1228,7 +1238,7 @@ subroutine initialise_embm(alon1,alat1,alon2,alat2,
 c     rsictscsf = ch_ice*u_tau_ice*rho0sea*cpo_ice
       rsictscsf = ch_ice*u_tau_ice*rho0*cpo_ice
       if (debug_init) print*,'rsictscsf = ',rsictscsf
-      rsictscsf = dsc*dz(kmax)*rho0*cpo_ice/(17.5*86400.0)
+      rsictscsf = dsc*dz(kmax)*rho0*cpo_ice/(17.5*sodaylen)
       if (debug_init) print*,'rsictscsf = ',rsictscsf
 c minimum average sea-ice thickness over a grid cell
       hmin = 0.01

genie-embm/src/fortran/inm_netcdf_embm.F

@@ -28,7 +28,7 @@ subroutine inm_netcdf_embm
       real    area
 
 c      timestep=24.0*60.0*60.0*yearlen/real(nyear)
-      timestep=24.0*60.0*60.0*yearlen/real(nyear*ndta)
+      timestep=sodaylen*yearlen/real(nyear*ndta)
 
 c AY (06/10/04) : append input filename to input directory name
 c      fnamein=trim(indir_name(1:lenin))//trim(filenetin)
@@ -66,7 +66,7 @@ subroutine inm_netcdf_embm
          day_rest=iday
          ioffset_rest=mod(ioffset_rest,nint(yearlen))
 
-      day_rest=day_rest+timestep/(24*60*60.)
+      day_rest=day_rest+timestep/sodaylen
 c     This bit so that we don't get too far out in our count....
 c     Anchor to a day if we start drifting.
 c     Means timestep can never be less than 1/1000 of a day!!!!

genie-embm/src/fortran/netcdf_embm.F

@@ -38,7 +38,7 @@ SUBROUTINE INI_NETCDF_EMBM(istep,imode)
 c
 c AY (17/03/04) : alon1, etc. now calculated in initialise_ocean.F
 c
-      day=istep*dtatm*tsc/86400.0
+      day=istep*dtatm*tsc/sodaylen
 c AY (17/03/04) : 365.25 days per GOLDSTEIN year (not 360.0)
 c AY (23/03/04) : small constant added to day calculation for round-off
 c                 reasons (i.e. 365.2499 vs. 365.25)

genie-goldstein/src/fortran/initialise_goldstein.F

@@ -204,7 +204,7 @@ subroutine initialise_goldstein(
       NAMELIST /ini_gold_nml/indir_name,outdir_name,rstdir_name
       NAMELIST /ini_gold_nml/igrid,world
       NAMELIST /ini_gold_nml/npstp,iwstp,itstp,ianav
-      NAMELIST /ini_gold_nml/conserv_per,ans,yearlen,nyear
+      NAMELIST /ini_gold_nml/conserv_per,ans,sodaylen,sidaylen,yearlen,nyear
       NAMELIST /ini_gold_nml/temp0,temp1,rel,scf,diff,adrag
       NAMELIST /ini_gold_nml/hosing,hosing_trend,nyears_hosing
       NAMELIST /ini_gold_nml/fwanomin,cmip_model,albocn
@@ -287,9 +287,9 @@ subroutine initialise_goldstein(
          close(56,iostat=ios)
          call check_iostat(ios,__LINE__,__FILE__)
       end if
-
+      
 c SG > syr re-defined here
-      syr = yearlen * 86400
+      syr = yearlen * sodaylen
       if (debug_init) print*, 'syr = ',syr
 c SG <
 
@@ -341,6 +341,16 @@ subroutine initialise_goldstein(
       if (debug_init) print*,'new or continuing run ?'
       if (debug_init) print*,ans
 
+c CL (01/15/24) : number of seconds per solar day
+c     sodaylen = 86400.0
+      if (debug_init) print*,'number of seconds per solar day'
+      if (debug_init) print*,sodaylen
+
+c CL (01/15/24) : number of seconds per sidereal day
+c     sidaylen = 86164.0
+      if (debug_init) print*,'number of seconds per sidereal day'
+      if (debug_init) print*,sidaylen
+
 c AY (05/05/04) : number of days per GOLDSTEIN year (usually 365.25)
 c     yearlen = 365.25
       if (debug_init) print*,'number of days per GOLDSTEIN year'
@@ -366,7 +376,13 @@ subroutine initialise_goldstein(
       rsc = 6.37e6
 c      dsc = 5e3
       dsc = par_dsc
-      fsc = 2*7.2921e-5
+c     CL (01/15/24) : used sidereal day length for Coriolis effect scaling factor
+c     AR (24/03/03) : check for sidereal and solar days being different, otherwise use original code for backwards compatabilty
+      if (abs(sodaylen-sidaylen).gt.0.001)then
+         fsc = 4*pi/sidaylen
+      else
+         fsc = 2*7.2921e-5
+      endif
       gsc = 9.81
       rh0sc = 1e3
       rhosc = rh0sc*fsc*usc*rsc/gsc/dsc
@@ -496,7 +512,7 @@ subroutine initialise_goldstein(
       if(nyear.gt.maxnyr)stop 'goldstein : nyear > maxnyr' 
       if (debug_init) print*,nyear
 
-      tv = 86400.0*yearlen/(nyear*tsc)
+      tv = sodaylen*yearlen/(nyear*tsc)
 
 c     dtatm = tv/ndta
 
@@ -507,7 +523,7 @@ subroutine initialise_goldstein(
          dzu(1,k) = 0
          dzu(2,k) = 0
       enddo
-      if (debug_init) print*,'dimensional ocean timestep',tv*tsc/86400
+      if (debug_init) print*,'dimensional ocean timestep',tv*tsc/sodaylen
 
 c set up grid
 c For variable (exponential) dz use ez0 > 0, else use ez0 < 0
@@ -1098,8 +1114,8 @@ subroutine initialise_goldstein(
 c drgf per gridpoint close to equator and in regions of
 c shallow water (k1>kmxdrg) ie land in the case kmxdrg=kmax
 c jeb = 1/2 width of equatorial region of maximum drag
-
-      adrag = 1.0/(adrag*86400*fsc)
+c CL (01/15/24) : used sidereal day length for inverse of frictional drag
+      adrag = 1.0/(adrag*sidaylen*fsc)
 c cross equator need * 4 if drag is constant ie if drgf=1
       drgf = 3.0
       kmxdrg = kmax/2
@@ -1363,7 +1379,7 @@ subroutine initialise_goldstein(
       if (debug_init) print*,'density variation scale',rhosc,' kg/m**3'
       if (debug_init) print*,
      & 'vertical velocity scale',usc*dsc/rsc,' m/s'
-      if (debug_init) print*,'time scale',tsc/86400/yearlen,' yrs'
+      if (debug_init) print*,'time scale',tsc/sodaylen/yearlen,' yrs'
       if (debug_init) print*,'overturning scale',dsc*usc*rsc*1e-6,' Sv'
       if (debug_init) print*,
      & 'vertical heat flux scale',dsc*usc*rh0sc*cpsc/rsc,' W/m**2'
@@ -1460,7 +1476,7 @@ subroutine initialise_goldstein(
       rhoice = 913.
 c useful constant proportional to inverse timscale for surface freezing
 c     rsictscsf = ch_ice*u_tau_ice*rho0*cpo_ice
-      rsictscsf = dsc*dz(kmax)*rho0*cpo_ice/(17.5*86400.0)
+      rsictscsf = dsc*dz(kmax)*rho0*cpo_ice/(17.5*sodaylen)
 c minimum average sea-ice thickness over a grid cell
       hmin = 0.01
       rhmin = 1.0/hmin

genie-goldstein/src/fortran/inm_netcdf.F

@@ -31,7 +31,7 @@ subroutine inm_netcdf(lrestart_genie)
 
       logical lrestart_genie
 
-      timestep=24.0*60.0*60.0*yearlen/real(nyear)
+      timestep=sodaylen*yearlen/real(nyear)
 
       fnamein=trim(filenetin)
 
@@ -77,7 +77,7 @@ subroutine inm_netcdf(lrestart_genie)
          day_rest=iday
          ioffset_rest=mod(ioffset_rest,nint(yearlen))
 
-      day_rest=day_rest+timestep/(24*60*60.)
+      day_rest=day_rest+timestep/sodaylen
 c     This bit so that we don't get too far out in our count....
 c     Anchor to a day if we start drifting.
 c     Means timestep can never be less than 1/1000 of a day!!!!

genie-goldstein/src/fortran/netcdf.F

@@ -43,7 +43,7 @@ SUBROUTINE INI_NETCDF_OCN(istep,imode)
 c
 c AY (17/03/04) : alon1, etc. now calculated in initialise_ocean.F
 c
-      day=(istep*dt(1)*tsc/86400.0)
+      day=(istep*dt(1)*tsc/sodaylen)
 c AY (17/03/04) : 365.25 days per GOLDSTEIN year (not 360.0)
 c AY (23/03/04) : small constant added to day calculation for round-off
 c                 reasons (i.e. 365.2499 vs. 365.25)

genie-goldstein/src/fortran/ocean.cmn

@@ -127,6 +127,14 @@ c AY (29/11/04) : ASURF grid cell area
       real asurf(maxj)
       common /ocn_asurf/asurf
 
+c CL (01/15/24) : seconds per solar day (necessary for GENIE)
+      real sodaylen
+      common /ocn_sodaylen/sodaylen
+
+c CL (01/15/24) : seconds per sidereal day (necessary for GENIE)
+      real sidaylen
+      common /ocn_sidaylen/sidaylen
+
 c AY (08/04/04) : days per year (necessary for GENIE)
       real yearlen
       common /ocn_yearlen/yearlen

genie-goldstein/src/fortran/surf_ocn_sic.F

@@ -214,9 +214,9 @@ subroutine surf_ocn_sic(istep,
             endif
          enddo
       enddo
-      diff_latent=diff_latent*3600.0*24.0*yearlen*rsc*rsc*4.0*pi/
+      diff_latent=diff_latent*sodaylen*yearlen*rsc*rsc*4.0*pi/
      :               real(nyear)
-      diff_sensible=diff_sensible*3600.0*24.0*yearlen*rsc*rsc*4.0*pi/
+      diff_sensible=diff_sensible*sodaylen*yearlen*rsc*rsc*4.0*pi/
      :               real(nyear)
       endif
 c-DJL
@@ -776,9 +776,9 @@ subroutine surf_ocn_sic(istep,
             endif
          enddo
       enddo
-      diff_netsolar=diff_netsolar*3600.0*24.0*yearlen*rsc*rsc*4.0*pi/
+      diff_netsolar=diff_netsolar*sodaylen*yearlen*rsc*rsc*4.0*pi/
      :               real(nyear)
-      diff_netlong=diff_netlong*3600.0*24.0*yearlen*rsc*rsc*4.0*pi/
+      diff_netlong=diff_netlong*sodaylen*yearlen*rsc*rsc*4.0*pi/
      :               real(nyear)
       test_energy_seaice=real(
      :              test_energy_seaice-diff_latent-diff_sensible-

genie-goldsteinseaice/src/fortran/gold_seaice.F

@@ -313,7 +313,7 @@ subroutine gold_seaice(istep,
 c THE ENERGY BIT (HERE AND IN SURF_OCN_SIC) IS WELL-DODGY AND
 c IS A COMPLETE FIX.  SEE DJL FOR MORE INFO.
        if (mod(istep,conserv_per).eq.0) then
-       tv = 86400.0*yearlen/(nyear*tsc)
+       tv = sodaylen*yearlen/(nyear*tsc)
        vsc = dphi*rsc*rsc
        tot_energy=0.0
        tot_water=0.0

genie-goldsteinseaice/src/fortran/initialise_seaice.F

@@ -115,7 +115,7 @@ subroutine initialise_seaice(
       NAMELIST /ini_sic_nml/indir_name,outdir_name,rstdir_name
       NAMELIST /ini_sic_nml/igrid,world
       NAMELIST /ini_sic_nml/npstp,iwstp,itstp,ianav,conserv_per
-      NAMELIST /ini_sic_nml/ans,yearlen,nyear,diffsic,lout
+      NAMELIST /ini_sic_nml/ans,sodaylen,sidaylen,yearlen,nyear,diffsic,lout
       NAMELIST /ini_sic_nml/netin,netout,ascout,filenetin
       NAMELIST /ini_sic_nml/dirnetout,lin
       NAMELIST /ini_sic_nml/dosc,impsic,debug_init,debug_end,debug_loop
@@ -217,6 +217,16 @@ subroutine initialise_seaice(
       if (debug_init) print*,'new or continuing run ?'
       if (debug_init) print*,ans
 
+c CL (01/15/24) : number of seconds per solar day
+c     sodaylen = 86400.0
+      if (debug_init) print*,'number of seconds per solar day'
+      if (debug_init) print*,sodaylen
+
+c CL (01/15/24) : number of seconds per sidereal day
+c     sidaylen = 86164.0
+      if (debug_init) print*,'number of seconds per sidereal day'
+      if (debug_init) print*,sidaylen
+
 c AY (05/05/04) : number of days per GOLDSTEIN sea-ice year (usually 365.25)
 c     yearlen = 365.25
       if (debug_init) print*,'number of days per GOLDSTEIN sea-ice year'
@@ -235,6 +245,13 @@ subroutine initialise_seaice(
       usc = 0.05
       rsc = 6.37e6
       dsc = 5e3
+c CL (01/15/24) : used sidereal day length for Coriolis effect scaling factor
+c     AR (24/03/03) : check for sidereal and solar days being different, otherwise use original code for backwards compatabilty
+      if (abs(sodaylen-sidaylen).gt.0.001)then
+         fsc = 4*pi/sidaylen
+      else
+         fsc = 2*7.2921e-5
+      endif
       fsc = 2*7.2921e-5
       gsc = 9.81
       rh0sc = 1e3
@@ -334,11 +351,11 @@ subroutine initialise_seaice(
 c v2 seasonality
       if (debug_init) print*,'timesteps per year'
       if (debug_init) print*,nyear
-      tv = 86400.0*yearlen/(nyear*tsc)
+      tv = sodaylen*yearlen/(nyear*tsc)
 
 c AY (13/01/04) : oops!  left this out
       dtsic = tv  
-      if (debug_init) print*,'sea-ice timestep in days',tv*tsc/86400
+      if (debug_init) print*,'sea-ice timestep in days',tv*tsc/sodaylen
 
       rdtdim = 1.0/(tsc*dtsic)
       if (debug_init) print*,'rdtdim = ',rdtdim
@@ -425,7 +442,7 @@ subroutine initialise_seaice(
 c read some scalings
 
 c AY (11/12/03) : year length inconsistent here - corrected to 365.25
-      ryear = 1.0/(yearlen*86400)
+      ryear = 1.0/(yearlen*sodaylen)
 
 c more constants
 

genie-goldsteinseaice/src/fortran/inm_netcdf_sic.F

@@ -28,7 +28,7 @@ subroutine inm_netcdf_sic
       integer i,j,l,icell
       real    tmp_val(4)
 
-      timestep=24.0*60.0*60.0*yearlen/real(nyear)
+      timestep=sodaylen*yearlen/real(nyear)
 
       fnamein=trim(filenetin)
 
@@ -68,7 +68,7 @@ subroutine inm_netcdf_sic
          day_rest=iday
          ioffset_rest=mod(ioffset_rest,nint(yearlen))
 
-      day_rest=day_rest+timestep/(24*60*60.)
+      day_rest=day_rest+timestep/sodaylen
 c     This bit so that we don't get too far out in our count....
 c     Anchor to a day if we start drifting.
 c     Means timestep can never be less than 1/1000 of a day!!!!

genie-goldsteinseaice/src/fortran/netcdf_sic.F

@@ -40,7 +40,7 @@ SUBROUTINE INI_NETCDF_SIC(istep,imode)
 c
 c AY (17/03/04) : alon1, etc. now calculated in initialise_ocean.F
 c
-      day=istep*dtsic*tsc/86400.0
+      day=istep*dtsic*tsc/sodaylen
 c AY (17/03/04) : 365.25 days per GOLDSTEIN year (not 360.0)
 c AY (23/03/04) : small constant added to day calculation for round-off
 c                 reasons (i.e. 365.2499 vs. 365.25)

genie-goldsteinseaice/src/fortran/seaice.cmn

@@ -58,6 +58,14 @@ c AY (23/09/04) : sea-ice albedo added (for restarts)
      3 ,varicedy,variceth                                                 ! yka101004
      4 ,maxalbice,scale_oheat
 
+c CL (01/15/24) : seconds per solar day (necessary for GENIE)
+      real sodaylen
+      common /sic_sodaylen/sodaylen
+
+c CL (01/15/24) : seconds per sidereal day (necessary for GENIE)
+      real sidaylen
+      common /sic_sidaylen/sidaylen
+
 c AY (08/04/04) : days per year (necessary for GENIE)
       real yearlen
       common /sic_yearlen/yearlen