Introduce shelf-sea water residence time tracer acc. to Liu et al 2019 #426

jmaerz · 2024-10-28T18:37:09Z

Hi @JorgSchwinger and @TomasTorsvik ,
I introduced a shelf-sea-water residence time tracer according to Liu et al. 2019: Simulating Water Residence Time in the Coastal Ocean: A Global Perspective. Since there might be a broader interest in using such a tracer, I thought it would be good to have that included in iHAMOCC in an optional way. It can be switched on via the xml-switch: HAMOCC_SHELFSEA_RES_TIME. The shelf-sea mask can be either set via a file (1. where ocean shelves, 0. else, variable name shelfmask) or the model creates the mask internally as a fall-back solution (current default 200m isobath).

@matsbn , I include you in the loop since it may be of interest for evaluating the ocean physics as well.

…. 2019

jmaerz · 2024-10-28T21:04:52Z

I currently still get some strange pattern. While the value range seem to be ok-ish, the pattern is off - not linked to bathymetry. Will look into it with fresh eyes tmw. If you've an immediate idea, what goes wrong, let me know. The run was performed using cgs and the CMIP6 settings as described in https://github.com/NorESMhub/BLOM/wiki/new-BLOM-with-CMIP6-settings

Outcome at some depth	mask file for shelf sea regions

jmaerz · 2024-10-29T12:27:10Z

