Add script to compute vertical coordinate for initial condition#495
Conversation
|
@vanroekel thanks for providing the script. I only altered formatting, but not content. |
|
I added a little plotting routine to the script. Parameters are in bottom right. |
mark-petersen
force-pushed
the
ocean/vertical_coordinate_script
branch
2 times, most recently
from
593322b to
36deeba
Compare
|
Here are the current settings for the 64-layer grid:
|
|
@mark-petersen: Andrew Roberts (@proteanplanet) has generated a 64 level profile. maybe using the same code, but you should make sure you guys are on the same page. |
|
@proteanplanet, please check the layer distribution above. Can you check against or post your 64-layer grid? Did you use the script from @vanroekel or something else? |
|
I generated the layers for @proteanplanet -- here are the layer thicknesses from that invocation of the script 1.98312923 2.22519441 2.49685395 2.80181737 3.1442291 |
|
sorry, accidentally clicked closed. I was going to say that I don't have the parameters I used anymore. but your graph looks close to what I post here. |
There was a problem hiding this comment.
@xylar, I've only tested the QU240 right now, but this file changes all of them. Is it best to call the python script as an command-line executable here? That's what I'm most familiar with.
@vanroekel thanks for the values. This new one looks extremely close: |
|
@vanroekel I'll run this through my scripts later today and post the results here and on Confluence. |
There was a problem hiding this comment.
These kinds of subprocess calls are not working well with the MPI versions of the COMPASS environment. It would be better to make this a function and call it that way. In general, calling python scripts from other python scripts with subprocess is not a good practice.
There was a problem hiding this comment.
@mark-petersen, do you want me to make these changes or do you want to do it?
There was a problem hiding this comment.
So it should look like this instead?
include make_vertical_grid
make_vertical_grid(arg1, arg2, ...)
There was a problem hiding this comment.
OK, that one is done.
There was a problem hiding this comment.
It should look like this:
#!/usr/bin/env python
from make_vertical_grid import create_vertical_grid
# Standard 64 layer vertical grid
create_vertical_grid(bottom_depth=5500, nz=64, dz1_in=2.0, dz2_in=200.0,
plot_vertical_grid=True, maxit=1000, epsilon=1e-2,
outfile='vertical_grid.nc')
|
OK, I've tested with QU240, QU240wISC, and EC60to30. I am confident the COMPASS a scripts are all working. |
testing_and_setup/compass/ocean/scripts/vertical_grid/make_vertical_grid.py
Outdated
Show resolved
Hide resolved
testing_and_setup/compass/ocean/scripts/vertical_grid/make_vertical_grid.py
Outdated
Show resolved
Hide resolved
testing_and_setup/compass/ocean/scripts/vertical_grid/make_vertical_grid.py
Outdated
Show resolved
Hide resolved
testing_and_setup/compass/ocean/scripts/vertical_grid/define_vertical_grid_64_layers.py
Outdated
Show resolved
Hide resolved
testing_and_setup/compass/ocean/scripts/vertical_grid/define_vertical_grid_64_layers.py
Outdated
Show resolved
Hide resolved
|
@vanroekel - An assorted set of comparisons between 60-, 64- and 80-levels: |
There was a problem hiding this comment.
I think links like this are needed for a lot more test cases, right?
- CUSP12
- CUSP8
- QU60
- SO60to10wISC
- SOQU60to15 (is this even a configuration we want to keep?)
xylar
left a comment
xylar
left a comment
There was a problem hiding this comment.
@mark-petersen and @vanroekel, I've spent time today trying to understand the algorithm by which layer thicknesses are selected and I really can't understand it. I would appreciate going through it carefully with @vanroekel if we want to keep the algorithm as it is but I would suggest switching to something simpler and easier to control, where we just use the root finder from scipy. I will post a suggested alternative algorithm in a gist below.
|
BTW, I will still address your remaining comments above on this PR. |
|
Thought a phone call would be easier... My concern about merging this is whether it makes it difficult to (temporarily) go back to the 60-level mesh for the EC60to30wISC initial condition that's currently our highest priority. I guess it shouldn't be hard to do that, though, right? |
|
During our meeting today, @maltrud suggested setting the maximum depth to 6000 or 6500 m. Here is a comparison of 5500 (black) with 6500 (red).
The lower right panel is just a zoom in on the upper right. |
This is a more straight-forward root finding approach using scipy. It finds a characteristic lengh scale `delta` over which layer thickenesses vary from `dz1` a the surface to `dz2` at `z --> -infinity` such that the sum of the layer thicknesses hits the target bottom depth in the given number of layers. Docstrings and formatting have been improved.
remove link for config_initial_state.xml copy ocean/develop version of config_initial_state.xml in
mark-petersen
force-pushed
the
ocean/vertical_coordinate_script
branch
from
5039d7c to
80630e9
Compare
|
Rebased on ocean/develop. Old 60 layer case is still available, but default is now the 64 layer case. I changed the test QU240 to 16 layers for speed. |
| bottom_depth=6000, | ||
| min_layer_thickness=2., | ||
| max_layer_thickness=210., | ||
| plot_vertical_grid=True, |
There was a problem hiding this comment.
Based on the comments at the meeting today I set the 64-layer bottom depth to 6000m, varying from 2m to 210m.
|
Here is the auto-generated plot: |
|
EC60to30 with 64 levels. This is the new plot that is produced at the end of the initial_state step: |
Previously, the vertical coordinate (layer depths) was read from a file when generating initial conditions in COMPASS. This adds a script with input parameters and adds it to the COMPASS init sequence. This will be used to create a 64-layer domain for the EC60to30 and others right now.
This step inadvertently got replaced with the non-ice-shelf-cavity version in #495.
This change is needed because of a recent modification (MPAS-Dev#495) to `build_mesh.py`
…evelop Previously, the vertical coordinate (layer depths) was read from a file when generating initial conditions in COMPASS. This adds a script with input parameters and adds it to the COMPASS init sequence. This will be used to create a 64-layer domain for the EC60to30 and others right now.
This change is needed because of a recent modification (MPAS-Dev#495) to `build_mesh.py`
6 participants
Previously, the vertical coordinate (layer depths) was read from a file when generating initial conditions in COMPASS. This adds a script with input parameters and adds it to the COMPASS init sequence. This will be used to create a 64-layer domain for the EC60to30 and others right now.