Identification of Flexural Rigidity in a Kirchhoff Plates Model Using a Convex Objective and Continuous Newton Method
This work provides a detailed theoretical and numerical study of the inverse problem of identifying flexural rigidity in Kirchhoff plate models. From a mathematical standpoint, this inverse problem requires estimating a variable coefficient in a fourth-order boundary value problem. This inverse problem and related estimation problems associated with general plates and shell models have been investigated by numerous researchers through an optimization framework using the output least-squares (OLSs) formulation. OLS yields a nonconvex framework and hence it is suitable for investigating only the local behavior of the solution. In this work, we propose a new convex framework for the inverse problem of identifying a variable parameter in a fourth-order inverse problem. Existence results, optimality conditions, and discretization issues are discussed in detail. The discrete inverse problem is solved by using a continuous Newton method. Numerical results show the feasibility of the proposed framework.
Let be a nonempty bounded domain in with a sufficiently smooth boundary . We consider the following fourth-order elliptic boundary value problem (BVP):where is the usual normal derivative.
The above BVP models the pure bending (without twisting) of a Kirchhoff plate occupying the region . Here, the parameter is the flexural rigidity and the solution is the lateral deflection under the force per unit area (see  for details). The parameter is connected to the thickness of the plate, Young’s modulus of elasticity, and the Poisson ratio of the material. Boundary conditions (1b) and (1c) are the so-called clamped boundary conditions. However, all the results given below can easily be transferred to the pinned boundary conditions given by
Moreover, the proposed framework can be extended to the following more general plate model :
Given the force and the flexural rigidity , the direct problem is to find the lateral deflection . In this work, our objective is to identify the flexural rigidity so that the corresponding deflection is the nearest to a deflection measurement in a certain sense. This inverse problem and its dynamic analogue have been explored mostly through the output least-squares formulation. In the present context, the output least-squares (OLSs) approach seeks to minimize the functional:defined by an appropriate norm. Here, is the data (a measurement of ) and is the unique solution of the weak form of (1a), (1b), and (1c) that corresponds to the flexural rigidity .
One of the major difficulties associated with the OLS approach is the fact that the coefficient-to-solution map is nonlinear and hence, in general, the OLS functional is nonconvex. Consequently, the OLS-based approach only allows the study of local properties of the inverse problem. Moreover, the numerical algorithms designed through the necessary optimality conditions can only ensure convergence to a local optimizer and such conditions are not sufficient optimality conditions.
In the following, we will give a quick overview of some related works. In an interesting paper, Lesnic et al.  investigate the coefficient identification problem in a fourth-order differential equation that governs the deflection of the Euler-Bernoulli beam and give useful stability estimates. An optimization-based approach for the same problem was recently given in . In a recent paper, Marinov and Marinova  developed a numerical method by using the so-called method of variational embedding for solving an inverse problem for bending stiffness estimation in a Kirchhoff-Love plate from overdetermined data. Ewing et al.  study a nonlinear analogue of this problem using the mixed finite element method.
Manservisi and Gunzburger  studied a control problem associated with a plate model and give interesting applications to a simplified model for the “sag bending process” in the manufacturing of automobile windscreens. Inspired by this inverse problem and its useful application, Engl and Kügler  and Kügler  presented a regularized Landweber method for its numerical solution. Salazar and Westbrook  studied the identification of the flexural rigidity of a nonhomogeneous plate of uniform thickness under the assumption that Poisson’s ratio is constant and only the Young modulus varies pointwise. Lopes et al.  used identification techniques to detect holes in plates. In a related work, Kim and Kreider  conducted a numerical study of the inverse problem of parameter identification in nonlinear elastic and viscoelastic plates. They presented numerical results for a variety of constitutive models. In an intriguing paper by Yuan and Yamamoto , they studied Lipschitz stability estimates for both inverse problems of determining spatially varying Lamé coefficients and determining the mass density by a finite number of boundary observations. This inverse problem for the case of interior measurement was studied in [14–16]. A systematic study of inverse problems of parameter identification in plates and shells models was conducted in a series of excellent papers by White [2, 17, 18] where the author investigates various aspects of the inverse problem. An overview of various approaches for parameter identification problems in inverse problems can be found in .
In this work, our primary objective is to employ a new convex functional for the identification of flexural rigidity . The convexity of the new functional circumvents one of the major drawbacks, nonconvexity, of the OLS functional. Furthermore, the proposed convex objective functional characterizes global optimality, and the associated variational inequality is a necessary as well as a sufficient optimality condition. The new objective functional also has a computational advantage in that to compute its gradient one does not need to compute the derivative of the coefficient-to-solution map.
This paper is divided into six sections. Section 2 begins with the essential background material and continues with the introduction of the related minimization problem which is subsequently shown to have a solution. The problem is discretized using finite elements in Section 3 and it is shown that the continuous minimization problem can be approximated by the discrete one. A description of the continuous Newton method is given in Section 4 and the details of numerical experiments and their results are provided in Section 5. The paper concludes with some remarks concerning the outlined approach and other related issues.
2. Modified Output Least-Squares (MOLSs)
The variational formulation of (1a), (1b), and (1c) plays an important role in formulating the MOLS approach. Assuming that the feasible parameters are in , , and by taking , the weak form of (1a), (1b), and (1c) reads: find such that
In the following, we take to be the set of admissible coefficients. A variety of other conditions can be imposed on the feasible coefficients so that the following two inequalities hold for some constants and :
Motivated by [20, 21], we propose the following objective functional to identify the flexural rigidity:where is the measured data and solves the weak form (5).
To give optimality conditions for the above minimization problem using (8) and to design some derivative-based numerical schemes, we need to compute the derivatives of . Since depends on , the derivatives of the solution map are instrumental to evaluate the derivatives of the functional . For this, we give the following result.
Lemma 1. For each in the interior of , is infinitely differentiable at Given , the first derivative is the unique solution of the variational equation:and the second derivative is the unique solution of the variational equation:
Proof. A direct proof relying on (6) and (7) can be extracted from . A proof based on the implicit function theorem can be found in .
The following result gives derivative formulas for (8).
Theorem 2. Let be an arbitrary element in the interior of . Then, consider the following:(1)The first derivative of is given by(2)The second derivative of is given by
Proof. The first derivative of is a direct consequence of the chain rule:We note that, from (9), we haveand henceIt now follows thatwhere in the last step we applied (9) again.
The following result gives a useful feature of MOLS.
Theorem 3. The modified output least-squares functional defined in (8) is convex in the interior of the set
Proof. In view of the form of the second derivative of , it follows from (7) that the following inequality holds for all in the interior of :which proves the convexity of the functional .
The inverse problem of identifying parameters is ill-posed and hence we need to regularize the MOLS functional. In the following, we assume that is the set of admissible coefficients which we further assume to be closed and convex subset of We also assume that is in the interior of the domain of the parameter-to-solution map so that we can use the convexity of MOLS. We now consider the following regularized minimization problem: find by solvingwhere is the regularization parameter, is the data, and is a regularization map. Throughout this work, we take . The associated inner product is denoted by
We have the following useful information concerning a minimizer of the MOLS functional.
Theorem 4. The regularized minimization problem (18) is uniquely solvable. Moreover, is a minimizer if and only if the following variational inequality holds:
Proof. Since , for every , we can find a minimizing sequence such thatThe inequalityconfirms that the sequence is bounded in . By the compact embedding of into , there exists a subsequence that converges weakly in and strongly in . By using the same notation for the subsequences as well, we obtain that in and in . Let be the corresponding sequence of the solutions of the weak form (5). It follows from (6) and (7) that is uniformly bounded in and consequently there exists a weakly convergent subsequence. Once again by using the same notion for the subsequences, we assume that . By the definition of , we havewhich, in view of the convergence properties of and , after passing to the limit implies thatwhere is arbitrary. Since the solution to the above variational problem is unique, we deduce that . In fact, converges to strongly in . From the identities,After a rearrangement of the terms, we deduce thatand the strong convergence follows at once from (6) and (7).
We next claim thatFor this, we notice that the following identities hold due to the definitions of , , , and :For (26), we note that and .
Finally, we havewhere we used the weak-lower semicontinuity of The uniqueness follows from the fact that the regularizer is strictly convex. Finally, since is a convex functional, the variational inequality is a necessary and sufficient optimality condition. The proof is complete.
The continuous problem (18) has to be discretized to obtain a numerical solution. In this work, we use the finite element discretization on a triangulation of . We choose to be the finite-dimensional space of the parameter space relative to Analogously, we assume that is a sequence of finite-dimensional subspaces of and, for each , is a projection operator. We suppose that and have the property that, for every ,
We define and assume that
To keep the structure of finite-dimensional set of admissible sets quite general and also to allow other regularizers, we impose approximation properties given by the following assumption: for any , there exists a sequence such that
The above property is an analogue of (29) for the set of feasible constraints.
We now define as the following discrete analogue of (8):where is the solution of the following discrete variational problem:
Here, the data is replaced by the projected data for computational convenience. Since converges to strongly in , all of the following results hold irrespective of whether or is used in the definition of . We then define the following regularized analogue of the above discrete MOLSand solve
The solvability of the above regularized problems can easily be proved by arguments used above.
We can now prove the following convergence result.
Theorem 5. Assume that, for each , is a solution of (35). Then, any weak limit point of the sequence is a solution of (18).
Proof. By (30), there exists a constant such that , for all . Therefore, the sequence is bounded in and hence it has a subsequence that converges weakly in and strongly in For simplicity of notation, we assume that in and in We will show that is a solution of (18). Let be the corresponding solution of (33). Evidently, is bounded and hence there is a subsequence that weakly converges to some in . We claim that In view of the definition of , we haveLet be arbitrary. We set in the above identity to getwhich after a rearrangement of terms implies that and after passing to the limit , we obtain that and because was chosen arbitrarily, we obtain that solves (5). That is,
By using the arguments given in Theorem 4, we can also show thatwhich, in view of the lower semicontinuity of , implies thatNow consider an arbitrary . By (31a), (31b), and (31c), there exists a sequence such that , for all , in , and . This, however, implies thatand hence, by the optimality of , we haveSince was arbitrary, this shows that is a solution of (8), which completes the proof.
3. Discrete Formulas
To represent the discrete weak form in a computable way, we proceed as follows. We represent bases for and by and , respectively. The space is then isomorphic to , and for any , we define by , , where is a nodal basis corresponding to the nodes . Conversely, each corresponds to defined by . Similarly, will correspond to , where , and . Here, are the nodes of the mesh defining .
Then, the stiffness matrix and the load vector are given by
We also define the matrix by the following condition:
Using the above notation, the following formulas can be verified:
4. Continuous Newton Method
After the discretization, the finite-dimensional minimization problem can be solved by any known optimization solver. In this paper, we will employ the continuous Newton method recently proposed by Zhang et al. . Continuous methods for optimization problems solve an optimization problem by following the solution trajectory of a system of ordinary differential equations. Although numerous authors have focused on continuous optimization methods for solving unconstrained and constrained optimization problems, most of these works have been of a theoretical nature. However, recently, many researchers have advocated the usefulness of continuous methods for solving large scale optimization problems and have given significant numerical evidence to show the competitiveness of these methods with the conventional optimization methods (see [23–26]).
In this work, we will use the continuous Newton method  to solve the finite-dimensional optimization problems given by the regularized MOLS and OLS. For simplicity, the side constraints will be ignored. Zhang et al.  propose following the trajectory given by the following initial value problem:where is defined by
Here, corresponds to either the regularized MOLS or the regularized OLS, the minimum eigenvalue of the Hessian of is given by , and . Moreover,
The idea of the proposed scheme is to combine the steepest descent and Newton’s direction through a convex combination. It involves a simple yet efficient strategy to get around the potential singularities of the Hessian. To be specific, when is between and , if its value is closer to , then will be larger, in turn putting more emphasis on the Newton direction. On the other hand, if is closer to , then will be larger, putting more weight on the gradient direction. The approach is particularly appealing as one of the main challenges in the nonlinear inverse problems is the ill-posedness of the Hessian. Moreover, since the continuous Newton method takes into account the behaviour of the Hessian by examining its minimum eigenvalue, we also gain insight into the relative convexity (or nonconvexity) of the OLS or MOLS functionals. Finally, this approach places less emphasis on computing an optimal regularization parameter, a difficult task in inverse problems.
5. Numerical Experiments
In this section, we present the results of numerical experiments designed to gauge the preliminary effectiveness of the MOLS approach when coupled with the continuous Newton method and applied to the flexural rigidity inverse problem.
5.1. Experiment Setup
We consider first an example boundary value problem (Example 1) derived from (1a), (1b), and (1c):where the solution and the exact parameter are given byand where is subsequently defined by (50). The domain is taken as the unit square, , with the boundary as the square’s outside edges.
Additionally, we considered a similar example (Example 2) with the recovery of the same parameter as Example 1, but on a simply supported circular plate (with appropriate changes to the forward problem applied).
Finally, we considered an example (Example 3) again on the unit square, but with the recovery of a different parameter and with more complicated mixed boundary conditions (i.e., clamped along the top and bottom edges, simply supported along the left edge, and with the right edge free). The parameter recovered in this case was given by
Discretization of the solution was performed using cubic Hermite elements for all examples. For Examples 1 and 3, the mesh used was a regular 25 × 25 triangular grid consisting of 1,352 triangles resulting in 3,529 total degrees of freedom for the solution. For the circular domain in Example 2, a mesh was generated using a Delaunay triangulation yielding 848 triangles with corresponding 2,219 degrees of freedom.
5.2. Algorithm Details
The optimization was performed using the continuous Newton method outlined in the previous section applied to both the discretized MOLS and OLS functionals. Additionally, for means of comparison, the optimization was also performed using a conjugate gradient trust-region-reflective method (MATLAB’s fminunc) using similar stopping criteria. In all examples, seminorm regularization was used for computational simplicity with a fixed regularization parameter used per example (see Table 1 for values).
In the continuous Newton algorithm, the Lanczos iteration was used to compute/estimate the minimum eigenvalue of and in those cases where the Lanczos method failed to converge in the given (limited) number of iterations, the method defaulted to the steepest descent trajectory. The resulting ordinary differential equation was solved using MATLAB’s ode113 solver with default parameters. The ODE was solved over the time domain with for MOLS and for OLS. The additional continuous Newton algorithmic parameters were fixed for each example (again, see Table 1 for values).
5.3. Numerical Results
Figures 1–3, respectively, show the computed parameter at several algorithm steps for the OLS and MOLS approaches using both the continuous Newton and fminunc algorithms. These figures also show the final output of the algorithm and the error between the computed and exact parameters.
In Figure 4, we compare the minimum eigenvalues of the Hessian of both the MOLS and OLS functional over a number of continuous Newton algorithm iterations applied to Example 1. Here, we see that, for the OLS method, the algorithm transitions to both the “mixed” step direction and, when the minimum eigenvalue becomes exceedingly small or negative, the steepest descent direction. We also note that, for the same regularization parameter, the MOLS method uses the full Newton direction. The minimum and maximum of along with the error and the total number of algorithm iterations for all examples are provided in Table 1.
The results of the numerical experiments suggest the following:(i)The MOLS functional converges faster than the OLS functional (Table 1) for the flexural rigidity inverse problem under similar conditions (i.e., similar regularization parameters, optimization algorithm).(ii)In at least one instance for the MOLS approach (see Example 2 in Table 1), the continuous Newton method shows accelerated convergence over the fminunc algorithm.(iii)The eigenvalues of the Hessian of the MOLS functional remain strictly positive compared with the Hessian of the OLS functional, whose minimum eigenvalues can become negative after a number of algorithm iterations (see both Table 1 and Figure 4).(iv)The more complicated domains and boundary conditions in Examples 2 and 3 have a negative impact in terms of algorithm performance and accuracy in recovery of the flexural rigidity parameter. As can be seen in Table 1, this is particularly true in the case of the OLS approach.(v)Overall, for either the MOLS or OLS approach, the error in the recovery of the coefficient is most prominent along the boundaries with significantly better reconstruction in the domain’s interior. Again, in particular, the OLS approach seems to have the most difficulty in reconstructing the flexural rigidity along the boundaries.
7. Concluding Remarks
Further directions for research are as follows:(i)Thorough comparison of the performance and stability differences between the OLS and MOLS approaches.(ii)Thorough comparison of differing approaches to solve the ODE system arising from the continuous Newton method including a comparison of the relative performance and stability of calculating/estimating the eigenvalues of the Hessian of .(iii)Comparison/exploration of differing regularization approaches.(iv)The use of mixed methods for solving the variational problem.
Conflict of Interests
The authors declare that there is no conflict of interests regarding the publication of this paper.
P. G. Ciarlet and P. Rabier, Les équations de von Kármán, vol. 826 of Lecture Notes in Mathematics, Springer, Berlin, Germany, 1980.
L. W. White, “Estimation of elastic parameters in a nonlinear elliptic model of a plate,” Applied Mathematics and Computation, vol. 42, no. 2, pp. 139–187, 1991.View at: Publisher Site | Google Scholar | MathSciNet
D. Lesnic, L. Elliott, and D. B. Ingham, “Analysis of coefficient identification problems associated to the inverse Euler-Bernoulli beam theory,” IMA Journal of Applied Mathematics, vol. 62, no. 2, pp. 101–116, 1999.View at: Publisher Site | Google Scholar | Zentralblatt MATH | MathSciNet
T. T. Marinov and R. S. Marinova, “Coefficient identification in Euler-Bernoulli equation from over-posed data,” Journal of Computational and Applied Mathematics, vol. 235, no. 2, pp. 450–459, 2010.View at: Publisher Site | Google Scholar | Zentralblatt MATH | MathSciNet
T. T. Marinov and R. S. Marinova, “An inverse problem for estimation of bending stiffness in Kirchhoff-Love plates,” Computers & Mathematics with Applications, vol. 65, no. 3, pp. 512–519, 2013.View at: Publisher Site | Google Scholar | MathSciNet
R. E. Ewing, T. Lin, and Y. Lin, “A mixed least-squares method for an inverse problem of a nonlinear beam equation,” Inverse Problems, vol. 15, no. 1, pp. 19–32, 1999.View at: Publisher Site | Google Scholar | MathSciNet
S. Manservisi and M. Gunzburger, “A variational inequality formulation of an inverse elasticity problem,” Applied Numerical Mathematics, vol. 34, no. 1, pp. 99–126, 2000.View at: Publisher Site | Google Scholar | MathSciNet
H. W. Engl and P. Kügler, “The influence of the equation type on iterative parameter identification problems which are elliptic or hyperbolic in the parameter,” European Journal of Applied Mathematics, vol. 14, no. 2, pp. 129–163, 2003.View at: Publisher Site | Google Scholar | MathSciNet
P. Kügler, “A parameter identification problem of mixed type related to the manufacture of car windshields,” SIAM Journal on Applied Mathematics, vol. 64, no. 3, pp. 858–877, 2004.View at: Publisher Site | Google Scholar | MathSciNet
D. Salazar and R. Westbrook, “Inverse problems of mixed type in linear plate theory,” European Journal of Applied Mathematics, vol. 15, no. 2, pp. 129–146, 2004.View at: Publisher Site | Google Scholar | MathSciNet
P. d. Lopes, A. B. Jorge, and J. Sebastião Cunha, “Detection of holes in a plate using global optimization and parameter identification techniques,” Inverse Problems in Science and Engineering, vol. 18, no. 4, pp. 439–463, 2010.View at: Publisher Site | Google Scholar | MathSciNet
S. Kim and K. L. Kreider, “Parameter identification for nonlinear elastic and viscoelastic plates,” Applied Numerical Mathematics, vol. 56, no. 12, pp. 1538–1554, 2006.View at: Publisher Site | Google Scholar | MathSciNet
G. Yuan and M. Yamamoto, “Lipschitz stability in inverse problems for a Kirchhoff plate equation,” Asymptotic Analysis, vol. 53, no. 1-2, pp. 29–60, 2007.View at: Google Scholar | MathSciNet
M. S. Gockenbach, B. Jadamba, and A. A. Khan, “Equation error approach for elliptic inverse problems with an application to the identification of Lamé parameters,” Inverse Problems in Science and Engineering, vol. 16, no. 3, pp. 349–367, 2008.View at: Publisher Site | Google Scholar | MathSciNet
M. S. Gockenbach and A. A. Khan, “Identification of Lamé parameters in linear elasticity: a fixed point approach,” Journal of Industrial and Management Optimization, vol. 1, no. 4, pp. 487–497, 2005.View at: Publisher Site | Google Scholar | MathSciNet
M. S. Gockenbach and A. A. Khan, “An abstract framework for elliptic inverse problems. II. An augmented Lagrangian approach,” Mathematics and Mechanics of Solids, vol. 14, no. 6, pp. 517–539, 2009.View at: Publisher Site | Google Scholar | MathSciNet
L. W. White, “Estimation of elastic parameters in beams and certain plates: regularization,” Journal of Optimization Theory and Applications, vol. 60, no. 2, pp. 305–326, 1989.View at: Publisher Site | Google Scholar | MathSciNet
L. W. White, “Identification of flexural rigidity in a dynamic plate model,” Journal of Mathematical Analysis and Applications, vol. 144, no. 1, pp. 275–303, 1989.View at: Publisher Site | Google Scholar | Zentralblatt MATH | MathSciNet
B. Jadamba, A. A. Khan, and M. Sama, “Inverse problems of parameter identification in partial differential equations,” in Mathematics in Science and Technology, pp. 228–258, World Science Publisher, Hackensack, NJ, USA, 2011.View at: Publisher Site | Google Scholar | MathSciNet
M. S. Gockenbach and A. A. Khan, “An abstract framework for elliptic inverse problems. I. An output least-squares approach,” Mathematics and Mechanics of Solids, vol. 12, no. 3, pp. 259–276, 2007.View at: Publisher Site | Google Scholar | Zentralblatt MATH | MathSciNet
B. Jadamba, A. A. Khan, G. Rus, M. Sama, and B. Winkler, “A new convex inversion framework for parameter identification in saddle point problems with an application to the elasticity imaging inverse problem of predicting tumor location,” SIAM Journal on Applied Mathematics, vol. 74, no. 5, pp. 1486–1510, 2014.View at: Publisher Site | Google Scholar | MathSciNet
L.-H. Zhang, C. T. Kelley, and L.-Z. Liao, “A continuous Newton-type method for unconstrained optimization,” Pacific Journal of Optimization, vol. 4, no. 2, pp. 259–277, 2008.View at: Google Scholar
H. Attouch and M. Teboulle, “Regularized Lotka-Volterra dynamical system as continuous proximal-like method in optimization,” Journal of Optimization Theory and Applications, vol. 121, no. 3, pp. 541–570, 2004.View at: Publisher Site | Google Scholar | MathSciNet
A. A. Brown and M. C. Bartholomew-Biggs, “ODE versus SQP methods for constrained optimization,” Journal of Optimization Theory and Applications, vol. 62, no. 3, pp. 371–386, 1989.View at: Publisher Site | Google Scholar | MathSciNet
A. A. Brown and M. C. Bartholomew-Biggs, “Some effective methods for unconstrained optimization based on the solution of systems of ordinary differential equations,” Journal of Optimization Theory and Applications, vol. 62, no. 2, pp. 211–224, 1989.View at: Publisher Site | Google Scholar | Zentralblatt MATH | MathSciNet
L.-Z. Liao, L. Qi, and H. W. Tam, “A gradient-based continuous method for large-scale optimization problems,” Journal of Global Optimization, vol. 31, no. 2, pp. 271–286, 2005.View at: Publisher Site | Google Scholar | MathSciNet