[BUG/ISSUE] Typo in julday_mod.F which may lead to failure when year of simulation spans across 2 centuries

[joeylamcy]

Hi all,

TL;DR

A typo in calculating Julian date leads to an incorrect simulation length and raises an error for simulation that spans across 2 centuries, e.g. from 2099 to 2100. Fortunately, the case from 1999 to 2000 is unaffected.

I discovered the typo in julday_mod.F upon checking with the original algorithm in the literature.

Description

My colleague encounters this error when he runs a simulation from 01-Dec-2099 to 01-Jan-2100:

GEOS-Chem ERROR: No diagnostic output will be created for collection: 
"Restart"!  Make sure that the length of the simulation as specified in 
"input.geos" (check the start and end times) is not shorter than the frequency 
setting in HISTORY.rc!  For example, if the frequency is 010000 (1 hour) but 
the simulation is set up to run for only 20 minutes, then this error will occur.
 -> at History_ReadCollectionData (in module History/history_mod.F90)

 -> ERROR occurred at (or near) line     67 of the HISTORY.rc file

Julian date algorithm

Upon investigation, I checked the book "Practical Astronomy With Your Calculator", Third Edition, by Peter Duffett-Smith (1992). Here is the relevant algorithm used in calculating Julian date:

Typo in calculating Julian date

The bug is in line 110 of GeosUtil/julday_mod.F:

geos-chem/GeosUtil/julday_mod.F

Lines 100 to 111 in 075d55e

	! Compute YEAR and MONTH1
	IF ( ( MM == 1 ) .OR. ( MM == 2 ) ) THEN
	YEAR1 = YYYY - 1
	MONTH1 = MM + 12
	ELSE
	YEAR1 = YYYY
	MONTH1 = MM
	ENDIF
	! Compute the "A" term.
	X1 = DBLE( YYYY ) / 100.0d0
	A = MINT( X1 )

YYYY (equivalent to y in the book) in Line 110 should instead be YEAR1 (equivalent to y' in the book), i.e.

      X1 = DBLE( YEAR1 ) / 100.0d0

This affects the calculation of SimLengthSec in History/history_mod.F90:

geos-chem/History/history_mod.F90

Lines 634 to 641 in 075d55e

	! Compute the Astronomical Julian Date corresponding to the yyyymmdd
	! and hhmmss values at the start and end of the simulation, which are
	! needed below. This can be done outside of the DO loop below.
	CALL Compute_Julian_Date( yyyymmdd, hhmmss, JulianDate )
	CALL Compute_Julian_Date( yyyymmdd_end, hhmmss_end, JulianDateEnd )
	! Compute the length of the simulation, in elapsed seconds
	SimLengthSec = NINT( ( JulianDateEnd - JulianDate ) * SECONDS_PER_DAY )

and, in our case, gives 2592000 seconds (=30 days) instead of 2678400 seconds (= 31 days). Meanwhile the FileWriteAlarm remains correct (31 days), so it gives the following error message:

geos-chem/History/history_mod.F90

Lines 1241 to 1258 in 075d55e

	!-----------------------------------------------------------------
	! ERROR CHECK: Make sure that the length of the simulation is
	! not shorter than the requested "File Write" interval. This
	! will prevent simulations without diagnostic output.
	!-----------------------------------------------------------------
	IF ( SimLengthSec < Container%FileWriteAlarm ) THEN
	! Construct error message
	ErrMsg = &
	'No diagnostic output will be created for collection: "' // &
	TRIM( CollectionName(C) ) // '"! Make sure that the ' // &
	'length of the simulation as specified in "input.geos" ' // &
	'(check the start and end times) is not shorter than ' // &
	'the frequency setting in HISTORY.rc! For example, if ' // &
	'the frequency is 010000 (1 hour) but the simulation ' // &
	'is set up to run for only 20 minutes, then this error ' // &
	'will occur.'

When will this bug lead to simulation failure?

I think this bug only affects simulation that spans across different centuries, e.g. from 1899 to 1900, 2099 to 2100, 2199 to 2200. But fortunately, the case 1999 to 2000 is unaffected because of the calculation of "B" term: B = 2 - A + INT( A/4 ) is exactly the same as if the bug was not there. Therefore, most simulation would not encounter this issue.

Relevant configurations

GEOS-Chem version 12.2.0
Excerpt from HISTORY.rc:

#==============================================================================
# %%%%% THE Restart COLLECTION %%%%%
#
# GEOS-Chem restart file fields
#
# Available for all simulations
#==============================================================================
  Restart.template:           '%y4%m2%d2_%h2%n2z.nc4',
  Restart.format:             'CFIO',
  Restart.frequency:          'End',
  Restart.duration:           'End',
  Restart.mode:               'instantaneous'
  Restart.fields:             'SpeciesRst_?ALL?               ', 'GIGCchem',
                              'Chem_H2O2AfterChem             ', 'GIGCchem',
                              'Chem_SO2AfterChem              ', 'GIGCchem',
                              'Chem_DryDepNitrogen            ', 'GIGCchem',
                              'Chem_WetDepNitrogen            ', 'GIGCchem',
                              'Chem_KPPHvalue                 ', 'GIGCchem',
                              'Met_DELPDRY                    ', 'GIGCchem',
                              'Met_PS1WET                     ', 'GIGCchem',
                              'Met_PS1DRY                     ', 'GIGCchem',
                              'Met_SPHU1                      ', 'GIGCchem',
                              'Met_TMPU1                      ', 'GIGCchem',
::

Excerpt from input.geos:

GEOS-CHEM UNIT TEST SIMULATION: merra2_2x25_tropchem
------------------------+------------------------------------------------------
%%% SIMULATION MENU %%% :
Start YYYYMMDD, hhmmss  : 20991201 000000
End   YYYYMMDD, hhmmss  : 21000101 000000
Run directory           : ./
Root data directory     : /project/TGABI/GEOS-Chem/ExtData
Global offsets I0, J0   : 0 0
------------------------+------------------------------------------------------

Best regards,
Joey
02 Sep 2019

[yantosca]

Thanks Joey for bringing this to our attention. I have pushed commit cf35e9b, which contains your recommended fix, to the dev/12.6.0 branch.

This will go live into GEOS-Chem 12.6.0, as GEOS-Chem 12.5.0 is already in the benchmarking process.

This commit has passed a geosfp_4x5_benchmark difference test.