User david ketcheson - MathOverflow

Answer by David Ketcheson for Global Error Analysis of Euler's Method

2012-10-12T11:54:43Z

Any analysis of global error must include information about how local errors are amplified in subsequent steps. So your statement

I know that the local error at each step of Euler's method is O(t^2), where t is the time step. And since there are (b-a)/t steps, the order of the global error is O(t).

isn't accurate without some assumption of stable error propagation.

You are correct that the "derivation of the global error" given does not say anything about global error.

Answer by David Ketcheson for Is is preferable to use a difference formula of higher order of accuracy for spatial derivatives to solve this IVP problem ?

2012-10-10T18:17:11Z

You could use the method of lines to solve this PDE. If you use an explicit finite difference method, you will need to take a rather small time step (${\mathcal O(\Delta x^3))}$ due to the $u_{xxx}$ term.

Given that you don't specify any boundary conditions, I will assume that you are solving the Cauchy problem. In that case (or in case of periodic boundary conditions), there is a much more efficient approach. Since your PDE is linear and the coefficients don't vary in space, you can solve in terms of Fourier modes. Decompose the initial data as usual:

$$u(x,0) = \sum_k c_k\exp{ikx}.$$

Then use the ansatz

$$u(x,0) = \sum_k g_k(t)\exp{ikx}$$

with $g_k(0) = c_k$. Substituting this in your PDE gives

$$g_k'(t) = ik a(t) - ik^3b + c.$$

This can be solved very easily by numerical quadrature (or exactly, if $a(t)$ can be integrated symbolically).

Answer by David Ketcheson for "Must read" papers in numerical analysis

2012-07-17T08:23:02Z

S.K. Godunov, A difference method for numerical calculation of discontinuous solutions of the equations of hydrodynamics, Matematicheskii Sbornik (1959). (in Russian)

Virtually all modern methods for nonlinear hyperbolic PDEs are based on this.

Differentiability of minimax objective function with respect to a decision variable

2012-06-29T13:02:37Z

I have the following optimization problem:

$$\text{find } x= \min_{a} \max_{\lambda\in\Lambda} |R(h\lambda)|$$

where $\Lambda$ is some finite, fixed set of complex numbers and $R(z)$ is a polynomial of degree at most $s$ of the form

$$R(z) = 1 + z + \sum_{j=2}^{s}a_j z^j$$

where $h>0$ is fixed and the real coefficients $a_j$ are the decision variables. Can it be shown that $x$ is a differentiable function of $h$? Are there general tools that might be useful for proving this type of property?

Answer by David Ketcheson for sum of maxima vs the maximum of the sum

2012-07-11T09:44:17Z

In some cases there exists no ordering of $\Gamma$ such that solving the problems in $\Gamma$ sequentially gives the optimal value of the original problem.

Here is a simple counterexample. Let $i,j$ range from 1 to 2 and set $c_j=1, f(i,j)=1$ for all $i,j$. Take $U_{i,j}=1$ for $i=j$ and $U_{i,j}=2$ for $i\ne j$. Then the original problem is solved by taking $x_{1,1}=x_{2,2}=0$ and $x_{1,2}=x_{2,1}=1$, giving an objective function value of 4. But solving $\Gamma$ in either order leads to an objective function value of 3.

Let's go through it step by step.

The problem is $$\max x_{11} + x_{22} + 2x_{12} + 2x_{21}$$ subject to $$x_{11}+x_{21}\le 1$$ $$x_{12}+x_{22}\le 1$$

For the other two problems, we have

$$\max x_{11} + 2x_{12}$$ subject to $$x_{11}+x_{21}\le 1$$ $$x_{12}+x_{22}\le 1$$ which leads to $x_{12}=1$ AND $x_{11}=1$, in order to achieve a maximum value of 3. Notice that the solution to this problem in the comment below is not correct.

Then the second problem becomes $$\max 2x_{21} + x_{22}$$ subject to $$x_{21}\le 0$$ $$x_{22}\le 0.$$

Answer by David Ketcheson for Basic software libraries for numerical analysis using modern programming languages?

2012-07-04T18:43:37Z

There are several good answers already, but none of them so far deal with high performance computing, which nowadays means computing on tens or hundreds of thousands of cores. I am aware of only three Python codes that have scaled to this level:

PyClaw (my code -- you can read about it and my take on HPC software development in Python in this paper)
GPAW
FiPy

Of course, you can't write a performant numerical code entirely in Python, so typically these use a compiled language for performance-critical parts (in PyClaw, we use both Fortran and C). One reason for choosing Python is that automated tools like f2py and Cython make it easy to incorporate compiled code.

