Add ISSM HISTORY, Restarts, and BCs

CHANGELOG.md

@@ -13,6 +13,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 
 - Added support for river routing.
 - Added optional NetCDF4 output mode for ObsFcstAna, including NetCDF metadata and runtime context. Changed namelist variable "out_ObsFcstAna" from logical to integer.
+- Added support for running ISSM (Ice-Sheet and Sea-level System Model). This includes a new collection in GEOSldas_HIST.rc (tavg24_1d_issm_Nt) as well as additions to the landice gridded collection (tavg24_2d_glac_Nx), and small additions to lenkf_j_template.py and ldas.py for handling ISSM boundary conditions and restarts.  
 
 ### Changed
 

GEOSldas_App/GEOSldas_HIST.rc

@@ -24,6 +24,7 @@ COLLECTIONS:
 #                'const_2d_lnd_Nx'
 #                'tavg24_2d_glc_Nx'
 #                'tavg24_1d_glc_Nt'
+#                'tavg24_1d_issm_Nt'
 #                'tavg24_1d_route'
            ::
 
@@ -531,6 +532,10 @@ EASEv2_M36.LM:           1
                                 'ALBNF'       , 'LANDICE' ,  'ALBNF_GL'           ,
                                 'ASNOW_GL'    , 'LANDICE' ,
                                 'DELTS'       , 'LANDICE' ,
+                                'ICESMB'      , 'LANDICE' ,
+>>>HIST_ISSM<<<                 'ICESURF'    ,  'LANDICE' ,
+>>>HIST_ISSM<<<                 'ICETHICK'   ,  'LANDICE' ,
+>>>HIST_ISSM<<<                 'ICEVEL'     ,  'LANDICE' ,
                                 'DNICFLX'     , 'LANDICE' ,
                                 'EVAPOUT'     , 'LANDICE' ,
                                 'QH'          , 'LANDICE' ,
@@ -597,4 +602,18 @@ EASEv2_M36.LM:           1
                                 'QINFLOW'     , 'ROUTE' ,
                                 'QRUNOFF_CAT' , 'ROUTE' ,
                                  ::                                 
+
+  tavg24_1d_issm_Nt.descr:         'Tile-space,Daily,Time-Averaged,Single-level,ISSM Diagnostics',
+  tavg24_1d_issm_Nt.nbits:         12,
+  tavg24_1d_issm_Nt.template:      '%y4%m2%d2_%h2%n2z.nc4' ,
+  tavg24_1d_issm_Nt.mode:          'time-averaged' ,
+  tavg24_1d_issm_Nt.format:        'CFIO',
+  tavg24_1d_issm_Nt.frequency:      240000 ,
+  tavg24_1d_issm_Nt.ref_time:       000000 ,
+  tavg24_1d_issm_Nt.fields:      'ICESURF'    , 'ISSM' ,
+                                 'ICETHICK'   , 'ISSM' ,
+                                 'ICEVX'      , 'ISSM' ,
+                                 'ICEVY'      , 'ISSM' ,
+                                 ::
+
 # ========================== EOF ==============================================================

GEOSldas_App/ldas.py

@@ -101,10 +101,12 @@ def __init__(self, cmdLineArgs):
         self.bcs_dir_land       = ''
         self.bcs_dir_geom       = ''
         self.bcs_dir_landshared = ''
+        self.bcs_dir_landiceshared = ''
         self.tile_types         = ''
         self.with_land          = False
         self.with_landice       = False
         self.run_route          = 0
+        self.with_issm          = False
         self.adas_expdir        = ''
 
         # assert necessary optional arguments in command line if exeinp file does not exsit
@@ -234,6 +236,10 @@ def __init__(self, cmdLineArgs):
           assert int(self.ExeInputs['LSM_CHOICE']) <= 2, "\nLSM_CHOICE=3 (Catchment-CN4.5) is no longer supported. Please set LSM_CHOICE to 1 (Catchment) or 2 (Catchment-CN4.0)"
         if  "20" in self.tile_types :
           self.with_landice = True
