Hamocc hybrid coord2

[TomasTorsvik]

For hybrid cooridinates, define depth of surface mixed layer in HAMOCC by the ocean boundary layer (OBL) variable provided by CVmix. Default value kml(i,j) = 2 used for isopycnic variables for consistency with the bulk mixed layer formulation.

[jmaerz]

I am a bit uncertain, but CVMix_kpp_compute_kopb_depth() i) seem to return a real value and ii) depending on which index one wants, it potentially requires to subtract 1 (not sure about the vertical grid layout in BLOM) iii) why calculating it again - and not using hOBL? - from CVMix:

! !DESCRIPTION:
!  Computes the index of the level and interface above OBL\_depth. The index is
!  stored as a real number, and the integer index can be solved for in the
!  following way:\\
!    \verb|kt| = index of cell center above OBL\_depth = \verb|nint(kOBL_depth)-1|
!    \verb|kw| = index of interface above OBL\_depth = \verb|floor(kOBL_depth)|

In any case, to be clean, I guess it would be good to explicitly convert from real to integer for the k-index. This applies again below.

[TomasTorsvik]

Thanks @jmaerz !

@matsbn @milicak I would like to hear from you about this issue.

I did not look into the CVMix routine itself when I made these code changes, but I agree that this looks suspicious. For iHAMOCC we probably want the index of the cell center above OBL_depth, not the index of the interface.

hOBL is computed twice, first by the subroutine CVMix_kpp_compute_OBL_depth and later by function CVMix_kpp_compute_kOBL_depth. The function CVMix_kpp_compute_kOBL_depth is used internally in CVMix_kpp_compute_OBL_depth, so there should be no need to do this again. Then kOBL is computed (twice) by the function CVMix_kpp_compute_kOBL_depth.

kOBL is used in the calculation of buoyancy flux acting on the OBL, so in this case the index may need to reference the interface. If so, we may want to define a new index OBL_cell_center = nint(hOBL) - 1 for use in iHAMOCC.

[TomasTorsvik]

Will update this PR depending on feedback for PR #182

[jmaerz]

Thanks for keeping that going!

[JorgSchwinger]

I haven't dived into the kOBL-issue, but from the HAMOCC side this looks fine to me.

[TomasTorsvik]

Hi,

I have revived this PR after changing the dependency between BLOM and iHAMOCC, so now I'm porting the real value surface mixed layer depth hOBL, corresponding to the fractional depth in the index space, instead of the integer index value kOBL.

@matsbn , @milicak
There is only a minor change in mod_difest.F where I have moved the hOBL variable initiation from the difest_vertical_hyb to the main module initiation section, and made it public so that it is accessible to iHAMOCC. This should not change anything substantial in the physics part.

@jmaerz , @JorgSchwinger
The iHAMOCC implementation is mostly as before, only that now the kmle variable is defined from hOBL instead of kOBL. I had a problem with the default value spval = 1.e33, which doesn't fit a regular integer type, hence I included a test to exclude hOBL instances with this value before attempting to do nint(hOBL), and I set kmle to the default value 2 if this test fails. I'm not sure this is the best solution to this problem, please let me know if you have suggestions how to improve this.

[matsbn]

In routine inivar_difest, you could initialise hOBL(:,:) with the following code block at the end of the routine:

c$OMP PARALLEL DO PRIVATE(l,i)
      do j=1,jj
        do l=1,isp(j)
        do i=max(1,ifp(j,l)),min(ii,ilp(j,l))
          hOBL(i,j)=3.
        enddo
        enddo
      enddo
c$OMP END PARALLEL DO

This should give kmle = 2 always at ocean grid cells, except when vcoord_type_tag == cntiso_hybrid. Maybe then you do not have to access the mod_vcoord variables in iHAMOCC for now.

[jmaerz]

Hi @TomasTorsvik , just to understand the problem with the spval value: do I see it right that the moment, the whole water column is mixed, you experience the problem that hOBL is still equal to spval? - or when does the problem happen? Is spval needed for initialization for hOBL or could you also use e.g. 9999 as initialization (expecting to never reach more than 9998 layers in BLOM in the future). In any case, we don't want to end up with kmle > number of layers.

As minor thing: are the line breaks (the re-formatting) in mo_apply_rivin needed?

Apart from that the changes look fine to me.

[TomasTorsvik]

As minor thing: are the line breaks (the re-formatting) in mo_apply_rivin needed?

No, I changed it back to the the original form.

[TomasTorsvik]

Hi @TomasTorsvik , just to understand the problem with the spval value: do I see it right that the moment, the whole water column is mixed, you experience the problem that hOBL is still equal to spval? - or when does the problem happen? Is spval needed for initialization for hOBL or could you also use e.g. 9999 as initialization (expecting to never reach more than 9998 layers in BLOM in the future). In any case, we don't want to end up with kmle > number of layers.

Hi @jmaerz ,
the update of OBLdepth and hOBL in mod_difest.F is only done for "wet" cells, so the spaval value remained for land cells. If I had done the looping differently, i.e.

   do j=1,jj
      do l=1,isp(j)
         do i=max(1,ifp(j,l)),min(ii,ilp(j,l))
            kmle(i,j) = nint(hOBL(i,j))-1
         enddo
      enddo
   enddo

the loop should not encounter any spval values. However, if for some reason an spval value was encountered, there is the potential for having kmle > number of layers, as you say.

I have changed the initiation for hOBL so that it gets a default value of 3. (defined by hOBL_static), so that kmle gets a default value of 2. In this way we don't need to introduce an explicit dependence on vcoord_type_tag in iHAMOCC, at least for now.

[jmaerz]

Ok, I see, so it was more of a theoretical problem and you haven't encountered it really. Good to know about it and thanks for the information! From my side, you've my approval to merge.

[jmaerz]

Approved.

[matsbn]

It should be OK now. I would prefer that hOBL is first set to spval everywhere, also including the halo, then initialised to hOBL_static only for wet points. The purpose of first initialising to a very large value, spval, is that the code would typically crash if using the variable where one shouldn't (e.g at land points). This can be corrected later.

[TomasTorsvik]

It should be OK now. I would prefer that hOBL is first set to spval everywhere, also including the halo, then initialised to hOBL_static only for wet points. The purpose of first initialising to a very large value, spval, is that the code would typically crash if using the variable where one shouldn't (e.g at land points). This can be corrected later.

Thanks @matsbn , I see that iHAMOCC is referencing the wet points indexing at some places, but I'm not sure if this is done consistently. I therefore prefer to merge this version of the PR, and make a separate issue on the use of wet points indexing in iHAMOCC.