Merge pull request mom-ocean#1031 from marshallward/tc4_fix

Tc4 fix

.testing/Makefile

@@ -251,6 +251,7 @@ results/%/ocean.stats.$(1): ../build/$(2)/MOM6
 	if [ $(3) ]; then find ../build/$(2) -name *.gcda -exec rm -f '{}' \; ; fi
 	mkdir -p work/$$*/$(1)
 	cp -rL $$*/* work/$$*/$(1)
+	cd work/$$*/$(1) && if [ -f Makefile ]; then make; fi
 	mkdir -p work/$$*/$(1)/RESTART
 	echo $(4) > work/$$*/$(1)/MOM_override
 	cd work/$$*/$(1) && $$(call MPIRUN_CMD,$(5)) -n $(6) ../../../$$< 2> debug.out > std.out \
@@ -285,6 +286,7 @@ results/%/ocean.stats.restart: ../build/symmetric/MOM6
 	rm -rf work/$*/restart
 	mkdir -p work/$*/restart
 	cp -rL $*/* work/$*/restart
+	cd work/$*/restart && if [ -f Makefile ]; then make; fi
 	mkdir -p work/$*/restart/RESTART
 	# Generate the half-period input namelist
 	# TODO: Assumes runtime set by DAYMAX, will fail if set by input.nml
@@ -294,20 +296,19 @@ results/%/ocean.stats.restart: ../build/symmetric/MOM6
 		&& if [ -z "$${timeunit}" ]; then timeunit="8.64e4"; fi \
 		&& printf -v timeunit_int "%.f" "$${timeunit}" \
 		&& halfperiod=$$(printf "%.f" $$(bc <<< "scale=10; 0.5 * $${daymax} * $${timeunit_int}")) \
-		&& printf "\n&ocean_solo_nml\n    seconds = $${halfperiod}\n/\n" >> input.nml \
-		&& echo $${daymax} $${timeunit}
+		&& printf "\n&ocean_solo_nml\n    seconds = $${halfperiod}\n/\n" >> input.nml
 	# Run the first half-period
-	cd work/$*/restart && $(MPIRUN) -n 1 ../../../$< 2> debug.out > std.out \
-		|| ! sed 's/^/$*.restart1: /' std.out debug.out \
-		&& sed 's/^/$*.restart1: /' std.out
+	cd work/$*/restart && $(MPIRUN) -n 1 ../../../$< 2> debug1.out > std1.out \
+		|| ! sed 's/^/$*.restart1: /' std1.out debug1.out \
+		&& sed 's/^/$*.restart1: /' std1.out
 	# Setup the next inputs
 	cd work/$*/restart && rm -rf INPUT && mv RESTART INPUT
 	mkdir work/$*/restart/RESTART
 	cd work/$*/restart && sed -i -e "s/input_filename *= *'n'/input_filename = 'r'/g" input.nml
 	# Run the second half-period
-	cd work/$*/restart && $(MPIRUN) -n 1 ../../../$< 2> debug.out > std.out \
-		|| ! sed 's/^/$*.restart2: /' std.out debug.out \
-		&& sed 's/^/$*.restart2: /' std.out
+	cd work/$*/restart && $(MPIRUN) -n 1 ../../../$< 2> debug2.out > std2.out \
+		|| ! sed 's/^/$*.restart2: /' std2.out debug2.out \
+		&& sed 's/^/$*.restart2: /' std2.out
 	# Archive the results and cleanup
 	mkdir -p $(@D)
 	cp work/$*/restart/ocean.stats $@

.testing/_tc4/MOM_input → .testing/tc4/MOM_input

@@ -7,10 +7,15 @@
 USE_REGRIDDING = True           !   [Boolean] default = False
                                 ! If True, use the ALE algorithm (regridding/remapping). If False, use the
                                 ! layered isopycnal algorithm.
-DT = 300.0                      !   [s]
+DT = 1200.0                      !   [s]
                                 ! The (baroclinic) dynamics time step.  The time-step that is actually used will
                                 ! be an integer fraction of the forcing time-step (DT_FORCING in ocean-only mode
                                 ! or the coupling timestep in coupled mode.)
+DT_THERM = 3600.0                !   [s] default = 300.0
+                                ! The thermodynamic and tracer advection time step. Ideally DT_THERM should be
+                                ! an integer multiple of DT and less than the forcing or coupling time-step,
+                                ! unless THERMO_SPANS_COUPLING is true, in which case DT_THERM can be an integer
+                                ! multiple of the coupling timestep.  By default DT_THERM is set to DT.				
 C_P = 3925.0                    !   [J kg-1 K-1] default = 3991.86795711963
                                 ! The heat capacity of sea water, approximated as a constant. This is only used
                                 ! if ENABLE_THERMODYNAMICS is true. The default value is from the TEOS-10
