Augmented Lagrangian Methods for Numerical Solutions to Higher Order Differential Equations

A large number of problems in engineering can be formulated as the optimization of certain functionals. In this paper, we present an algorithm that uses the augmented Lagrangian methods for finding numerical solutions to engineering problems. These engineering problems are described by differential equations with boundary values and are formulated as optimization of some functionals. The algorithm achieves its simplicity and versatility by choosing linear equality relations recursively for the augmented Lagrangian associated with an optimization problem. We demonstrate the formulation of an optimization functional for a 4th order nonlinear differential equation with boundary values. We also derive the associated augmented Lagrangian for this 4th order differential equation. Numerical test results are included that match up with well-established experimental outcomes. These numerical results indicate that the new algorithm is fully capable of producing accurate and stable solutions to differential equations.


Introduction
Many problems in engineering are described by boundary value problems of differential equations of order 2k in the form of ( ) , , , , , 0, , , However, a strong solution to (1) may not exist for some problems.Even when a strong solution to (1) exists, such a solution may be too costly to calculate numerically in practice.
Therefore, weak formulations of (1) are preferred.Assume that (1) admits a weak formulation in the following form. ( where  is usually a subset of the Hilbert space ( ) subject to some boundary conditions, 0  is a subspace of ( ) satisfying some homogeneous boundary conditions, and ( ) In other words, a solution to (1) also satisfies (2).
A solution to (2) is referred to as a weak solution to (1) because of ( ) We also assume that for y ∈  , there exists a certain functional ( ,  such that a stationary point of functional ( ) J y satisfies (2).A stationary point y of (3) satisfies ( ) where η is referred to as a test function.Equivalently, a stationary point y of (3) satisfies Refer to [1] for more details.
We see that a weak solution of (1) corresponds to an optimization problem of functional ( ) J y defined in (3).Such a weak solution is simply a stationary point of (3).
As pointed out in [1], a boundary value problem of inhomogeneous boundary conditions can be effectively treated as a problem with homogeneous boundary conditions by introducing a special function satisfying the inhomogeneous boundary conditions.For simplicity, we will assume homogeneous boundary conditions for the rest of this paper.Under this assumption, Hestenes [2] and Powell [3] introduced the augmented Lagrangian methods or the method of multipliers in 1969 for the study of optimization problems.
Comprehensive applications of augmented Lagrangian methods for optimization and boundary value problems were investigated by Bertsekas [4], Fortin and
Even though most research on the augmented Lagrangian methods have been focusing on their direct applications in optimization, there have been consistent interests in using the augmented Lagrangian methods in other fields over the years.Recent applications of the augmented Lagrangian methods include image processing and computer vision [6], numerical solutions to Laplace equation with various boundary values [7], mechanics [7] [8] [9], medical imaging and computational biology [10], geophysics [11], ontology regarding big data representation and storage [12], and elastica theory [13], to name just a few.
In this paper, we'll investigate the applications of the augmented Lagrangian methods in boundary value problems of differential equations associated with problems from engineering.These are usually higher order differential equations (4th or higher).The feature of higher order of these problems is exploited to allow the development of simple and stable methods for their numerical solutions using finite elements.
The rest of the paper is organized in the following way.In section 2, we briefly review the augmented Lagrangian methods for optimization problems, and the finite element methods for solving weak formulation (2).In section 3, we'll introduce a special formulation for the construction of an augmented Lagrangian for finding numerical solutions to higher order differential equations using finite elements.We show a sample boundary value problem of a higher order nonlinear differential equation, its weak formulation, and the associated optimization problem in section 4. We'll present numerical tests for this sample differential equation in section 5.The numerical results demonstrate that this new version of the augmented Lagrangian methods is accurate, stable, versatile, and easy to implement.

Augmented Lagrangian Methods and Finite Element Solutions
Generally speaking, an analytic solution to an optimization problem is not available.We therefore focus on numerical solutions using augmented Lagrangian methods and finite elements.

Augmented Lagrangian Methods
For a constrained optimization problem ( ) where K is a subset of the set of all natural numbers, and  is some functional space, an augmented Lagrangian method for (7) where for each j K ∈ , j µ is a Lagrangian multiplier, and j r is a pre-chosen constant serving as a penalty when ( ) j c y is not sufficiently small.
Notice that the function spaces in (7) and ( 9) are identical.In particular, both (7) and (9) impose the same regularity requirements on y through  .
When an optimization problem is associated with a boundary value problem of a differential equation, we have the leeway for constructing the associated augmented Lagrangian that imposes weaker regularity requirements than those of the original optimization problem.We'll present the formulation later in a separate section.

Finite Element Methods
Assume that N  is a N -dimensional subspace that approximates ( ) The Galerkin method of finite elements amounts to finding a function where ( ) Because (10) is true for all N η ∈  , it must be true for every ,1 In fact, (12) forms a system of N equations for N unknown coefficients . These equations are linear or nonlinear depending on ( ) , a ⋅ ⋅ being linear or nonlinear.Finding an approximate solution to (2) corresponds to solving the system of Equations ( 12) for coefficients { } 1 Such regularity requirement may become too demanding on the basis for N  when k is big.Using weaker regularity requirements generally results in the following advantages.
1) We may use simpler test function η that improves overall numerical effi- ciency.
Our goal is to find new approaches that use the least regularity requirement possible, to be introduced next.

Augmented Lagrangian Methods for Differential Equations
In this section, we propose a new approach for the construction of augmented Lagrangian associated with boundary value problems of differential equations.
This new approach exploits the feature of higher order of the differential equations to arrive at the least regularity requirements on weak solutions to boundary value problems.The approach coincides with those discussed in [5] for differential equations of second order or less.
Recall that to find a weak solution y to (1) means to solve for y in (5).
That means all test functions must be from ( ) . We now propose the following augmented Lagrangian associated with (3).
where ,0 < that is associated with some of the homogeneous boundary conditions of (1) (see assumption for We immediately recognize the major difference between functional ( ) J y of (3) and the augmented Lagrangian of (13).Whereas the regularity requirement for ( ) the regularity requirements for augmented Lagrangian of ( 13) are independent of k , the order of the differential equation.
Using results in (25) along with (32), we have proved that ( , ) That completes the proof.We therefore have three approaches for finding approximate solutions to (1).

1) Solve (1) for an approximate strong solution in 2 ( , )
k a b  directly, e.g., by a certain finite difference method.This approach demands the regularity requirement where functions involved must belong to 2 ( , ) k a b  .
2) Solve ( 1) for an approximate weak solution in ( , ) k a b  via ( 2), ( 10) and ( 12).This approach demands that functions involved belong to ( , ) k a b  .This has been the standard approach found in most work related to finite element methods.
3) Solve ( 1) for an approximate weak solution in 1 ( , ) a b  by finding stationary point of the augmented Lagrangian ( 13).This approach demands that functions involved belong to 1 ( , ) a b  , regardless of the order 2k of differential Equation (1).This approach achieves the minimum regularity requirement possible for any differential equation of order second or higher.
The method that uses the augmented Lagrangian (13) has the least regularity requirement among the three approaches.
Regardless of the order of differential Equation (1), finding approximate solution using the augmented Lagrangian (13) with finite elements is simple and standard as illustrated below.
• Because of ( 16), we can always use Hermite cubic polynomials to approximate ,0 j p j k ≤ < , and use piecewise linear functions to approximate k p and ,1 • The bases for the corresponding finite dimensional spaces for finite elements are the Hermite cubic shape functions and linear shape functions, respectively.
• The approximate solutions are obtained from (12) by substituting Hermite cubic shape functions and linear shape functions for { } j φ , respectively.We'll implement the augmented Lagrangian methods associated with Lagrangian (13) for finding approximate solution to a boundary value problem of a 4th order nonlinear differential equation in the next section.

