<p><b>mhoffman@lanl.gov</b> 2013-03-11 11:37:48 -0600 (Mon, 11 Mar 2013)</p><p>BRANCH COMMIT - landice-ocean coupling<br>
<br>
The script now succesfully runs and restarts the land ice model. Added the framework to allow either the land ice or ocean model to be data (useful for testing, may also simplify running ISOMIP with the ocean model). I am also committing a land ice namelist file with comments explaining how certain settings need to be set for using the coupler.<br>
<br>
Next steps are a) getting the script to successfully run and restart the ocean model; and b) adding actual calculations made by the coupler (boundary layer physics and ocean surface pressure).<br>
</p><hr noshade><pre><font color="gray">Modified: branches/landice_ocean_coupler/landice-ocean-coupler.py
===================================================================
--- branches/landice_ocean_coupler/landice-ocean-coupler.py 2013-03-11 14:55:04 UTC (rev 2580)
+++ branches/landice_ocean_coupler/landice-ocean-coupler.py 2013-03-11 17:37:48 UTC (rev 2581)
@@ -25,19 +25,28 @@
# The script assumes you have setup a directory structure:
# ocean model is setup in ./ocean
# land ice model is setup in ./land_ice
+# and each has a subdirectory called 'output-files' where all the output is accumulated.
+
# =======================================
+
+# Set if you are using data or running a model
+runLICE = True
+runOCN = False
+
+# Setup needed executables
+runLICEcmd = 'mpirun -np 4 land_ice_model.exe'
+runOCNcmd = 'mpirun -np 4 ocean_model.exe'
+
# Setup needed files
+liceInput = 'land_ice_grid.nc'
+liceOutput = 'output.nc'
ocnInput = 'grid.nc'
-liceInput = 'land_ice_grid.nc'
+ocnOutput = 'output.nc'
-# Setup needed executables
-runLICE = 'mpirun -np 4 land_ice_model.exe'
-runOCN = 'mpirun -np 4 ocean_model.exe'
-
# Define number of coupling intervals (assumes you have set dt for each model accordingly, and ocean dt is evenly divisble in land ice dt)
# (We may want to do more sophisticated time-keeping of the coupling if there is a need.)
-numCoupleIntervals = 2
+numCoupleIntervals = 8
# Some parameters
rhoi = 900.0 # land ice model's ice density (kg m^-3)
@@ -45,59 +54,101 @@
# =======================================
+
+
# clean up previous runs
-os.system('rm -rf land_ice/output*')
-os.system('rm -rf land_ice/restart*')
-os.system('rm -rf ocean/output*')
-os.system('rm -rf ocean/restart*')
+os.system('rm -rf land_ice/output.*')
+os.system('rm -rf land_ice/restart.*')
+os.system('rm -rf land_ice/log.*')
+os.system('rm -rf land_ice/output-files/*')
+os.system('rm -rf ocean/output.*')
+os.system('rm -rf ocean/restart.*')
+os.system('rm -rf ocean/output-files/*')
for t in range(numCoupleIntervals):
print '============================================================'
print 'Starting coupling number ', t
- # Get the appropriate source and destination files for the input needed by the coupler
+ # =======================================
+ # 1. Get the appropriate files and variables for the coupler
+ # =======================================
+
+ # Get the appropriate source and destination files for the input/output of the coupler
if t==0:
- # Use I.C. from input files at initial time
- liceSourceFile = NetCDFFile('./land_ice/'+liceInput,'r')
-# ocnSourceFile = NetCDFFile('./ocean/'+ocnInput,'r')
+ if runLICE:
+ # Use I.C. from input files as Source at initial time
+ liceSourceFile = NetCDFFile('./land_ice/'+liceInput,'r')
+ # At initial time, we want to write to the input files
+ liceDestFile = NetCDFFile('./land_ice/'+liceInput,'r+')
+ else:
+ # TODO Read from a file
+ pass
- # At initial time, we want to write to the input files
- liceDestFile = NetCDFFile('./land_ice/'+liceInput,'r+')
-# ocnDestFile = NetCDFFile('./ocean/'+ocnInput,'r+')
+ if runOCN:
+ # TODO - set this up properly for the ocean model
+ # Use I.C. from input files as Source at initial time
+ ocnSourceFile = NetCDFFile('./ocean/'+ocnInput,'r')
+ # At initial time, we want to write to the input files
+ ocnDestFile = NetCDFFile('./ocean/'+ocnInput,'r+')
+ else:
+ # TODO Read from a file
+ pass
else:
- # Use values from output files at other times.
- p = os.popen('ls -1t ./land_ice/output.*.nc | head -n 1') # get the most recent output file
- liceOutput = p.read().rstrip() # get stdout
- liceSourceFile = NetCDFFile(liceOutput, 'r')
-# p = os.popen('ls -1t ./ocean/output.*.nc | head -n 1') # get the most recent output file
-# ocnOutput = p.read().rstrip() # get stdout
-# ocnSourceFile = NetCDFFile(ocnOutput, 'r')
+ if runLICE:
+ # Use values from output files as Source at other times.
+ liceSourceFile = NetCDFFile('./land_ice/' + liceOutput, 'r')
- # At other times we want to write to the restart files
- # Need to choose the appropriate restart file - there may be smarter ways to do this, but an easy way is to take the most recent
- # os.listdir does not have a sort ability, so use a shell command (unpythonic but works)
- p = os.popen('ls -1t ./land_ice/restart.*.nc | head -n 1')
- liceRestart = p.read().rstrip() # get stdout
- liceDestFile = NetCDFFile(liceRestart,'r+')
-# p = os.popen('ls -1t ./ocean/restart.*.nc | head -n 1')
-# ocnRestart = p.read().rstrip() # get stdout
-# ocnDestFile = NetCDFFile(ocnRestart,'r+')
-
+ # At other times we want to write to the restart files as the Destination
+ # Need to choose the appropriate restart file - there may be smarter ways to do this, but an easy way is to take the most recent
+ # os.listdir does not have a sort ability, so use a shell command (unpythonic but works)
+ p = os.popen('ls -1t ./land_ice/restart.*.nc | head -n 1')
+ liceRestart = p.read().rstrip() # get stdout
+ print 'For land ice model restart data, using file: ' + liceRestart
+ liceDestFile = NetCDFFile(liceRestart,'r+')
+ else:
+ # TODO Read from a file
+ pass
- # Get the source fields needed by the coupler
+ if runOCN:
+ # Use values from output files as Source at other times.
+ # TODO - set this up properly for the ocean model
+ pass
+
+ else:
+ # TODO Read from a file
+ pass
+
+ # Get the Source fields needed by the coupler
# (Get the -1 time level for all variables - this will be the last whether this is from input or output
- liceThickness = liceSourceFile.variables['thickness'][-1,:]
- liceBasalTemp = liceSourceFile.variables['temperature'][-1,:,-1]
-# ocnSurfTemp = ocnSourceFile.variables['temperature'][-1,:,0]
-# ocnSurfSalin = ocnSourceFile.variables['salinity'][-1,:,0]
+ try:
+ liceThickness = liceSourceFile.variables['thickness'][-1,:]
+ liceBasalTemp = liceSourceFile.variables['temperature'][-1,:,-1]
+ except:
+ print "Problem getting needed Source fields from the land ice files"
+ try:
+ ocnSurfTemp = ocnSourceFile.variables['temperature'][-1,:,0]
+ ocnSurfSalin = ocnSourceFile.variables['salinity'][-1,:,0]
+ except:
+ print "Problem getting needed Source fields from the ocean files"
# Get the destination fields needed by the coupler
- liceBMB = liceDestFile.variables['marineBasalMassBalTimeSeries'][:,0] # There should only be on time level
- #liceBHF = liceDestFile.variables['basalHeatFluxTimeSeries'][:,0] # NEED TO RESTRICT THIS TO CELLS WITH ICE
-# ocnSurfPressure = ocnDestFile.variables['ssh'][0,:] # TODO WHAT IS THE PROPER VARIABLE HERE?
- # TODO WHAT ARE THE PROPER VARIABLE FOR MASS, HEAT, AND SALINITY FLUXES?
+ try:
+ liceBMB = liceDestFile.variables['marineBasalMassBalTimeSeries'][:,0] # There should only be on time level
+ #liceBHF = liceDestFile.variables['basalHeatFluxTimeSeries'][:,0] # NEED TO RESTRICT THIS TO CELLS WITH ICE
+ except:
+ print "Problem getting needed Destination fields from the land ice files"
+ try:
+ ocnSurfPressure = ocnDestFile.variables['ssh'][0,:] # TODO WHAT IS THE PROPER VARIABLE HERE?
+ # TODO WHAT ARE THE PROPER VARIABLE FOR MASS, HEAT, AND SALINITY FLUXES?
+ except:
+ print "Problem getting needed Destination fields from the ocean files"
+
+ # =======================================
+ # 2. Do the coupler calculations
+ # =======================================
+
# Calculate the ocean surface pressure - this could be a function call
# TODO WHAT TO DO ABOUT NON-ICE SHELF OCEAN CELLS?
# ocnSurfPressure = rhoi * grav * liceThickness # TODO ARE THESE THE PROPER UNITS?
@@ -105,86 +156,124 @@
# Calculate boundary layer fluxes - this could be an internal or external function call
#calculateBdyLyrFluxes(liceBasalTemp, ocnSurfTemp, ocnSurfSalin, liceBHF, ocnMassFlux, ocnHeatFlux, ocnSalinFlux) #TODO
+
+ # =======================================
+ # 3. Write the new forcing data to the files
+ # =======================================
+
# Write new quantities to the destination files
# TODO Does this need to be done explicitly? I always forget how netCDF i/o works in python precisely.
+
# =======================================
- # Run the models # TODO Do we want to do this in parallel?
-# try:
-# os.chdir('ocean')
-# # Set the namelist file properly for a restart or no
-# print 'Setting up ocean namelist file'
-# if t==0:
-# subprocess.check_call('sed -i.SEDBACKUP "s/^.*config_do_restart.*/ config_do_restart = .false./" namelist.input' , shell=True, executable='/bin/bash') # no restart
-# subprocess.check_call('sed -i.SEDBACKUP "s/^.*config_write_output_on_startup.*/ config_write_output_on_startup = .true./" namelist.input' , shell=True, executable='/bin/bash') # write the initial time
-# subprocess.check_call('sed -i.SEDBACKUP "s/^.*config_start_time.*/ config_start_time = \'0000-01-01_00:00:00\'/" namelist.input' , shell=True, executable='/bin/bash') # set the start time
-# else:
-# subprocess.check_call('sed -i.SEDBACKUP "s/^.*config_do_restart.*/ config_do_restart = .true./" namelist.input' , shell=True, executable='/bin/bash') # do a restart
-# subprocess.check_call('sed -i.SEDBACKUP "s/^.*config_write_output_on_startup.*/ config_write_output_on_startup = .false./" namelist.input' , shell=True, executable='/bin/bash') # don't write the initial time
-# # Get the old stop time so we can set it to be the new start time.
-# p = os.popen('ls -1t output*nc | head -n 1')
-# oldoutputfile=p.read().rstrip()
-# p = os.popen('ncks -v xtime -H -s "%c" ' + oldoutputfile + ' | tail -r -c 65')
-# oldstoptime=p.read().rstrip()
-# subprocess.check_call('sed -i.SEDBACKUP "s/^.*config_start_time.*/ config_start_time = ' + oldstoptime + '/" namelist.input ' , shell=True, executable='/bin/bash') # set the start time to the old stop time
+ # 4. Run the models
+ # =======================================
+ # TODO Do we want to do this in parallel? Probably no reason to do so.
-# print 'Starting ocean model.'
-# subprocess.check_call(runOCN, shell=True, executable='/bin/bash')
-# print 'Ocean model completed.'
-# # Don't need to copy output file if the namelist is setup properly
-# ### Copy output file so we don't lose it when the model restarts - could use Python's shutil module. TODO Also may want to change format of the string conversion of t
-# ##subprocess.check_call('cp ' + ocnOutput + ' output.' + str(t) + '.nc', shell=True, executable='/bin/bash')
-# os.chdir('..')
-# except:
-# sys.exit('Ocean model failed!')
-
+ # === LAND ICE MODEL ===
try:
- os.chdir('land_ice')
- # Set the namelist file properly for a restart or no and set the start and stop times
- print 'Setting up land ice namelist file'
- if t==0:
- subprocess.check_call('sed -i.SEDBACKUP "s/^.*config_do_restart.*/ config_do_restart = .false./" namelist.input' , shell=True, executable='/bin/bash') # no restart
- print 'set no restart'
- subprocess.check_call('sed -i.SEDBACKUP "s/^.*config_write_output_on_startup.*/ config_write_output_on_startup = .true./" namelist.input' , shell=True, executable='/bin/bash') # write the initial time
- print 'set write initial time'
- subprocess.check_call('sed -i.SEDBACKUP "s/^.*config_start_time.*/ config_start_time = \'0000-01-01_00:00:00\'/" namelist.input' , shell=True, executable='/bin/bash') # set the start time
- print 'set start time'
- else:
- subprocess.check_call('sed -i.SEDBACKUP "s/^.*config_do_restart.*/ config_do_restart = .true./" namelist.input' , shell=True, executable='/bin/bash') # do a restart
- print 'set restart'
- subprocess.check_call('sed -i.SEDBACKUP "s/^.*config_write_output_on_startup.*/ config_write_output_on_startup = .true./" namelist.input' , shell=True, executable='/bin/bash') # don't write the initial time
- print 'disabled initial write'
- # Get the old stop time so we can set it to be the new start time. - it seems more reliable to search for the most recent restart file rather than the most recent output file.
- p = os.popen('ls -1t restart*nc | head -n 1')
- lastrestartfile=p.read().rstrip()
- p = os.popen('ncks -v xtime -H -s "%c" ' + lastrestartfile + ' | tail -r -c 65')
- oldstoptime=p.read().rstrip()
- subprocess.check_call('sed -i.SEDBACKUP "s/^.*config_start_time.*/ config_start_time = ' + oldstoptime + '/" namelist.input ' , shell=True, executable='/bin/bash') # set the start time to the old stop time
- print 'set start time'
+ if runLICE:
+ print 'hello'
+ os.chdir('land_ice')
+ # Set the namelist file properly for a restart or no and set the start and stop times
+ print 'Setting up land ice namelist file'
+ if t==0:
+ subprocess.check_call('sed -i.SEDBACKUP "s/^.*config_do_restart.*/ config_do_restart = .false./" namelist.input' , shell=True, executable='/bin/bash') # no restart
+ print 'set no restart'
+ subprocess.check_call('sed -i.SEDBACKUP "s/^.*config_write_output_on_startup.*/ config_write_output_on_startup = .true./" namelist.input' , shell=True, executable='/bin/bash') # write the initial time
+ print 'set write initial time'
+ subprocess.check_call('sed -i.SEDBACKUP "s/^.*config_start_time.*/ config_start_time = \'0000-01-01_00:00:00\'/" namelist.input' , shell=True, executable='/bin/bash') # set the start time
+ print 'set start time'
+ else:
+ subprocess.check_call('sed -i.SEDBACKUP "s/^.*config_do_restart.*/ config_do_restart = .true./" namelist.input' , shell=True, executable='/bin/bash') # do a restart
+ print 'set restart'
+ subprocess.check_call('sed -i.SEDBACKUP "s/^.*config_write_output_on_startup.*/ config_write_output_on_startup = .true./" namelist.input' , shell=True, executable='/bin/bash') # don't write the initial time
+ print 'disabled initial write'
+ # Get the old stop time so we can set it to be the new start time. - it seems more reliable to search for the most recent restart file rather than the most recent output file.
+ p = os.popen('ls -1t restart*nc | head -n 1')
+ lastrestartfilename=p.read().rstrip()
+ p = os.popen('ncks -v xtime -H -s "%c" ' + lastrestartfilename + ' | tail -r -c 65')
+ oldstoptime=p.read().rstrip()
+ print 'For land ice, getting the restart time from file: ' + lastrestartfilename + ' which gives a time of ' + oldstoptime
+ subprocess.check_call('sed -i.SEDBACKUP "s/^.*config_start_time.*/ config_start_time = ' + oldstoptime + '/" namelist.input ' , shell=True, executable='/bin/bash') # set the start time to the old stop time
+ print 'set start time'
- print 'Starting land ice model.'
- subprocess.check_call(runLICE, shell=True, executable='/bin/bash')
- print 'Land ice model completed.'
- # Don't need to copy output file if the namelist is setup properly
- ### Copy output file so we don't lose it when the model restarts - could use Python's shutil module. TODO Also may want to change format of the string conversion of t
- ##subprocess.check_call('cp ' + liceOutput + ' output.' + str(t) + '.nc', shell=True, executable='/bin/bash')
- os.chdir('..')
+ print 'Starting land ice model.'
+ subprocess.check_call(runLICEcmd, shell=True, executable='/bin/bash')
+ print 'Land ice model completed.'
+ # Don't need to copy output file if the namelist is setup properly
+ ### Copy output file so we don't lose it when the model restarts - could use Python's shutil module. TODO Also may want to change format of the string conversion of t
+ subprocess.check_call('cp ' + liceOutput + ' ./output-files/output.' + '{0:04d}'.format(t+1) + '.nc', shell=True, executable='/bin/bash') # the number is the ending time level (if there are more than one)
+ subprocess.check_call('cp log.0000.out ./output-files/log.0000.' + '{0:04d}'.format(t+1) + '.out', shell=True, executable='/bin/bash') # the number is the ending time level (if there are more than one)
+ subprocess.check_call('cp log.0000.err ./output-files/log.0000.' + '{0:04d}'.format(t+1) + '.err', shell=True, executable='/bin/bash') # the number is the ending time level (if there are more than one)
+
+ os.chdir('..')
except:
sys.exit('Land Ice model failed!')
- #time.sleep(62)
+ # === OCEAN MODEL ===
+ try:
+ if runOCN:
+ os.chdir('ocean')
+ # Set the namelist file properly for a restart or no
+ print 'Setting up ocean namelist file'
+ if t==0:
+ subprocess.check_call('sed -i.SEDBACKUP "s/^.*config_do_restart.*/ config_do_restart = .false./" namelist.input' , shell=True, executable='/bin/bash') # no restart
+ subprocess.check_call('sed -i.SEDBACKUP "s/^.*config_write_output_on_startup.*/ config_write_output_on_startup = .true./" namelist.input' , shell=True, executable='/bin/bash') # write the initial time
+ subprocess.check_call('sed -i.SEDBACKUP "s/^.*config_start_time.*/ config_start_time = \'0000-01-01_00:00:00\'/" namelist.input' , shell=True, executable='/bin/bash') # set the start time
+ else:
+ subprocess.check_call('sed -i.SEDBACKUP "s/^.*config_do_restart.*/ config_do_restart = .true./" namelist.input' , shell=True, executable='/bin/bash') # do a restart
+ subprocess.check_call('sed -i.SEDBACKUP "s/^.*config_write_output_on_startup.*/ config_write_output_on_startup = .false./" namelist.input' , shell=True, executable='/bin/bash') # don't write the initial time
+ # Get the old stop time so we can set it to be the new start time.
+ p = os.popen('ls -1t output*nc | head -n 1')
+ oldoutputfile=p.read().rstrip()
+ p = os.popen('ncks -v xtime -H -s "%c" ' + oldoutputfile + ' | tail -r -c 65')
+ oldstoptime=p.read().rstrip()
+ subprocess.check_call('sed -i.SEDBACKUP "s/^.*config_start_time.*/ config_start_time = ' + oldstoptime + '/" namelist.input ' , shell=True, executable='/bin/bash') # set the start time to the old stop time
+
+
+ print 'Starting ocean model.'
+ subprocess.check_call(runOCNcmd, shell=True, executable='/bin/bash')
+ print 'Ocean model completed.'
+ # Don't need to copy output file if the namelist is setup properly
+ ### Copy output file so we don't lose it when the model restarts - could use Python's shutil module. TODO Also may want to change format of the string conversion of t
+ ##subprocess.check_call('cp ' + ocnOutput + ' output.' + str(t) + '.nc', shell=True, executable='/bin/bash')
+ os.chdir('..')
+ except:
+ sys.exit('Ocean model failed!')
+
+
# =======================================
print 'Coupled model run complete!'
-# After the coupled run is completed, we may want to go in and use NCO to patch together all of the output files for each model.
+# =======================================
+# 5. Post-processing
+# =======================================
+# After the coupled run is completed, we may want to go in and use NCO to patch together all of the output files for each model.
+if runLICE:
+ try:
+ os.chdir('land_ice/output-files')
+ subprocess.check_call('ncrcat output.*.nc output-all.nc', shell=True, executable='/bin/bash')
+ os.chdir('../..')
+ except:
+ sys.exit('Failed combining land ice output files!')
+# TODO Does the ocean want to do this too? If not, delete these lines.
+if runOCN:
+ try:
+ os.chdir('ocean/output-files')
+ subprocess.check_call('ncrcat output.*.nc output-all.nc', shell=True, executable='/bin/bash')
+ os.chdir('../..')
+ except:
+ sys.exit('Failed combining ocean output files!')
+print 'Combined output files.'
+print 'Successful competion.'
-
Added: branches/landice_ocean_coupler/namelist.input.land_ice
===================================================================
--- branches/landice_ocean_coupler/namelist.input.land_ice (rev 0)
+++ branches/landice_ocean_coupler/namelist.input.land_ice 2013-03-11 17:37:48 UTC (rev 2581)
@@ -0,0 +1,55 @@
+&land_ice_model
+ config_time_integration = 'ForwardEuler'
+ config_thickness_advection = 'FO-Upwind'
+ config_tracer_advection = 'None'
+ config_stats_interval = 1
+ config_dycore = 'SIA'
+ config_forcing_frequency = 'Annual'
+ config_allow_additional_advance = true
+ config_marine_advance = 'CFBC'
+! config_calving_law = 'critical-thicknes'
+! config_calving_critical_thickness = 300.0
+ config_dynamic_thickness = 0.0
+
+! use the same calendar type in both models
+ config_calendar_type = 'gregorian_noleap'
+! set the dt of the land ice model to be evenly divisble by the dt of the ocean model
+ config_dt_years = 1.0
+! the start time does not matter because the script will set it. But it needs to be present.
+ config_start_time = 0007-01-01_00:00:00
+! the script uses run_duration instead of stop_time. Be sure run_duration is present and stop_time is not. The run duration of both the land ice and ocean models should be the same
+ config_run_duration = '365_00:00:00'
+!!! config_stop_time = '0001-01-01_00:00:00'
+! write_output_on_startup must be present. value does not matter because script will set it.
+ config_write_output_on_startup = .true.
+/
+
+&advection
+ config_vert_tracer_adv_order = 2
+ config_horiz_tracer_adv_order = 3
+ config_coef_3rd_order = 0.25
+ config_monotonic = true
+ config_check_monotonicity = true
+/
+
+&io
+! you can use any name for input, but script assumes restart_name is restart.nc and output_name is output.nc
+ config_input_name = 'land_ice_grid.nc'
+ config_output_name = 'output.nc'
+ config_restart_name = 'restart.nc'
+! make sure you get at least one output per run
+ config_output_interval = '00_00:01:00'
+! frames_per_outfile needs to be 0. The script assumes that the output from each run is named 'output.nc',
+! which will be the case when 0 is used. If a positive integer is used, then MPAS adds a timestamp to each file name,
+! but this mode causes problems because it also currently writes an empty file for the initial time of output write_output_on_startup=false
+ config_frames_per_outfile = 0
+/
+
+&restart
+! You want to make sure a restart file is written with the last time step of the run (if not more frequent).
+! For the ice model we will probably not do more than a single time step per run, so set this to write a restart every time step.
+ config_restart_interval = '001_00:00:00'
+! do_restart needs to be present but the value doesn't matter because script will set it.
+ config_do_restart = .true.
+
+/
</font>
</pre>