@@ -377,10 +382,15 @@ WIND_CONFIG = "zero"            !
 ! === module MOM_restart ===
 
 ! === module MOM_main (MOM_driver) ===
-DAYMAX = 1.0                    !   [days]
+DAYMAX = 0.25                    !   [days]
                                 ! The final time of the whole simulation, in units of TIMEUNIT seconds.  This
                                 ! also sets the potential end time of the present run segment if the end time is
                                 ! not set via ocean_solo_nml in input.nml.
+
+ENERGYSAVEDAYS = 0.125          !   [days] default = 1.44E+04
+                                ! The interval in units of TIMEUNIT between saves of the
+                                ! energies of the run and other globally summed diagnostics.
+
 RESTART_CONTROL = 3             ! default = 1
                                 ! An integer whose bits encode which restart files are written. Add 2 (bit 1)
                                 ! for a time-stamped file, and odd (bit 0) for a non-time-stamped file. A
@@ -397,3 +407,6 @@ MAXCPU = 2.88E+04               !   [wall-clock seconds] default = -1.0
                                 ! processors used.
 
 ! === module MOM_file_parser ===
+
+DIAG_AS_CHKSUM = True
+DEBUG = True

.testing/tc4/Makefile

@@ -0,0 +1,3 @@
+all:
+	python build_grid.py
+	python build_data.py

.testing/tc4/build_data.py

@@ -0,0 +1,80 @@
+import netCDF4 as nc
+import numpy as np
+
+x = nc.Dataset('ocean_hgrid.nc').variables['x'][1::2, 1::2]
+y = nc.Dataset('ocean_hgrid.nc').variables['y'][1::2, 1::2]
+zbot = nc.Dataset('topog.nc').variables['depth'][:]
+zbot0 = zbot.max()
+
+
+def t_fc(x, y, z, radius=5.0, tmag=1.0):
+    """a radially symmetric anomaly in the center of the domain.
+    units are meters and degC.
+    """
+    ny, nx = x.shape
+    nz = z.shape[0]
+
+    x0 = x[int(ny/2), int(nx/2)]
+    y0 = y[int(ny/2), int(nx/2)]
+
+    tl = np.zeros((nz, ny, nx))
+    zb = z[-1]
+    if len(z) > 1:
+        zd = z / zb
+    else:
+        zd = [0.]
+    for k in np.arange(len(zd)):
+        r = np.sqrt((x - x0)**2 + (y - y0)**2)
+        tl[k, :] += (1.0 - np.minimum(r / radius, 1.0)) * tmag * (1.0 - zd[k])
+    return tl
+
+
+ny, nx = x.shape
+nz = 3
+z = (np.arange(nz) * zbot0) / nz
+
+temp = t_fc(x, y, z)
+salt = np.zeros(temp.shape)+35.0
+fl = nc.Dataset('temp_salt_ic.nc', 'w', format='NETCDF3_CLASSIC')
+fl.createDimension('lon', nx)
+fl.createDimension('lat', ny)
+fl.createDimension('depth', nz)
+fl.createDimension('Time', None)
+zv = fl.createVariable('depth', 'f8', ('depth'))
+lonv = fl.createVariable('lon', 'f8', ('lon'))
+latv = fl.createVariable('lat', 'f8', ('lat'))
+timev = fl.createVariable('Time', 'f8', ('Time'))
+timev.calendar = 'noleap'
+timev.units = 'days since 0001-01-01 00:00:00.0'
+timev.modulo = ' '
+tv = fl.createVariable('ptemp', 'f8', ('Time', 'depth', 'lat', 'lon'),
+                       fill_value=-1.e20)
+sv = fl.createVariable('salt', 'f8', ('Time', 'depth', 'lat', 'lon'),
+                       fill_value=-1.e20)
+tv[:] = temp[np.newaxis, :]
+sv[:] = salt[np.newaxis, :]
+zv[:] = z
+lonv[:] = x[0, :]
+latv[:] = y[:, 0]
+timev[0] = 0.
+fl.sync()
+fl.close()
+
+
+# Make Sponge forcing file
+dampTime = 20.0     # days
+secDays = 8.64e4
+fl = nc.Dataset('sponge.nc', 'w', format='NETCDF3_CLASSIC')
+fl.createDimension('lon', nx)
+fl.createDimension('lat', ny)
+lonv = fl.createVariable('lon', 'f8', ('lon'))
+latv = fl.createVariable('lat', 'f8', ('lat'))
+spv = fl.createVariable('Idamp', 'f8', ('lat', 'lon'), fill_value=-1.e20)
+Idamp = np.zeros((ny, nx))
+if dampTime > 0.:
+    Idamp = 0.0 + 1.0 / (dampTime * secDays)
+spv[:] = Idamp
+lonv[:] = x[0, :]
+latv[:] = y[:, 0]
+fl.sync()
+fl.close()