A Sample Differential Equation
Many of the problems in engineering are described by boundary value problems of differential equations of orders two or higher.An example is the following nonlinear Euler-Bernoulli beam equation [14] [15] [16], with various boundary conditions (natural boundary conditions are not enforced by standard arguments), such as: ( ) ( ) for a rectangular beam of length L, width w and thickness t .Here, v is the deflection of the beam, E Young's modulus, I moment of inertia of cross-section of the beam, r T represents residual force which is independent of v , represents the axial force in the beam, and ( ) f v is the intensity of external force exerted on the beam which is assumed to be continuous in v .We see that (33) along with (34) is just a special case of (1) for 2 k = .A strong solution v of (33) requires that 4 (0, ) v L ∈  .
Subspace  associated with (33) becomes Because boundary conditions in (34) are themselves homogeneous, the subspace 0  of all test functions is identical to  , as has been assumed back in (6).
When multiplying any test function 0 η ∈  to (33), integrating the result over [0, ] L , and applying integration by parts, we arrive at the following.
To find a solution to (33) becomes to find a function v ∈  such that (37) is true.A solution to (37) is therefore called a weak solution to (33) because the regularity requirement has become We also call (37) the weak formulation of (33).
Furthermore, for any v ∈  , if we define ( ) ( ) , we can verify with some algebraic manipulations that a stationary point of ( ) J v in (39) is a weak solution to (33).Details are omitted here for brevity.
We've shown that (33) indeed admits a weak formulation (37).A weak solution of (33) corresponds to a stationary point of (39).Such a stationary point is a solution to an optimization problem of functional ( ) J v of (39).
Based on functional ( ) J v in (39), we can define the associated augmented Lagrangian as shown below.
Notice that boundary conditions (34) become ( ) ( ) and ( ) Finding an exact solution to a stationary point from (43) is highly unlikely, if not impossible.Fortin and Glowinski [5] suggested an iterative algorithm (called ALG1) for finding approximate solution to (43) as shown below.
Steps of ALG1.
3) Check for convergence.Repeat step 2 when needed.Convergence is reached when the following conditions are met, for a pre-chosen relative tolerance r ε .
It seems that ( 48) is just as difficult to solve as (43) is.Fortin and Glowinski [5] suggested yet another iterative algorithm (called ALG2) for finding approximate solution to (48) as shown below.
3) Repeat step 2 for a certain number of times, or until convergence.Then In particular, we solve (51) and (52) using the standard finite element approach where we approximate ( , )    We'll present numerical results from solving (33) using ALG1 and ALG2 in the next section.

