WIP: New speed_limit parameter to avoid too many dt reductions

[rjleveque]

A new parameter rundata.geo_data.speed_limit can be set in setrun.py and is used in b4step2.f90 and getmaxspeed.f90 to scale down the fluid speed sqrt(u**2 + v**2) if it is larger than this value (preserving the flow direction).

This helps avoid crashing the code with "too many dt reductions" in cases when flow off very steep topography with delta B larger than the fluid depth gives big speeds in the Riemann solution (shallow water equations aren't valid for flow off a cliff).

If not specified in setrun.py, the default value speed_limit = 50 m/s (approx 100 knots) is set in data.py. This is much larger than physically reasonable speeds for most problems, so this should give consistent results with past version of GeoClaw for runs that don't die along the way.

In practice it is better to use a smaller value for most problems. Setting speed_limit = 20 has been found to work well for several tsunami inundation problems, eliminating the tiny time steps that are sometimes forced by unphysical huge speeds in one cell, and allowing jobs to run to completion that otherwise fail in sevaral tests. This speed is probably larger than fluid speeds expected in most tsunami, storm surge, or flooding problems.

But more experiments are needed to try to optimize this. And for some extreme problems, such as asteroid-impact mega-tsunamis or catastrophic dam failures, larger values of speed_limit may be appropriate.

[mjberger]

why did you set the default to 50 m/s, instead of the 20 we've tested? — Marsha

…

On Feb 17, 2025, at 4:31 PM, Randall J. LeVeque ***@***.***> wrote: Assigned #637 <#637> to @mjberger <https://github.com/mjberger>. — Reply to this email directly, view it on GitHub <#637 (comment)>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/AAGUGCZ7CGQAJQ532KZ6YOD2QJIMZAVCNFSM6AAAAABXKGE6NCVHI2DSMVQWIX3LMV45UABCJFZXG5LFIV3GK3TUJZXXI2LGNFRWC5DJN5XDWMJWGMZTSMJRGQZDONY>. You are receiving this because you were assigned.

[rjleveque]

I set it to 50 so it would be backwards compatible in case anyone is running a case where speeds are actually greater than 20 and doesn't know about this new limit. I'm not sure what the best value is in general, but I recommend setting it to 20 in setrun.py for now. I need to dig into past work more to see if there were cases where speed > 20 m/s was correct some places.

At any rate, for the few cases I've looked at so far, we got "too many dt reductions" in cells where the speed was much larger than 50 at times (supersonic in some cases), so setting it to 50 might be sufficient to keep it from crashing. Recall that it dies only if it's trying to reduce the timestep by more than a factor of 100. On the other hand, if a smaller speed_limit keeps the Courant number smaller then larger timesteps can be taken, which is good for both efficiency and accuracy. But I need to do more experimenting and I hope others will too with this WIP PR.

[kbarnhart]

@rjleveque this looks very useful. I'm testing it in some D-Claw simulations and will let you know if I have any feedback.

[mandli]

I imagine that 20 m/s would be plenty high to avoid anything in storm surge but I can check our examples.

Couple general thoughts:

1. Do we want to print this to the terminal, the log file, or both?
2. It may be nice to also have an option to abort. Even better would be to write a restart file and abort but that's probably not possible where this is.

[kbarnhart]

@rjleveque A question/thought: Why not put this in stepgrid right before src2 is called? I think this may be a better location because if there are elements of src2 that depend on hu and hv (like manning friction) it may be better to have this adjustment made immediately after step2 but before src2.

[rjleveque]

@mandli, Maybe we should just write to the log file fort.amr, eventually at least. Regarding aborting, the main point of this modification is to keep the code from aborting with "too many dt reductions" (at which point it outputs the levels where the error occurred to help with debugging, but it's not possible to do a checkpoint mid-cycle). I am hoping these huge speeds only occur in isolated cells on cliffs as in the the cases I've run so far. But I suppose there may be cases where the user wants to abort if this happens. Maybe just leave it to the user to put a stop in the Fortran code at this point if they really need to do this, since I think it would mostly be a debugging thing?

@kbarnhart: Interesting idea to do it in stepgrid before calling src2. I'll have to dig into it some more to see if that might work instead of (or perhaps in addition to) doing it in b4step2. The reason it is also in getmaxspeed has to do with ghost cells when regridding, I think. Without this addition, some cases ran fine but others still died if there was a grid patch edge right at the troublesome cell on the cliff. But @mjberger suggested we should also revisit the need for this routine and where it is called, at some point.

[kbarnhart]

manning-explore.zip
@rjleveque Thanks for the thoughts. Inspired by this conversation, I did a little work up of how manning friction (with a variety of h, u, coeffmanning, and dt values) reduces the velocity. For the types of fall-off-cliff velocities that I've seen, manning friction doesn't do much. Code an figures attached.

[mandli]

That makes sense to me. I think I was thinking about a "too many dt reductions" type thing mostly. If this is being implemented but probably never reached without explicit knowledge it's probably something for future thought.

@mandli, Maybe we should just write to the log file fort.amr, eventually at least. Regarding aborting, the main point of this modification is to keep the code from aborting with "too many dt reductions" (at which point it outputs the levels where the error occurred to help with debugging, but it's not possible to do a checkpoint mid-cycle). I am hoping these huge speeds only occur in isolated cells on cliffs as in the the cases I've run so far. But I suppose there may be cases where the user wants to abort if this happens. Maybe just leave it to the user to put a stop in the Fortran code at this point if they really need to do this, since I think it would mostly be a debugging thing?

[rjleveque]

@kbarnhart: I tried imposing the speed limit in stepgrid after the call to step2 and just before the call to src2, and I found that it was invoked in many cells in the inundation zone where the depth was very small, h < 0.01, and the speed was very high after taking a step. But if left uncorrected, all of these were damped out by the friction term with the value 0.025 we normally use, and even with values as small as 0.005. (With Manning coefficient 0.001 or no friction, the code does die.)

So I don't think it's necessary to reduce the speeds in these cells before applying friction, and in fact doesn't seem necessary because of friction (unless trying to do a frictionless run).

Moreover, I checked that doing so does not eliminate the need to also do it before the call to step2 (as now done in b4step2) or in getmaxspeed.