Skip to content
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

Processes at the boundary #17

Open
MatthisAuger opened this issue Jan 16, 2023 · 63 comments
Open

Processes at the boundary #17

MatthisAuger opened this issue Jan 16, 2023 · 63 comments
Assignees

Comments

@MatthisAuger
Copy link

Last week we discussed the processes at the model frontier and how the model was behaving with boundary conditions. It seems like they are not the candidates for the large jump and trends between the models anymore, but it might be still interesting to look at what is happening there.

I made a few videos, looking at the bottom layer of Temperature, Salinity, U, and V outputs from the hycom1 model and boundary conditions. This is focused on the bottom layer as initially the interest was on the lower overturning cell, but I can make these videos in any other layer or region if needed.

Here are all the parameters along the boundary over longitudes 50W-10W. These are the outputs over the first few months of the model, boundary conditions, and differences between outputs and boundary conditions.

TempSaltUVdiff_lon50w-10w-4800_3000_panant-hycom1.mp4

Same thing but across the boundary this time, at longitude 70E if I remember correctly.

TempSaltUV_accross_lat39s-37s-4800_3000_panant-hycom1.mp4

For both videos, T and S do not vary a lot as it is the bottom layer, but U outputs seem to go a bit crazy compared to boundary conditions.
As I said just tell me if you are interested in similar plots for other regions or layers. I can make those closer to the surface if you want to see how the model reacts to stronger changes in boundary conditions temperature or salinity.

@MatthisAuger
Copy link
Author

Same thing for the surface layer:

Along the boundary:
https://user-images.githubusercontent.com/68908668/213604623-a53fcd56-c315-4064-a4e9-12e0b7fc95a5.mp4
And location of the transect
along_pos (1)

Across the boundary
https://user-images.githubusercontent.com/68908668/213604621-acaa47b4-d1f3-4aa0-8237-36f9e76d21c7.mp4
And location of the transect
accross_pos (1)

@adele-morrison
Copy link
Contributor

Update from our meeting 20/1/23:

The movies above show that there are particularly large u velocities at the northern open boundary. To investigate this more, we can:

  1. Plot SSH transects at the boundary. Are the large u being driven by SSH gradients at the boundary?
  2. Do some short (maybe 1 month is sufficient?) runs with panan-01 changing the reservoir volume and the nudging timescale. See how the rms of u near the boundary changes with these parameters.

@MatthisAuger and @PaulSpence and @pedrocol volunteered to look into this.

@PaulSpence
Copy link

@adele157 @AndyHoggANU
Can you please point me to the latest mom6 panan config that we should use to test the boundary condition settings? We aim to start testing them today. Thanks!

@AndyHoggANU
Copy link
Contributor

Probably this branch, for 0.1° zstar:
https://github.com/COSIMA/mom6-panan/tree/panan-01-zstar-run

@MatthisAuger
Copy link
Author

MatthisAuger commented Jan 25, 2023

Just a quick plot of sea level at the boundary.
Those are sea level variations over one year for the 1st year of simulation (red) and the boundary conditions (black) at one particular longitude.

sea_level_boundary

If both variables are supposed to be the sea level above the geoid then it explains the strong U along the boundary.
Looking at the Mean Sea Surface plots at 37S it looks like the model has the good order of magnitude of sea level and not the boundary conditions files (the plot is at longitude 70E)

image

@adele-morrison
Copy link
Contributor

Wow that's a big difference! Might be useful to compare lat/lon maps of SSH in the panantarctic and ACCESS-OM2-01-RYF. Is the RYF just really biased or is there something wrong with the boundary forcing?

@AndyHoggANU
Copy link
Contributor

Yes, it is a big difference. But the offset surprises me -- I had thought U oscillated a bit and was positive at some times, negative at others.

Maybe it's just a definition of sea level? If yes, SSH gradient perpendicular to the boundary would tell us?

@PaulSpence
Copy link

Hi Folks,

We are running 1 month sims varying these params in MOM_input:
OBC_SEGMENT_001_VELOCITY_NUDGING_TIMESCALES = .3, 360.0 ! inflow and outflow timescales
OBC_TRACER_RESERVOIR_LENGTH_SCALE_OUT = 30000
OBC_TRACER_RESERVOIR_LENGTH_SCALE_IN = 3000