+
+        if self.with_landice == True:
+            if int(self.ExeInputs.get('DO_ISSM'))==1:
+                self.with_issm = True
 
         self.run_route = int(self.ExeInputs.get('RUN_ROUTE',0))
 
@@ -261,6 +267,8 @@ def __init__(self, cmdLineArgs):
         self.bcs_dir_land       = self.ExeInputs['BCS_PATH']+ '/land/'       + self.ExeInputs['BCS_RESOLUTION']+'/'
         self.bcs_dir_geom       = self.ExeInputs['BCS_PATH']+ '/geometry/'   + self.ExeInputs['BCS_RESOLUTION']+'/'
         self.bcs_dir_landshared = self.ExeInputs['BCS_PATH']+ '/land/shared/'
+        self.bcs_dir_landiceshared = self.ExeInputs['BCS_PATH']+ '/landice/shared/'
+
 
         # make sure MET_PATH and RESTART_PATH have trailing '/'
         if self.ExeInputs['MET_PATH'][-1] != '/':
@@ -448,6 +456,11 @@ def __init__(self, cmdLineArgs):
               routeRstFile=self.in_rstdir+'/'+tmpFile
               assert os.path.isfile(routeRstFile), 'route_internal_rst file [%s] does not exist!' %(routeRstFile)
 
+           if self.with_issm:
+              tmpFile=self.ExeInputs['RESTART_ID']+'.issm_internal_rst.'+y4m2d2_h2m2
+              issmRstFile=self.in_rstdir+'/'+tmpFile
+              assert os.path.isfile(issmRstFile), 'issm_internal_rst file [%s] does not exist!' %(issmRstFile)
+
         # DEAL WITH mwRTM input from exec
         self.assim = True if self.ExeInputs.get('LAND_ASSIM', 'NO').upper() == 'YES' and self.with_land else False
         # verify mwrtm file
@@ -767,6 +780,9 @@ def createLnRstBc(self) :
         # if running routing, make sure domain is global
         if self.run_route>0 and self.isZoomIn:
             exit( "Must have global domain to run routing model, RUN_ROUTE=" + self.run_route )
+        # if running ISSM, make sure domain is global
+        if self.with_issm and self.isZoomIn:
+            exit( "Must have global domain to run ISSM (DO_ISSM: 1)")    
 
         # update tile domain
         if self.isZoomIn:
@@ -830,7 +846,13 @@ def createLnRstBc(self) :
            if ("catchcn" in self.catch):
               os.symlink(self.bcs_dir_landshared + 'CO2_MonthlyMean_DiurnalCycle.nc4', \
                           self.inpdir+'/CO2_MonthlyMean_DiurnalCycle.nc4')
-
+
+        if self.with_issm:
+           for pattern in ('*.bin', '*.toolkits', '*.nc*'):
+              for issmbc in glob.glob(self.bcs_dir_landiceshared + pattern):
+                 myISSMBC = self.inpdir+'/'+os.path.basename(issmbc)
+                 os.symlink(issmbc, myISSMBC)        
+
         # create and link restart
         print ("Creating and linking restart...")
         _start = self.begDates[0]
@@ -910,12 +932,15 @@ def createLnRstBc(self) :
         catchRstFile0  = ''
         vegdynRstFile0 = ''
         landiceRstFile0 = ''
+        issmRstFile0   = ''
+
         for iens in range(self.nens) :
             ensdir   = self.ensdirs[iens]
             ensid    = self.ensids[iens]
-            myCatchRst   = myRstDir+'/'+ self.catch +ensid +'_internal_rst'
-            myLandiceRst = myRstDir+'/'+ 'landice'  +ensid +'_internal_rst'
-            myVegRst     = myRstDir+'/'+ 'vegdyn'   +ensid +'_internal_rst'
+            myCatchRst   = myRstDir+'/'+self.catch +ensid +'_internal_rst'
+            myLandiceRst = myRstDir+'/'+ 'landice' +ensid +'_internal_rst'
+            myIssmRst = myRstDir+'/'+ 'issm' +ensid +'_internal_rst'
+            myVegRst     = myRstDir+'/'+ 'vegdyn'+ensid +'_internal_rst'
             myPertRst    = myRstDir+'/'+ 'landpert' +ensid +'_internal_rst'
             myRouteRst   = myRstDir+'/'+ 'route'    +ensid +'_internal_rst'
 