@TomasTorsvik and @JorgSchwinger , I tested a bit further and fixed some issues, with no avail. I manipulated the mask file for debugging purposes (introduced a big square, since I wasn't sure if the mask file reading introduced some turns, etc. of the mask). This resulted in a strange pattern (I would have expected a square with some blurred edges as result). I will further check tmw. My suspicion is some not properly initialized field, since the switch use_shelfsea_res_time seem to be overridden with time.

Outcome at some depth	mask file for shelf sea regions

…metry in case maskfile file is missing

jmaerz · 2024-10-30T19:31:31Z

Hi @JorgSchwinger , after implementing a fall-back initialization of the mask initialization via internal bathymetry values, it seems that indeed the mask reading isn't functioning properly. Since the reading of the input file follows the implementation of ocean alkalinization, my suspicion is that the format of the mask file isn't right, or that my settings are wrong for reading the file: could you provide the ncdump -h for one of your ocean alkalinization files? - potentially just via email?

One month calculation with fall-back solution - 200m bathymetry line:

jmaerz · 2024-10-31T14:20:18Z

Hi @TomasTorsvik and @JorgSchwinger , this PR is now ready. I eventually got it working in both ways - either via input file or via internally calculated mask (200m isobath generated from BLOM variable depths). Using the file as input can be useful to e.g. define the Antarctic shelf different from the rest of the global shelf regions (i.e. 400m isobath).

This time with using the shelf-sea mask file:

TomasTorsvik

@jmaerz - thanks, I think in principle this is fine.
Will the masking routine work if the ocean grid file is changed? The 200m depth limit should work, but I can't judge about the input file. I assume this is the case (hence I have approved the PR). If not, there should be some testing for the grid configuration as well, presumably before reading the mask file.

TomasTorsvik · 2024-11-01T13:21:56Z

hamocc/mo_read_shelfmask.F90

+    integer :: i,j,errstat,ncid,ncstat
+    real,allocatable  :: mask(:,:)
+
+    ! Check, if we are goiung to run with shelf-sea water residence time tracers


Suggested change

! Check, if we are goiung to run with shelf-sea water residence time tracers

! Check, if we are going to run with shelf-sea water residence time tracers

fix typo

TomasTorsvik · 2024-11-01T13:42:38Z

hamocc/mo_read_shelfmask.F90

+
+  public :: read_shelfmask
+
+  character(len=512),               public :: shelfsea_maskfile =''


Suggested change

character(len=512), public :: shelfsea_maskfile =''

character(len=:), allocatable, public :: shelfsea_maskfile

Is it possible to use the Fortran 2003 convention with allocatable character array here? I'm not sure how this works when reading from namelist file. I think it's better if we don't restrict string lengths.
Just a suggestion, leave it out if inconvenient to implement.

I tried this now locally and then the namelist reading fails to correctly initialize the file name (it silently ignores the entry and leaves the file name as empty - and just runs further). Maybe we should introduce a string-length parameter to be able to set those things at one place... - up for later.

OK, no problem, just leave it as it is.

@TomasTorsvik , not sure, if this can be extended to a full namelist: https://stackoverflow.com/questions/14886787/fortran-character-input-at-undefined-length - for now, I keep the old style.

jmaerz · 2024-11-01T14:32:04Z

Hi @TomasTorsvik : what is your concern with regards to the maskfile reading? - I fear, I simply went for the wrong netcdf-reading routine initially, which caused the trouble... For myself, I prepared two files, tnx2v1 and tnx1v4 - I can set up a test run for tnx1v4 (thus far tested with tnx2v1).

TomasTorsvik · 2024-11-01T14:59:55Z

Hi @TomasTorsvik : what is your concern with regards to the maskfile reading? - I fear, I simply went for the wrong netcdf-reading routine initially, which caused the trouble... For myself, I prepared two files, tnx2v1 and tnx1v4 - I can set up a test run for tnx1v4 (thus far tested with tnx2v1).

I'm guessing the mask file that the user provides should correspond to a specific grid, is that right?
If so, is there some way to do a sanity check on the input mask file? Or to put it in another way, what happens for instance if you provide a tnx1v4 maskfile for a tnx2v1 run?

jmaerz · 2024-11-01T15:05:01Z

@TomasTorsvik : Correct - as it is for e.g. oafx. Not sure, if it is a real issue for those who will switch on the optional tracer.

JorgSchwinger · 2024-11-01T16:35:25Z

In general, I'm not sure I understand why a file is needed. If the definition of shelf-sea is just based on the bathymetry on the model grid, this can be easily done online without the need to read a file. I understand that you would like to have a different definition for Antarctica, but that is also very easy to implement based on model bathymetry (where lat < -60 ... elsewhere ...). Unless there are not more complicated definitions to implement, I'm not sure I see the need to implement reading of a file?

jmaerz · 2024-11-01T16:40:28Z

I guess, in my opinion, I prefer to tinker around with the file than inside the code for such use cases... - but if you want me to delete it, fine - I think, the reasoning of Tomas would also apply to e.g. the oafx use case (and other examples), if I am not wrong.

JorgSchwinger · 2024-11-01T16:53:53Z

Mhm, I was assuming here that we would want to use a fixed definition for "shelf-sea", whereas for the OAFX, you can have many different cases (e.g. OAE is done in the EEZs of the US, EU, and China (all could have different rates), or alternatively, is done in the form of electrochemical weathering at the location of desalination plants, or whatever)

JorgSchwinger · 2024-11-01T16:56:30Z

... and yes, otherwise the reasoning of Tomas also applies to the oafx-files - a user could do stupid things here (but we assume that someone using this is an expert user)

jmaerz · 2024-11-01T16:58:09Z

So you would implement other switches for Antarctic and potentially other shelf regions (e.g. South America) to enable to investigate those separately? I feel this as easier with an input file.... But if the reading of a file is of real concern, I delete the section (while I am not going to implement the Antarctic shelf region hard-coded now - this is much easier with one cdo-command with a file). Well, the shelf tracer will/should be optional - so expert user as well ;-)

jmaerz · 2024-11-01T17:02:26Z

Maybe we need to invest time in such things - I see the trouble not only for this mask - but also for other readings (while with nuopc, this phases out now). If there is a template how to do it, it's easier to apply.

JorgSchwinger · 2024-11-01T17:04:20Z

I think my concern is that this would be another file that could not be automatically generated for a new grid-resolution, could it? So in view of Mariana's NUOPC-magic, I would think it is good to introduce such files only if really necessary (and putting a few if-blocks in the code that generates a shelf-mask based on model bathymetry seems to me quite doable)

JorgSchwinger · 2024-11-01T17:06:12Z

So you would implement other switches for Antarctic and potentially other shelf regions

I'm not sure I understand what you mean. The definition of "shelf" is fixed, isn't it?

