Analysis of an Adaptive Safeguarded
Newton-Anderson Algorithm of Depth One with Applications to
Fluid Problems
Matt Dallas
111Corresponding author: Matt Dallas. MD and SP were
supported in part by the US NSF
Project DMS 2011519 (PI: Pollock). LR was supported by
US NSF Project DMS 2011490 (PI: Rebholz).
Keywords:
Anderson acceleration, Newton’s Method, safeguarding,
singular points, bifurcation.
2020 Mathematics Subject Classification: 65J15,1
Sara Pollock2
Leo G. Rebholz3 1Department of Mathematics, University of Dallas, USA
(mdallas@udallas.edu)
2Department of Mathematics, University of Florida,
USA (s.pollock@ufl.edu)
3Department of Mathematics, Clemson University, USA (rebholz@clemson.edu)
(August 21, 2024)
Abstract
The purpose of this paper is to develop a
practical strategy to accelerate Newton’s method in
the vicinity of singular points.
We present an adaptive safeguarding scheme
with a tunable parameter, which we call
adaptive -safeguarding, that
one can use in tandem with Anderson acceleration
to improve the performance of Newton’s
method when solving
problems at or near singular points. The key features of adaptive
-safeguarding are that it converges locally
for singular problems, and it can detect nonsingular problems
automatically,
in which case the Newton-Anderson
iterates are scaled towards a standard Newton step.
The result is a flexible algorithm that performs well for singular and nonsingular problems, and can recover convergence from both standard Newton and Newton-Anderson with the
right parameter choice.
This leads to faster local convergence compared to both
Newton’s method, and Newton-Anderson without safeguarding,
with effectively no additional computational cost.
We demonstrate three strategies
one can use when implementing Newton-Anderson and
-safeguarded Newton-Anderson to solve parameter-dependent
problems near singular points. For our benchmark problems, we
take two parameter-dependent incompressible flow
systems: flow
in a channel and Rayleigh-Bénard convection.
1 Introduction
Nonlinear systems of equations
of the form , with ,
arise frequently in applications.
Many times the solutions depend on parameters
(see [6, 7, 33, 39, 40, 41, 53]
and references therein) that can have a significant effect
on the solution.
Of particular interest are bifurcation points, which
are characterized by the breakdown of
local uniqueness
of a solution for a particular parameter, and
correspond to a qualitative change in the
solution’s behavior [6, 32, 39].
Sets of solutions of similar qualitative
behavior are called branches [39].
Studying these branches provides a
more complete understanding of the solutions,
and has applications in a wide variety of fields such as
echocardiography [41], economics [56], physics [40], and engineering
[57].
A necessary condition for bifurcation comes from
the Implicit Function Theorem, which says that
the Jacobian at a solution , , is
necessarily singular
if is a bifurcation point [32, p. 8].
Bifurcation points are thus examples of singular points, i.e.,
points for which is singular.
Techniques for computing solution branches such
as continuation [53] and deflation
[7, 17] require many solves of
as the parameter is varied, and these problems become
singular or nearly singular near bifurcation points.
A popular method for solving nonlinear equations
is Newton’s method defined in Algorithm 1.
1: Choose .
2:for k=1,2,… do
3:
4:
5:endfor
Algorithm 1 Newton
When is Lipschitz continuous and is
nonsingular,
Newton’s method exhibits local quadratic convergence in a ball
centered at . This is essentially the celebrated
Newton-Kantorovich theorem [35].
If we remove the nonsingular assumption and let
be singular, then the convergence
behavior changes dramatically.
Rather than local quadratic convergence from any in a
sufficiently small ball around , we see local linear convergence
in a starlike domain of convergence
around [11, 12, 22, 48, 49]. Since bifurcation points are necessarily
singular points, this means
that continuation or deflation algorithms using Newton’s method
may converge slowly or fail to converge at or near a bifurcation
point. This challenge has motivated the study of modifications
[11, 13, 18, 19, 23, 22, 27, 31]
or alternatives
[3, 4, 5, 16, 30, 50] to Newton’s
method that can improve convergence behavior at singular points.
Among the modifications, Richardson extrapolation and overrelaxation
are popular and perform well. They can achieve superlinear and
arbitrarily fast linear convergence respectively under certain
conditions [23] at the cost of additional
function evaluations and some knowledge of
the order of the singularity [22, 23]. The order may be inferred from
monitoring the singular values of the Jacobian as the solve
progresses. A popular alternative to Newton’s
method for singular problems is the Levenberg-Marquardt method
[3, 4, 5, 16, 30]. Under standard
assumptions and the
local error bound, or local Lipschitzian error bound, the
Levenberg-Marquardt method can achieve local quadratic
convergence [30]. The local error bound is known to be much weaker
than the standard nonsingularity assumption. Indeed, it can hold
even for singular problems [4]. However, without
the local error bound or nonsingularity,
Levenberg-Marquardt is not guaranteed
such success. In the absence of the local error bound, one
may insist that the function is 2-regular
[18, 19, 28, 29], in which case
Levenberg-Marquardt converges locally linearly in a starlike domain
much like Newton’s method [28].
The method of interest
in this paper, Anderson acceleration, has a long track record of
accelerating linearly converging fixed-point methods, and
has been applied in many different fields
[1, 34, 37, 38, 44, 51, 55].
Further, when applied as a modification
to Newton’s method at singular points,
in contrast to Richardson extrapolation or overrelaxation
discussed in [23], this success requires
no knowledge of the order of the root or additional function
evaluations. It was also shown recently that under a condition
equivalent to 2-regularity (discussed further in Section
2), Anderson accelerates Newton’s
method when applied to singular problems [10],
and therefore outperforms Levenberg-Marquardt in
the absence of the local error bound. This was demonstrated
numerically in [10].
Also known as
Anderson extrapolation or Anderson mixing [54], Anderson
acceleration was first introduced in 1965 by
D.G. Anderson in [2] to improve the convergence of fixed-point iterations applied to
integral equations. The algorithm combines the previous
iterates and update steps into a
new iterate at each step of the solve. The number is
commonly known as the
algorithmic depth. The combination of the iterates often
involves solving a least-squares problem,
but since is typically small, the computational cost of this step is in general orders of magnitude less than that of a single linear solve.
There are problems for which taking much larger can be beneficial [42, 54], and later in this paper we will see numerically
that increasing
can improve or recover convergence near bifurcation points.
With greater depths, the least-squares problem
may suffer from ill-conditioning if proper care
is not taken in the implementation [43].
In [10], the authors developed a convergence and
acceleration theory for Anderson accelerated
Newton’s method with depth , defined
in Algorithm 2,
applied to singular problems. This is a special case of depth
given in Algorithm 3.
1: Choose . Set , and .
2:for k=1,2,… do
3:
4:
5:
6:endfor
Algorithm 2 Newton-Anderson(1)
1: Choose and .
Set , and .
2:for k=1,2,… do
3:
4:
5:
6:
7:
8:
9:endfor
Algorithm 3 Newton-Anderson(m)
A major challenge when proving
convergence of Newton-like methods near singular points
is ensuring that the iterates remain well-defined.
The authors in [10] introduced a novel
safeguarding scheme called
-safeguarding, defined in Algorithm 4
in Section 2,
to deal with this problem.
The result was a convergence proof
for -safeguarded Newton-Anderson,
and
it was
observed numerically to
perform better or no worse
compared to standard
Newton-Anderson, particularly when applied to nonsingular
problems.
The purpose of this paper is to extend these
ideas of
[10] by developing
an adaptive version of
-safeguarding
that automatically
detects nonsingular problems, and to demonstrate the effectiveness
of Newton-Anderson and adaptive
-safeguarded Newton-Anderson at solving parameter-dependent
PDEs near bifurcation points in fluid problems.
This new adaptive scheme is
proven to be locally convergent under the
same conditions as the non-adaptive scheme,
and the automatic detection of nonsingular problems
enables local quadratic
convergence
when applied to nonsingular problems.
Such a property is desirable
when solving nonlinear problems near bifurcation points, because even if the problem
itself is not singular, convergence can still be affected if it is close to a singular
problem [14].
One would like to enjoy the benefits of
Newton-Anderson in the preasymptotic regime such as
a larger domain of convergence [45], but not
lose quadratic convergence in the
asymptotic regime if the problem is nonsingular. It
is not always known a priori if a problem is
singular or nonsingular, and
Anderson acceleration can reduce the order of convergence
when applied to superlinearly converging iterations
such as Newton’s method applied to a nonsingular
problem [47]. This problem is solved with adaptive
-safeguarded Newton-Anderson
with effectively no additional computational cost.
We also show numerically
that increasing the algorithmic depth
of the Newton-Anderson algorithm can recover convergence
when Newton fails
near bifurcation
points, but only for specific choices of .
The algorithms of interest in this paper are
Newton, defined in Algorithm 1, Newton-Anderson
with algorithmic depth (NA())
defined in Algorithm 2,
Newton-Anderson with algorithmic depth (NA())
defined in Algorithm 3,
-safeguarded Newton-Anderson (NA()), defined in Algorithm 5, and adaptive -safeguarded Newton-Anderson (), defined in
Algorithm 7.
We implement NA() and by replacing line
5 in Algorithm 2 with, respectively, -safeguarding (Algorithm 4) and
adaptive -safeguarding (Algorithm 6).
The norm in the algorithms is the Euclidean norm.
There should be no confusion
between NA() and NA() or since -safeguarding is currently
only developed for depth . The rest of the paper
is organized as follows. In Section 2, we
review the original -safeguarding algorithm and its
role in the convergence theory developed in [10].
In Section 3,
we introduce the new adaptive -safeguarding algorithm
and prove that can recover local quadratic
convergence when applied to nonsingular problems
in Corollary 3.1.
We conclude in Section
4 by applying NA and to two
parameter-dependent incompressible flow systems, and discussing
various strategies one can use when employing -safeguarding.
2 The -safeguarding algorithm
For the theory discussed in Section 2 and Section 3,
unless stated otherwise, we take
to be a function,
, , , , and and
to be the
orthogonal projections onto and respectively.
The assumption that is common in the
singular Newton literature
[11, 21], and no generality
is lost in finite dimensions.
Indeed, Newton’s method is essentially
invariant under nonsingular affine transformations of the
domain and nonsingular linear transformations of the range [15]. Thus to determine the convergence behavior
of Newton’s method applied to a general
function , it suffices to study
that of , where and come from the
single value decomposition . Since
, it follows that
.
Lastly, let denote the Euclidean 2-norm, denote a ball of radius centered at ,
, ,
and
(1)
where is computed via Algorithm 2.
The term
is known as the optimization gain,
and is key to determining when Anderson acceleration
is successful both in the singular and nonsingular
cases [10, 42].
Algorithm 4 -safeguarding
1: Given , , , , , and . Set .
2:
3:ifor then
4:
5:elseifthen
6:
7:endif
8:
1: Choose and . Set , and
2:for k=1,2,… do
3:
4:
5:
6:
7:ifor then
8:
9:elseifthen
10:
11:endif
12:
13:endfor
Algorithm 5 -Safeguarded Newton-Anderson
To implement
NA, one replaces line 5 in Algorithm 2
with Algorithm 4 giving Algorithm
5.
The resulting steps, , can be viewed as
NA steps scaled by towards a Newton step based on a
user-chosen parameter, , which is set at the start of the solve. This parameter determines
how strongly the NA steps are scaled towards a Newton step,
i.e., how close a NA step is to a
Newton step.
When ,
the NA iterates will be heavily scaled towards the latest
Newton step, and when , the NA
iterates will behave more like
standard NA. The idea behind
-safeguarding is to take advantage of the particular
convergence behavior of Newton’s method near singular points, which
we now describe.
Since , is spanned by
some nonzero . If the linear operator
is nonsingular as a map from to ,
then there exists
and such that
is nonsingular for all and Newton’s method
converges linearly to from any
[11].
The assumption that is nonsingular as a linear map on
is equivalent to the assumption that is 2-regular at
in the direction
[19]. A stronger result due to Griewank
[22, Theorem 6.1] says that when is
nonsingular, Newton’s method converges from every
in a starlike region with density 1 with respect to ,
and the iterates lead into provided
is sufficiently small. Thus for our purposes in this
paper studying local convergence of NA, it suffices to study
the behavior of iterates in . Though the focus of this work is the case, we note that Griewank’s Theorem [22, Theorem 6.1] holds for
, and numerical experiments from [10] on
small-scale problems indicate that NA and are effective when . While the current theoretical results
for require , extensions of these results for will be studied in future work.
The main challenge in accelerating Newton’s method is ensuring
the iterates remain in , which requires
to remain bounded. In
other words, the iterates can’t be accelerated “too much”
along the null space.
Evidently, Anderson may accelerate Newton significantly,
especially in when .
This is demonstrated by the following proposition. We will use
the notation . Similarly, , , and
. Here
denotes a linear map defined in the proof of Proposition 2.1. We also note that
is nonsingular if and only if is nonsingular for all
with .
Proposition 2.1.
Let , , nonsingular, and so that
is well-defined. If and ,
then for sufficiently small and
there is a constant such that
(2)
where .
Proof.
First note that .
Since by Algorithm 2, and ,
we have that for any ,
(3)
Taking gives .
By Proposition 3.1 in [10], we can write
(4)
where
is a linear
map whose range lies in [10], and for
a constant determined by . Hence
We also have from Proposition 3.1 in [10] that
. Let . This expansion of
combined with
Equation (4) gives
(5)
where we define . One can show that as and
tend to zero [10]. Thus for sufficiently small and we have
. This completes the proof.
∎
Note that can be very large when ,
which may correspond to the terminal phase of the solve. So this bound is meaningful in
the asymptotic regime when there is still a significant decrease in the residual at each step.
Such a bound is good for a single step,
but this dramatic acceleration
of could place outside the
domain of invertibility , i.e., may
stray too far from .
This is where -safeguarding is useful.
It ensures
that the NA iterates remain within
by taking advantage of the way Newton steps
are attracted to and scaling NA steps towards
Newton steps when the conditions of Algorithm 4
are met. This is also the key to the convergence proof of
NA given in [10].
However, given the results of [47], if the problem is
nonsingular one should use Newton’s method without
Anderson, but it is not always obvious a priori if the
problem at hand is singular or nonsingular.
In
the next section, we develop an adaptive version
of NA
that enjoys guaranteed local convergence when
applied to singular problems, but can also detect nonsingular
problems automatically, at no additional computational
cost, and “turn off” NA in response.
This leads to local quadratic convergence if the problem is
nonsingular.
3 Adaptive -safeguarding
It was observed in [10] that
NA performed competitively
with standard NA(1). For example, NA(0.5)
could outperform NA(1) when
applied to certain nonsingular problems. This is not surprising given
the results of [47], and it had previously been
observed numerically in [45] that
NA does not necessarily improve convergence when applied to nonsingular problems.
Thus, by setting the
-safeguarding parameter , thereby scaling the
iterates
closer to Newton steps, we see less of the effect of a full NA step on the
order of convergence.
Even if it is known that the problem is nonsingular, one
may still wish to use NA to take advantage of the larger
domain of convergence [45].
In such a scenario one could set close to zero, but this means
the effect of NA is
never completely eliminated which could lead to a smaller order of convergence
relative to Newton for nonsingular problems.
This, and the fact that often it is not known a priori
if the problem is singular or nonsingular, motivates the development of an
adaptive form of NA() that can automatically
detect nonsingular problems, and scale an NA step accordingly,
without sacrificing local convergence and
acceleration for singular
problems. We will denote this
adaptive choice of from Algorithm 4 by , with denoting the iteration count.
There are three criteria that the choice of
should satisfy, which we record in Criteria 3.1.
Criteria 3.1.
An adaptive -safeguarding tolerance should satisfy the following.
0.
if ;
0.
if ; and
0.
if is nonsingular.
Criterion 3.1. ‣ 3.1 says that if is very small, then
we want to scale the NA step generated from and heavily
towards . In the singular case, this will (locally) keep
within the domain of invertibility. Alternatively, if the problem is nonsingular, but close to a singular problem,
then scaling NA towards a Newton step is also preferred when
is small since then we do not
slow Newton’s fast local quadratic convergence.
Criterion 3.1. ‣ 3.1 says
that if
is close to one, then the error is not decreasing
significantly, and we want to allow
NA to act on and without significant scaling of from safeguarding.
Lastly, Criterion 3.1. ‣ 3.1 is important because if is nonsingular,
Newton’s method will converge quadratically in a neighborhood of . We therefore
want to “turn off” NA near , and insisting that
asymptotically achieves this.
Observing Criteria 3.1. ‣ 3.1- ‣ 3.1, one may note that we essentially want to
behave like
within the domain of convergence .
Of course, we can not compute , but
Equation (4) says that in ,
. So
if we take , we can expect
Criteria 3.1. ‣ 3.1- ‣ 3.1 to
be enforced locally. For Criterion 3.1. ‣ 3.1, if is nonsingular, then locally we will
have . Thus satisfies the
three criteria within the domain of convergence. With
this choice of , we have adaptive -safeguarding and .
Algorithm 6 Adaptive -safeguarding
1: Given , , , , , and , set
, , and .
2:
3:ifor then
4:
5:elseifthen
6:
7:endif
8:
1: Choose and . Set and
2:for k=1,2,… do
3:
4:
5:
6:
7:
8:
9:ifor then
10:
11:elseifthen
12:
13:endif
14:
15:endfor
Algorithm 7 Adaptive -Safeguarded Newton-Anderson
Adaptive
-safeguarding differs from
Algorithm 4 only in line 2. In Algorithm 4, whereas in Algorithm
6 with .
This one
change can have a significant impact on convergence as demonstrated
in Section 4, and can enable locally quadratic
convergence when applied to nonsingular problems (see Corollary
3.1). Similar to NA,
one implements by replacing line 5 in
Algorithm 2 with Algorithm 6
and setting at the start of the solve. The
result is Algorithm 7.
The choice of here sets the weakest safeguarding the user wants to impose.
Thus we are always safeguarding at least as strictly as standard -safeguarding
with . Stated concisely, we have .
Local convergence of then follows from
Theorem in [10]. To state this precisely, let
be the value of computed by Algorithm 6 at step ,
, , and
.
Then we have the following theorem.
Theorem 3.1.
Let , and let be invertible as a map on .
Let . If is
chosen so that and , for
sufficiently small and , , and
for ,
then for all and .
That is,
remains well-defined and converges to . Furthermore,
there exist constants and such that
(6)
(7)
for all .
Under the assumptions of Theorem 3.1, Griewank’s Theorem [22, Theorem 6.1] says that Newton’s method almost surely leads into the domain of
convergence provided is sufficiently close to . This effectively means that if the sequence generated by approaches we will have almost sure convergence eventually, and this convergence will be faster than Newton. The precise improvement is determined by the asymptotic behavior of . Globalization techniques such as linesearch methods may be used to bring the iterates
closer to . In particular, with an Armijo linesearch was shown to be effective in [10].
Our choice of in Algorithm 6 is partially motivated by Criterion 3.1. ‣ 3.1.
That is, in the case of a nonsingular
problem, we prefer to use Newton asymptotically rather
than NA. Hence we want to tend to
zero as our solver converges, thereby scaling the
iterates heavily towards pure Newton steps in the
asymptotic regime and enjoying quadratic
convergence locally. The remainder of
this section is dedicated to quantifying how close
a iterate
is to a standard Newton
iterate in the asymptotic regime when
is nonsingular.
The main result is
Theorem 3.2
below which bounds locally, where
, and . This notation is introduced to
emphasize the Newton and Newton-Anderson iterates in the
comparison.
We will also define
.
The following lemma will be used in the proof of
Theorem 3.2.
Lemma 3.1
bounds , the scaled
returned by at iteration , in terms of
and .
The proof consists of walking through
the cases in Algorithm 6, and is therefore
left to the interested reader.
Lemma 3.1.
Let and .
Define
and
as in
Algorithm 6.
Let be the value computed by
Algorithm 6 at iteration . If , then returned by
Algorithm 6 satisfies
(8)
when , and
(9)
when .
With Lemma 3.1, we can bound
in terms of ,
, and .
Theorem 3.2.
If is nonsingular, then there exists a
and a constant depending only on
such that for and
in and ,
(10)
when , and
(11)
when .
Proof.
Using our notation from the discussion preceding
Lemma 3.1, an iterate generated by
takes the form . Therefore
.
Since is nonsingular, we can take
sufficiently small to ensure that , where is a constant determined by
, when for . Hence, upon
adding and subtracting we obtain
.
To complete the proof, we write
and apply Lemma 3.1 to bound
for the cases
and .
∎
We conclude this section with Corollary 3.1,
which proves that can recover local quadratic convergence
from NA when applied to nonsingular problems.
Corollary 3.1.
If is nonsingular, , and , then there exists a and constants and depending only on such that
(12)
for .
Proof.
Adding and subtracting to
gives
By Theorem 3.2, , and
for , . Since , we have . Moreover, when is nonsingular, Taylor expansion shows that
(13)
for . Thus , and therefore
(14)
∎
4 Numerics
In this section, we demonstrate the effectiveness of NA and near bifurcation points by applying
these algorithms to the following parameter-dependent PDEs. All
computations are performed on an M1 MacBook with GNU Octave 8.2.0.
4.1 Test Problems
1.
Navier-Stokes Flow in a Channel
(15)
2.
Rayleigh-Bénard Convection
(16)
In both models, denotes the
fluid velocity, the pressure, and
the outward normal. In Model (15), denotes
the viscosity parameter. In Model (16), denotes
the temperature of the fluid, and denotes the
Richardson number, which is the parameter of interest
for Model (16). We set .
These parameters values are chosen
so as to replicate the results seen in [20] as ranges from to .
For the flow in a channel, Model (15), we use Taylor-Hood elements [8, p. 164]. The channel,
shown in Figure 1, is arranged such that the left most boundary lies at , the right most boundary lies
at , and
the boundary components are given by
, , and
.
For the Rayleigh-Béndard Model, we use
Scott-Vogelius elements [24], where denotes piece-wise
linear discontinuous elements.
The elements are stable
on Alfeld-split, also known as barycenter-split, triangulations
[46, p. 77].
The meshes used for
each model are shown below.
Figure 1: Meshes used for benchmark problems. Top: Mesh used for Rayleigh-Bénard
model. Bottom: mesh used for flow in a channel.
For Model (15), it is known [39] that there exists a critical
viscosity at which a bifurcation occurs. For ,
the stable velocity solution is symmetric about
the center horizontal () as seen in the top plot
of Figure 2. For ,
there is still a symmetric solution, but it is unstable. Stability is inherited by
two asymmetric solutions seen in the bottom two plots of
Figure 2.
Figure 2: Solutions to channel flow
problem (15) for different . Top: Representative
symmetric solution for . Middle: Representative
asymmetric solution with positive vertical velocity
upon exiting the narrow channel for .
Bottom: Representative
asymmetric solution with negative vertical velocity
upon exiting the narrow channel for .
The parameter region of interest for Model (16) is . In this range,
the flow appears to be in transition from a single eddy
in the center of the domain to
two eddies as seen in Figure 3 below.
Figure 3: Velocity streamlines for velocity
from Model (16) showing transition from
one eddy to two eddies. Top Left: Ri = 3.0. Top Right: Ri = 3.1.
Bottom Left: Ri = 3.2. Bottom Right: Ri = 3.4.
For Model (15), we take the zero vector
as our initial iterate. For Model (16), we initialize
our iterate by applying a Picard step and then applying
. The seminorm is used in our implementations [8].
Before we discuss the numerical results, we recall that
NA, Algorithm 5, is Algorithm 2
with line 5 replaced by Algorithm 4.
•
, Algorithm 7, is Algorithm 2 with line
5 replaced by Algorithm 6.
For NA and , we set and
respectively to a fixed quantity for all iterations.
4.2 General Discussion of Results
The following experiments demonstrate three strategies for solving nonlinear problems
near bifurcations using NA and . The first two are
asymptotic safeguarding and preasymptotic
safeguarding. With asymptotic safeguarding,
we run NA until
the residual is smaller
than some user-chosen threshold, and
we use for all subsequent iterates. Hence the
solve will behave like NA until the last few iterations when
is applied. For the
experiments in Section 4.3, we activated
when .
We chose since
we want to activate -safeguarding
as early as possible.
This allows for earlier detection of a nonsingular problem, and thus
faster convergence. We found that
activating -safeguarding for will not
necessarily break convergence,
but if the problem is nonsingular, this will
not be detected as early.
Activating
when can break convergence,
though this
seems to be problem-dependent.
Convergence for Model (15) was virtually unaffected
with threshold , but for Model (16)
with ,
setting the activation threshold to caused
to diverge when it had converged with threshold .
With activation threshold , asymptotic safeguarding
is shown to be effective
close to the bifurcation point, where Newton’s method
can fail to converge.
Preasymptotic safeguarding, on the other hand, applies
at each step of the solve. Interestingly, we observe
that applied preasymptotically can outperform NA when applied to Model (16), and even recover
convergence when both Newton and diverge (see
Figure 11). A theoretical explanation for
this requires a better understanding of these methods in the
preasymptotic regime. In particular, a better understanding
of the descent properties of Anderson acceleration would be
of great value. It is known [52] that for
singular problems in the
preasymptotic regime, where is not small, the
Newton update step ,
or in our notation, can be large
and nearly orthogonal to the gradient of .
With , our update step takes the form
.
It is clear that is a descent direction for
sufficiently small , since in this case it is
nearly . It is possible that for certain values of
, is a stronger descent direction
than .
Another possible explanation as to why can outperform
in some cases is its resemblance to restarted Anderson
acceleration methods. Restarted versions of Anderson
acceleration are often applied in various
forms for depth , and
have been shown to effective in practice [9, 25, 26, 36]. In the special case of Newton-Anderson with depth , every odd iterate is simply
a Newton step, rather than
a combination of the previous two Newton steps. Hence
the algorithm is “restarted” every other step. Explicitly,
we have for ,
In other words, for all . With
, is not necessarily set to zero, but
it is scaled towards zero, significantly so depending on
. In this way, one may think of as
a quasi-restarted Anderson scheme when the depth
. At the moment, -safeguarding is not
developed for , but this interpretation of
as a quasi-restarted method could lead to such
a development.
Presently, these are only heuristics, but they provide
interesting questions for future projects.
The third technique we demonstrate to solve these problems
near bifurcation points is increasing the depth
. Evidently, the right choice of can significantly
improve convergence by reducing the number of iterations
to convergence by half, and increase the domain of
convergence with respect to the parameter. We found, however, that such
performance was very sensitive to the choice of . So while
the results suggest this could be developed into a viable
strategy, more work is required to achieve this.
In all experiments we take
since this is the range in which local convergence of
is guaranteed by Theorem 3.1. In
practice, neither Algorithm 5 nor Algorithm
7 breaks down if one sets or
respectively to zero, one, or a value greater than one.
Setting or to zero reduces the iteration to
Newton, and choosing one leads to a more NA like iteration in
the preasymptotic regime. A systematic
study of with has not been performed,
but experiments performed
thus far show no significant advantage over .
What the best choice of is remains an open
question. Numerical experiments suggest the best choice
depends on the initial guess . For example, when applied to the channel flow problem preasymptotically (see Section
4.4), setting results in faster convergnece than Newton if is the zero vector.
If we perturb this (discussed in Section 4.4.1), then the choice of leads to converging slower than Newton, while outperforms NA and Newton. This phenomena, that the
best choice of in varies with , is
seen with other choices of as well. Elucidating this
dependence is the subject of ongoing work.
4.3 Asymptotic Safeguarding
Under the assumptions of Theorem 3.1
or Theorem 3.2, is guaranteed to converge locally. This motivates
the strategy of this subsection.
As discussed above, asymptotic safeguarding is when we only
apply once the residual is smaller than a set
threshold.
This allows one to take full advantage of NA in the
preasymptotic regime, and
ensures fast quadratic convergence for nonsingular problems.
In practice, this means fast local convergence provided
NA reaches the domain of convergence.
We activate when
in the experiments below. When ,
we run NA.
4.3.1 Results for Channel Flow Model
The results when applied
to Model
(15) are shown below in
Figures 4, 5, and 6.
The takeaway is that when applied asymptotically, convergence
of is not as sensitive to the choice of as it is when applied in the preasymptotic regime (see
Section 4.4), and
Algorithm 6
is working as intended
by detecting that the problem is
nonsingular.
This is seen in the plot on the right of
Figures 4, 5, and 6.
Recall that is the adaptive parameter in Algorithm 6
that determines how close a
step is to a Newton step.
When ,
and the criteria in Algorithm 6 is met,
the
step will be close to a Newton step. When
, this scaling will be much less severe,
and the step will
be close to an NA step. From our discussion
in Section 3, we want when the problem is
nonsingular in order to enjoy local quadratic convergence. This is not guaranteed
with NA [47]. Observing the plots in Figures 4,
5, and 6, one notes that as the solver converges. Since , this is equivalent to as the solve converges. Thus by Theorem 3.2, the iterates converge to
Newton iterates asymptotically. This is precisely what was designed to do: detect nonsingular problems, and respond by scaling the iterates towards Newton
asymptotically. Since
is only activated when in these examples of asymptotic safeguarding, is only computed, and plotted, for the last two iterations. In the next section on preasymptotic safeguarding, a more interesting
history is seen. Our methods, including
, failed to converge for . If one wishes
to solve a problem at a particular parameter, and a direct
solve fails like we see here for , one could still
employ or to solve the problem directly for
a parameter value close to the desired one to obtain an
initial guess for continuation. The benefit here is that
the continuation would be required in a smaller parameter
range, thus reducing the total number of solves required.
We will see in Section 4.5 that
increasing can lead to
convergence for a wider range of parameters.
Figure 4: Comparison of
applied asymptotically with Newton and NA applied
to Model (15) with .
Figure 5: Comparison of
applied asymptotically with Newton and NA applied
to Model (15) with .
Figure 6: Comparison of
applied asymptotically with Newton and NA applied
to (15) with .
4.3.2 Results for Rayleigh-Bénard Model
The results of asymptotic safeguarding with activation
threshold applied to Model (16) are similar to
those of Model (15) seen in the previous section.
In this case, however, NA diverged for and
. For , Newton’s method converges, but
all methods diverged for . In Section
4.4, we are able to recover convergence
for with preasymptotic safeguarding.
For ,
, , and , Newton diverged
while NA and converged. Like with Model (15), we see as the solve converges. The results for
, shown below, are representative of the others for which NA converged.
Figure 7: Comparison of
applied asymptotically with Newton and NA applied
to Model (16) with .
4.4 Preasymptotic Safeguarding
In this section, we demonstrate the preasymptotic safeguarding strategy, where is activated starting at iteration , the first iteration where can be applied. Compared to asymptotic safeguarding, preasymptotic safeguarding is more sensitive to the choice of , but, with the right choice of , it can recover convergence when both Newton and fail.
4.4.1 Results for Channel Flow Model
Here we present the results from applying the preasymptotic
safeguarding strategy to the channel flow Model
(15). Evidently, Newton’s method can still converge quickly with the right initial guess when
. Figure 8 demonstrates that in this case, and also
perform well. The right-most plot in Figure 8 demonstrates that the is
working as intended. That is, , and therefore
is detecting that Newton is converging quickly,
and responds by scaling its update steps
towards a pure Newton iteration.
With , we are closer to the bifurcation point, and observing the left plot in Figure
9,
we see that
Newton takes many more iterations to converge. However, after a long preasymptotic
phase, Newton does eventually converge quickly, which suggests that the problem is not
truly singular for . This is again detected by .
Even though the algorithms take a few more iterations to converge than NA with
,
the terminal order of convergence of is greater than that of NA. Approximating
the rate by at each step , and letting
denote the terminal order, we found that
for NA, for , and
for .
In the right most plot in Figure
9, we take as our initial iterate the
zero vector, but with the fifth entry set to 50. From this
perturbed initial guess, Newton’s method is seen to
perform better than NA. Further, and
outperform both Newton and NA, with
converging in about half as many
iterations as NA. This demonstrates the flexibility
offered by . That is, whether Newton
or NA is the best choice for a particular problem
and initial guess, is more agnostic to
these choices, and can perform well in either
case. There is still, however, sensitivity
to .
We found that failed to converge for
, , , , and . converged with
, but only after 74 iterations. We observed that
NA(0.5), the non-adaptive version of -safeguarding, managed to
converge. The reason for this variation is likely due to the complex behavior of Newton and NA
in the preasymptotic regime, leading to sensitivity to the
choice of . With preasymptotic
safeguarding, our methods failed for , hence
continuation may still be required in some cases,
but NA and can be used to efficiently solve the
problem closer to the bifurcation point compared to
Newton, thereby reducing the overall computational
cost.
The point is that applying NA and in the
preasymptotic regime can still lead to faster convergence
than Newton, but this convergence is again
sensitive to the choice of , and further
work is required to understand this sensitivity.
Figure 8:
Results of preasymptotic safeguarding applied to
Model (15) for . Left: Convergence history for Newton, NA, and
for . Right: The value of
at each iteration for the adaptive methods.
Figure 9: Results of preasymptotic safeguarding applied to
Model (15) for . Left: Convergence history for Newton, NA, and
for with
zero vector initial guess.
Right: Convergence history for Newton, NA,
and for
with perturbed initial guess.
4.4.2 Results for Rayleigh-Bénard Model
With preasymptotic safeguarding applied to
Model (16), we again see more varied
behavior since we have -safeguarding activated from
the beginning of the solve. The results are shown in Figures 10, 11, and
12. The plots are only
shown for those values of for which converged.
We found that for and ,
failed to converge with , , or
. Convergence is recovered with and
respectively. We also ran with these
values for , , , and . The results were similar for
, , and . Hence
we only show the results for
in Figure 12.
The theme demonstrated in Figures 10, 11, and
12 is that
when preasymptotic
safeguarding is employed, it is possible to converge faster
than standard NA. Moreover, as seen in Figure
11, can converge when both Newton and NA diverge. One may note that in the right-most
plot in Figure 10,
prior to the asymptotic regime where ,
we observe only twice. Using the restarted Anderson
interpretation discussed in Section 4.2,
we could say that there are two quasi-restarts prior
to the asymptotic regime. Evidently, these two
quasi-restarts are essential, since we found that
non-adaptive -safeguarding with ,
Algorithm 4, diverges. Similar
behavior of is seen in Figures
11 and 12.
It remains unclear
precisely how the choice of affects convergence,
e.g., in Figure 10, why does
converge, but diverges for
, , , and ? As
previously discussed, a better
understanding of these methods in the preasymptotic regime
would help answer questions like these, and this will be
the focus of future studies.
Figure 10: Comparison of
applied preasymptotically with Newton and NA applied
to Model (16) with .
Figure 11: Comparison of
applied preasymptotically with Newton and NA applied
to Model (16) with .
Figure 12: Comparison of
applied preasymptotically with Newton and NA applied
to Model (16) with .
4.5 Increasing Anderson Depth
The strategy employed in this section is to increase the
depth of the NA from used in previous sections.
4.5.1 Results for Channel Flow Model
The plots in Figure 13 demonstrate the
effectiveness of increasing the algorithmic depth
to solve Model (15) near the bifurcation point. We also experimented with applying asymptotically. When , we set and activated
with .
These results are seen as the
dashed lines in the left-most plot in Figure
13.
The philosophy is similar to that of
preasymptotic safeguarding. We use NA(3) to reach the
asymptotic regime, and then allow adaptive -safeguarding
to detect if the problem is nonsingular.
We set because,
for the present, -safeguarding is only designed
for . The left-most plot in Figure 13
shows how, from the same initial iterate, we are able
to solve Model (15) for a wider range of
values, including in the regime where Newton,
, and failed to converge.
The right-most
plot in Figure 13 focuses on the
results of applying to Model (15)
with . The point here is that even in the
regime where converges, converges
in about half as many iterations. Thus increasing
the algorithmic depth can lead to faster convergence.
The improved convergence seen with increasing the depth suggests that
a generalization of -safeguarding for greater depths could be
useful as a generalization of the strategies presented in previous sections.
However, for this strategy to be effective in general,
further study is needed on the proper choice of depth .
In the chosen parameter regime,
was the only value of observed to improve convergence
when NA(1) failed.
Figure 13: Results from increasing Anderson
depth . Left: Convergence
history of NA(3) for various choices of near . Solid lines are NA, and
dashed lines are NA with asymptotic .
Right: Convergence history
for NA(1) and NA(3) with .
4.5.2 Results for Rayleigh-Bénard Model
The results of increasing to solve Model (16)
are shown below in
Figure 14.
As with Model (15), we experimented with reducing
and activating -safeguarding asymptotically. The dashed lines in Figure 14 are the results of these experiments.
For this problem, activation occurred when
since we found that activating
with when broke convergence like it did
in Section 4.3.2.
We once again
observe that the right choice of
can improve convergence
significantly. This can be seen in Figure 14
for . Moreover, Figure 14 demonstrates that increasing can recover convergence when
fails to converge for and .
However, for 3.1, 3.4, and 3.5, there was no significant improvement
gained from increasing . The results for
shown below in Figure
14 are representative of the results for and .
The most significant difference between
without asymptotic , and with asymptotic , is
seen with .
It also appears that is more sensitive to the
activation threshold than smaller choices of . This is
seen in Figure 14 for
and . These results are promising, and
motivate further investigation to fully understand
how the choice of affects convergence near singularities,
and in particular, near bifurcation points.
Figure 14: Convergence history
of for 1,2,3,4,5, and 10 for 3.0, 3.2, 3.3, and 3.4. Dashed lines denote with asymptotic
with activated when .
5 Conclusion
We have presented a modification of Anderson accelerated
Newton’s method for solving nonlinear equations near
bifurcation points. We proved that, locally, this modified
scheme can detect nonsingular problems and scale the
iterates towards a pure Newton step, which leads to faster
local convergence compared to standard NA.
We numerically demonstrated two strategies one can
employ when using our modified NA scheme
to solve nonlinear problems near bifurcation points, with
our test problems being two Navier-Stokes type parameter-dependent PDEs. Asymptotic safeguarding was shown to recover
local quadratic convergence when the problem is nonsingular,
and it shows virtually no sensitivity to the choice of
parameter . It can, however, be sensitive to the
choice of activation threshold. Preasymptotic safeguarding
is shown to significantly improve convergence, and
can recover convergence when both Newton and fail.
There is strong sensitivity to the choice of though,
and future work will clarify this dependence.
We also demonstrated that increasing the Anderson depth
can improve convergence, and increase the domain of
convergence with respect to the problem parameter. Future
projects will further study how the choice of impacts
convergence, and work towards developing -safeguarding
for greater algorithmic depths.
6 Acknowledgements
MD and SP are supported in part by the National Science
Foundation under Grant No. DMS-2011519. LR is supported in
part by the National Science Foundation under Grant No.
DMS-2011490. This material is based upon work supported by
the National Science Foundation under Grant No. DMS-1929284
while the authors were in residence at the Institute for
Computational and Experimental Research in Mathematics in
Providence, RI, during the Acceleration and Extrapolation
Methods (MD, SP and LR), and the Numerical PDEs: Analysis,
Algorithms and Data Challenges (SP and LR), programs.
References
[1]
H. An, X. Jia, and H.F. Walker.
Anderson acceleration and application to the three-temperature energy
equations.
J. Comput. Phys., 347:1–19, 2017.
[2]
D.G. Anderson.
Iterative procedures for nonlinear integral equations.
J. Assoc. Comput. Mach., 12(4):547–560, 1965.
[3]
R. Behling and A. Fischer.
A unified local convergence analysis of inexact constrained
Levenberg-Marquardt methods.
Optim. Lett., 6:927–940, 2012.
[4]
R. Behling, A. Fischer, M. Herrich, A. Iusem, and Y. Ye.
A Levenberg-Marquardt method with approximate projections.
Comput. Optim. Appl., 59:5–26, 2014.
[5]
S. Bellavia and B. Morini.
Strong local convergence properties of adaptive regularized methods
for nonlinear least squares.
IMA J. Numer. Anal., 35, 2015.
[6]
R.G. Bettiol and P. Piccione.
Instability and bifurcation.
Notices Am. Math. Soc., 67(11):1679–1691, 2020.
[7]
N. Boullé, V. Dallas, and P.E. Farrell.
Bifurcation analysis of two-dimensional Rayleigh-Bénard
convection using deflation.
Phys. Rev. E, 105, 2022.
[8]
D. Braess.
Finite Elements: Theory, Fast Solvers, and Applications in Solid
Mechanics.
Cambridge University Press, 2007.
[9]
M. Chupin, M. Dupuy, G. Legendre, and E. Séré.
Convergence analysis of adaptive DIIS algorithms with application
to electronic ground state calculations.
ESAIM Math. Model. Numer. Anal., 55(6):2785–2825, 2021.
[10]
M. Dallas and S. Pollock.
Newton-Anderson at singular points.
Int. J. Numer. Anal. and Mod., 20(5):667–692, 2023.
[11]
D.W. Decker, H.B. Keller, and C.T. Kelley.
Convergence rates for Newton’s method at singular points.
SIAM J. Numer. Anal., 20(2):296–314, 1983.
[12]
D.W. Decker and C.T. Kelley.
Newton’s method at singular points I.
SIAM J. Numer. Anal., 17(1):66–70, 1980.
[13]
D.W. Decker and C.T. Kelley.
Convergence acceleration for Newton’s method at singular points.
SIAM J. Numer. Anal., 19(1):219–229, 1982.
[14]
D.W. Decker and C.T. Kelley.
Expanded convergence domains for Newton’s method at nearly
singular roots.
SIAM J. Sci. Stat. Comput., 6(4), 1985.
[15]
P. Deuflhard.
Newton Methods for Nonlinear Problems.
Springer Series in Computational Mathematics. Springer, 2005.
[16]
J. Fan and J. Zeng.
A Levenberg–Marquardt algorithm with correction for singular
system of nonlinear equations.
Appl. Math. Comput., 219(17):9438–9446, 2013.
[17]
P.E. Farrell, A. Birkisson, and S.W. Funke.
Deflation techniques for finding distinct solutions of nonlinear
partial differential equations.
SIAM J. Sci. Comput., 37, 2015.
[18]
A. Fischer, A.F. Izmailov, and M.V. Solodov.
Accelerating convergence of the globalized Newton method to
critical solutions of nonlinear equations.
Comput. Optim. Appl., 78:273–286, 2021.
[19]
A. Fischer, A.F. Izmailov, and M.V. Solodov.
Unit stepsize for the Newton method close to critical solutions.
Math. Program., 187:697–721, 2021.
[20]
J.K. Galvin, A. Linke, L.G. Rebholz, and N. Wilson.
Stabilizing poor mass conservation in incompressible flow problems
with large irrotational forcing and application to thermal convection.
Comput. Methods Appl. Mech. Engrg., pages 166–176, 2012.
[21]
A.O. Griewank.
Analysis and Modification of Newton’s Method at
Singularities.
PhD thesis, 1980.
[22]
A.O. Griewank.
Starlike domains of convergence for Newton’s method at
singularities.
Numer. Math., 35:95–111, 1980.
[23]
A.O. Griewank.
On solving nonlinear equations with simple singularities or nearly
singular solutions.
SIAM Review, 27(4):537–563, 1985.
[24]
J. Guzḿan and L.R. Scott.
The Scott-Vogelius finite elements revisited.
Math. Comp., 88(316):515–529, 2019.
[25]
H. He, S. Zhao, Y. Xi, J.C. Ho, and Y. Saad.
Solve minimax optimization by Anderson acceleration.
International Conference on Learning Representations, 2022.
[26]
N.C. Henderson and R. Varadhan.
Damped Anderson acceleration with restarts and monotonicity control
for accelerating EM and EM-like algorithms.
Journal of computational and graphical statistics,
28(4):834–846, 2019.
[27]
J.L. Hueso, E.M., and J.R. Torregrosa.
Modified Newton’s method for systems of nonlinear equations with
singular Jacobian.
J. Comput. Appl. Math., 224(1):77–83, 2008.
[28]
A.F. Izmailov, A.S. Kurennoy, and M.V. Solodov.
Critical solutions of nonlinear equations: local attraction for
Newton-type methods.
Math. Progam., 167:355–379, 2018.
[29]
A.F. Izmailov, A.S. Kurennoy, and M.V. Solodov.
Critical solutions of nonlinear equations: stability issues.
Math. Progam., 168:475–507, 2018.
[30]
C. Kanzow, N. Yamashita, and M. Fukushima.
Levenberg–Marquardt methods with strong local convergence
properties for solving nonlinear equations with convex constraints.
J. Comput. Appl. Math, 172:375–397, 2004.
[31]
C.T. Kelley and R. Suresh.
A new acceleration method for Newton’s method at singular points.
SIAM J. Numer. Anal, 20(5):1001–1009, 1983.
[32]
H. Kielhöfer.
Bifurcation Theory: An Introduction with Applications to Partial
Differential Equations.
Springer Science+Business Media, 2012.
[33]
C. Kuehn.
PDE Dynamics: An Introduction.
SIAM, 2019.
[34]
P.A. Lott, H.F. Walker, C.S. Woodward, and U.M. Yang.
An accelerated Picard method for nonlinear systems related to
variably saturated flow.
Adv. Water Resour., 38:92–101, 2012.
[36]
M.L. Pasini.
Convergence analysis of Anderson-type acceleration of
Richardson’s iteration.
Numer. Linear Algebra Appl., 26, 2019.
[37]
M.L. Pasini, J. Yin, V. Reshniak, and M.K. Stoyanov.
Anderson acceleration for distributed training of deep learning
models.
In SoutheastCon 2022, pages 289–295, 2022.
[38]
Y. Peng, B. Deng, J. Zhang, F. Geng, W. Qin, and L. Liu.
Anderson acceleration for geometry optimization and physics
simulation.
ACM Trans. Graph., 37(4), 2018.
[39]
F. Pichi.
Reduced order models for parametric bifurcation problems in
nonlinear PDEs.
PhD thesis, 2020.
[40]
F. Pichi, A. Quaini, and G. Rozza.
A reduced order modeling technique to study bifurcating phenomena:
Application to the Gross–Pitaevskii equation.
SIAM J. Sci. Comput., 42(5), 2020.
[41]
F. Pichi, M. Strazzullo, F. Ballarin, and G. Rozza.
Driving bifurcating parametrized nonlinear pdes by optimal control
strategies: application to Navier–Stokes equations with model order
reduction.
ESAIM: M2AN, 56(4):1361–1400, 2022.
[42]
S. Pollock and L.G. Rebholz.
Anderson acceleration for contractive and noncontractive operators.
IMA J. Numer. Anal., 41(4):2841–2872, 2021.
[43]
S. Pollock and L.G. Rebholz.
Filtering for Anderson acceleration.
SIAM J. Sci. Comput., 45(4):A1571–A1590, 2023.
[44]
S. Pollock, L.G. Rebholz, and M. Xiao.
Anderson-accelerated convergence of Picard iterations for
incompressible Navier-Stokes equations.
SIAM J. Numer. Anal., 57(2):615–637, 2019.
[45]
S. Pollock and H. Schwartz.
Benchmarking results for the Newton–Anderson method.
Results Appl. Math., 8:100095, 2020.
[46]
J. Qin.
On the convergence of some low order mixed finite elements for
incompressible fluids.
PhD thesis, 1994.
[47]
L.G. Rebholz and M. Xiao.
The effect of Anderson acceleration on superlinear and sublinear
convergence.
J. Sci. Comput., 96(2), 2023.
[48]
G.W. Reddien.
On Newton’s method for singular problems.
SIAM J. Numer. Anal., 15(5):993–996, 1978.
[49]
G.W. Reddien.
Newton’s method and high order singularities.
Comput. Math. Appl., 5(2):79–86, 1979.
[50]
R.B. Schnabel and P.D. Frank.
Tensor methods for nonlinear equations.
SIAM J. Numer. Anal., 21(5):815–843, 1984.
[51]
R.B. Thompson, K.O. Rasmussen, and T. Lookman.
Improved convergence in block copolymer self-consistent field theory
by Anderson mixing.
J. Chem. Phys., 120(1):31–34, 2004.
[52]
R.S. Tuminaro, H.F. Walker, and J. N. Shadid.
On backtracking failure in Newton-GMRES methods with
demonstration for the Navier-Stokes equations.
J. Comput. Phys., 180:549–558, 2002.
[53]
H. Uecker.
Continuation and bifurcation in nonlinear pdes – algorithms,
applications, and experiments.
Jahresbericht der Deutschen Mathematiker-Vereinigung,
124:43–80, 2021.
[54]
H.F. Walker and P. Ni.
Anderson acceleration for fixed-point iterations.
SIAM J. Numer. Anal., 49(4):1715–1735, 2011.
[55]
D. Wang, Y. He, and H. De Sterck.
On the asymptotic linear convergence speed of Anderson acceleration
applied to ADMM.
J. Sci. Comput., 88(2):38, 2021.
[56]
Q. Yuan, Y. Sun, and J. Ren.
How interest rate influences a business cycle model.
Discrete Contin. Dyn. Syst., 13(11):3231–3251, 2020.
[57]
H. Zhou and Y. Qian.
Double Hopf bifurcation analysis for coupled van der
Pol–Rayleigh system with time delay.
J. Vib. Eng. Technol., 12:6075–6087, 2024.