-
Notifications
You must be signed in to change notification settings - Fork 37
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
Port the MISMIP+ spinup testcase from the old version of compass #737
Conversation
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
@andrewdnolan, this is looking good! I made some comments regarding compass conventions and whatnot, but this has come a long way in a short time. However, I'm getting an error when running this on Perlmutter. Here's what I take to be the important points from the traceback:
File "/global/cfs/cdirs/fanssie/users/trhille/compass/compass/landice/tests/mismipplus/setup_mesh.py", line 265, in _setup_MISMPPlus_IC
src['dirichletVelocityMask'].loc[:] = xr.where(mask, 1, 0)
...
ValueError: could not broadcast input array from shape (1,1,3456) into shape (1,3456,11)
It looks like it created the landice_grid.nc
file successfully, but is throwing an error when trying to populate the mesh.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
@andrewdnolan , this looks great! MISMIP+ meshes will be very useful to have available, and recently we've just been avoiding using this test for lack of meshes, so this will be great to have. Overall, this PR does everything it sets out to do, and the code is very clean and easy to follow. I made a number of mostly minor comments/suggestions, plus raised a couple questions (mostly for @trhille ).
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
@andrewdnolan , I did a quick review now to get you quick feedback. Overall, this looks fantastic. All of my comments are nitpicky things. I appreciate the level of detail that went in to getting a robust product here, and despite what I said initially about this mesh probably not getting used much, we keep finding new reasons we'll need it in the coming months, so I think this was worth doing right. I'll review again once the docs are added, and I will then actually test it by running it.
<var name="allowableDtACFL"/> | ||
<var name="allowableDtDCFL"/> | ||
<var name="calvingCFLdt"/> | ||
<var name="dtCalvingCFLratio"/> |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
@andrewdnolan , on further reflection, let's not remove these calving related variables here. They are all scalars so don't really hurt leaving in. That way if someone modifies this test to include calving, they don't have to remember to add these in.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Okay! Works for me. Are they only needed in the globalStatsOutput
stream? Or would you like them in the restart
and/or output
streams?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Just in the globalStats file is good.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Added!
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Another quick pass with a few changes. I'll follow up with actually running this in the next round!
<var name="allowableDtACFL"/> | ||
<var name="allowableDtDCFL"/> | ||
<var name="calvingCFLdt"/> | ||
<var name="dtCalvingCFLratio"/> |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Just in the globalStats file is good.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
@andrewdnolan, this is excellent stuff! Thanks for your hard work. I have a handful of requested changes, but nothing major.
# Number of vertical levels | ||
levels = 10 | ||
|
||
# Initial ice thickness (m) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Needs another sentence of explanation to this comment to tell the user why we need an initial ice thickness, where on the domain it is applied, and under what circumstances a user might want to change this.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Added a bit more info to the comment in commit 7ef36a1
. Tried to explain where/why this is needed.
# create the calvingMask data array and add Time dimension along axis 0 | ||
calvingMask = xr.where(mask, 1, 0).expand_dims({"Time": 1}, axis=0) | ||
# assign data array to dataset and ensure it's a 32 bit int field | ||
src['calvingMask'] = calvingMask.astype('int32') |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Does this not need .loc[:]
as well?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
If we included .loc[:]
it would throw an error, because calvingMask
doesn't exist as a variable in src
until line 452.
Those other variables, where .loc[:]
is used, already exist in src
since they were created in the create_landice_grid_from_generic_MPAS_grid.py
subprocess call. calvingMask
isn't initialized within that script, so there are no existing dimension to broadcast/align with by using .loc[:]
.
You'll notice here there a little bit of extra stuff, specifically (...).expand_dims({"Time": 1}, axis=0)
in line 450. Since calvingMask
doesn't exist in the dataset, I have to manually add the Time
dimension to the dataarray, before placing it in the dataset.
(I could manually add the missing dimensions to all the local variables through the expand_dims
function, but I find that a lot more cumbersome than the .loc[:]
syntax. Plus it's more error prone since you need to specify the axis to add the missing dimension on, whereas .loc[:]
will automatically determine the correct location).
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Okay, thanks for the explanation. No need to change how you have it here.
@andrewdnolan , when I try to set up the mismip+ smoke test, I get a cfg error:
Looks like we need tweak the cfg hierarchy a bit. |
# check if the resolution has been changed since the `compass setup` | ||
# command was run | ||
if self.resolution != resolution: | ||
raise Exception(f'Resolution was set at {self.resolution:2d}m' |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I got this exception to be raised with the -f
config testing I mentioned above. In so doing, I got this value error for this line:
File "/usr/projects/climate/mhoffman/compass/mismip+_spinup/compass/landice/tests/mismipplus/setup_mesh.py", line 94, in run
raise Exception(f'Resolution was set at {self.resolution:2d}m'
ValueError: Unknown format code 'd' for object of type 'float'
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Thanks for testing this! I wasn't really aware of the -f
option.
I've revamped the way the directory structure is created (see response to your original comment). With that revamp, this bug should be fixed and configuration options passed from the command line should be properly parsed. (Also the error you got in the SmokeTest
should be fixed).
All that was pushed in 4b0f35a.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
@andrewdnolan, I just have a handful of minor comments. I tested this by creating and running 1-, 2-, 4-, 8-, and 16-km meshes, and everything seems to be functioning as expected!
@@ -0,0 +1,5 @@ | |||
config_run_duration = '20000-00-00_00:00:00' |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I'd prefer using config_stop_time
rather than config_run_duration
because in practice the latter often causes overshoots of the desired stop time when you do a restart.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Gotcha, I've switched to config_stop_time
for both the SmokeTest
and the SpinUp
.
@@ -56,17 +72,20 @@ def __init__(self, test_case, name='run_model', subdir=None, | |||
if suffixes is None: | |||
suffixes = ['landice'] | |||
self.suffixes = suffixes | |||
|
|||
if resolution is None: | |||
resolution = 0.0 |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
How does this work?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Yeah this is definitely strange but has to do with the very different way resolution
is set between SmokeTest
and SpinUp
test cases.
SmokeTest
So, this last round of commits removes the option to specify the SmokeTest
resolution in a config file. The problem with the former approach is a user could request a resolution not supported by test case (currently, any value other than 2000m
). With the current code, the supported resolutions are added as individual test cases at the test group level:
https://github.com/andrewdnolan/compass/blob/67b26525f73724c1d6ba0b2f5d7c724427a2c9e0/compass/landice/tests/mismipplus/__init__.py#L19-L21
(This also has the added benefit of the supported resolutions being shown in compass list
). So, for the SmokeTest
a resolution is always passed to the initialization and in the code above resolution
will never be None
.
SpinUp
Okay, so now we get to why a "dummy" value of 0.0
is specified. For the SpinUp
we want to support arbitrary resolution, which should be set in the users config file. The big hiccup here is that we want to use the requested resolution in the directory structure, but the config file which isn't fully parsed until after the step is initialized. Ideally we'd be able to pass resolution
into the initialization of the SpinUp
testcase, but that doesn't work. (We can parse the resolution value in the compass repo's config before initialization, but can not parse a user specified config file (e.g. compass setup -f ....
) prior to the testcases initialization). So, instead we have this dummy value of 0.0
and instead set the value of the resolution
attribute in the SpinUp.configure
method so that we have access to all the config options:
https://github.com/andrewdnolan/compass/blob/67b26525f73724c1d6ba0b2f5d7c724427a2c9e0/compass/landice/tests/mismipplus/spin_up/__init__.py#L49-L83
So, in short the SmokeTest
knows the value of the resolution
attribute at the time of initialization, whereas the SpinUp
needs to parse resolution
from the config options, so a temporary dummy value is used until the resolution
attribute is set by the SpinUp.configure
method.
Now, you might ask: Why set the dummy value of 0.0
at all? Why not just leave the resolution
attribute unset for the SpinUp
testcase until the configure method is called? That's definitely an option, but in talking to Xylar that seems like poor python practice. I guess it's best to initialize all attributes within the __init__
method.
Happy to ditch this if need be in order to have a bit clearer but less "proper" python code.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Okay, this is fine. I just wanted to make sure it was actually doing something intentional.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I added a comment in f916f60 which summarizes the message above to provide a brief explanation to new/future users.
The minimum number of tasks | ||
""" | ||
|
||
# read MPI related parameters from configuration file |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Can we just use config.getint()
instead of getting floats and then converting to floats below?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I've switched to parsing all the parameters as int
s and then casting them as floats where need be. This was done as of f916f60. Thanks for catching this!
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Do we need to cast these as floats at all? It does no harm, but it will give exactly the same answer as using ints because integer division yields a float anyway.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Your right, in python casting as floats is kind of overkill. Especially since all of the floating point operations with int
s we are doing are division, where we are explicitly using the floating point operator (e.g. /
vs. //
).
I was just trying to follow best practice, but seems a little unnecessary here and at the expense of readable code. I'll remove the casting as floats to make things a little more concise.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Thanks. I wasn't sure what was considered best practice here, but I think it is more readable without casting them as floats.
Moved the SmokeTest to it's own case directory. Package imports also reorganized (by isort hook) in order to be alphabetical.
Added a mock testcase which I'm using to debug and confirm the mesh creation step is workin as planned.
Added a step which creates the mesh and sets the initial conditions, following the example of the previous compass version. Also added some experiment related parameters to the config file.
Fixed a problem where array dimensions were overwritten when the initial condition was set. Previously the only the dimensions present in the variables used to set the IC were broadcast, which means `Time` was almost always dropped as a dimension of the variables with initial condition. Now I use `xr.DataArray.loc[:] = ...` to set the initial conditions, which ensures the dimensions are broadcast properly.
Enforce the Y-axis to centered about the trough in the MISMIP+ domain. Whereas previously, the origin was placed in the lower left hand corner, I now ensure the Y-axis is centered about the bed trough. While this is fairly minor, I thik it will be usefull for easy visulaization. Also cleaned up the number of output files and their naming conventions to match more closely with out landice compass tests.
Renamed the coordinate translation function in order to reflect that we now center Y-axis about the bed trough. ALso added some more comments to explain in a little more detail exactly what and how we are doing the translation.
Create subdirecotry sturcutre based on the resolution as defined in the `.cfg` file at the time `compass setup` is run. Directory naming based on resolution changed after `compass setup` is run, is not currently supported.
Since the directory tree structure can not be changed at runtime, any changes to the resolution within the cfg file at runtime will produce a warning. The resolution as it was set at the time of `compass setup` wil be used, and wanring will be printed. If you want to change the mesh resolution, by editing the `.cfg` file, then you'll need to re-run `compass setup` for each resolution of interest.
Apply Trevor suggestion which fix code style, spelling, and clean up comments. Co-authored-by: Trevor Hillebrand <[email protected]>
This allows the resolution to be read from the config file and the `ntasks` to be dynamically calculated. While only one resolution is supported now, this provided the infrastructure to support multiple resolutions in the future.
Switched `clobber` modes to `overwrite` and increased the number of digits for the years postiion in filenames.
Updated the user guide, developers guide, and the API.
- `smoke_test` step was being passed `resoluion` as a string when it should have been a float. This has been fixed, and directory structure now use resolution value. - Error message raised by `setup_mesh` when value of resolution is changed at runtime had incorrect units (i.e. km). Switched to units of "m" to properly reflect the config file. - Switched the `clobber_mode` for the `output` stream from `truncate` to `overwrite` since the filename interval is equal to the simulation length, so there will only be one output file.
Spelling and grammatical errors caught by Matt. Co-authored-by: Matt Hoffman <[email protected]>
- Enforce the MISMIP+ rate factor value of 6.338e-25 - Ensure the sliding exponent is set to 1/3 per MISMIP+
- Add the density value from the config file to the namelist so that IC and simulation have a consistent density value - Density wasn't being parsed in setting of initial conditions This has been fixed so that config density is used consistently through spin_up test case.
Correct more typos caught by Matt. Co-authored-by: Matt Hoffman <[email protected]>
- ensured the `approx_cell_count` function returns a (correct) value and that it is called the the parent function. - Added claving timestep realted variables back the globalstats stream file per Matt's request.
Co-authored-by: Trevor Hillebrand <[email protected]>
- Updated namelist and streams with additional options Trevor suggested. - Updated comment about `init_thickness`, reorganized config a little, and fixed a couple typos in config. - Updated directory to look for smoke_test output in (for validation), based on the `work_dir` of `run_model` step. - Removed heuristic for resolution scaling from test group config. Now have it hard coded within the `tasks.py` file where its actually used. - Updated scaling of heuristic to not fail when no entry for `gutter_length` is found for the `SmokeTest`. The `try`/`except` solution is probably only temporary.
- reconfigured how/when `compass.step.subdir` attribute is set such that config options pass via `compass setup -f ...` are properly parsed. As part of this reconfiguration, I've created a common function for creating subdirectory structure, so that there more code sharing and consitent direcotry strucutres across steps and testcases (within this test group). - streamlined the way `ntasks` is calculated, which better handles the different config options present between the test cases. This should fix the error Matt got when he tried to run the `SmokeTest` with the version of code run during the review.
- Previously mixed camel-case and snake-case for variables names. Have tried to use consitent snake-case for all varibales, while keeping camel-case for classes.
- Removed the resolution option, and the configuration file as a whole for the `SmokeTest`. Instead supported resolutions are added as different test cases, with a different `subdir` for each resolution. Per Xylar, this is more standard practice for actual "test". (c.f. SpinUp which is more meant to set up files for longer run outside of `compass`, where arbitary resolution is fine. - Added comment about equilateral trinagles, per Trevor request. - Minor stylistic changes based on a limited review with Xylar.
Co-authored-by: Trevor Hillebrand <[email protected]>
Also switched the reference/start time to start at year one instead of year zero.
- Added comment about why resolution attribute of the `RunModel` step has a strange default value. - Make sure core/node info is parsed as integers but recast as floating point numbers, where floating point arithmetic is needed. - Changed the way intial conditions are set, so the friction coefficent from the MISMIP+ protocol is put in the `muFriction` field (c.f. `effectivePressure`) and now `effectivePresure` field is set to a sptially constant value of 1.0. - Updated the `albany_input.yaml` so that both Mu and effective pressure are "Fields", which will enable easier integration with coupled hydrology runs.
ebf6b1f
to
6d5391f
Compare
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
@andrewdnolan, amazing work! I have two tiny comments, but this is ready to go.
The minimum number of tasks | ||
""" | ||
|
||
# read MPI related parameters from configuration file |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Do we need to cast these as floats at all? It does no harm, but it will give exactly the same answer as using ints because integer division yields a float anyway.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
@andrewdnolan , thanks for all the revisions on this PR! I'm happy to merge it. Before I do that, can you update the PR description describing the final product? You're also welcome to make the small changes Trevor suggested if you like. Once you're happy with things, let me know and I will merge.
Co-authored-by: Trevor Hillebrand <[email protected]>
Was casting integer value as floats in calculation of `ntasks` and `min_task`, but as Trevor pointed out because we are using the floating point operator for division (`/`), this isn't really necessary. Removed the float casting to make the code a bit more legible.
@trhille @matthewhoffman Thanks for all the help reviewing / work shopping this PR! I've updated the PR description to reflect the final product. I've also implemented those two requested changes from Trevor. So, I'm happy with where things are if you guys are. |
Description:
Currently the only test case within the MISMIP+ test group is a
SmokeTest
, which uses an archived "spun up" initial condition. This PR adds aSpinUp
test case that generates the initial conditions and other input files (e.g. namelist and streams) needed to run MALI for a ~20,000 year spin up. TheSpinUp
testcase supports an arbitrary resolution, which needs to be specified in the config file at the time of runningcompass setup
. Furthermore, the user can specify agutterLength
in the config file to extend the eastern boundary of the MISMIP+ domain, which was added to support simulations using a dynamic claving law that will results in an irregular grounding line. To this end, we have added aSetupMesh
step, used by theSpinUp
testcase that creates a MALI mesh (and initial conditions) following the MISMIP+ requirements (Asay-Davis et al. 2016).We've also slightly altered the
SmokeTest
to explicit list the supported resolutions (e.g.compass list | grep mismipplus
). Currently this is only2000m
, but in the future theSpinUp
testcase will be used to generated the spun up conditions needed to support additional resolutions of theSmokeTest
.Checklist:
api.rst
) has any new or modified class, method and/or functions listedTesting
in this PR) any testing that was used to verify the changes