jmaerz · 2024-11-01T17:15:29Z

I can easily provide such script.... - (in fact, I have it already - it's a one-liner for constant 200m isobath). Well, I perceive the definition kind of fixed, while as written, the Antarctic Shelf is usually differently suggested. And I don't know about the Patagonian shelf region, potentially such deep narrow deep trenches like the Norwegian trench could also be of interest for some users, etc. Having a clear control and hand on an input file (easily by eye) feels safer than having to hard-code such grid-dependent features... (I would end up doing things first in cdo and then implement it in iHAMOCC :-D)

jmaerz · 2024-11-01T17:22:54Z

Technically, one could use the mask e.g. also whole basins, etc. - like the Arctic Ocean - I see a number of use cases, where I would want to use a mask different from the 200/400m isobath.

JorgSchwinger · 2024-11-01T17:26:52Z

Mhm, my perception would be that a hard coded definition would be the easiest and safest option... Particularly with regards to the upcoming NUOPC based interpolation to different grids. But never mind, I am also fine with a file-solution

jmaerz · 2024-11-01T17:35:15Z

Well, the bathymetry grid hopefully doesn't change that frequent - will it be really read by nuopc? I mean currently, there is a default option: if there is no file provided by the user, the 'classical' 200m isobar is used (so no need for the file reading - maybe one should just reformulate the log-message - hinting the user that there is a potential to use a file?). Sorry, I guess, I was bit perplexed about the discussion...

jmaerz · 2024-11-04T10:27:09Z

Hi @JorgSchwinger and @TomasTorsvik , just to summarize the status again on this PR:
The current default setting is that the code isn't requiring an input file for the shelf-sea mask. It then produces a mask internally for the 200m isobath, based on the depths provided through BLOM. This is fully in-line with the publication by Liu et al. 2019, who also used the 200m isobath globally (while they didn't look at the Antarctic shelf, though).
With respect to reading the file: I understand that the dimensions can be an issue - and could be worth checking at some point (not sure, how applicable it is, but the pio-package seem to provide this functionality - used e.g. in MOSART/CESM through shared libraries - likely requires adding a few lines in BLOM as well to provide this information for the bathymetry grid). However, I would assume that the usage of a shelf-sea mask file is for 'advanced' applications of the shelf-sea residence time tracers. The latter can go beyond pure shelf sea regions (e.g. including the 400m isobath for the Antarctic shelf region, considering the age of whole basin regions, etc.). Having the capability to read a mask file makes the water age tracers bound to regions more flexible to my understanding (while again, the default setup without reading a mask file is completely scientifically meaning- and useful).
If there is a deeper interest to study the Antarctic shelf residence time, I feel the user can be expected to either set up a shelf mask properly or to include it hard-wired and less flexible in the code. In either case, the right usage of the mask can be checked within a time step of model run time - hence not too computationally costly for setting up a proper mask (while I feel setting up a mask for dedicated regions is easier and more flexible via an input file).

In summary, I would currently go with the code as is - rewriting the info message so that the mask reading is more regarded as further option (and not regarded as default).

TomasTorsvik · 2024-11-04T10:56:08Z

@jmaerz - I think this is fine, at least with the default depth limiter.
In principle, the mask doesn't need to be shelf seas, it could also be e.g. a sub-basin (for instance, a mask of the Arctic Ocean), if you introduce it through a mask file. At least, I don't see why the mask would have to be confined to shelf seas with this implementation. Optionally, this could be implemented as a more general "region residence time" with a default to shelf seas. I'm not sure how much added value there would be in generalizing this routine, but maybe something to think about now that we introduce non-global tracers.

jmaerz · 2024-11-04T12:10:35Z

Hi @TomasTorsvik , yes, indeed, it could be the Arctic Ocean or any other region - which makes running it with the mask more flexible. I was debating, whether one should rename the tracer (region age - reg_age?), but for now remained with the present name (also due to the default use as a shelf-sea water residence time tracer). I am not entirely sure, what you mean with 'generalizing' the routine further.

TomasTorsvik · 2024-11-04T12:20:56Z