See relevant code here: https://github.com/mom-ocean/MOM6/blob/7467a63efea7025ceb9118448d593709dc1cdf47/src/core/MOM_open_boundary.F90

Since we don't understand this code, the plan is to run the following *2 and *0.5 sensitivity tests:

  1. OBC_SEGMENT_001_VELOCITY_NUDGING_TIMESCALES = .6, 720.0 ! inflow and outflow timescales *2
  2. OBC_SEGMENT_001_VELOCITY_NUDGING_TIMESCALES = .15, 180.0 ! inflow and outflow timescales *.5
  3. OBC_TRACER_RESERVOIR_LENGTH_SCALE_OUT = 60000
    OBC_TRACER_RESERVOIR_LENGTH_SCALE_IN = 6000
  4. OBC_TRACER_RESERVOIR_LENGTH_SCALE_OUT = 15000
    OBC_TRACER_RESERVOIR_LENGTH_SCALE_IN = 1500

Trial by brushfire :)

@angus-g
Copy link
Collaborator

angus-g commented Jan 27, 2023

I think one of the configs that inspired our choice of parameters was https://github.com/ESMG/ESMG-configs/tree/dev/esmg/CCS2, although the reservoir length scales are different to that one: they're swapped. I think we adopted that because the other configuration we looked at (https://github.com/ESMG/Arctic6) doesn't override them at all, but we were getting issues leaving them at the defaults.

@aekiss
Copy link

aekiss commented Jan 30, 2023

Thanks @MatthisAuger.

Does the panan model directly use the ACCESS-OM2-01-RYF SSH at the boundary value when calculating its pressure gradient?

It doesn't appear to, at least from your plot. If I'm reading it correctly, it shows SSH increases to the south, implying a westward geostrophic current, but I think your movies show it's mostly eastward. Also a ~0.7m difference over 0.08°=8.9km implies a zonal geostrophic current of ~8.8 m/s at 37S, which is pretty extreme and seems more than what your movies were showing. I also share Andy's surprise that the SSH offset is nearly constant in time, in contrast to the U anomalies.

So all up, I don't think the SSH difference is consistent with the anomalous velocity, but I agree with Adele that lat/lon maps (ideally movies) of SSH difference between the panantarctic and ACCESS-OM2-01-RYF would be helpful, to get a fuller picture of the spatial pattern.

Although I'm not sure whether this SSH difference is relevant to the U anomaly, so this is probably a rabbithole, I'm still wondering how the SSH difference can be maintained.

  • Is SSH one of the fields handled by the open boundary conditions?
  • If so, how is it applied?
  • why aren't the SSH values brought closer?
  • Is there some condition of zero net volume flux through the boundaries?
  • Does the panan take its SSH initial condition from ACCESS-OM2-01-RYF? edit: Angus said it was uniformly zero, so it seems to me that would explain the SSH offset

@MatthisAuger
Copy link
Author

MatthisAuger commented Feb 2, 2023

Here is the zonal mean of the SSH at 37S for the cases listed by @PaulSpence

sea_level_boundary_sensitivity_zonalmean

Nudging timescale seems seem to have a strong impact, with a SSH difference of 1cm between the *0.5 and *2 cases, which might be a lot considering that this is a zonal mean over all the longitudes.

