<p>The usual panacea for numeric stability problems is to decrease the timestep. It doesn't always help (and, indeed, in some cases can even make the instability worse), but it's generally the first thing to try.</p>

<p>Of course, a shorter timestep costs more computer time.  Another approach to try, if reducing the timestep doesn't help, or leads to unacceptably poor performance, is to either refine or, somewhat paradoxically, coarsen the mesh in the problematic regions.  (Which one of these is more likely to help depends on the exact cause of the instability.)  If this alone doesn't help either, try adjusting the mesh <em>and</em> decreasing the timestep.</p>

<p>Ultimately, the <em>best</em> way to fix such problems is usually to switch to a more stable timestepping method.  As a general rule of thumb, higher-order methods are more stable than lower-order methods (up to some point, anyway), and implicit methods are more stable than explicit.  In both cases, stability comes at the cost of extra computational effort, but it's still generally cheaper than just decreasing the timestep enough to achieve comparable stability.</p>

<p>Unfortunately, I'm not personally familiar with the MIKE3 software, and have no idea what timestepping algorithm(s) it uses, or whether there's any way for the end user to adjust the algorithm.</p>

CC BY-SA 3.0
2014-04-19T00:02:24.037Z
Ps. Sorry for giving such a vague answer, but I figured that, in case no-one else comes up with a better answer, it might at least be better than nothing.
Ilmari Karonen
CC BY-SA 3.0
2014-04-19T05:46:34.937Z
I am not sure if higher-order methods are generally cheaper than smaller time steps.
BHF
CC BY-SA 3.0
2014-04-19T05:52:03.247Z
@BHF: Sometimes they are, sometimes they aren't. But if they weren't, often enough, more efficient than just decreasing the timestep,  there wouldn't  be much point in developing them, would there?
Ilmari Karonen
CC BY-SA 3.0
2014-04-19T06:10:01.157Z
I agree. By the way, stability without worrying about the time step could be such a reason, independently of computational costs.
BHF
CC BY-SA 3.0
2014-04-19T00:00:45.507Z
2014-04-19T00:00:45.507Z
Ilmari Karonen
<p>Almost all finite methods that use forward time models adhere to the Courant-Friedrichs-Lewy law which calculates a courant number and compares it to a $C_{max}$, which is what determines stability, so for 2-D:</p>

<p>$$C = \frac {V_xdt}{dx}  + \frac{V_ydt}{dy} \ge C_{max}$$</p>

<p>Where $C$ is the courant number, $V_i$ is the velocity in subscripted direction, $dx$ or $dy$ is the length interval in specified direction and $dt$ is the specified time step. If you are getting numerical instability, increasing the mesh spacing will lower the courant number, but then you might be missing the dynamics you want to capture. If you decrease the time step, it will take longer for the dynamics to develop computationally. </p>

<p>One way to solve this computational problem is to code in a dynamic time step: Keep calculating the $C$ and $C_{max}$. If $C &gt; C_{max}$, then have the model lower the step by half. Check again, then reset for the next step so you can have larger timesteps for "less dynamic" times. I have implemented this many times in other models and it has worked well.</p>

<p>Another possibility might be that you have high aspect ratios in your mesh: if your $x$ direction is 10 km and your $y$ is 1 km, a $10:1$ aspect ratio is very high. In your case, since the model package you use are triangles, you want the angles in the triangle grid to be as acute as possible. Obtuse angles would mean high aspect ration, in general. </p>

<p>Like llmari, I am not familiar with hydrodynamic models, so I do not do specifics. My expertise are in both mantle convection and plate flexure models. I hope this helps.</p>

CC BY-SA 3.0
2014-04-20T07:23:42.950Z
I think you have acute/obtuse backwards :-) Actually, the best triangular mesh will tend to have equilateral triangles. Both are good answers, but accepting this one since it's higher-voted :-)
Semidiurnal Simon
CC BY-SA 3.0
2014-04-19T00:18:55.133Z
2014-04-24T15:00:07.403Z
2014-04-24T15:00:07.403Z
Neo
<p>I am working with a coastal hydrodynamic model in the MIKE3 modelling package. MIKE3 uses the finite volume method to solve the shallow water equations with explicit timestepping, over a triangular mesh.</p>

