22 votes
Accepted

Puzzling remark about stability region of fifth-order Runge-Kutta method

van der Houwen's statement is correct, but it is not a statement about all fifth-order Runge-Kutta methods. The "Taylor polynomials" he is referring to are (as you seem to know) just the polynomials ...
user avatar
13 votes
Accepted

Numerically stable approach for calculating x in Ax=b

If the solution of $Ax=b$ is unstable, the matrix is very ill-conditioned (i.e., has a very large condition number), and (paraphrasing Lanczos) no amount of mathematical trickery can make it stable. ...
user avatar
13 votes
Accepted

Forcing an ODE solver to preserve the norm

The best approach is to use an ODE solver that is guaranteed to conserve the norm of the initial condition, i.e., for which $\|y_n\| = \|y_0\|$ for all $n\in\mathbb{N}$. Such solvers exist, and are ...
user avatar
13 votes

How much regularization to add to make SVD stable?

The singular value decomposition for a symmetric matrix $A=A^{T}$ is one and the same as its canonical eigendecomposition (i.e. with an orthonormal matrix-of-eigenvectors), while the same thing for a ...
user avatar
12 votes
Accepted

Clenshaw-type recurrence for derivative of Chebyshev series

You can just take Clenshaw's recurrence $$ u_k(x) = 2xu_{k+1}(x)-u_{k+2}(x)+\color{red}{a_k},\\ f(x) = x u_1(x)-u_2(x)+\color{red}{a_0} $$ and differentiate it directly: $$ u_k'(x) = 2xu_{k+1}'(x)-u_{...
user avatar
  • 11.4k
11 votes
Accepted

Stability of hyperbolic PDE and DG-FEM

Stability does indeed mean that small changes in the data lead to small changes in the solution. This can be shown for the linear advection equation through the energy method; in proving the stability ...
user avatar
8 votes
Accepted

Stability analysis of Heun's method

Notice that $\hat C_*=(1-r F(h))\hat C_n$, but the sign in front of $r$ is lost when you use $\hat C_*$ inside $\hat C_{n+1}$. I'm looking at p.149 of Numerical Solution of Time-Dependent Advection-...
user avatar
  • 11.4k
8 votes
Accepted

Order of operations, numerical algorithms

Let's denote by $\otimes,\oplus,\ominus$ (I was lazy trying to get circled version of division operator) the floating-point analogs of exact multiplication ($\times$), addition ($+$), and subtraction (...
user avatar
  • 8,362
8 votes
Accepted

More Smearing with decreasing timestep in advection problems