I do not believe there are any high performance computing codes written in the languages you mention. I suspect modern supercomputers don't support those languages.

Answer by David Ketcheson for what method can I employ to solve this optimization problem which involves \min?

2012-07-04T09:06:39Z

A common transformation when faced with a problem of this type:

$${\rm maximize} \min (f(x), g(x))$$

is to instead solve the equivalent problem

$${\rm maximize} \ \ \ z $$ subject to $$ z\le f(x); z\le g(x).$$

This can be helpful, for instance, in making the problem more tractable for some numerical optimization methods (if they are better at handling inequality constraints rather than complicated objective).

Many other transformations are possible; for instance, you can replace the inequality constraints above with equality constraints via slack variables. For a discussion of transformations of optimization problems, I recommend Section 4.1.3 of Boyd and Vandenberghe (free PDF here). Actually, I recommend the whole book.

Answer by David Ketcheson for convergence of finite difference method for boundary value ODE

2012-06-22T22:06:54Z

A very straightforward explanation is given in R.J. LeVeque's text. In Chapter 2 there are simple explanations of how to show convergence for the linear problem in both the maximum norm and the Euclidean norm. There is also a discussion of the nonlinear problem. Numerical solution of nonlinear BVPs typically requires a nonlinear iterative solver, such as Newton's method, and the numerical discretization inherits the properties of Newton's method (i.e., it is only guaranteed to converge if you start sufficiently close to a solution, and the solution it converges to depends on the initial guess.

The nonlinear example considered in the text is $y''(x)=\sin(x)$, and it is demonstrated that the numerical solution is non-unique (as is the analytic solution).

Answer by David Ketcheson for numerically track spectrum curves of a parameter dependent linear operator

2012-01-13T06:43:09Z

I faced this problem a few years ago. In that case, I obtained a satisfactory approach along the lines of one of your suggestions. Specifically, I found that the eigenvectors changed relatively slowly with $t$. So I could associate the corresponding eigenvalue curves after an intersection with the right ones before an intersection by considering the similarity of the corresponding eigenvectors. I found that the inner-product $X^T_i Y_j$, where $X_i$ is one of the eigenvectors after the intersection and $Y_j$ is one of the eigenvectors after the intersection was a good measure of their similarity.

By the way, you might get more answers by posting this on http://scicomp.stackexchange.com/

Comment by David Ketcheson

2013-05-07T15:25:54Z

The links are dead.

Comment by David Ketcheson

2012-12-06T08:35:56Z

What does "*" mean here? If just multiplication, it would be clearer to omit the symbol.

Comment by David Ketcheson

2012-11-14T16:57:58Z

This is neither numerical analysis nor functional analysis (unless there is some application I'm missing).

Comment by David Ketcheson

2012-11-14T16:56:26Z

You might ask on scicomp.stackexchange.com.

Comment by David Ketcheson

2012-11-09T12:10:30Z

In the literature, such matrices (without the sparseness requirement) are referred to as $L$-matrices. If they are diagonally dominant, they are referred to as $M$-matrices. Many things are known about such matrices, so those keywords may help.

Comment by David Ketcheson

2012-10-22T04:02:14Z

You might try asking on scicomp.stackexchange.com.

Comment by David Ketcheson

2012-10-22T04:00:08Z

This certainly won't give all solutions.

Comment by David Ketcheson

2012-10-18T10:28:17Z

@Pat M: sorry, I missed the "all solutions" part. For general nonlinear functions, I believe there is no such algorithm.

Comment by David Ketcheson

2012-10-12T11:55:13Z

This question would be more appropriate on Math.SE, as it pertains to undergraduate numerical analysis.

Comment by David Ketcheson

2012-10-10T18:19:55Z

I recommend asking this question on scicomp.stackexchange.com.

Comment by David Ketcheson

2012-09-13T10:15:22Z

You may want to ask this question at scicomp.stackexchange.com.

Comment by David Ketcheson

2012-07-14T07:16:40Z

@john, Your solution of this problem (in comments) is not correct. I have made it completely explicit for you now.

Comment by David Ketcheson

2012-07-12T10:16:00Z

I suggest asking on scicomp.stackexchange.com.

Comment by David Ketcheson

2012-07-12T10:14:35Z

I suggest asking on scicomp.stackexchange.com. But you'll want to clarify the question.

Comment by David Ketcheson

2012-07-04T18:32:51Z

Supercomputers don't have Java compilers. Of course, one may argue about whether this is a cause or an effect.