@@ -972,11 +997,22 @@ def createLnRstBc(self) :
                    vegdynRstFile = vegdynRstFile0
 
             landiceRstFile = ''
+            issmRstFile   = ''
             if self.with_landice :
                if RESTART_str in ['1', '3'] :
                   landiceRstFile = rstpath+ensdir +'/'+ y4m2+'/'+self.ExeInputs['RESTART_ID']+'.'+'landice_internal_rst.'+y4m2d2_h2m2
+                  if self.with_issm:
+                     issmRstFile = rstpath+ensdir +'/'+ y4m2+'/'+self.ExeInputs['RESTART_ID']+'.'+'issm_internal_rst.'+y4m2d2_h2m2
+
+
+
                if RESTART_str in ['2', 'M']:
                   landiceRstFile = glob.glob(self.exphome+'/'+exp_id+'/mk_restarts/*'+'landice_internal_rst.'+YYYYMMDD+'*')[0]
+                  if self.with_issm and RESTART_str in ['2']:
+                     exit( "Restart '2' not supported when running ISSM (DO_ISSM: 1)") 
+                  if self.with_issm and RESTART_str in ['M']:
+                     print("issm_internal_rst will be bootstrapped for Restart 'M' ")
+
 
                if os.path.isfile(landiceRstFile) :
                   landiceLocal = self.rstdir+ensdir +'/'+ y4m2+'/'+self.ExeInputs['EXP_ID']+'.landice_internal_rst.'+y4m2d2_h2m2
@@ -995,6 +1031,17 @@ def createLnRstBc(self) :
                else :
                    landiceRstFile = landiceRstFile0
 
+               if self.with_issm:
+                  if os.path.isfile(issmRstFile) :
+                     issmLocal = self.rstdir+ensdir +'/'+ y4m2+'/'+self.ExeInputs['EXP_ID']+'.issm_internal_rst.'+y4m2d2_h2m2
+                     shutil.copy(issmRstFile,issmLocal)
+                     issmRstFile = issmLocal
+
+                     if '0000' in ensdir :
+                        issmRstFile0 = issmRstFile
+                  else :
+                     issmRstFile = issmRstFile0
+
             routeRstFile = ''
             if self.run_route > 0 :
                if RESTART_str in ['1', '2'] :
@@ -1037,6 +1084,11 @@ def createLnRstBc(self) :
                print("link route restart: " + myRouteRst)
                os.symlink(routeRstFile, myRouteRst)
 
+            if self.with_issm:
+               if RESTART_str in ['1', '3']:
+                  print("link issm restart: " + myIssmRst)
+                  os.symlink(issmRstFile, myIssmRst)
+
             if ( self.has_geos_pert and  self.perturb == 1 ):
                os.symlink(pertRstFile,    myPertRst)
 
@@ -1189,7 +1241,8 @@ def createRCFiles(self):
                         + str(self.ExeInputs['LSM_CHOICE'])         + ' ' \
                         + str(self.ExeInputs['AEROSOL_DEPOSITION']) + ' ' \
                         + str(self.ExeInputs['RUN_IRRIG'])          + ' ' \
-                        + str(self.nens)
+                        + str(self.nens)                            + ' ' \
+                        + str(self.ExeInputs['DO_ISSM'])   
                     print(cmd)
                     #os.system(cmd)
                     sp.call(shlex.split(cmd))
@@ -1298,6 +1351,9 @@ def createRCFiles(self):
                 if self.with_landice:
                   rstkey.append('LANDICE')
                   rstval.append('landice')