Numerical Tests
In this section, we study a microbeam switch that is electrostatically actuated by an applied voltage.Hu, Chang and Huang [17] first studied such a problem in 2004.The structure [17] can be schematically described as a thin metal beam hanging over a substrate separated by some insulator, where one end of the microbeam is fixed and the other is free (fixed-free beam).The beam is pulled (deflected) towards the substrate when a voltage is applied between the beam and the substrate.The action of such a structure is modelled using (33) subject to boundary conditions (34).We are interested in the gap between the free end of the microbeam and the substrate subject to different applied voltages.
The specifics of the microbeam structure are listed below.• Beam length, width and thickness are 20 mm, 5 mm and 57 µm, respectively.
• Initial gap between microbeam and substrate is 92 µm.
• Poisson's ratio of 0.06 is used because the microbeam is considered wide, i.e., width is significantly greater than thickness.
We present in Table 1 numerical results from solving (33) subject to boundary conditions (34).Our numerical results indicated that the newly introduced augmented Lagrangian methods (ALM) are fully capable of producing quality solutions that match up well with experiment data.The number of finite elements used is denoted by ne in Table 1.
The relative errors of results obtained from ALM are comparable to those in [17] and [16], for all cases of finite element approximations.In fact, it appears

Conclusion
An algorithm is developed based on the augmented Lagrangian methods and the finite element, exploiting the order of the differential equation it solves.Independent of the order of the differential equation, we are always able to use only Hermite cubic and linear finite elements to approximate variables involved.As a result, this algorithm is easy to implement, and is capable of producing accurate and stable solutions to engineering problems that admit weak formulations associated with optimization of some functionals.Extensions of this algorithm for solving engineering problems described by higher order partial differential equations are being investigated by the author.Results will be submitted for publication in the near future.
to certain boundary conditions, where x is the independent variable in interval ( ) , a b ⊂  , k is a positive integer, to(1).Such a solution is also referred to as a strong solution because of the requirement that

are
to be linked to y by the recursive linear equality relations Lagrangian multipliers, which are also functions of x , and { } 1 k j r are pre-chosen constants serving as penalty in case any of the linear equality relations in (14) is not met satisfactorily.

Table 1 .
Comparison of results, ρ 1 = ρ 2 = 0.0208, r = s = 330, and ε r = 0.0001.the augmented Lagrangian methods have already converged with relatively small number (i.e., 12) of finite elements because increasing the number of finite elements used does not cause much of a change in accuracy. that