The SOR-Newton method

Directly to the input form

• The SOR-Newton method is meant for the solution of large nonlinear systems of equations F(x)=0.
• It is a combination of the ideas of Newton's method and the SOR method for linear systems, in the form first relax then linearize.
• Given an approximation x_k for the zero x^* in a cycle of n substeps new coordinates are computed as
  x_k+1,i = x_k,i -ω*F_i(x_k+1,1,..,x_k+1,i-1,x_k,i,..,x_k,n)/((d/d x_i) F_i(x_k+1,1,..,x_k+1,i-1,x_k,i,..,x_k,n)), i=1,...,n .
  That means we take a onedimensional Newton's step with stepsize variation for the i-th equation, considered as a function of the i-th variable only.
• The relaxation parameter ω has the same role as in the linear case.
• If F(x)=Ax-b, then this is identical to the SOR method for linear equations.
• Proofs of local convergence are possible under the same assumptions on the Jacobian of F at x^* as used in the linear case. Global convergence can be shown if F is the gradient of an uniformly convex function on Rⁿ (this means that the Jabobian is symmetric and its eigenvalues have a globally valid upper and lower bound > 0). In addition ω needs to be ''sufficiently small'' or otherwise one needs to use a test for sufficient decrease of the principal function of F. Even if this principal function is available its computation (for each k,i) is as costly as a complete cycle in this iteration. Hence this is never realized in practice, rather one tries a ''small'' omega first, tries to check the ''optimal'' one and continues. Here we work with a fixed user defined ω.
• The advantage of this method is the fact that it is completely ''matrix free'', the required diagonal of the Jacobian might well be approximated by finite differences.
• The convergence conditions however are very restrictive. Local convergence can be shown under the condition that the Jacobian of F at the solution has properties which are sufficient for proving convergence in the case of a linear system. For example for an irreducible and positive definite symmetrix tridiagonal matrix one knows that an unique optimal relaxation parameter exists which is larger than one and approximately
  ω_opt = 2/(1+2/√(cond(J_F(x^*))))
  with the convergence rate
  ρ = ω_opt -1 = 1 - O(1/√(cond(J_F(x^*)))).

Input

• For the nonlinear system
  F(x) = (F₁(x₁,...,x_n) ,..., F_n(x₁,...,x_n)) = 0
  there is a predefined case with chooseable dimension (in {60,1200}), namely the standard discretization of the nonlinear two point boundary problem
  y''(t) = γ sin(π y(t)) + t(1-t), 0 <= t <= 1 , y(0) = 1 , y(1) = 2
  by finite differences of second order consistency. This is well defined for γ <= 1/2 but for larger values it becomes illdefined. You may set γ as you wish, but don't complain if it does not work! If you know a little bit more about the theory: think about the eigenvalues of the tridiagonal matrix with the typical entries
  -1 , 2 + γ π cos(π y(j))/(n+1)² , -1 .
  You have also the possibility to provide your own testcase by code for computing the i-th component of F yourself. The partial derivatives (d/dx_i)F_i are computed by a high precision numerical differentiation (sixth order Richardson extrapolation), hence you need not to bother about that.
• In the case of a selfdefined problem input is the dimension n, four indices I in the range {1,...,n} for which you will see the graph of x_i, i in I over the stepnumber k and the initial point (x_0,1,...,x_0,n).
• in both cases:
• The relaxation parameter ω. Important: 0 < ω < 2 ! This does not necessarily yield convergence!
• The required precision ε in the norm of F. We require
  ||F||₂ <= ε ||diag(J_F)||
  for the current point for successful termination.
• The maximum number of iterations allowed.
• A bound C for the allowed groth of ||F|| (there is normally no monotonic decrease!). We diagnose ''divergence'' and give up if
  ||F(x_k)|| > C ||F(x₀)|| .
• You might want to see in addition to the graphical output a printed table of some components of x and ||F||. Since the step numbers might be quite large, you can decide to see only every k-th step.

Output

• Three 2D-Plots showing the sequence x_I for the selected indices and ||F||, and x_i=1,...,n over i.
• A table of results with this sequence if you did require this.
• A table of the final solution.

Questions?!

• When can you observe monotonic convergence in the euclidean sense?
• Which influence has the relaxation parameter ?
• In which cases you observe convergence, in which divergence?

To the input form