You solve the 1D-advection equation with $c$ a constant velocity: $$\frac{\partial u}{\partial t} + c\frac{\partial u}{\partial x}=0~~~~~~~~(1)$$ When you discretize this equation (with an explicit ...
user avatar
  • 649
8 votes

Why is my simulation of a first-order wave equation not stable?

Numerical solution of the advection equation with centered differences in space and forward Euler in time is unconditionally unstable. So the behavior you are seeing is expected. Here is a nice ...
user avatar
  • 5,744
7 votes
Accepted

Method to solve linear, first order ODE of generalized matrix matrix form

There are a few ways you can do this. Even though RK4 fails unless stepsizes are really small, one thing that can happen a lot is that the model can start out more stiff than it really is in the full ...
user avatar
7 votes
Accepted

What is the flaw in my stability analysis?

But clearly, this is not the case as my programs do come up with (an approximate) solution though. I believe you did not continue the integration until you see that your integration is not convergent ...
user avatar
6 votes
Accepted

What is the meaning of stability in numerical analysis? How to deterimne the stability of a numerical method?

You can define a method/algorithm stable if during the various steps it not amplifies excessively the errors on the data. So you introduce some bounded estimate. Note that this kind of errors, ...
user avatar
6 votes
Accepted

Von Neumann stability analysis with a constant term

From your link, consider the definition of the round off error and the statement "Since the exact solution must satisfy the discretized equation exactly, the error must also satisfy the discretized ...
user avatar
  • 2,199
6 votes
Accepted

Why is subtraction a stable operation?

(My apologies for the terrible typesetting) Let me start with the definition of backward stability of an algorithm: An algorithm $A$ to solve a problem $P$ is called backward stable, if for all $x$ in ...
user avatar
6 votes
Accepted

How to interpret if $\displaystyle \sum_{j = 0}^{n} \frac {1}{j!}$ is a stable algorithm for computing $e$?

As hardmath mentioned in the comments above, if data space $D$ is empty, then the first half of the condition "for all $X\in D$ there exist $Y\in D$ such that $||X-Y||/||X||\approx\epsilon_m$ and ...
user avatar
5 votes

what do zero real parts of eigenvalues mean? Any good references?

Eigenvalues with zero eigenvalue correspond to purely oscillatory modes. You can see it by diagonalising the system. Your matrix $A$ can be written as \begin{equation} A = P \Sigma P^{-1} \end{...
user avatar
  • 1,198
5 votes
Accepted

stabilizing advection-diffusion with multi-grid?

There are two different flavors of smoothing and stability. The spurious oscillations in convection-diffusion problems are not an artifact of the linear solver but an inevitable artifact of the ...
user avatar
  • 10.8k
5 votes

Finite Difference Method Stability

Analytical problem : what you are expecting is positive diffusion : you want the $T_i$ values to spread over your domain as time passes to eventually reach $T_i(t\rightarrow \infty) = cte$ if $\...
user avatar
5 votes
Accepted

Improve numeric stability of subtraction in C++

This is comment, but too long for comment so I write here. I think is better you edit the post with more information, now I explain why. Here you can find some general advice to improve precision. ...
user avatar
5 votes

What is the origin of the spurious oscillations in the Crank-Nicolson scheme?

There is something very basic that you should know about hyperbolic problems. Consider the most basic example $\partial_tu+a\partial_xu=0$ with a numerical marching scheme of the form $$u_j^{n+1}=\...
user avatar
5 votes

Lagrange multipliers space is too rich in a mathematical view

Each Lagrange multiplier corresponds to a constraint. So if the space of Lagrange multipliers is too large, then you have too many constraints that can no longer be all satisfied at the same time ...
user avatar
5 votes
Accepted

Lagrange multipliers space is too rich in a mathematical view

The saddle point matrix you must solve takes the following form: $$\begin{bmatrix} A & B^T \\ B \end{bmatrix},$$ where $A$ is the unconstrained matrix and $B$ is the matrix that observes the ...
user avatar
  • 3,013
5 votes

Numerical stability in the product of many matrices

Orthogonal matrices are about as well-conditioned as you can get, but numerical errors still occur. One common error is loss of orthogonality. A fix for this could be to re-orthogonalize your columns ...
user avatar
  • 1,097
5 votes
Accepted

Advantage of diagonal "jitter" for numerical stability?

So, you want to invert your matrix $A=\Phi^T\Phi$. For $A$ to be invertible it must not have zero eigenvalues. We can show that $A$ is positive semi-definite as follows. Positive semi-definite means ...
user avatar
  • 428
5 votes

stability of a numercial scheme for a hyperbolic system?

It is worth making some additional points. What you set out is just one version of the Lax-Wendroff method. That scheme is unique in one space dimension but has several free parameters in two or three ...
user avatar
4 votes

Finite Difference Method Stability

For stability you need $\omega \leq 1$, as stated here. If your diffusivity $\alpha$ is positive (that should be based on a physical basis), $\omega$ should also be positive since $\Delta x$ and $\...
user avatar
  • 8,006
4 votes

How to avoid negative values of numerical solution of transport equation using FEM scheme?

Typically you would use a slope limiter (or artificial diffusion and just cross your fingers) which detects where the solution has gone negative and modifies the solution to restore positivity (often ...
user avatar
  • 2,961
4 votes

Avoid arithmetic overflow in matrix multiplication

This isn't a direct answer to your question but rather an alternative approach. It looks like you are solving a least squares (LS) problem using the normal equations. The normal equations are known ...
user avatar
  • 5,744

Only top scored, non community-wiki answers of a minimum length are eligible