diff --git a/CHANGELOG.md b/CHANGELOG.md index 2b7fde4838b..35000718a4c 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -9,14 +9,22 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 ### Added -- Added subroutine `MAPL_SunGetLocalSolarHourAngle()` in `base/MAPL_sun_uc.F90`. This provides a -  convenient local solar hour angle diagnostic which will be used to detect local -  solar noon via the EXAMPLE OF USE in the subroutine header. See DESCRIPTION in code for more -  details. Provides the TRUE local solar hour angle (i.e., with equation of time included), but -  can also provide the MEAN value (without EOT) via FORCE_MLSHA=.TRUE. optional argument. +- Added subroutine `MAPL_SunGetLocalSolarHourAngle()` in `base/MAPL_sun_uc.F90`. This + provides a convenient local solar hour angle diagnostic which will be used to detect local + solar noon via the `EXAMPLE OF USE` in the subroutine header. See `DESCRIPTION` in code + for more details. Provides the TRUE local solar hour angle (i.e., with equation of time + included), but can also provide the MEAN value (without EOT) via `FORCE_MLSHA=.TRUE.` + optional argument. ### Changed +- Changed call to `MAPL_SunOrbitCreate()` inside `MAPL_Generic.F90` to call to new function + `MAPL_SunOrbitCreateFromConfig()`, the latter which get the orbital parameters from the MAPL + state's Config. In this way no default orbital parameter values need appear in `MAPL_Generic.F90`. + Rather, these default values are encapsulated where they belong in `Sun_Mod` in `base/MAPL_sun_uc.F90` + and are now explicitly named and commented on at the head of the module. This is a structural + zero-diff change. + ### Fixed - Added the correct values to halo corner of LatLon grid diff --git a/base/MAPL_sun_uc.F90 b/base/MAPL_sun_uc.F90 index cd4b4306516..3ebd68166a1 100644 --- a/base/MAPL_sun_uc.F90 +++ b/base/MAPL_sun_uc.F90 @@ -34,6 +34,7 @@ module MAPL_SunMod ! !PUBLIC MEMBER FUNCTIONS: public MAPL_SunOrbitCreate + public MAPL_SunOrbitCreateFromConfig public MAPL_SunOrbitCreated public MAPL_SunOrbitDestroy public MAPL_SunOrbitQuery @@ -56,6 +57,53 @@ module MAPL_SunMod integer, public, parameter :: MAPL_SunDailyMean = 5 integer, public, parameter :: MAPL_SunAnnualMean = 6 +! Default solar orbital system parameters (private). +! Dont change these unless you know what you are doing. +! They are appropriate for the current modern epoch circa 2000. +! ------------------------------------------------------------- + + ! Parameters of old orbital system (tabularized intercalation cycle) + ! ------------------------------------------------------------------ + real, parameter :: DEFAULT_ORBIT_ECCENTRICITY = 0.0167 + real, parameter :: DEFAULT_ORBIT_OBLIQUITY = 23.45 ! degrees + real, parameter :: DEFAULT_ORBIT_PERIHELION = 102.0 ! degrees + integer, parameter :: DEFAULT_ORBIT_EQUINOX = 80 ! days + + ! Parameters of new orbital system (analytic two-body), which allows some + ! time-varying behavior, namely, linear variation in LAMBDAP, ECC, and OBQ. + ! ------------------------------------------------------------------------- + + ! Fixed anomalistic year length in mean solar days + real, parameter :: DEFAULT_ORB2B_YEARLEN = 365.2596 + + ! Reference date and time for orbital parameters + ! (defaults to J2000 = 01Jan2000 12:00:00 TT = 11:58:56 UTC) + integer, parameter :: DEFAULT_ORB2B_REF_YYYYMMDD = 20000101 + integer, parameter :: DEFAULT_ORB2B_REF_HHMMSS = 115856 + + ! Orbital eccentricity at reference date + real, parameter :: DEFAULT_ORB2B_ECC_REF = 0.016710 + ! Rate of change of orbital eccentricity per Julian century + real, parameter :: DEFAULT_ORB2B_ECC_RATE = -4.2e-5 + + ! Earth's obliquity (axial tilt) at reference date [degrees] + real, parameter :: DEFAULT_ORB2B_OBQ_REF = 23.44 + ! Rate of change of obliquity [degrees per Julian century] + real, parameter :: DEFAULT_ORB2B_OBQ_RATE = -1.3e-2 + + ! Longitude of perihelion at reference date [degrees] + ! (from March equinox to perihelion in direction of earth's motion) + real, parameter :: DEFAULT_ORB2B_LAMBDAP_REF = 282.947 + ! Rate of change of LAMBDAP [degrees per Julian century] + ! (Combines both equatorial and ecliptic precession) + real, parameter :: DEFAULT_ORB2B_LAMBDAP_RATE = 1.7195 + + ! March Equinox date and time + ! (defaults to March 20, 2000 at 07:35:00 UTC) + integer, parameter :: DEFAULT_ORB2B_EQUINOX_YYYYMMDD = 20000320 + integer, parameter :: DEFAULT_ORB2B_EQUINOX_HHMMSS = 73500 + +! ------------------------------------------------------------- interface MAPL_SunGetInsolation module procedure SOLAR_1D @@ -759,6 +807,160 @@ end function MAPL_SunOrbitCreate !========================================================================== +!BOPI + +! !IROUTINE: MAPL_SunOrbitCreateFromConfig + +! !DESCRIPTION: + +! Like MAPL_SunOrbitCreate() but gets orbital parameters from Config CF. + +! !INTERFACE: + + function MAPL_SunOrbitCreateFromConfig ( & + CF, CLOCK, FIX_SUN, RC) result (ORBIT) + +! !ARGUMENTS: + + type (ESMF_Config), intent(INOUT) :: CF + type (ESMF_Clock), intent(IN ) :: CLOCK + logical, intent(IN ) :: FIX_SUN + integer, optional, intent(OUT ) :: RC + + type (MAPL_SunOrbit) :: ORBIT + +!EOPI + + character(len=ESMF_MAXSTR), parameter :: IAm = "SunOrbitCreateFromConfig" + integer :: STATUS + + real :: ECC, OB, PER + integer :: EQNX + + logical :: EOT, ORBIT_ANAL2B + integer :: ORB2B_REF_YYYYMMDD, ORB2B_REF_HHMMSS, & + ORB2B_EQUINOX_YYYYMMDD, ORB2B_EQUINOX_HHMMSS + real :: ORB2B_YEARLEN, & + ORB2B_ECC_REF, ORB2B_ECC_RATE, & + ORB2B_OBQ_REF, ORB2B_OBQ_RATE, & + ORB2B_LAMBDAP_REF, ORB2B_LAMBDAP_RATE + + ! pmn: There is one orbit is per STATE, so, for example, the MAPL states of the + ! solar and land gridded components can potentially have independent solar orbits. + ! Usually these "independent orbits" will be IDENTICAL because the configuration + ! resources such as "ECCENTRICITY:" or "EOT:" will not be qualified by the name + ! of the gridded component. But for example, if the resource file specifies + ! "EOT: .FALSE." + ! but + ! "SOLAR_EOT: .TRUE." + ! then only SOLAR will have an EOT correction. The same goes for the new orbital + ! system choice ORBIT_ANAL2B. + ! A state's orbit is actually created in this routine by requesting the ORBIT + ! object. If its not already created then it will be made below. GridComps that + ! don't needed an orbit and dont request one will not have one. + + ! Parameters of standard orbital system (tabularized intercalation cycle) + ! ----------------------------------------------------------------------- + call ESMF_ConfigGetAttribute (CF, & + ECC, label="ECCENTRICITY:", & + default=DEFAULT_ORBIT_ECCENTRICITY, _RC) + + call ESMF_ConfigGetAttribute (CF, & + OB, label="OBLIQUITY:", & + default=DEFAULT_ORBIT_OBLIQUITY, _RC) + + call ESMF_ConfigGetAttribute (CF, & + PER, label="PERIHELION:", & + default=DEFAULT_ORBIT_PERIHELION, _RC) + + call ESMF_ConfigGetAttribute (CF, & + EQNX, label="EQUINOX:", & + default=DEFAULT_ORBIT_EQUINOX, _RC) + + ! Apply Equation of Time correction? + ! ---------------------------------- + call ESMF_ConfigGetAttribute (CF, & + EOT, label="EOT:", & + default=.FALSE., _RC) + + ! New orbital system (analytic two-body) allows some time-varying + ! behavior, namely, linear variation in LAMBDAP, ECC, and OBQ. + ! --------------------------------------------------------------- + + call ESMF_ConfigGetAttribute (CF, & + ORBIT_ANAL2B, label="ORBIT_ANAL2B:", & + default=.FALSE., _RC) + + ! Fixed anomalistic year length in mean solar days + call ESMF_ConfigGetAttribute (CF, & + ORB2B_YEARLEN, label="ORB2B_YEARLEN:", & + default=DEFAULT_ORB2B_YEARLEN, _RC) + + ! Reference date and time for orbital parameters + call ESMF_ConfigGetAttribute (CF, & + ORB2B_REF_YYYYMMDD, label="ORB2B_REF_YYYYMMDD:", & + default=DEFAULT_ORB2B_REF_YYYYMMDD, _RC) + call ESMF_ConfigGetAttribute (CF, & + ORB2B_REF_HHMMSS, label="ORB2B_REF_HHMMSS:", & + default=DEFAULT_ORB2B_REF_HHMMSS, _RC) + + ! Orbital eccentricity at reference date + call ESMF_ConfigGetAttribute (CF, & + ORB2B_ECC_REF, label="ORB2B_ECC_REF:", & + default=DEFAULT_ORB2B_ECC_REF, _RC) + + ! Rate of change of orbital eccentricity per Julian century + call ESMF_ConfigGetAttribute (CF, & + ORB2B_ECC_RATE, label="ORB2B_ECC_RATE:", & + default=DEFAULT_ORB2B_ECC_RATE, _RC) + + ! Earth's obliquity (axial tilt) at reference date [degrees] + call ESMF_ConfigGetAttribute (CF, & + ORB2B_OBQ_REF, label="ORB2B_OBQ_REF:", & + default=DEFAULT_ORB2B_OBQ_REF, _RC) + + ! Rate of change of obliquity [degrees per Julian century] + call ESMF_ConfigGetAttribute (CF, & + ORB2B_OBQ_RATE, label="ORB2B_OBQ_RATE:", & + default=DEFAULT_ORB2B_OBQ_RATE, _RC) + + ! Longitude of perihelion at reference date [degrees] + ! (from March equinox to perihelion in direction of earth's motion) + call ESMF_ConfigGetAttribute (CF, & + ORB2B_LAMBDAP_REF, label="ORB2B_LAMBDAP_REF:", & + default=DEFAULT_ORB2B_LAMBDAP_REF, _RC) + + ! Rate of change of LAMBDAP [degrees per Julian century] + ! (Combines both equatorial and ecliptic precession) + call ESMF_ConfigGetAttribute (CF, & + ORB2B_LAMBDAP_RATE, label="ORB2B_LAMBDAP_RATE:", & + default=DEFAULT_ORB2B_LAMBDAP_RATE, _RC) + + ! March Equinox date and time + call ESMF_ConfigGetAttribute (CF, & + ORB2B_EQUINOX_YYYYMMDD, label="ORB2B_EQUINOX_YYYYMMDD:", & + default=DEFAULT_ORB2B_EQUINOX_YYYYMMDD, _RC) + call ESMF_ConfigGetAttribute (CF, & + ORB2B_EQUINOX_HHMMSS, label="ORB2B_EQUINOX_HHMMSS:", & + default=DEFAULT_ORB2B_EQUINOX_HHMMSS, _RC) + + ! create the orbit object + ORBIT = MAPL_SunOrbitCreate ( & + CLOCK, ECC, OB, PER, EQNX, & + EOT, ORBIT_ANAL2B, ORB2B_YEARLEN, & + ORB2B_REF_YYYYMMDD, ORB2B_REF_HHMMSS, & + ORB2B_ECC_REF, ORB2B_ECC_RATE, & + ORB2B_OBQ_REF, ORB2B_OBQ_RATE, & + ORB2B_LAMBDAP_REF, ORB2B_LAMBDAP_RATE, & + ORB2B_EQUINOX_YYYYMMDD, ORB2B_EQUINOX_HHMMSS, & + FIX_SUN=FIX_SUN,_RC) + + _RETURN(ESMF_SUCCESS) + + end function MAPL_SunOrbitCreateFromConfig + +!========================================================================== + !BOP ! !IROUTINE: MAPL_SunOrbitDestroy diff --git a/generic/MAPL_Generic.F90 b/generic/MAPL_Generic.F90 index ed77a88b2b3..d29db2dcce8 100644 --- a/generic/MAPL_Generic.F90 +++ b/generic/MAPL_Generic.F90 @@ -4031,20 +4031,9 @@ subroutine MAPL_GenericStateGet (STATE, IM, JM, LM, VERTDIM, & character(len=ESMF_MAXSTR), parameter :: IAm = "MAPL_GenericStateGet" integer :: status - real :: ECC - real :: OB - real :: PER - integer :: EQNX logical :: FIX_SUN character(len=ESMF_MAXSTR) :: gname - logical :: EOT, ORBIT_ANAL2B - integer :: ORB2B_REF_YYYYMMDD, ORB2B_REF_HHMMSS, & - ORB2B_EQUINOX_YYYYMMDD, ORB2B_EQUINOX_HHMMSS - real :: ORB2B_YEARLEN, & - ORB2B_ECC_REF, ORB2B_ECC_RATE, & - ORB2B_OBQ_REF, ORB2B_OBQ_RATE, & - ORB2B_LAMBDAP_REF, ORB2B_LAMBDAP_RATE type(MaplGrid), pointer :: temp_grid if(present(IM)) then @@ -4088,142 +4077,15 @@ subroutine MAPL_GenericStateGet (STATE, IM, JM, LM, VERTDIM, & CF=STATE%CF endif - ! pmn: There is one orbit is per STATE, so, for example, the MAPL states of the - ! solar and land gridded components can potentially have independent solar orbits. - ! Usually these "independent orbits" will be IDENTICAL because the configuration - ! resources such as "ECCENTRICITY:" or "EOT:" will not be qualified by the name - ! of the gridded component. But for example, if the resource file specifies - ! "EOT: .FALSE." - ! but - ! "SOLAR_EOT: .TRUE." - ! then only SOLAR will have an EOT correction. The same goes for the new orbital - ! system choice ORBIT_ANAL2B. - ! A state's orbit is actually created in this routine by requesting the ORBIT - ! object. If its not already created then it will be made below. GridComps that - ! don't needed an orbit and dont request one will not have one. - if(present(ORBIT)) then if(.not.MAPL_SunOrbitCreated(STATE%ORBIT)) then - call ESMF_GridGet(STATE%GRID%ESMFGRID,name=gname,rc=status) - _VERIFY(status) - if (index(gname,"DP")>0) then - FIX_SUN=.true. - else - FIX_SUN=.false. - end if - - ! Fixed parameters of standard orbital system (tabularized intercalation cycle) - ! ----------------------------------------------------------------------------- - - call MAPL_GetResource(STATE, ECC, Label="ECCENTRICITY:", default=0.0167, & - RC=status) - _VERIFY(status) - - call MAPL_GetResource(STATE, OB, Label="OBLIQUITY:", default=23.45, & - RC=status) - _VERIFY(status) - - call MAPL_GetResource(STATE, PER, Label="PERIHELION:", default=102.0, & - RC=status) - _VERIFY(status) - - call MAPL_GetResource(STATE, EQNX, Label="EQUINOX:", default=80, & - RC=status) - _VERIFY(status) - - ! Apply Equation of Time correction? - ! ---------------------------------- - call MAPL_GetResource(STATE, EOT, Label="EOT:", default=.FALSE., & - RC=status) - _VERIFY(status) - - ! New orbital system (analytic two-body) allows some time-varying - ! behavior, namely, linear variation in LAMBDAP, ECC, and OBQ. - ! --------------------------------------------------------------- - - call MAPL_GetResource(STATE, & - ORBIT_ANAL2B, Label="ORBIT_ANAL2B:", default=.FALSE., & - RC=status) - _VERIFY(status) - - ! Fixed anomalistic year length in mean solar days - call MAPL_GetResource(STATE, & - ORB2B_YEARLEN, Label="ORB2B_YEARLEN:", default=365.2596, & - RC=status) - _VERIFY(status) - - ! Reference date and time for orbital parameters - ! (defaults to J2000 = 01Jan2000 12:00:00 TT = 11:58:56 UTC) - call MAPL_GetResource(STATE, & - ORB2B_REF_YYYYMMDD, Label="ORB2B_REF_YYYYMMDD:", default=20000101, & - RC=status) - _VERIFY(status) - call MAPL_GetResource(STATE, & - ORB2B_REF_HHMMSS, Label="ORB2B_REF_HHMMSS:", default=115856, & - RC=status) - _VERIFY(status) - - ! Orbital eccentricity at reference date - call MAPL_GetResource(STATE, & - ORB2B_ECC_REF, Label="ORB2B_ECC_REF:", default=0.016710, & - RC=status) - _VERIFY(status) - - ! Rate of change of orbital eccentricity per Julian century - call MAPL_GetResource(STATE, & - ORB2B_ECC_RATE, Label="ORB2B_ECC_RATE:", default=-4.2e-5, & - RC=status) - _VERIFY(status) - - ! Earth's obliquity (axial tilt) at reference date [degrees] - call MAPL_GetResource(STATE, & - ORB2B_OBQ_REF, Label="ORB2B_OBQ_REF:", default=23.44, & - RC=status) - _VERIFY(status) - - ! Rate of change of obliquity [degrees per Julian century] - call MAPL_GetResource(STATE, & - ORB2B_OBQ_RATE, Label="ORB2B_OBQ_RATE:", default=-1.3e-2, & - RC=status) - _VERIFY(status) - - ! Longitude of perihelion at reference date [degrees] - ! (from March equinox to perihelion in direction of earth's motion) - call MAPL_GetResource(STATE, & - ORB2B_LAMBDAP_REF, Label="ORB2B_LAMBDAP_REF:", default=282.947, & - RC=status) - _VERIFY(status) - - ! Rate of change of LAMBDAP [degrees per Julian century] - ! (Combines both equatorial and ecliptic precession) - call MAPL_GetResource(STATE, & - ORB2B_LAMBDAP_RATE, Label="ORB2B_LAMBDAP_RATE:", default=1.7195, & - RC=status) - _VERIFY(status) - - ! March Equinox date and time - ! (defaults to March 20, 2000 at 07:35:00 UTC) - call MAPL_GetResource(STATE, & - ORB2B_EQUINOX_YYYYMMDD, Label="ORB2B_EQUINOX_YYYYMMDD:", default=20000320, & - RC=status) - _VERIFY(status) - call MAPL_GetResource(STATE, & - ORB2B_EQUINOX_HHMMSS, Label="ORB2B_EQUINOX_HHMMSS:", default=073500, & - RC=status) - _VERIFY(status) + call ESMF_GridGet(STATE%GRID%ESMFGRID,name=gname,_RC) + FIX_SUN = (index(gname,"DP")>0) ! create the orbit object - STATE%ORBIT = MAPL_SunOrbitCreate(STATE%CLOCK, ECC, OB, PER, EQNX, & - EOT, ORBIT_ANAL2B, ORB2B_YEARLEN, & - ORB2B_REF_YYYYMMDD, ORB2B_REF_HHMMSS, & - ORB2B_ECC_REF, ORB2B_ECC_RATE, & - ORB2B_OBQ_REF, ORB2B_OBQ_RATE, & - ORB2B_LAMBDAP_REF, ORB2B_LAMBDAP_RATE, & - ORB2B_EQUINOX_YYYYMMDD, ORB2B_EQUINOX_HHMMSS, & - FIX_SUN=FIX_SUN,RC=status) - _VERIFY(status) + STATE%ORBIT = MAPL_SunOrbitCreateFromConfig (STATE%CF, STATE%CLOCK, FIX_SUN, _RC) end if ORBIT=STATE%ORBIT