Repository navigation
Read surface topography from Elmer-Xios output - #12
Conversation
BenjaminRodenberg
left a comment
There was a problem hiding this comment.
I merged main into this branch and cleaned up some formatting errors (will push in a moment). One part is still problematic (see comment below).
We introduced some checks in #13 and #14. These should help to have a uniform formatting and make sure all code is properly licensed. Please refer to https://github.com/wardvp/EBFM?tab=readme-ov-file#developer-notes for documentation.
|
I did another merge with Main change: The last commit before the merge was 29350ad (just in case the merge causes trouble and you need to go back). |
BenjaminRodenberg
left a comment
There was a problem hiding this comment.
I reviewed the design of this PR and provided some general comments.
I can also take care of the implementation of (some of) these parts if you want. Please just let me know in the comments.
|
|
||
|
|
||
| def main(OUTFILE, io, OUT, grid, t, time, C): | ||
| def main(OUTFILE, io, OUT, grid, t, time, C, gridtype): |
There was a problem hiding this comment.
Is it possible to use grid["input_type"] here? Introducing a new type is no problem in my opinion. Extending the interface is something I would only do if necessary.
There was a problem hiding this comment.
So this and everything that is related to changes in the save_netCDF_file(gridtype) function is something, we should probably discuss whether it belongs to this pull request. I started implementing writing, a netcdf output file, but did not finish the work. So the question now is whether, I scratch the code bits that belong to this and open a new PR once this is merged or whether I keep this "unfinished" work in there. I am thinking the former is the better option.
There was a problem hiding this comment.
I assume this is a problem originating from the general design of the code. The design is not very extensible yet.
Can you provide a link to your sketch? This helps me to better understand what you are trying to achieve with the writer and how we could update the design to make it a bot more extensible and avoid merge conflicts etc.
(maybe you noticed: I tried to follow a similar goal with #44 to decouple the configuration from the implementation in order to avoid merge conflicts and allow for easier parallel development)
| return True | ||
|
|
||
| def save_netCDF_file(): | ||
| def save_netCDF_file(gridtype): |
There was a problem hiding this comment.
This function is a bit problematic in my opinion due to the variables that are used from the outer scope (this makes the input/output hard to understand). Doing a big refactoring of the codebase is out of scope. But I think we could start here (since we are touching this function) by turning this function from a nested function (inside main) into an independent function on the root level of this file. I.e.:
Not
def main(...):
...
elif io["output_type"] == 2:
save_netCDF_file(gridtype=gridtype)
...
def save_netCDF_file(gridtype):
...
but
def save_netCDF_file(gridtype, io, ...):
...
def main(...):
...
elif io["output_type"] == 2:
save_netCDF_file(gridtype=gridtype, io=io, ...)
...
There was a problem hiding this comment.
After some refactoring this problem still exists but it's actually not really the scope of this PR (it's already a problem of the original code). I would ignore this problem for now and do the refactoring at some later point when it fits better.
@ClemensSchannwell feel free to resolve this discussion if you agree.
| IN["h"] = data_from_elmer["h"] | ||
| grid["z"] = IN["h"][0].ravel() | ||
| # TODO add gradient field later | ||
| # IN['h'] = data_from_elmer['h'] | ||
| # grid['z'] = IN['h'][0].ravel() |
There was a problem hiding this comment.
I assume we need to revert this in order to not break existing functionality?
Why did you have to comment it out in the first place? Is this something that is easy to fix or should we add a guard around this part? I.e.
if grid["grid_type"] == GridInputType.ELMERXIOS:
logger.warning(f"Grid type {grid["grid_type"]} does not support outputting IN["h"] and grid["z"]. Ignoring these fields.") # todo: find a fix if necessary
else:
# usual behaviorThere was a problem hiding this comment.
Yes, we can revert this. It does actually not break anything. Since I am passing 42 from Elmer, I simply wanted to see whether it would do something more sensible if I did not update the topography (h in this case). Generally, a little further down the line it might actually not bad to also allow a flag --with-topography-update or similar.
| grid["slope_gamma"][(grid["slope_x"] > 0) & (grid["slope_y"] == 0)] = np.pi / 2 | ||
| grid["slope_gamma"][(grid["slope_x"] < 0) & (grid["slope_y"] == 0)] = -np.pi / 2 | ||
| grid["slope_gamma"] = -grid["slope_gamma"] | ||
| elif grid["input_type"] is GridInputType.ELMER_XIOS: |
There was a problem hiding this comment.
Does this work?
It's generally a good idea to make sure that at least one of the cases in this chain of if-elif-statements is entered. I always append a final else: raise Exception(f"Cannot handle unexpected {grid["input_type"]=}"). This will make sure that the program exits.
| io, OUTFILE = LOOP_write_to_file.main(OUTFILE, io, OUT, grid, t, time2, C) | ||
| pass | ||
| if grid["input_type"] is GridInputType.MATLAB: | ||
| gridtype = "structured" |
There was a problem hiding this comment.
To make sure: gridtype = "unstructured" is not tested/supported yet? Or do you manually set gridtype = "unstructured" during development here?
| grid["slope_gamma"][(grid["slope_x"] < 0) & (grid["slope_y"] == 0)] = -np.pi / 2 | ||
| grid["slope_gamma"] = -grid["slope_gamma"] | ||
| elif grid["input_type"] is GridInputType.ELMER_XIOS: | ||
| grid = read_elmer_xios_grid(grid=grid, gridfile=args.elmer_xios_mesh) |
There was a problem hiding this comment.
| grid = read_elmer_xios_grid(grid=grid, gridfile=args.elmer_xios_mesh) | |
| grid = read_elmer_xios_grid(grid=grid, gridfile=args.netcdf_mesh_unstructured) |
args.elmer_xios_mesh does not exist. Can you clarify?
|
I performed another merge of |
BenjaminRodenberg
left a comment
There was a problem hiding this comment.
Some minor comments on the design. I can take care of this as soon as you are done. Please let me know when I can take over here again.
| help="End time of the simulation in format 'DD-Mon-YYYY HH:MM', e.g., '02-Jan-1979 09:00'.", | ||
| default="2-Jan-1979 09:00", | ||
| help="End time of the simulation in format 'DD-Mon-YYYY HH:MM'", | ||
| default="2-Jan-1979 00:00", |
There was a problem hiding this comment.
I did another merge of main into this PR. The default values are quite arbitrary. So I'm fine to just change the default to the value you are using if this saves you some typing.
You can also use the new options to provide different values.
b8e6e52 to
db2c652
Compare
BenjaminRodenberg
left a comment
There was a problem hiding this comment.
Looks good to me.
|
Yes, I have tested the new version this morning and see no regression. So also from my side it looks good to go! |
In order to read the surface topography from an actual Elmer/Ice simulation, I included the option to
--netcdf-mesh-unstructured. Since this surface topography uses the libraryXIOSto write the file, it lives on the same grid as Elmer. This simplifies the reading of the file as the surface topography comes in the same order as the Elmer mesh.I think this feature will be helpful in the case of initialisation, because we can read the same surface topography in with which Elmer is initialised and also if we decide that we do not need to update the surface topography within a EBFM run, we can simply read it it at the start of the next simulation.
I also wanted to write some model output for easy inspection using
ncviewor similar, but while it works in principle, this is still in its early stage. This does not hold up further development of theEBFM-dummybut would be good to have to check whether results from EBFM are physically plausible.I made some test with the
feature_greenlandbranch and all previous documented running commands of the model still work as before (as far as I can tell at least).To test the new model on levante, you can run:
The
elmer-xiosfile should be accessible for all project members on levante.