@jmaerz - I was only thinking that once you start with sub-regions, there are multiple options/configurations one could imagine, e.g. regional limits, depth range limits come to mind. But I think we don't need to implement a more general method unless there is a need.

jmaerz · 2024-11-04T12:26:30Z

Mmmh - except for partial parts of the water column, all this is already included via the capability to set a mask, I would say?! Of course, one is currently limited to one tracer, but that's a totally different question again and could be elaborated, if really needed. Thus far, my use case is for better understanding BLOM capabilities in representing shelves with a priority on understanding riverine impacts.

JorgSchwinger

Sorry that it took quite a bit of time... Looks good to me, just a few minor comments

JorgSchwinger · 2024-11-06T10:23:55Z

hamocc/mo_ncout_hamocc.F90

@@ -792,6 +803,9 @@ subroutine ncwrt_bgc(iogrp)
      call wrtlvl(jlvlprefalk(iogrp),  LVL_PREFALK(iogrp),  rnacc*1e3,      0.,cmpflg,'p_talklvl')
      call wrtlvl(jlvlprefdic(iogrp),  LVL_PREFDIC(iogrp),  rnacc*1e3,      0.,cmpflg,'p_diclvl')
    endif
+    if (use_shelfsea_res_time) then
+      call wrtlvl(jlvlshelfage(iogrp),  LVL_SHELFAGE(iogrp),  rnacc,      0.,cmpflg,'shelfage_lvl')


I think we should stick to the convention that level-output is named <variable_name>lvl (without underscore), where <variable_name> is the name of the layer-variable. My plotting stuff relies on this convention.

JorgSchwinger · 2024-11-06T10:26:00Z

hamocc/mo_ncout_hamocc.F90

@@ -1892,6 +1916,10 @@ subroutine hamoccvardef(iogrp,timeunits,calendar,cmpflg)
      call ncdefvar3d(LVL_PREFDIC(iogrp),cmpflg,'p',                              &
           &   'p_diclvl','Preformed DIC',' ','mol C m-3',2)
    endif
+    if (use_shelfsea_res_time) then
+      call ncdefvar3d(LVL_SHELFAGE(iogrp),cmpflg,'p',                           &
+           &   'shelfage_lvl','Residence time of shelf water',' ','d',2)


See my comment above, I think this should be 'shelfagelvl'

JorgSchwinger · 2024-11-06T10:38:25Z

hamocc/mo_read_shelfmask.F90

+      return
+    endif
+
+    ! Check if shelfsea mask file exists. If not, abort.


Remove the "If not, abort" since this is not done.

JorgSchwinger · 2024-11-06T10:40:14Z

hamocc/mo_shelfsea_restime.F90

@@ -0,0 +1,73 @@
+! Copyright (C) 2021-2022  J. Schwinger


Update this?

You mean to 2021-2024? (- or do you mean your name?)

I did now, but maybe we can have a short conversation on this at some point - thus far, I never updated those entries in any of the code, where I was working on in iHAMOCC (at least not in what was present) - I see that in BLOM, this is somewhat handled differently. I guess, I to this date see the PRs as main documentation of this (and not being sure, how the copyright holding note should be handled as long as it is one codebase). Let's have a chat together, when all being back.

I also always forget to update this for changes made to existing files. For new files, I think I have usually updated it. Here I just found it weird to have my name in a routine that was written by someone else... I'm happy to have a chat

JorgSchwinger · 2024-11-06T10:48:02Z

hamocc/mo_hamocc_init.F90

@@ -203,6 +204,7 @@ subroutine hamocc_init(read_rest,rstfnm_hamocc)
       call ini_read_ndep(idm,jdm)
    end if
    call ini_read_rivin(idm,jdm,omask)
+    call read_shelfmask(idm,jdm,nbdy,depths,omask)


I would prefer the name ini_shelfmask, since this is what the routine does (yes, it might read a file for doing this, but not necessarily) and for consistency with the other names.

renamed to ini_read_shelfmask

jmaerz · 2024-11-06T14:02:53Z

