Set flag transition ocean continent

[anne-glerum]

I noticed that the amount of sediments in the marine domain increases for increasing the marine diffusivity, even for zero continental sediment input. I also noticed that the flag that defines whether a node is on the boundary between the continental and the marine domain is used, but not set to something other than its default value 0.

After discussing this with Xiaoping, he suggested making this modification to Marine.f90 that does set the flag to 1 at the boundary. Xiaoping has used this fix in his version of the code as well. I verified that the flag values change in space over time following the boundary.

The snapshot below shows the surface of two simulations without continental erosion and a very small continent at a certain moment in time. One run does not change the flag from 0 (i.e., current version, red interface), while the other sets the flag to 1 at the transition from continent to ocean (green interface). The initially vertical margin does not change over time in the green version, as there is no continental input. In the red simulation artificial sediment is built up as the margin diffuses.

I will add more examples demonstrating the difference later.

[benbovy]

Thanks @anne-glerum.

I'm not very familiar with the marine component but it seems to make more sense delineating the land/ocean interface using elevation and sea-level than using a flux.

Some quick thoughts / questions (I haven't looked in detail yet):

• Is it good enough using the single flow direction graph here? Or should the flag be computed from the multi-direction flow graph?
• flag is used only in SiltSandCouplingDiffusion. Would it be possible to remove this array and do the checks directly in the latter routine instead?

[xpyuan1]

Dear Anne,

I have checked your statement, and agree with it. As we have discussed, this flag setting is important to make a more realistic boundary between ocean and continent. Otherwise, there will be ‘automatic' diffusion along this boundary which will cause your unrealistic red-interface simulation.

For Benoit’s questions,

1. You already use the sea-level elevation to set the flag (if (h(ij).ge.sealevel.and.h(ijr).le.sealevel) flag(ijr)=1). So it is correct.
2. I think single-flow direction is enough here, not necessary to complicate things (to use, e.g., multi-direction flow algorithm). The results will look very similar. I also use single-flow algorithm in the continental area, no one questioned me. ;-)
3. flag is only used in this marine-diffusion code, and it is fine to keep it here.

All the best,
Xiaoping

[benbovy]

Thanks @xpyuan1 for your input.

I think single-flow direction is enough here, not necessary to complicate things

OK, good. Maybe it would be nice to add a small comment explaining that just above the lines computing the flag values.

flag is only used in this marine-diffusion code, and it is fine to keep it here.

OK. Not sure to really see what is the point to create outside of a function an internal variable that is used only within that function, but if it is more readable to define the continent / ocean boundary where it is currently defined lets keep it there then. flag is a generic name and not very meaningful here, so it would be nice to choose a better name for that variable.

If these two tiny concerns could be addressed, that would be great and I'd be happy to merge this PR promptly!

[anne-glerum]

Hi @benbovy, apologies for the late reply.

I'm not very familiar with the marine component but it seems to make more sense delineating the land/ocean interface using elevation and sea-level than using a flux.
Note that the flag was not set to anything but the default 0 before, the line using the flux was commented out. So the interface was not taken into account.

If these two tiny concerns could be addressed, that would be great and I'd be happy to merge this PR promptly!

I've added a small comment about how the flag is computed and I've changed the flag's name to COTflag.

flag is used only in SiltSandCouplingDiffusion. Would it be possible to remove this array and do the checks directly in the latter routine instead?

I agree that either the COTflag array could be filled directly in SiltSandCouplingDiffusion or only the if-statement itself could be done in SiltSandCouplingDiffusion. However, for simplicity, I've kept the filling of the array in the location is was before.

Let me know if you want anything else to be changed.

[anne-glerum]

Do I need to worry about the one failing test? The error messages don't seem to be related to my changes.

[anne-glerum]

Anything else I need to address @benbovy ?

[benbovy]

Thanks @anne-glerum, I think it is all good. The failing test is not related to this PR, I'll investigate this in #59 next week and merge this as well.

[benbovy]

#59 is merged and this should be good to go in now. Thanks again @anne-glerum and sorry for the delay.