search

edoddridge / aronnax

An idealised isopycnal model that can be run either with n+1/2 layers, or with n layers and variable bathymetry.
http://aronnax.readthedocs.io/en/latest/
MIT License
24 stars 6 forks source link

Physically motivated tests #66

Open edoddridge opened 7 years ago

edoddridge commented 7 years ago

The current test suite checks for consistency between the new outputs and a number of outputs designated as "good".

We should augment these checks with some physics motivated tests. These new physics motivated tests should include:

axch commented 7 years ago

Quick question: are not the boundary conditions enforced? If they are, is there any content to checking that they are respected? Or am I missing something?

edoddridge commented 7 years ago

They are enforced. I guess the point of this would be to check that they've been done correctly.

axch commented 7 years ago

More generically: these all look great. I suggest the following form for the tests: As one stage, measure the quantity in question, ideally in a fashion that's amenable to aggregation and plotting; and then as a second stage, if the answer is apparent, make an assertion about the result of the measurement.

For example, the current test suite does and does not follow that advice. It does in the sense that it aggregates array-wise relative error for each dumped output file, instead of asserting element-wise relative error at each point. It doesn't in the sense that it asserts a relative error bound for each dumped output file separately, instead of aggregating a relative error profile for the whole run.

The reason for this advice is that all the above metrics seem to be much easier to specify as measurements than to know what the "right" tolerances to permit actually are.

edoddridge commented 7 years ago

I've been working on this in the background, and have found that in n+1/2 layer mode the model receives 1% more momentum than I expect it to. I've gone through the code and am still not sure what is causing the discrepancy.

f_plane_momentum_test Momentum evolution through time

percent_error Percentage error between the simulated momentum and the expected momentum.

Obviously 1% is pretty small, but I'd like to work out where it is coming from. Provided I keep the wind forcing smaller than 2e-3 the error is 1%. With stronger wind forcings things get a bit weird - probably because I'm running these tests with zero viscosity, and the numerics get unhappy.

axch commented 7 years ago

Re: "numerics get unhappy": Perhaps invest the time in some diagnostic plots that will conclusively detect various forms of numerical unhappiness? A priori, I would expect an actual instability to be quite easy to detect.

Another generic piece of advice: What's your development cycle time for those two plots? Can we get it down small enough to just fiddle, and see what parameters make that 1% bigger or smaller? (Resolution, size of the pool, depth(s), etc)

Also, is this test code on some branch(es) somewhere?

edoddridge commented 7 years ago

I agree - I don't think this is a numerical instability. For larger forcings the momentum increases more than it should, but also sporadically decreases back towards the expected value.

The plots are very quick to generate, and the code is on the pysics-tests branch. Apparently correct spelling costs extra...

edoddridge commented 7 years ago

One thing I've noticed while doing this is that the convergence criterion used for the linear solve makes quite a difference to the solution. This may be what caused the large discrepancy you found between the results with the SOR solver and the external solver. Since the external one converges much faster, it likely overshot the criterion used by the internal SOR solver.

If you run the test script in the pysics-test branch you'll see that it also looks at the n-layer configuration, where the momentum is ~4 times larger than expected. I haven't had a chance to look at that closely yet, but I think I may have done something silly with the layer thicknesses used for the expected momentum calculation.

axch commented 7 years ago

It occurs to me that if the model dumps output on every time-step, it may be possible to compute from the output how much volume or momentum is lost (or gained) due to sponging (especially if we add an auxiliary dump before sponging is applied). With that, it may be possible to do conservation studies on an infinite domain, where the infinity is implemented by a sponge region around the border of the simulation. This may be better than trying to predict edge effects, and may be faster than running a simulation on a sufficiently large finite domain that the edges are far enough away not to matter.

edoddridge commented 7 years ago

That's a good thought. We can definitely use sponges to mop up edge effects.

edoddridge commented 7 years ago

I have some code on the pysics-test branch that I think closes #99 and provides more detail for #100. But, because of my earlier, rather cavalier, approach to merging the external-matrix-solver branch into it, any pull request on the branch with the tests is currently blocked on #89. The behaviour of the error with resolution in the reduced gravity mode is not quite what I expected, but does suggest truncation error. The graph is included below.

error_by_resolution_semilogx

Depending on where you're at with reviewing that branch, I'm happy to split the relevant files off on to a separate branch - I'm not sure if I can do that in a history preserving way, but I can manually recreate it - and then put that in as a PR.

axch commented 7 years ago

Sent next round of feedback on #89. My inclination for now is to still try to get that merged first, so no history surgery becomes necessary.