Newton.f90

$ \renewcommand{\vec}[1]{ {\bf #1} } \newcommand{\bnabla}{ \vec{\nabla} } \newcommand{\Rey}{Re} \def\vechat#1{ \hat{ \vec{#1} } } \def\mat#1{#1} $

This page describes usage of the Newton solver for pipe flow. See also the generic routine File:NewtonHook.f90.

Inputs

This is the contents of a typical params.in file:

 0
 150 100 1 5
 3d-6
 -1d0 1d-7 1d+7
 1d-3 1d-6

You are unlikely to need to change these settings, except perhaps the value 5, which is the number of leading eigenvalues that are calculated if newton converges. 150 is the max number of GMRES iterations. If it is taking this many iterations (per Newton iteration), convergence is unlikely. 100 is the max number of Newton iterations.

A typical shifts.in file:

  10000
0.  -1.  100.

Respectively:

• ndts - Number of timesteps to take in period T,
• sth - shift in $\theta$,
• sz - shift in $z$,
• T - period.

A value -1. for sz informs the code to search for a suitable initial guess for the shift.

The timestep size follows from the inputs, T/ndts (=0.01 in the example) which is stretched slightly if T changes. The shifts are shifts backwards, i.e. if a travelling wave has phase speed c (>0), then after time T we expect sz=c*T (>0).

The value T is arbitrary for a travelling wave, but T=10 usually works well. If the number of GMRES iterations is observed to be large (more than 50), try doubling T. If it is small (less than 10), halve T.

The initial guess should be in state.cdf.in.

Outputs

The output file of most interest is newton.dat:

       14679          150       279939
 0     0  0.21542E-02  0.21542E-03  0.25636E-01  0.84033E-01
 1    18  0.46324E-03  0.21542E-03  0.55202E-02  0.83918E-01
 2    19  0.12734E-03  0.83324E-04  0.15176E-02  0.83905E-01
 3    23  0.20499E-04  0.83324E-04  0.24427E-03  0.83920E-01

Along the top: ndts, max GMRES iterations, degrees of freedom.

Then columns: Newton iteration, GMRES iterations taken, $||g(x(T))-x(0)||$, $\delta$, $||g(x(T))-x(0)||\,/\,||x(0)||$ and $||x(0)||$, where $g(x)$ denotes the spatial shift of $x$ (backwards by sz).

The penultimate column, relative error, is the most important one. If it drops below 3d-6, this is usually considered sufficient. $\delta$ is the step size being taken (trust-region approach).

Other outputs:

• shift.dat records the shift at each iteration. Each line is like shifts.in preceded by the iteration number.

• state0???.cdf.dat is the state at each iteration.

• state1000.cdf.dat is the converged state.

• state1001.cdf.dat, state1002.cdf.dat... are eigenvectors.

Eigenvalues are found in arnoldi.dat. If the descending eigenvalues are e.g.

  real, complex, real, real, complex, real,...

then the respective state numbers for the eigenvectors are

  1001, 1002+3, 1004, 1005, 1006+7, 1008,...

where 1002+3 are real and imaginary components respectively.

Continuation of solutions

If the file params_ct.in is present, the code will attempt continuation from the input state.cdf.in. A typical params_ct.in is:

 2 3 1
  -50.

On the first line, 2 and 3 are the minimum and maximum number of newton iterations per continuation step, and 1 indicates the parameter to be continued, 1==$Re$, 2==$\alpha$. On the second line -50 is the size of the first step in the chosen parameter from the (assumed to be converged) input state.cdf.in.

The outputs from each continuation step are moved to directories 0001, 0002, 0003... and a brief summary of the data in each is written to newton.dat_ct. Its columns are

number of continuation step (= directory number), 
number of newton iterations required, 
$Re$, $\alpha$, wavespeed $c$.

Continuation restarting from two previous steps is possible. If state10.cdf.in is present, then for example, params_ct.in

 2 3 1
   1.5

asks the code to start from a linear extrapolation: if state.cdf.in=$u_0$ and state10.cdf.in=$u_1$, then it would restart with initial guess $u=u_1+ 1.5*(u_1-u_0)$. The corresponding extrapolation for the initial shifts should be supplied in shifts.in.

Note: actual shifts should be given for the starting point. Do not use the -1. flag in shifts.in.