-
Notifications
You must be signed in to change notification settings - Fork 25
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
Merge extended nitrogen cycle, M4AGO scheme and preformed Si tracer into master #269
Conversation
Introducing a static 3D sediment porosity field that can be optionally read in with effects on molecular pore water diffusion and shifting.
…way split function in nitrification
Bug fixes for denitrification stoichiometry and C-isotopes in the sediment
…fixes as well as better oxygen limitation in ocprod
@jmaerz - Made a quick attempt at regression test from your branch, building on noresm2.1.1. I get a failure at build time:
Need to look further into this tomorrow, but thought it best to notify you right away. |
Briefly discussed this with @TomasTorsvik - setting up a NorESM2.1 requires to set the NorESM
(and then use
|
@jmaerz - I ran the regression test again today. Now it passes through with comparison to
|
I'm looking at the way iHAMOCC is set to interact with M4AGO. I see the |
Before merging this PR, we should have a look at how this will impact noresm2.5 development. @mvertens is preparing a code sandbox for running regression tests based on the development branches. |
Hi @TomasTorsvik , I am open for discussions on this - and how to improve the interface. I perceive the |
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.
@jmaerz - thanks for the great contribution! I have looked through the files for a technical review. I don't see any red flags, so I approve it now. However, I suggest to put the PR merge on hold until we have had a look at the impact for noresm2.5 development.
hamocc/mo_hamocc_init.F90
Outdated
@@ -63,6 +64,10 @@ subroutine hamocc_init(read_rest,rstfnm_hamocc) | |||
sedlay2,powtra2,burial2,blom2hamocc,atm2 | |||
use mo_ini_fields, only: ini_fields_ocean,ini_fields_atm | |||
use mo_aufr_bgc, only: aufr_bgc | |||
use mo_extNsediment,only: alloc_mem_extNsediment_diag | |||
! use mo_m4ago, only: init_m4ago_nml_params, init_m4ago_params |
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.
Remove comment to unused module
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.
Done.
</entry> | ||
|
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.
trailing whitespace
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.
Done (and some further places).
|
||
<entry id="flx_ndepnoy"> |
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.
leading whitespace
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.
Done.
cime_config/ocn_in.readme
Outdated
@@ -372,6 +374,8 @@ | |||
! OXY, NO3, SIL, D13C, and D14C | |||
! WITH_DMSPH : Logical switch to activate DMS calculation as function of pH | |||
! PI_PH_FILE : File name (incl. full path) for surface PI pH input data. | |||
! LM4AGO : Switch for M4AGO settling scheme |
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 would prefer the switch to be named "with_m4ago" instead of "lm4ago", as a more descriptive variable name.
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.
@TomasTorsvik and @JorgSchwinger , wrt naming the switches: we currently don't have a real naming convention (when looking into mo_control_bgc
) which could be changed (on the longer term). I am happy to change the switch name for M4AGO - but maybe better use_m4ago
(as complementary to use_wlin
)?!
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.
sure, that's fine for me.
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 would also vote for use_m4ago
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.
Changed switch name from lm4ago
to use_M4AGO
.
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 a lot for all the work that has gone into this major upgrade of HAMOCC! I have left a couple of comments - nothing really critical. Some things I might not have understood correctly. Most of these things can also be considered later (e.g. I would suggest to abandon the reading of old ndep files).
cesm/mod_cesm.F90
Outdated
|
||
! If SSS restoring is requested, read climatological sea surface salinity. | ||
! If SSS restoring is requested, read climatological sea surface salinity. | ||
if (srxday > 0._r8) then | ||
call rdcsss | ||
end if |
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.
Looks to me these two changes in indentation should not be there?
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.
Came in before through master
- fixed.
cime_config/ocn_in.readme
Outdated
@@ -372,6 +374,8 @@ | |||
! OXY, NO3, SIL, D13C, and D14C | |||
! WITH_DMSPH : Logical switch to activate DMS calculation as function of pH | |||
! PI_PH_FILE : File name (incl. full path) for surface PI pH input data. | |||
! LM4AGO : Switch for M4AGO settling scheme |
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 would also vote for use_m4ago
cime_config/ocn_in.readme
Outdated
@@ -363,7 +363,9 @@ | |||
! DO_RIVINPT : Logical switch to activate riverine input | |||
! RIVINFILE : File name (incl. full path) for riverine input data | |||
! DO_NDEP : Logical switch to activate N-deposition | |||
! DO_NDEP_COUPLED: Logical to apply N-deposition fluxes received from the atmosphere (true=atm, false=clim) |
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.
Should be clarified: is DO_NDEP=false, DO_NDEP_COUPLED=true a valid setting?
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.
Nope. Clarified now in the description. Further clarity when the code restructuring happened for the atmosphere fluxes.
@@ -186,7 +225,7 @@ subroutine accfields(kpie,kpje,kpke,pdlxp,pdlyp,pddpo,omask) | |||
bgct2d(i,j,jpodisi) = bgct2d(i,j,jpodisi) + sedfluxo(i,j,ipowasi)/2.0 | |||
endif | |||
! N-deposition, ocean alkalinization, and riverine input fluxes | |||
bgct2d(i,j,jndep) = bgct2d(i,j,jndep) + ndepflx(i,j)/2.0 | |||
bgct2d(i,j,jndepnoy) = bgct2d(i,j,jndepnoy) + ndepnoyflx(i,j)/2.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.
Just to understand this: If climatological n-deposition is used, the n-deposition from file is assigned to be noy (as a pragmatic solution). The data for NOX and NOY deposition is actually available, just not carried through into the inputfile for HAMOCC. I would suggest to do this (as a future task) then we would always have NOX and NOY deposition (and lump it together when use_extNcycle=false).
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 am not entirely sure, if I understand your comment. In any case, NOy is deposited. Only when the extended N cycle is switched on, NHx is also deposited as individual tracer (see a few lines above, l. 211). Currently, in the default version, ndepnoy=NHx+NOy
(lumped together) while in the extended nitrogen cycle, ndepnoy=NOy
and ndepnhx=NHx
. But I'll see, if I will change the reading of the input file - depends a bit on whether I have write access to the input data...
if (use_extNcycle) then | ||
ndepnhxflx(i,j) = ndep(i,j,idepnhx)*dtb/365. | ||
ocetra(i,j,1,ianh4) = ocetra(i,j,1,ianh4) + ndepnhxflx(i,j)/pddpo(i,j,1) | ||
ocetra(i,j,1,ialkali) = ocetra(i,j,1,ialkali) + ndepnhxflx(i,j)/pddpo(i,j,1) | ||
endif | ||
endif | ||
enddo |
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.
We should have this as a task to make use of the CMIP input4MIPs NOX/NOY fields instead of lumping them together. Would simplify things(?) Maybe something for release 2.5?
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.
Forget about this comment I have seen (and remembered) that you have already produced new input files. I would suggest to abandon the reading of old files (see below)
atn2ov = atm(i,j,iatmn2o) | ||
else | ||
atn2ov = atm_n2o | ||
endif |
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.
One could follow the treatment of o2 and n2 also for n2o, by writing the n2o field into atm(i,j,iatmn2o) in ini_fields_atm.
It's of course an unnecessary overhead in terms of having a 2d array instead of a scalar, but that's what we were doing for o2 and n2, too.
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.
Done.
if (omask(i,j).gt.0.5) then | ||
do k=1,kmle(i,j) | ||
if (omask(i,j) > 0.5) then | ||
do k=1,merge(kwrbioz(i,j),kmle(i,j),leuphotic_cya) ! if leuphotic_cya=.true., do bluefix only in euphotic zone |
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 really want to keep the option to do n-fixation in the mixed layer? I would be fine to replace the old mixed layer treatment with the euphotic zone treatment.
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.
It was thus far a pragmatic choice to enable investigating its effect when switching between isopycnal and hybrid coordinates. Maybe we can/should keep it until the hybrid coordinates become default. We can make the euphotic zone default, though. Any thoughts @TomasTorsvik ?
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.
NorESM2.3 is not ready for scientific experiments with hybrid coordinates, so I think it makes sense to keep the mixed layer formulation here, and aim for a replacement in NorESM2.5.
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.
@JorgSchwinger , for now, I would keep the option and make euphotic zone only for the cyanobacteria (bluefix) the default (which makes scientifically more sense - probably more in line with the aim of NorESM2.3 making it a scientifically used release). I'll do this as last step after the regression testing since it is answer changing.
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 don't think the vertical grid is a strong argument. I wouldn't expect a big difference anyway: Mixed layers deeper than 100m are only found when it is cold enough and then there is anyway not much N-fixation going on. In very warm waters, the use of euphotic zone might increase N-fixation, but given the low PP we have in these waters that would actually not be so bad.
In summary, I still don't see a good reason to keep the old "mixed layer based N-fixation", but up to you of course.
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.
@JorgSchwinger , I gained the impression that it would be worth keeping it until one can (again) compare the hybrid coordinates with the isopycnic coordinates (at least from discussions with @TomasTorsvik, since he wanted to test this). I don't have strong feelings about it, but it also doesn't cost much to leave it in or taking it out. We can put this also into the #340 issue for iHAMOCC towards NorESM2.5 to collect changes to be (potentially) done - to not forget about the to do list. @TomasTorsvik: any further opinion on this switch?
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.
Ok, good plan.
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 this task to #340.
hamocc/mo_hamocc4bcm.F90
Outdated
real, intent(in) :: patmnh3(1-kbnd:kpie+kbnd,1-kbnd:kpje+kbnd) ! atmospheric ammonia concentration [ppt] used in fully coupled mode | ||
real, intent(out) :: pflxnh3(1-kbnd:kpie+kbnd,1-kbnd:kpje+kbnd) ! Ammonia flux [kg NH3 m-2 s-1]. | ||
real, intent(in) :: patmnhxdep(1-kbnd:kpie+kbnd,1-kbnd:kpje+kbnd) ! Atmospheric NHx deposition [kgN m-2 s-1] | ||
real, intent(in) :: patmnoydep(1-kbnd:kpie+kbnd,1-kbnd:kpje+kbnd) ! Atmospheric NOy deposition [kgN m-2 s-1] |
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 would somehow find it less confusing if ndep
is also used in the case that the deposition comes from the atmosphere. In this case the code below copying the interactive deposition would be in blom2hamocc
(and read_ndep
would return if do_ndep_coupled=.true.
). Or is there a specific reason why it has been done like 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.
I can only follow partially here. I was mainly following the bromoform implementation (which might not have been the best choice in this particular case). Indeed ndep
is used (the fields patmnhxdep
and patmnoydep
are/were written to ndep
in hamocc4bcm
, before ndep
is applied). In case that do_ndep_coupled=true
get_ndep
wasn't entered at all (through if statement in mo_hamocc_step
- I felt that clearer than having to open the file mo_read_ndep
just to see that it immediately returns - but it's a matter of taste). However, I do see the point that writing the atmosphere fluxes to ndep
could be closer to where ndep
is filled in the case that input files are used. Hence, the new interpretation of the module name mo_read_ndep
now also includes 'module reading from file or stream' and I put the writing of the atmosphere-delivered fluxes to ndep
into get_ndep
to not clutter mo_hamocc_step
with this information. Note that I left the gaseous fluxes (N2O and NH3) in hamocc4bcm
since all other atmosphere gaseous fluxes are treated the same way. Before pushing, I'll do a coupled test run tmw.
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.
My point was just that I would suggest to avoid two (three) different variables for the same thing. Core HAMOCC does not care wether the N-deposition is read from file or coming from the atmosphere...
I would have put the copying of atmosphere n-deposition into ndep
into blom2hamocc
(since it is a field that is passed through BLOM, but I'm also fine with your solution.
Regarding whether the read_*
routines are not entered at all or return immediately, I would vote for a consistent treatment (this is the same problem for all of these routines where an input can be switched off). I was preferring the solution with entering and returning, because then the output field can be set to zero if needed.
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.
Hi @JorgSchwinger , I guess, I am a bit confused about what you mean by 2-3 different variables for the same thing... - but maybe we can quickly talk about this tmw. Wrt blom2hamocc
- then it would make sense to do the very same for all sorts of fields coming from BLOM - e.g. also the gaseous fluxes, etc. I do see the mixture of approaches that we use, but then I would prefer to restructure all of it. Can we keep the status as is for now and merge it in (after some further testing, potentially still re-writing the reading of the N-deposition files - dependent of progress with the provision of input files)? - and we can discuss further re-organization afterward. As you might have seen, I reverted back to going in and out of the read_ndep
(your preferred way of handling this).
hamocc/mo_read_ndep.F90
Outdated
ndep(i,j,idepnoy) = noydepread(i,j) | ||
ndep(i,j,idepnhx) = nhxdepread(i,j) | ||
else | ||
ndep(i,j,idepnoy) = ndepread(i,j) |
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 would favor abandoning the reading of the old ndep-files and always read the new files here. Only in apply_ndep
I would lump together the fields in the extended N-cycle is not activated.
@@ -77,7 +79,7 @@ subroutine sedshi(kpie,kpje,omask) | |||
& + oplfa*sedlay(i,j,k,issssil) & | |||
& + clafa*sedlay(i,j,k,issster) | |||
! "full sediment has sedlo=1 | |||
wsed(i,j)=max(0.,(sedlo-1.)/(sedlo+1.e-10)) | |||
wsed(i,j)=max(0.,(sedlo-1.)/(abs(sedlo)+1.e-10)) |
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.
Was this a bug?
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.
It doesn't change anything in the 3D, but it fixed running the sediment in the 1D setup.
Hi @TomasTorsvik and @JorgSchwinger , thanks for your comments thus far. I'll tackle them tomorrow with a fresh mind. @JorgSchwinger , since we briefly spoke about the compsets - you basically suggested to delete those parts
in the |
write(io_stdo_bgc,*) ' ' | ||
write(io_stdo_bgc,*) 'AUFR_BGC info: preformed silica not in restart file ' | ||
write(io_stdo_bgc,*) 'Initialising preformed tracer from scratch' | ||
endif |
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.
A post-review comment: I think the option to initialize preformed Si from scratch while the other preformed tracers are read from restart file is potentially dangerous if users are not aware of this feature.
I do see that being able to start from older restart files is useful for now, but I would suggest to remove this option again for release NorESM 2.5
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.
Hi @JorgSchwinger , agreed, so keeping it for now. For NorESM2.5, I will check, if there are climatologies available that can be used (I vaguely remember to have seen some). I add this to #340.
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.
But wouldn't it be important to have the preformed tracers consistent - either all preformed tracers initialized from scratch (set to zero) or all preformed tracers initialized to a climatological estimate?
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 guess, my thought went along the lines of receiving a better spunup preformed silicate tracer field with the new tests and then provide them with the next release - since the longer the spinup time, the better it represents the model, if I am not mistaken?
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.
Maybe this could be implemented as an option - having a climatology of all preformed tracers to be used as an initial condition. I would have this optional, because I don't think it will make sense for higher resolution runs (unless we interpolate from lower resolutions).
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.
For this PR, I put this to the to-do-list #340.
@jmaerz , @mvertens - I ran the regression test based on
|
Hi @mvertens , @JorgSchwinger , and @TomasTorsvik , |
@jmaerz - if you think this can be done within a few days, I agree that it's worth trying to get into |
Latest regression on the
|
…h recent master Before this last commit leuphotic_cya=true, it was bfb with default settings with master
Final merge note: I perform the merge of the extended nitrogen cycle (++, see merge description above) with |
This merge request aims at bringing in i) the extended nitrogen cycle developments carried out in the ESM2025 project, ii) the preformed silicate tracer, iii) optionally bluefix (cyanobacteria) only in the euphotic zone, and iv) the M4AGO sinking scheme (Maerz et al. 2020) as submodule into iHAMOCC/BLOM/NorESM into
master
.The merge is carried out to enable a better overview of potentially needed adaptations wrt e.g. the nuopc cap. Interactive coupling of the extended nitrogen cycle to the atmosphere is currently only available through the MCT coupler and only with a certain CAM branch (so non-default for the extended nitrogen cycle, while it is anyways only run optionally).
Note 1: This is a branch still in (minor) progress (mainly concerns tuning and thus model parameter values). At current stage and in preparation for NorESM2.3, the extended nitrogen cycle and the M4AGO scheme are optional and can be switched on via
xmlchange
:where the last two are in place to enable the interactive coupling to the atmosphere (non-default - CAM branch under development).
Note 2: The M4AGO submodule is providing the M4AGO sinking scheme as expected, while the source code is up to further changes. Note that the ocean biogeochemistry currently needs about 1-2 year(s) until the sinking scheme gives reasonable sinking velocities and fluxes, which is due to changes in POM inventories.
About the extended nitrogen cycle:
The extended nitrogen cycle incorporates two new tracers: ammonium and nitrite.
The extended nitrogen cycle features an explicit representation of the following processes:
are represented (if applicable, in both the water column and the sediment).
The merge shall be carried out in squash mode
Closes #147