About the SSH difference between ACCESS OM and panant, the fields look very different. (I am even wondering if I'm getting the right variables?)
Sorry for the not convenient image formats (edit: 2nd was the wrong figure)

2 snapshots on the same day (notice the different colorbar)
zos_snap
access_snap (1)

mean difference over 1 year

sla_diffmean

And movie of the difference (mean difference removed):

sla_diff.mp4

@adele-morrison
Copy link
Contributor

Discussion from today's meeting:
The SSH offset may be arising because we are setting a 0m initial condition for sea level. Averaged south of 37S, the sea level in obs and ACCESS-OM2 is way less than 0m (perhaps equal to the 0.8m offset we are seeing?).

Suggestions from today's meeting:

  • Check if the boundary condition is mass conserving (check mom6.out and plot a time series of total mass).
  • Can we set a non-zero initial condition for sea level. @angus-g will look into this.

@angus-g
Copy link
Collaborator

angus-g commented Feb 2, 2023

Check if the boundary condition is mass conserving (check mom6.out and plot a time series of total mass).

There's also ocean.stats and ocean.stats.nc to avoid as much manual parsing.

@adele-morrison
Copy link
Contributor

@willaguiar do you feel like looking into whether the regional model and boundary condition are mass conserving? (see last two comments above)

@willaguiar
Copy link
Collaborator

@adele157 Yes! - I can do that and post here the results.

@willaguiar do you feel like looking into whether the regional model and boundary condition are mass conserving? (see last two comments above)

Yes! - I will do it. (I'll post the update/figs here)

@aekiss
Copy link

aekiss commented Feb 2, 2023

But does the SSH difference at the boundary play any role in the U anomalies at the boundary? It doesn't look like it to me - see above

@adele-morrison
Copy link
Contributor

adele-morrison commented Feb 3, 2023 via email

@willaguiar
Copy link
Collaborator

The total mass time series is plotted below.
It seems to me that there is a slight positive trend, but the variability would be about ~0.1% only. So maybe not strictly mass conserving?

panant-01-zstar-v13_Totalmass

@aekiss
Copy link

aekiss commented Feb 3, 2023

We have a ~1m SSH mismatch with ~4000m mean depth, so it would require a mass change of only about 0.025% for the SSH mismatch to disappear.

Your plot shows an initial very rapid adjustment of 0.018% which is remarkably close to the estimated required value. The subsequent variability and drift is about 0.003%.

So it looks like it's trying to make the SSH match, but for some reason it hasn't adjusted remove the SSH mismatch completely.

Can you please calculate the spatial-mean SSH anomaly? That will let us quantify this more accurately.

@aekiss
Copy link

aekiss commented Feb 3, 2023

The spatial-mean SSH anomaly timeseries would be even better...

@willaguiar
Copy link
Collaborator

Here it is the timeseries for the spatial-mean SSH anomaly (Anomaly referent to the first output). @aekiss
panant-01-zstar-v13_aSSH

@adele-morrison
Copy link
Contributor

adele-morrison commented Feb 8, 2023

@MatthisAuger Perhaps it would also be useful to compare the global-01-v2 simulation with the panan-01 simulation on the time series of SSH at the boundary. That might tell us if it's just a mismatch of the physics between the mom6 model and the boundary forcing, or whether the boundary forcing is doing something weird.

(The first couple of years of global-01-v2 output are missing, but you could compare for a later time period.)

@aekiss
Copy link

aekiss commented Feb 9, 2023

thanks @willaguiar

  • there's an initial very rapid drop of 0.2m, so that suggests a 0.2m/4000m = 0.005% mass change, which is about a third of the initial 0.018% mass loss you show above
  • the subsequent drift of about 0.03m, ie 0.00075% is also smaller than the mass drift

So that's pretty confusing - it's hard to see how it's possible for the mass to change by more than the SSH would allow.
Is the sea ice growing, and does it depress the ocean surface? Perhaps that might account for some of it.

@MatthisAuger
Copy link
Author

@adele157 Here is the plot of year 1994 SSH zonal mean time series at the boundary, for boundary forcing, mom6 and panan.

boundary_globalv2_panan_SSH

@aekiss
Copy link

aekiss commented Feb 13, 2023

@angus-g said the sea ice initial condition is no ice. So the model would definitely increase the sea ice mass, which presumably reduces the ocean mass. And if sea ice depresses the surface, that would give us a smaller mean SSH change than we'd expect from the ocean mass change, which is the same sign as the inconsistency we're trying to explain.

However if the sea ice is at all realistic, it won't quantitatively explain this little mystery. We need ~4e16 kg of new sea ice to account for the remaining 2/3 in this graph, but the Antarctic sea ice mass in ACCESS-OM2-01 varies seasonally between about 0.3e16 and 1.7e16kg, so that's too little mass, with too much variation, to explain the inconsistency between SSH and ocean mass changes.

Do we have a timeseries of total sea ice volume or mass?

@aekiss
Copy link

aekiss commented Feb 13, 2023

A timeseries of the total mass of ocean salt would also help distinguish between ocean mass loss due to the open boundaries vs sea ice formation & evap-precip-runoff.

@AndyHoggANU
Copy link
Contributor

seaice.stats should give you the sea ice stats, and I think includes mass.
A comparison with the global case (which also starts without sea ice in the first year) might be informative.

@MatthisAuger
Copy link
Author

Just plotted SSH spatial average time series over the first few years of simulation and got a different time series than this one. I find very small variations around 0 and no big changes during spin up. That might come from weighting the mean with the grid cell area?

image

In that case the SSH average would remain constant and not follow the total mass change?

Currently looking at mass budget and integrated in and out fluxes at the boundary so more figures to come.

@MatthisAuger
Copy link
Author

And same figure I presented last week with total ocean mass and sea ice mass, but with separated plots this time.

image

@PaulSpence
Copy link

Net surface water flux shows a small increasing trend. Need more diagnostics to see where this trend is coming from.
Screen Shot 2023-03-16 at 11 55 43 am

@PaulSpence
Copy link

PaulSpence commented Mar 16, 2023

WFO Anom relative to 2000-2060 mean shows positive trend, that is much smaller than the net input from rivers.

Screen Shot 2023-03-16 at 12 20 51 pm

Screen Shot 2023-03-16 at 12 20 37 pm

P-E should be also be positive (though we are missing the diagnostics to check). The OBC are fluxing mass into the model as well. Everything is increasing the volume of the SO. Where does this water go? How does the model handle a constant increase in volume/mass over time?

@PaulSpence
Copy link

PaulSpence commented Mar 27, 2023

Can we close the mass balance in panan, e.g. total mass change = fw surface + boundary? is there a missing process or does the budget make sense? Where does the excess mass go? Does SSH just increase and drive an outflow (into boundary)?

Ask Chris Chapman - John Rielly

@PaulSpence
Copy link

Does the SSH increase with the mass? Why not?

@MatthisAuger
Copy link
Author

From ocean.stats on the panant-01-zstar-v13 experiment:
It confirms that the excess of mass goes into sea level. Not sure how this mean sea level is computed though as it is different from the two previous SSH plots obtained with the model output.

image

@adele-morrison
Copy link
Contributor

Did we figure out if sea ice levitates or not? That could make a difference in different sea level diagnostics.

@MatthisAuger
Copy link
Author

Order of magnitude of the total mass change budget:

From previous plot: #17 (comment)
Total mass trend after spin up ~2.56 * 10^6 kg/s

From Paul's comment: #17 (comment)
Net surface mass flux ~1.07 * 10^9 kg/s

From meridional transport into the model: #17 (comment)
Mean flux into the model at the boundary (after spin up)= 0.85Sv, ~ 0.87 * 10^9 kg/s

So right now there is way more mass going in than the total mass change. Where does all this mass go?

@MatthisAuger
Copy link
Author

@adele157 Still not sure if it levitates.
Looking at #17 (comment) I would say yes, as the ocean and sea ice mass seasonal cycles should compensate each other?

@PaulSpence
Copy link

@access-hive-bot
Copy link

This issue has been mentioned on ACCESS Hive Community Forum. There might be relevant details there:

https://forum.access-hive.org.au/t/mom6-regional-modelling-mini-hackathon-may/725/2

@aekiss
Copy link

aekiss commented May 3, 2023

re. #17 (comment), if the mass flux into the model is about 1000x more than the mass change, is there a units error somewhere in the diagnostics or our calculations?

@aekiss
Copy link

aekiss commented May 3, 2023

@PaulSpence does Net surface water flux (precip+melt+|runoff+ice calving-evap) in #17 (comment) take into account sea ice formation as well as melt, or only sea ice melt?
And does positive sign indicate a flux into the model?

@aekiss
Copy link

aekiss commented May 3, 2023

@MatthisAuger re. #17 (comment), let me see if I understand your argument for seeing the effects of levitation from these ocean & ice mass timeseries #17 (comment)

In a global model, levitation would affect how much the ocean SSH is affected by mass transfer into sea ice:

  • if levitating, increasing sea ice mass will reduce SSH. Also snow accumulation on sea ice (if SIS2 includes that) would further decrease SSH since it would be precip withheld from the ocean. So we'd see a strong seasonal cycle in SSH.
  • if not levitating, sea ice (and any accumulated snow) depresses the surface so the the SSH in ice-free regions would be nearly unchanged, and show little seasonal signal.

However things are complicated by the open boundary - in the non-levitating case, would the pressure gradient at the boundary due to a seasonal SSH cycle at the boundary (from to the seasonal sea ice) cause seasonal pumping of mass in/out of the model, reducing the SSH and ocean mass signals? Is that what causes the seasonal cycle in the bottom figure here #17 (comment)? It flows outward in the melt season and inward in the freeze season, which is what is expected from this scenario.

The SSH doesn't show a clear seasonal cycle here #17 (comment) and and even less in this other SSH diagnostic #17 (comment) (it's odd that these look so different - do we know why?)
and the seasonal cycle in ocean mass is much weaker than in sea ice mass #17 (comment) so that is consistent with mass transfer to sea ice being compensated by inflow at the open boundary. Presumably there is enough of a seasonal cycle in SSH at the open boundary to drive these compensating flows.

So if this scenario is correct, while the fairly non-seasonal SSH timeseries might look like evidence for non-levitating, it is actually also inconsistent with levitating ice with seasonal boundary inflow/outflow. Does that sound right?

@aekiss
Copy link

aekiss commented May 3, 2023

Of course to settle the levitation question, we should look at the model parameters. In MOM5 this is controlled by max_ice_thickness - does anyone know the equivalent in MOM6? There's DEPRESS_INITIAL_SURFACE and TRIM_IC_FOR_P_SURF, but from a quick glance they seem to be intended for setting up initial conditions with ice shelves, rather than adjustments with sea ice changes during a model run. But I could be wrong. They both default to false anyway.

@aekiss
Copy link

aekiss commented May 3, 2023

Are we using a Boussinesq formulation? If so, that will conserve volume rather than mass (but of course the difference will be way too small to explain our factor of 1000 mismatch).

@access-hive-bot
Copy link

This issue has been mentioned on ACCESS Hive Community Forum. There might be relevant details there:

https://forum.access-hive.org.au/t/mom6-regional-modelling-mini-hackathon-may/725/4

@adele-morrison
Copy link
Contributor

BOUSSINESQ = True

Not sure about whether the sea ice levitates - I can't see any parameters associated with this. I suspect it doesn't levitate, because the reason for doing this in MOM5 was due to ocean - sea ice instabilities (see page 38 here). Updating SIS to SIS2 and more tightly coupling with MOM6 may have removed the need to levitate.

@aekiss
Copy link

aekiss commented May 3, 2023

Are we able to close the mass budget for the ice-free regional models, e.g. https://github.com/COSIMA/mom6-eac and what @ashjbarnes is developing? Maybe we're missing a term that has nothing to do with sea ice?

@PaulSpence
Copy link

Real vs virtual FW fluxes ... P-E+R is treated as a virtual freshwater flux that impacts the salt budget (not the volume) in MOM
https://www.sciencedirect.com/science/article/pii/S2095927317300932#b0095

@MatthisAuger
Copy link
Author

Discussion with Stephen Griffie.

Transport at the boundary and P-E+R should change the ocean volume accordingly.

Use ssh_ga to have the global averaged sea level, not zos_ga.

Check .wfo files. (Where are they?)

Double check the volume flux in P-E+R and at the boundary.

@aekiss
Copy link

aekiss commented May 9, 2023

from #17 (comment)

Mean flux into the model at the boundary (after spin up)= 0.85Sv, ~ 0.87 * 10^9 kg/s

Order of magnitude, the open boundary area would be about 1e11m^2 (40,000km * cos(37) * 4000m), so 0.85Sv implies a mean velocity of about 8.5e-6m/s, which seems low...?

But if it was bigger the problem would be worse.

@aekiss
Copy link

aekiss commented May 9, 2023

More order-of-magnitude calcs:

from #17 (comment)

From Paul's comment: #17 (comment)
Net surface mass flux ~1.07 * 10^9 kg/s

From meridional transport into the model: #17 (comment)
Mean flux into the model at the boundary (after spin up)= 0.85Sv, ~ 0.87 * 10^9 kg/s

So total boundary flux into model is ~1.9e9 kg/s

Total ocean mass =3.4e20kg, so timescale to double ocean volume is 3.4e20/1.9e9 = 1.8e11s = 5700 yr. Does that seem short?

But ocean area ~volume/depth ~ 3.4e20kg/4000m/1000 = 8.5e13 m2 so net surface mass flux ~1.07 * 10^9 kg/s is 1.2e-5 kg/s/m2 ~ 1.2e-8 m/s = 1 mm/day which seems a reasonable average precipitation - evap rate.

@aekiss
Copy link

aekiss commented May 9, 2023

The total ocean mass also looks ok to me. Rough numbers 1000* 4pi *6000km^2 *4000m/(say) 4 = 4.5e20 kg, close to the 3.4e20kg above.

So why is the mass trend so much smaller than the net mass flux?

@MatthisAuger
Copy link
Author

MatthisAuger commented May 9, 2023

vmo2d across the boundary = 1.062e+09 kg/s out of the model NB: opposite sign to above, and now including eddy mass flux
Ocean mass trend = 1.074e+05 kg/s NB: 25x smaller than previous estimate 2.134e+6kg/s
Panan mean of wfo = 1.083e+09 kg/s into the model

image

image

@aekiss
Copy link

aekiss commented May 9, 2023

vmo2d across the boundary = 1.062e+09 kg/s out of the ocean NB: opposite sign to above, and now including eddy mass flux
Ocean mass trend = 1.074e+05 kg/s NB: 25x smaller than previous estimate 2.134e+6kg/s
Panan mean of wfo = 1.083e+09 kg/s into the ocean

So net boundary flux into ocean is 1.083e+09 kg/s - 1.062e+09 kg/s = 2.1e7 kg/s. 4.97e6 kg/s - see below

So 200x 10x 2.3x bigger than ocean mass trend.

We would need very close (1 part in 10,000) cancellation of the boundary fluxes to match the mass trend. So everything would need to be calculated very carefully and consistently.

@MatthisAuger
Copy link
Author

Very sorry for all the diagnostics that turned out to be wrong. Still not sure to understand how my mass flux at the boundary was so far (even opposite) from the new one computed from the vmo variable. Also spotted a mistake in the computation of Tuesday's total ocean mass trend so please refer to the new value of 2.134e+6kg/s.

I merged @PaulSpence notebook on surface fluxes and mines to compute time series of mass and fluxes. Here is the notebook used to make the following figures:

Total mass flux:
image

Mass budget: Total mass from diagnostics and recomputed total mass from monthly mass fluxes.
image

Increase of mass associated with mass fluxes is still twice higher as the total ocean mass increase. Could be because of residuals due to averaging of the surface net flux?
Surface net flux is a monthly average while the vmo output at the boundary is computed from the total mass that went through the grid cell?

@PaulSpence
Copy link

PaulSpence commented May 11, 2023 via email

@StephenGriffies
Copy link

Have you confirmed that your notebook shows that everything is self-consistent for a simulation with vmo=0 at the open boundaries? Just need to run for a few days.

@aekiss
Copy link

aekiss commented May 11, 2023

No need to apologise, thanks for your work on this. It seems pretty clear from your plot that the surface mass inflow drives a compensating outflow at the open boundary, so the mass accumulation involves a subtraction of two terms that are very close, so it will show up every little difference in how those numbers are calculated.

Nice to see that the discrepancy is getting smaller. Your green trend line of 4.97e6 kg/s is much less than the 2.1e7 kg/s I put in my (now updated) post above, so the mismatch is now a factor of 2.3.

Are the surface and open boundary fluxes both calculated from monthly averages, not snapshots? And are the varying month lengths taken into account for dt?

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

9 participants