Ok, I guess, I push this now - and whenever there are different needs/wishes, we can get back to it - for now, it totally satisfies a number of needs.

Introduce shelfsea water residence time tracer according to Liu et al…

16c901f

…. 2019

jmaerz added enhancement New feature or request iHAMOCC Issue mainly concerns the iHAMOCC code base labels Oct 28, 2024

jmaerz requested review from JorgSchwinger and TomasTorsvik October 28, 2024 18:37

jmaerz self-assigned this Oct 28, 2024

jmaerz added 4 commits October 28, 2024 20:04

fix namelist

2ad8ae7

fix buildnml

37152dc

fix buildnml

3e13739

fix buildnml

5fee9ae

jmaerz changed the title ~~Introduce shelfsea water residence time tracer acc. to Liu et al 2019~~ Introduce shelf-sea water residence time tracer acc. to Liu et al 2019 Oct 28, 2024

fix declaration, add switch to output

7cbe225

jmaerz added 2 commits October 28, 2024 23:48

fix ocetra index

d40a1c1

fix strange entry

076e568

jmaerz added 2 commits October 30, 2024 17:38

Introduce fall-back shelf-sea mask reconstruction from internal bathy…

fd2296d

…metry in case maskfile file is missing

remove allocation

dc7b2ad

jmaerz added 5 commits October 31, 2024 09:38

merge with master

4aafd09

revert to allocation of mask field, add description and header

4deebae

merge with master

e950b24

change netcdf read for shelf sea mask

46e706c

cleanup

264db36

TomasTorsvik approved these changes Nov 1, 2024

View reviewed changes

fix typo

f25b1a8

NorESMhub deleted a comment from tomas Nov 3, 2024

JorgSchwinger approved these changes Nov 6, 2024

View reviewed changes

jmaerz added 3 commits November 6, 2024 14:28

minor changes in response to review

e242ce9

Update copyright date

d3f10e4

Update copyright data

8e0a4d4

jmaerz merged commit 79081d8 into NorESMhub:master Nov 6, 2024
5 checks passed

jmaerz deleted the shelf-residence-time branch November 6, 2024 14:04

	! Check, if we are goiung to run with shelf-sea water residence time tracers
	! Check, if we are going to run with shelf-sea water residence time tracers


		public :: read_shelfmask

		character(len=512), public :: shelfsea_maskfile =''

	character(len=512), public :: shelfsea_maskfile =''
	character(len=:), allocatable, public :: shelfsea_maskfile

Introduce shelf-sea water residence time tracer acc. to Liu et al 2019 #426

Introduce shelf-sea water residence time tracer acc. to Liu et al 2019 #426

Conversation

jmaerz commented Oct 28, 2024 • edited Loading

jmaerz commented Oct 28, 2024

jmaerz commented Oct 29, 2024

jmaerz commented Oct 30, 2024

jmaerz commented Oct 31, 2024

TomasTorsvik left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

jmaerz commented Nov 1, 2024 • edited Loading

TomasTorsvik commented Nov 1, 2024 • edited Loading

jmaerz commented Nov 1, 2024

JorgSchwinger commented Nov 1, 2024

jmaerz commented Nov 1, 2024

JorgSchwinger commented Nov 1, 2024

JorgSchwinger commented Nov 1, 2024

jmaerz commented Nov 1, 2024

jmaerz commented Nov 1, 2024

JorgSchwinger commented Nov 1, 2024

JorgSchwinger commented Nov 1, 2024

jmaerz commented Nov 1, 2024

jmaerz commented Nov 1, 2024

JorgSchwinger commented Nov 1, 2024

jmaerz commented Nov 1, 2024

jmaerz commented Nov 4, 2024

TomasTorsvik commented Nov 4, 2024

jmaerz commented Nov 4, 2024

TomasTorsvik commented Nov 4, 2024

jmaerz commented Nov 4, 2024

JorgSchwinger left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

jmaerz commented Nov 6, 2024

jmaerz commented Oct 28, 2024 •

edited

Loading

jmaerz commented Nov 1, 2024 •

edited

Loading

TomasTorsvik commented Nov 1, 2024 •

edited

Loading