+                  if self.with_issm:
+                     rstkey.append('ISSM')
+                     rstval.append('issm')
 
                 if self.run_route > 0:
                   rstkey.append('ROUTE')
@@ -1336,6 +1392,12 @@ def createRCFiles(self):
                        valn = val + tmpl_ + '_internal_checkpoint'
                        ldasrcInp[keyn] = valn
 
+                    if self.with_issm:
+                       keyn = 'ISSM_INTERNAL_CHECKPOINT_FILE'
+                       valn = 'issm'+tmpl_+'_internal_checkpoint'
+                       ldasrcInp[keyn]= valn
+
+
                 # specify LANDPERT restart file
                 if (self.perturb == 1):
                     keyn = 'LANDPERT_INTERNAL_RESTART_FILE'

GEOSldas_App/lenkf_j_template.py

@@ -158,6 +158,12 @@
 /bin/cp -f  $HOMDIR/*.rc .
 /bin/cp -f  $HOMDIR/*.nml .
 
+# copy ISSM input files to scratch
+/bin/find $EXPDIR/input -maxdepth 1 -name "ISSM*" -exec /bin/cp {{}} . \\;
+
+# move ISSM auxiliary mesh file to rc_out
+if (-e $EXPDIR/input/ISSM_MESH.nc) /bin/mv $EXPDIR/input/ISSM_MESH.nc $EXPDIR/output/$EXPDOMAIN/rc_out
+
 set LSMCHOICE = `grep -n -m 1 "LSM_CHOICE" $HOMDIR/LDAS.rc | cut -d':' -f3`
 
 #######################################################################
@@ -721,7 +727,7 @@
        set THISDIR = $EXPDIR/output/$EXPDOMAIN/rs/$ENSDIR/Y${{eYEAR}}/M${{eMON}}/
        if (! -e $THISDIR            ) mkdir -p $THISDIR
 
-       set rstfs = (${{LANDMODEL}} 'landice' 'route')
+       set rstfs = (${{LANDMODEL}} 'landice' 'route' 'issm')
        foreach rstf ( $rstfs )
           if (-f ${{rstf}}${{ENSID}}_internal_checkpoint ) then
              set tmp_file = $EXPDIR/output/$EXPDOMAIN/rs/$ENSDIR/Y${{eYEAR}}/M${{eMON}}/${{EXPID}}.${{rstf}}_internal_rst.${{eYEAR}}${{eMON}}${{eDAY}}_${{eHour}}${{eMin}}
@@ -759,8 +765,9 @@
        set rstfiles2 = `ls landpert${{ENSID}}_internal_checkpoint.*`
        set rstfiles3 = `ls landassim_obspertrseed${{ENSID}}_checkpoint.*`
        set rstfiles4 = `ls landice${{ENSID}}_internal_checkpoint.*`
+       set rstfiles5 = `ls issm${{ENSID}}_internal_checkpoint.*`
 
-       foreach rfile ( $rstfiles1 $rstfiles4 ) 
+       foreach rfile ( $rstfiles1 $rstfiles4 $rstfiles5 ) 
           set ThisTime = `echo $rfile | rev | cut -d'.' -f2 | rev`
           set TY = `echo $ThisTime | cut -c1-4`
           set TM = `echo $ThisTime | cut -c5-6`

GEOSldas_App/process_hist.csh

@@ -17,6 +17,8 @@ setenv LSM_CHOICE         $4
 setenv AEROSOL_DEPOSITION $5
 setenv RUN_IRRIG          $6
 setenv NENS               $7
+setenv DO_ISSM            $8
+
 
 # -------------------------------------------------
 
@@ -67,6 +69,12 @@ else
    sed -i 's/>>>HIST_IRRIG<<</''/g'        $HISTRC
 endif
 
+if($DO_ISSM == 0) then
+   sed -i '/^>>>HIST_ISSM<<</d'            $HISTRC
+else
+   sed -i 's/>>>HIST_ISSM<<</''/g'         $HISTRC
+endif
+
 # for ensemble simulations, set 'GridComp' to ENSAVG 
 
 if($NENS > 1) then