<p>My model frequently "blows up" (numerically diverges) in a few cells around a headland, at times of peak flow. What steps can I take to improve the numerical stability, either through modifications to the mesh or other measures?</p>

<oceanography><models><coastal>
How can I improve the numerical stability of a coastal hydrodynamic model?
Semidiurnal Simon
# How can I improve the numerical stability of a coastal hydrodynamic model?

(Asked by: username_0 on 2014-04-18)

I am working with a coastal hydrodynamic model in the MIKE3 modelling package. MIKE3 uses the finite volume method to solve the shallow water equations with explicit timestepping, over a triangular mesh.

My model frequently "blows up" (numerically diverges) in a few cells around a headland, at times of peak flow. What steps can I take to improve the numerical stability, either through modifications to the mesh or other measures?

**Tags:** oceanography, models, coastal


---

## 2 answers

### Answer by username_3 on 2014-04-19. Score: 12 | Accepted answer

Almost all finite methods that use forward time models adhere to the Courant-Friedrichs-Lewy law which calculates a courant number and compares it to a $C\_{max}$, which is what determines stability, so for 2-D:

$$C = \frac {V\_xdt}{dx} + \frac{V\_ydt}{dy} \ge C\_{max}$$

Where $C$ is the courant number, $V\_i$ is the velocity in subscripted direction, $dx$ or $dy$ is the length interval in specified direction and $dt$ is the specified time step. If you are getting numerical instability, increasing the mesh spacing will lower the courant number, but then you might be missing the dynamics you want to capture. If you decrease the time step, it will take longer for the dynamics to develop computationally.

One way to solve this computational problem is to code in a dynamic time step: Keep calculating the $C$ and $C\_{max}$. If $C > C\_{max}$, then have the model lower the step by half. Check again, then reset for the next step so you can have larger timesteps for "less dynamic" times. I have implemented this many times in other models and it has worked well.

Another possibility might be that you have high aspect ratios in your mesh: if your $x$ direction is 10 km and your $y$ is 1 km, a $10:1$ aspect ratio is very high. In your case, since the model package you use are triangles, you want the angles in the triangle grid to be as acute as possible. Obtuse angles would mean high aspect ration, in general.

Like llmari, I am not familiar with hydrodynamic models, so I do not do specifics. My expertise are in both mantle convection and plate flexure models. I hope this helps.

> **Comments:**
> **username_0**: I think you have acute/obtuse backwards :-) Actually, the best triangular mesh will tend to have equilateral triangles. Both are good answers, but accepting this one since it's higher-voted :-)  

---

### Answer by username_1 on 2014-04-19. Score: 10

The usual panacea for numeric stability problems is to decrease the timestep. It doesn't always help (and, indeed, in some cases can even make the instability worse), but it's generally the first thing to try.

Of course, a shorter timestep costs more computer time. Another approach to try, if reducing the timestep doesn't help, or leads to unacceptably poor performance, is to either refine or, somewhat paradoxically, coarsen the mesh in the problematic regions. (Which one of these is more likely to help depends on the exact cause of the instability.) If this alone doesn't help either, try adjusting the mesh *and* decreasing the timestep.

Ultimately, the *best* way to fix such problems is usually to switch to a more stable timestepping method. As a general rule of thumb, higher-order methods are more stable than lower-order methods (up to some point, anyway), and implicit methods are more stable than explicit. In both cases, stability comes at the cost of extra computational effort, but it's still generally cheaper than just decreasing the timestep enough to achieve comparable stability.

Unfortunately, I'm not personally familiar with the MIKE3 software, and have no idea what timestepping algorithm(s) it uses, or whether there's any way for the end user to adjust the algorithm.

> **Comments:**
> **username_1**: Ps. Sorry for giving such a vague answer, but I figured that, in case no-one else comes up with a better answer, it might at least be better than nothing.  
> **username_2**: I am not sure if higher-order methods are generally cheaper than smaller time steps.  
> **username_1**: @username_2: Sometimes they are, sometimes they aren't. But if they weren't, often enough, more efficient than just decreasing the timestep,  there wouldn't  be much point in developing them, would there?  
> **username_2**: I agree. By the way, stability without worrying about the time step could be such a reason, independently of computational costs.  

---

