Loading [MathJax]/jax/output/SVG/jax.js
Research article Special Issues

A numerical method based on barycentric interpolation collocation for nonlinear convection-diffusion optimal control problems

  • This paper describes a study of the barycentric interpolation collocation method for the optimal control problem governed by a nonlinear convection-diffusion equation. Using Lagrangian multipliers, we obtain the continuous optimality system which is composed of state equations, adjoint equations and optimality conditions. Then, barycentric interpolation collocation methods are applied to discretize the optimality system and the nonlinear term is treated by Newton's iteration. Furthermore, the corresponding consistency analyses of discrete schemes are demonstrated. Finally, the validity of the proposed schemes is demonstrated through several numerical experiments. Compared with the classical finite difference method, collocation schemes can yield the higher-order accurate solutions with fewer nodes.

    Citation: Rong Huang, Zhifeng Weng. A numerical method based on barycentric interpolation collocation for nonlinear convection-diffusion optimal control problems[J]. Networks and Heterogeneous Media, 2023, 18(2): 562-580. doi: 10.3934/nhm.2023024

    Related Papers:

    [1] Jin Li, Jinzheng Qu, Xibo Duan, Xiaoning Su . An image encryption algorithm based on heat flow cryptosystems. Networks and Heterogeneous Media, 2023, 18(3): 1260-1287. doi: 10.3934/nhm.2023055
    [2] Xiongfa Mai, Ciwen Zhu, Libin Liu . An adaptive grid method for a singularly perturbed convection-diffusion equation with a discontinuous convection coefficient. Networks and Heterogeneous Media, 2023, 18(4): 1528-1538. doi: 10.3934/nhm.2023067
    [3] Claudio Canuto, Anna Cattani . The derivation of continuum limits of neuronal networks with gap-junction couplings. Networks and Heterogeneous Media, 2014, 9(1): 111-133. doi: 10.3934/nhm.2014.9.111
    [4] Liping Yu, Hans Kleppe, Terje Kaarstad, Svein M. Skjaeveland, Steinar Evje, Ingebret Fjelde . Modelling of wettability alteration processes in carbonate oil reservoirs. Networks and Heterogeneous Media, 2008, 3(1): 149-183. doi: 10.3934/nhm.2008.3.149
    [5] Stefan Berres, Ricardo Ruiz-Baier, Hartmut Schwandt, Elmer M. Tory . An adaptive finite-volume method for a model of two-phase pedestrian flow. Networks and Heterogeneous Media, 2011, 6(3): 401-423. doi: 10.3934/nhm.2011.6.401
    [6] Yaxin Hou, Cao Wen, Yang Liu, Hong Li . A two-grid ADI finite element approximation for a nonlinear distributed-order fractional sub-diffusion equation. Networks and Heterogeneous Media, 2023, 18(2): 855-876. doi: 10.3934/nhm.2023037
    [7] Jérôme Droniou . Remarks on discretizations of convection terms in Hybrid mimetic mixed methods. Networks and Heterogeneous Media, 2010, 5(3): 545-563. doi: 10.3934/nhm.2010.5.545
    [8] Iryna Pankratova, Andrey Piatnitski . Homogenization of convection-diffusion equation in infinite cylinder. Networks and Heterogeneous Media, 2011, 6(1): 111-126. doi: 10.3934/nhm.2011.6.111
    [9] Didier Georges . Infinite-dimensional nonlinear predictive control design for open-channel hydraulic systems. Networks and Heterogeneous Media, 2009, 4(2): 267-285. doi: 10.3934/nhm.2009.4.267
    [10] Jinhu Zhao . Natural convection flow and heat transfer of generalized Maxwell fluid with distributed order time fractional derivatives embedded in the porous medium. Networks and Heterogeneous Media, 2024, 19(2): 753-770. doi: 10.3934/nhm.2024034
  • This paper describes a study of the barycentric interpolation collocation method for the optimal control problem governed by a nonlinear convection-diffusion equation. Using Lagrangian multipliers, we obtain the continuous optimality system which is composed of state equations, adjoint equations and optimality conditions. Then, barycentric interpolation collocation methods are applied to discretize the optimality system and the nonlinear term is treated by Newton's iteration. Furthermore, the corresponding consistency analyses of discrete schemes are demonstrated. Finally, the validity of the proposed schemes is demonstrated through several numerical experiments. Compared with the classical finite difference method, collocation schemes can yield the higher-order accurate solutions with fewer nodes.



    Optimal control problems governed by convection-diffusion equations arise in many real-world problems such as wastewater treatment problems [1] and air pollution problems [2]. Especially it is also an important step toward optimal control problems of the fluid. Such problems have profound application prospects, which simulate some chemical or biological process between different species. It is difficult or impossible to obtain their analytical solutions, and the only way to solve them is to seek approximate solutions by using numerical methods. Thus, the extraction of their numerical algorithms is of great practical significance for scholars. In this paper, we will propose the barycentric interpolation collocation methods for optimal control problems governed by the nonlinear convection-diffusion equation, which can be written as

    minJ(v,u)=12Ω(vˆv)2dΩ+w2Ωu2dΩ,s.t αΔv+βv+γv+ϕ(v)=f+u, vΩ,v=0, vΩ. (1.1)

    where Δ is the Laplace operator, is the gradient operator and v and u are the state variable and control variable, respectively. ˆv is the desired variable, and ϕ(v) represents the nonlinear term. w is the so-called regularization parameter. J(v,u) stands for the cost function. Ω is a bounded domain with the boundary Ω.

    In recent years, a large number of numerical methods have been presented for optimal control problems.

    Weng et al. [3] presented a stabilized finite element method for optimal control problems governed by a convection dominated diffusion equation. Sandilya and Kumar [4] proposed a discontinuous interpolated finite volume approximation of semilinear elliptic optimal control problems. Y¨ucel et al. [5] applied a discontinuous Galerkin method for optimal control problems governed by a system of convection-diffusion PDEs with nonlinear reaction terms. Lu et al. [6] presented a priori error estimates of interpolation coefficients mixed finite element methods for semilinear Dirichlet boundary optimal control problems. Lin et al. [7] developed Galerkin spectral methods for optimal control problem governed by elliptic equations. Porwal and Shakya [8] presented a finite element method for an elliptic optimal control problem with integral state constraints. In [9], a bilinear pseudo-spectral method based on Chebyshev polynomials was used for convection-diffusion optimal control problems, and the coupled Sylvester system was solved by some direct and iterative methods. Wang et al. [10] investigated a spectral Galerkin approximation of an optimal control problems governed by a fractional advection-diffusion-reaction equation with an integral fractional Laplacian. Casanova et al. [11] adopted radial basis function techniques for optimal constrained optimization problems governed by linear convection-diffusion PDEs.

    Recently, many researchers have developed a variety of meshless methods [12,13,14] for convection-diffusion problems. In particular, the barycentric interpolation collocation is a novel meshless method, which is different from the traditional element-based numerical method. It can effectively avoid the cumulative error caused by the difference scheme. Additionally, it can handle irregular domains. This method has attracted the attention of scholars due to its high precision, fast speed, good stability and convenient program execution. It was extended to solve various PDEs including linear optimal control problems [15,16], Sine-Gordon equations [17], high-dimensional Fredholm integral equations [18], telegraph equations [19], viscoelastic wave equations [20], the Allen-Cahn equation [21], nonlinear parabolic partial differential equations [22], etc. Lately, Yi and Yao [23] presented error analysis of a barycentric Lagrange interpolation (BLI) collocation scheme.

    To the best of our knowledge, there are few studies about using the barycentric interpolation collocation method for the nonlinear optimal control problem described by Eq (1.1). Based on the above work, we proposed two numerical schemes for the nonlinear optimal control in two-dimensions by using BLI collocation or barycentric rational interpolation (BRI) collocation. Moreover, consistency analysis of discrete schemes is presented.

    The rest of this paper is organized as follows: First, we obtain the continuous optimality system using Lagrangian multipliers in Section 2. Then, barycentric interpolation collocation methods are briefly introduced in Section 3. Two barycentric interpolation collocation schemes are proposed in Section 4, and we also provide the consistency analysis of discrete schemes. In Section 5, numerical simulations are presented to verify the accuracy and efficiency of the proposed collocation methods. Finally, conclusions are drawn in Section 6.

    In this section, we consider the optimize-then-discretize approach based on the Lagrangian multiplier method for solving Eq (1.1). First, we define the Lagrangian multiplier p on Ω as follows

    L=12Ω(vˆv)2dΩ+w2Ωu2dΩ+Ω(αΔv+βv+γv+ϕ(v)fu)pdΩ+Ωvpds . (2.1)

    Then, the state equations are obtained by differentiating the derivative of L with respect to p

    {αΔv+βv+γv+ϕ(v)=f+u, vΩ,v=0, vΩ. (2.2)

    Similarly, the adjoint equations are gained by calculating the Frˊechet derivative of L with respect to v

    {αΔpβp+γp+ϕ(v)p=vˆv, pΩ,p=0, pΩ. (2.3)

    Finally, by differentiating with respect to u, the optimality equation is written as follows

    wup=0, onΩ. (2.4)

    Combining Eqs (2.2)–(2.4), Eq (1.1) can be transformed into a continuous optimality system

    {αΔv+βv+γv+ϕ(v)=f+u, vΩ,v=0, vΩ,αΔpβp+γp+ϕ(v)p=vˆv, pΩ,p=0, pΩ,wup=0,onΩ. (2.5)

    In this section, we mainly introduce BLI and BRI collocation.

    This part presents a novel meshless method named BLI to solve nonlinear optimal control problems. Suppose that m+1 distinct nodes xk(k=0,1,,m) are given together with a set of functional values hk at the discrete nodes xk correspondingly. Let q(x) denote the approximate interpolation polynomial of h(x), satisfying q(xk)=hk (k=0,1,,m).

    As q(x) is unique, we can rewrite it in the Lagrange interpolation polynomial form

    h(x)q(x)=mk=0Lk(x)hk, k=0,1,,m, (3.1)

    where Lk(x) is the Lagrange basis function, and Lk(x)=mik(xxi)mik(xkxi), i=0,1,,m.

    Define ˜L(x)=(xx0)(xx1)(xxm); thus, the basis function Lk(x) can be expressed as

    Lk(x)=˜L(x)ωkxxk, k=0,1,,m, (3.2)

    where ωk=1ki(xkxi) is the barycentric weight.

    Supposing that q(x)=1 in Eq (3.1) and inserting Eq (3.2) into Eq (3.1), the results are shown as

    {1=˜L(x)mk=0ωkxxk,(3.3a)q(x)=˜L(x)mk=0ωkxxkhk.(3.3b)

    By dividing Eq (3.3b) by Eq (3.3a) and canceling the common polynomial ˜L(x), the BLI formula for v(x) can be expressed as

    q(x)=mk=0ωkxxkhkMk=0ωkxxk:=mk=0ψk(x)hk. (3.4)

    As a matter of fact, the BLI formula has excellent numerical stability when Chebyshev points are selected, so the nodes selected in this article and the corresponding barycentric weights are

    xk=cos(kmπ),k=0,1,,m. (3.5)
    ωk=(1)kδk,δk={12,k=0orm,1,else. (3.6)

    In this subsection, we consider the derivation process of the BRI collocation method. Suppose that (x0,h0),(x1,h1),(xm,hm) are given m+1 interpolation nodes of the interval [a,b], and that the function h(x) satisfies h(xk)=hk,k=0,1,2m. Additionally, d(0dm) is a chosen integer; We define the following BRI formula:

    r(x)=mk=0wkxxkhkmk=0wkxxk=mk=0rk(x)hk, (3.7)

    where rk(x)=wkxxkmk=0wkxxk is the basic function and the interpolation weight wk=iJk(1)ii+dl=i,lk1xkxl, Jk={iJ|kdik}, J={0,1,,md}. More details about the BRI formula are available [24,25].

    Similarly, the BRI formula for the binary function h(x,y) with rM,N(x,y) can be obtained. Let us partition the domain Ω=[a,b]×[c,d] into Ωh=(xk,yj),k=0,1,M,j=0,1,2,N. Then d1andd2 are two integers, where 0d1M, 0d2N, h(xk,yj)=hkj. Thus, rM,N(x,y) can be expressed as

    h(x,y):≈rM,N(x,y)=Mk=0Nj=0rk(x)rj(y)hkj, (3.8)

    where rk(x) and rj(y) are basis functions, which can be expressed as

    rk(x)=ˉwkxxkMs=0ˉwsxxs,ˉwk=iJk(1)ii+d1l=i,lk1xkxl,Jk={k0,1,,Md1:kd1ik}, (3.9)
    rj(y)=ˉvjyxjNs=0ˉvsyys,ˉvj=iJj(1)ii+d2l=i,lj1yjyl,Jj={i0,1,,Nd2:jd2ij}. (3.10)

    In order to solve the optimal control problem, derivatives expressions of the barycentric interpolation formula should be derived. Obviously, the μ-order derivative of q(x) at the different nodes xi (i=0,1,,M) can be obtained as follows:

    q(μ)(xi)=Mk=0ψ(μ)k(xi)hk:=Mk=0D(μ)ikhk, μ=1,2,, (3.11)

    where D(μ)ik denotes the element of the differential matrix D(μ), the recursive formula of which [26,27] is as follows according to the principles of mathematical induction

    {D(1)ik=wiwk1(xixk),ik,D(1)ii=Mk=0,kiD(1)ik. (3.12)
    {D(μ)ik=μ(D(1)ikD(μ1)iiD(μ1)ikxixk),ik,D(μ)ii=Mk=0,kiD(μ)ik. (3.13)

    In this part, we will establish numerical schemes for Eq (1.1) by discretizing x and y directions with barycentric interpolation formulae. We consider the domain Ω=[0,1]×[0,1], and that the x direction is discretized by M+1 nodes with N+1 nodes along the y direction. Supposing that β=(β1,β2), the optimality conditions of Eq (1.1) are as follows

    {α2vx2α2vy2+β1vx+β2vy+γv+ϕ(v)1wp=f, vΩ,v(0,y)=0,v(1,y)=0, y[0,1],v(x,0)=0,v(x,1)=0, x[0,1]. (4.1)
    {α2px2α2py2β1pxβ2py+γp+ϕ(v)pv=ˆv, pΩ,p(0,y)=0,p(1,y)=0, y[0,1],p(x,0)=0,p(x,1)=0, x[0,1]. (4.2)

    According to the barycentric interpolation formula, v(x,y) and p(x,y) can be written as follows

    v(x,y)=Mk=0Nj=0Ck(x)Dj(y)vkj, p(x,y)=Mk=0Nj=0Ck(x)Dj(y)pkj,  (4.3)

    where Ck(x) and Dj(y) represent basis functions of the barycentric interpolation method along the x and y directions. And vkj=v(xk,yj)andpkj=p(xk,yj).

    Consider the λ+μ order partial derivative of v(x,y) with respect to the spatial variables x and y; we have

    λ+μv(x,y)xλyμ:=v(λ,μ)(x,y)=Mk=0Nj=0C(λ)k(x)D(μ)j(y)vkj, λ,μ=0,1,2,. (4.4)

    Thus, the approximate values of the partial derivatives of v(x,y) at the node (xq,ys) can be written as

    v(λ,μ)(xq,ys)=Mk=0Nj=0C(λ)k(xq)D(μ)j(ys)vkj, q=0,1,2,,M,s=0,1,2,,N. (4.5)

    Similarly, we consider the λ+μ order partial derivative of p(x,y) with respect to the variables x and y, and its approximate value at the node (xq,ys) is

    {p(λ,μ)(x,y)=Mk=0Nj=0C(λ)k(x)D(μ)j(y)pkj, λ,μ=0,1,2,,p(λ,μ)(xq,ys)=Mk=0Nj=0C(λ)k(xq)D(μ)j(ys)pkj, q=0,1,,M,s=0,1,,N. (4.6)

    Selecting a suitable λ and μ and then substituting Eqs (4.4)–(4.6) into Eqs (4.1) and (4.2), is implied that

    {αMk=0Nj=0C(2)k(xq)D(0)j(ys)vkjαMk=0Nj=0C(0)i(xq)D(2)j(ys)vkj+β1Mk=0Nj=0C(1)k(xq)D(0)j(ys)vkj+β2Mk=0Nj=0C(0)k(xq)D(1)j(ys)vkj+γvkj+ϕ(vkj)1wpkj=fkj,αMk=0Nj=0ϕC(2)k(xq)D(0)j(ys)pkjαMk=0Nj=0C(0)i(xq)D(2)j(ys)pkjβ1Mk=0Nj=0C(1)k(xq)D(0)j(ys)pkjβ2Mk=0Nj=0C(0)k(xq)D(1)j(ys)pkj+γpkj+ϕ(vkj)pkjvkj=ˆvkj, (4.7)

    with the following boundary conditions

    {v0j=v(x0,yj)=0,vMj=v(xM,yj)=0, j=1,N1,vk0=v(xk,y0)=0,vkN=v(xk,yN)=0, k=1,M1,p0j=p(x0,yj)=0,pMj=p(xM,yj)=0, j=1,N1,pk0=p(xk,y0)=0,pkN=p(xk,yN)=0, k=1,M1. (4.8)

    This system can be derived by applying the matrix form

    {diag(α)(C(2)IN+IMD(2))V+diag(β1)(C(1)IN)V+diag(β2)(IMD(1))V+diag(γ)V+ϕ(V)1wP=F,diag(α)(C(1)IN+IMD(2))Pdiag(β1)(C(1)IN)Pdiag(β2)(IMD(1))P+diag(γ)PV+ϕ(V)P=ˆV, (4.9)

    where stands for the Kronecker product of the matrix, C(i),D(i)(i=1,2) are differentiating matrices of order i and V,P,F,andˆV are column vectors, which can be respectively expressed as follows

    {V=[v00,,v0N,v10,,v1N,vM0,,vMN]T,P=[p00,,p0N,p10,,p1N,pM0,,pMN]T,F=[f00,,f0N,f10,,f1N,,fM0,,fMN]T,ˆV=[ˆv00,,ˆv0N,ˆv10,,ˆv1N,,ˆvM0,,ˆvMN]T. (4.10)

    In this section, we present consistency estimates of proposed schemes based on the approximation properties of the BLI or BRI. For the unknown function v(x,y), the corresponding barycentric Lagrange interpolation is vh(x,y), so the error ζ(x,y) can be defined as

    ζ(x,y)=v(x,y)vh(x,y). (4.11)

    The approximation properties of BLI have been presented [23] as follows

    Lemma 1. If ζ(x,y)C(n+1)(Ω), Ω=[a,b]×[c,d] is a smooth and bounded domain, it holds that

    {|ζ(x,y)|≤∥v(n+1)(cx(eLx2M)M+cy(eLy2N)N),|ζ(x,y)x|≤∥v(n+1)(˙cx(eLx2(M1))M1+cy(eLy2N)N),|ζ(x,y)y|≤∥v(n+1)(cx(eLx2M)M+˙cy(eLy2(N1))N1),|2ζ(x,y)x2|≤∥v(n+1)(¨cx(eLx2(M2))M2+cy(eLy2N)N),|2ζ(x,y)y2|≤∥v(n+1)(cx(eLx2M)M+¨cy(eLy2(N2))N2), (4.12)

    where e is the natural logarithm, cx,cy,˙cx,˙cy,¨cxand¨cy are positive constants, Lx=ba2 and Ly=dc2.

    Suppose that v(xq,ys) and p(xq,ys) are the corresponding BLIs of v(x,y) and p(x,y), and we define the differential operators G1 and G2 as follows

    G1v(x,y)=f(x,y),G2p(x,y)=ˆv(x,y), (4.13)

    and

    lim(M,N)+G1v(x,y)=f(x,y),lim(M,N)+G2p(x,y)=ˆv(x,y). (4.14)

    Next, the consistency analysis of the BLI scheme is presented.

    Theorem 1. If v(x,y), p(x,y) and u(x,y) C(n+1)(Ω), Ω=[a,b]×[c,d] and the nonlinear term ϕ(v)W2,(Ω) satisfies the Lipschitz condition, it holds that

    |v(x,y)v(xq,ys)|Cv(n+1)((eLx2(M2))M2+(eLy2(N2))N2), (4.15)
    |p(x,y)p(xq,ys)|ˆCp(n+1)((eLx2(M2))M2+(eLy2(N2))N2), (4.16)
    |u(x,y)u(xq,ys)|˜Cp(n+1)((eLx2(M2))M2+(eLy2(N2))N2). (4.17)

    where C,ˆCand˜C are positive constants and W2,(Ω)={vL(Ω)αvL(Ω),|α|2} is the Sobolev space.

    Proof. The definitions of the operators G1 and G2 lead to the following result

    G1v(x,y)+G2p(x,y)=α[vxx(x,y)+vyy(x,y)]+β1vx(x,y)+β2vy(x,y)+ϕ(v(x,y))+(γ1)v(x,y)1wp(x,y)α[pxx(x,y)+pyy(x,y)]β1px(x,y)β2py(x,y)+γp(x,y)+ϕ(v(x,y))p(x,y)=f(x,y)ˆv(x,y). (4.18)

    The discrete equation corresponding to Eq (4.18) is as follows

    G1v(xq,ys)+G2p(xq,ys)=α[vxx(xq,ys)+vyy(xq,ys)]+β1vx(xq,ys)+β2vy(xq,ys)+ϕ(v(xq,ys))+(γ1)v(xq,ys)1wp(xq,ys)α[pxx(xq,ys)+pyy(xq,ys)]β1px(xq,ys)β2py(xq,ys)+γp(xq,ys)+ϕ(v(xq,ys))p(xq,ys)=f(xq,ys)ˆv(xq,ys). (4.19)

    To simplify the analysis process, we first consider the terms related to the state variable v(x,y). Combining Eq (4.18) with Eq (4.19), we obtain

    G1[v(x,y)v(xq,ys)]=α[vxx(x,y)vxx(xq,ys)]α[vyy(x,y)vyy(xq,ys)]+β1[vx(x,y)vx(xq,ys)]+β2[vy(x,y)vy(xq,ys)]+ϕ(v(x,y))ϕ(v(xq,ys))+(γ1)[v(x,y)v(xq,ys)]:=ζ1+ζ2+ζ3+ζ4+ζ5+ζ6, (4.20)

    where ζ1,ζ2,ζ3,ζ4,ζ5andζ6 represent

    ζ1:=α[vxx(x,y)vxx(xq,ys)],ζ2:=α[vyy(x,y)vyy(xq,ys)],ζ3:=β1[vx(x,y)vx(xq,ys)], ζ4:=β2[vy(x,y)vy(xq,ys)],ζ5:=ϕ(v(x,y))ϕ(v(xq,ys)),ζ6:=(γ1)[v(x,y)v(xq,ys)]. (4.21)

    Similarly, we have

    G2[p(x,y)p(xq,ys)]=α[pxx(x,y)pxx(xq,ys)]α[pyy(x,y)pyy(xq,ys)]β1[px(x,y)px(xq,ys)]β2[py(x,y)py(xq,ys)]+(γ1w)[p(x,y)p(xq,ys)]+ϕ(v(x,y))p(x,y)ϕ(v(xq,ys))p(xq,ys):=R1+R2+R3+R4+R5+R6, (4.22)

    where R1,R2,R3,R4,R5andR6 represent

    R1:=α[pxx(x,y)pxx(xq,ys)],R2:=α[pyy(x,y)pyy(xq,ys)],R3:=β2[px(x,y)px(xq,ys)],R4:=β2[py(x,y)py(xq,ys)],R5:=(γ1w)[p(x,y)p(xq,ys)],R6:=ϕ(v(x,y))p(x,y)ϕ(v(xq,ys))p(xq,ys). (4.23)

    Applying Lemma 1, it holds that

    |ζ1|=α|vxx(x,y)vxx(xq,y)+vxx(xq,y)vxx(xq,ys)|α[|vxx(x,y)vxx(xq,y)|+|vxx(xq,y)vxx(xq,ys)|]˜C1αv(n+1)(¨cx(eLx2(M2))M2+cy(eLy2N)N)C1v(n+1)((eLx2(M2))M2+(eLy2N)N), (4.24)
    |ζ2|=α|vyy(x,y)vyy(x,ys)+vyy(x,ys)vyy(xq,ys)|C2v(n+1)((eLx2M)M+(eLy2(N2))N2), (4.25)
    |ζ3|=β1|vx(x,y)vx(xq,y)+vx(xq,y)vx(xq,ys)|C3v(n+1)((eLx2(M1))M1+(eLy2N)N), (4.26)
    |ζ4|=β2|vx(x,y)vx(x,ys)+vx(x,ys)vx(xq,ys)|C4v(n+1)((eLx2M)M+(eLy2(N1))N1). (4.27)

    Since ϕ(v) satisfies the Lipschitz condition, assuming a positive constant K, we obtain

    |ζ5|Kv(x,y)v(xq,ys)C5v(n+1)((eLx2M)M+(eLy2N)N). (4.28)

    As for the ζ6, we have

    |ζ6|≤∣γ1∣∣v(x,y)v(xq,ys)C6v(n+1)((eLx2M)M+(eLy2N)N). (4.29)

    Taking C=max{C1,C2,C3,C4,C5,C6}, and substituting Eqs (4.24)–(4.29) into Eq (4.20), we complete the proof of Eq (4.15).

    Similarly, the approximation errors of R1,R2,R3andR4 can be straightforwardly proved. Now we consider the derivation process of R5andR6. Applying the Lemma 1, we derive

    |R5|=|(γ1w)[p(x,y)p(xq,ys)]|˜C5p(n+1)((eLx2M)M+(eLy2N)N). (4.30)

    As ϕ(v(x,y))W2,(Ω), we have

    |R6|=|ϕ(v(x,y))p(x,y)ϕ(v(xq,ys))p(xq,ys)|=|ϕ(v(x,y))p(x,y)ϕ(v(xq,ys))p(x,y)+ϕ(v(xq,ys))p(x,y)ϕ(v(xq,ys))p(xq,ys)||[ϕ(v(x,y))ϕ(v(xq,ys))]p(x,y)|+|ϕ(v(xq,ys))[p(x,y)p(xq,ys)]||[ϕ(v(xγ,yξ))(v(x,y)v(xq,ys))]||p(x,y)|+˙C6[p(x,y)p(xq,ys)]|˜C6v(n+1)((eLx2M)M+(eLy2N)N). (4.31)

    Finally, substituting the error results of R1,R2,R3,R4,R5andR6 into Eq (4.22), the proof of Eq (4.16) can be completed.

    As for the variable u(x,y), it can be proved that

    |u(x,y)u(xq,ys)|=1w|p(x,y)p(xq,ys)|˜Cp(n+1)((eLx2(M2))M2+(eLy2(N2))N2). (4.32)

    This completes the proof of Eq (4.17).

    In the following part, we present the error analysis of Eq (1.1), as solved by the BRI scheme. The compatibility errors of two-dimensional nonlinear optimal control problems solved by the BRI method are as follows based on Theorem 1 in [19]

    Theorem 2. If v(x,y), p(x,y) and u(x,y) Cd1+4[a,b]×Cd1+4[c,d], v(xq,ys) and p(xq,ys) are corresponding numerical solutions solved by BRI. The nonlinear term ϕ(v)W2,(Ω) satisfies the Lipschitz condition, so we have

    |v(x,y)v(xq,ys)|C(hd11x+hd11y), (4.33)
    |p(x,y)p(xq,ys)|C(hd11x+hd11y), (4.34)
    |u(x,y)u(xq,ys)|C(hd11x+hd11y). (4.35)

    where C is the positive constant independent of the grid size hx and hy; hxandhy are step sizes of the discretized grids in the x and y directions, respectively.

    The nonlinear term ϕ(Vk+1) is expanded by using Taylor's formula at Vk, and ϕ(Vk+1)Pk+1 is treated with

    ϕ(Vk+1)ϕ(Vk)+ϕ(Vk)(VVK),ϕ(Vk+1)Pk+1ϕ(Vk)Pk+1. (4.36)

    Thus, we derive the numerical scheme of Eq (1.1) based on Eq (4.9) as follows

    {{diag(α)(C(2)IN+IMD(2))+diag(β1)(C(1)IN)+diag(β2)(IMD(1))}Vk+1+diag(γ+ϕ(Vk))Vk+11wPk+1=ϕ(Vk)+ϕ(Vk)Vk+Fk+1,{diag(α)(C(2)IN+IMD(2))diag(β1)(C(1)IN)diag(β2)(IMD(1))}Pk+1+diag(γ+ϕ(Vk))Pk+1Vk+1=ˆVk+1. (4.37)

    Rewrite Eq (4.37) in matrix form

    [A+B+Cdiag(1w)IAB+C][Vk+1Pk+1]=[ϕ(Vk)+ϕ(Vk)Vk+Fk+1ˆVk+1], (4.38)

    where A=diag(α)(C(2)IN+IMD(2), B=diag(β1)(C(1)IN)+diag(β2)(IMD(1))} and C=diag(γ+ϕ(Vk)).

    We will present two numerical examples to illustrate the accuracy and stability of our proposed methods in this section.

    Consider Ω=[0,1]×[0,1]; the exact solutions of the state and adjoint functions are as follows

    {v(x,y)=sin(πx)sin(πy),p(x,y)=π2sin(πx)sin(πy), (5.1)

    with

    f=(1+2απ2π2w)sin(πx)sin(πy)+β1πcos(πx)sin(πy)+β2πsin(πx)cos(πy)+(sin(πx)sin(πy))3,
    ˆv=(12απ4π2)sin(πx)sin(πy)+β1π3cos(πx)sin(πy)+β2π3sin(πx)cos(πy)3π2(sin(πx)sin(πy))3.

    Take the simulation parameters as follows: α=0.3, β=[1,1], γ=1, w=0.5, d1=[3M4], d2=[3N4], and vary the numbers of mesh nodes M,N, ϕ(v(x,y))=(v(x,y))3. The L norm errors of different variables v,p,u solved by a BLI, BRI or five-point finite difference (FD) scheme are shown in Tables 13. The computational results in Tables 1, 2 show that. Obviously, both BLI and BRI schemes for Eq (1.1) produce the higher-order accurate solutions when selecting fewer nodes. Furthermore, it is observed that the accuracy of the BLI method is slightly higher than that of the BRI method when taking the same mesh nodes. In Table 3, the obtained results illustrate that the convergence rates of FD scheme are almost second order, and thus are in line with the theoretical analysis. It can be seen that the barycentric interpolation collocation method achieves a high accuracy on the order of 107 with 8×8 mesh nodes, while the five-point FD approach can only achieve an accuracy on the order of 104 using 64×64 mesh nodes. This shows that the barycentric interpolation collocation method possesses higher accuracy than the FD approach.

    Table 1.  The L errors of different variables v,p,u solved by the BLI scheme.
    (N,M) E(v) E(p) E(u)
    (4, 4) 4.3856 ×103 9.0961×104 9.0961×104
    (8, 8) 1.7191×107 1.0134×107 1.0134×107
    (12, 12) 1.6189×1012 1.1602×1012 1.1602×1012
    (16, 16) 9.0836×1014 4.8873×1014 4.8873×1014

     | Show Table
    DownLoad: CSV
    Table 2.  The L errors of different variables v,p,u solved by the BRI scheme.
    (N,M) E(v) E(p) E(u)
    (4, 4) 4.3856 ×103 9.0691×104 9.0691×104
    (8, 8) 2.1043×106 1.1426×106 1.1426×106
    (12, 12) 1.4114×1010 9.5459×1011 9.5459×1011
    (16, 16) 2.0288×1013 1.0036×1013 1.0036×1013

     | Show Table
    DownLoad: CSV
    Table 3.  The L errors, convergence rates of different variables v,p,u solved by FD scheme.
    (N,M) E(v) Order E(p) Order E(u) Order
    (8, 8) 1.6718 ×102 6.6866×103 6.6866×103
    (16, 16) 4.2523×103 1.9751 1.7503×103 1.9336 1.7503×103 1.9336
    (32, 32) 1.0727×103 1.9869 4.3889×104 1.9957 4.3889×104 1.9957
    (64, 64) 2.6842×104 1.9988 1.0996×104 1.9969 1.0996×104 1.9969

     | Show Table
    DownLoad: CSV

    Select the nodes M=N=16 to obtain the exact solution images, numerical solution images and error images of v, p, as shown in Figures 1 and 2. We can conclude that numerical solution images solved by the BLI and BRI schemes are all consistent with analytical solution images, indicating that collocation schemes are stable. Figure 3 shows that BLI and BRI schemes have exponential convergence effect when choosing different nodes. In addition, the five-point FD scheme has a second-order convergence rate.

    Figure 1.  Exact solution images, numerical solution images and error images of v, p solved by BLI scheme for Example 1 (M=N=16).
    Figure 2.  Exact solution images, numerical solution images and error images of v, p solved by BRI scheme for Example 1 (M=N=16).
    Figure 3.  Comparison results of convergence rates in space solved by different schemes for Problem 1.

    In this simulation, we consider the following exact solutions

    v(x,y)=(x2x)(y2y)ex+y,
    p(x,y)=4π(x2x)(y2y)ex+y,

    with

    ϕ(v(x,y))=(v(x,y))3.

    We consider Eq (1.1) on the computational domain Ω=[0,1]×[0,1] with homogeneous Dirichlet boundary condition and set the parameters α=0.1,β=[1,1],γ=1,w=0.3 respectively, d1=[3M4], d2=[3N4]. The simulation has confirmed the high accuracy and efficiency of the proposed schemes through the calculated errors in Tables 4, 5. Furthermore, the convergence rate of the FD method is almost second order in Table 6.

    Table 4.  The L errors of different variables v,p,u solved by BLI scheme.
    (N,M) E(v) E(p) E(u)
    (4, 4) 3.9671 ×102 3.7965×103 3.7965×103
    (8, 8) 2.6812×108 7.9556×109 7.9556×109
    (12, 12) 1.3786×1013 3.9897×1014 3.9897×1014
    (16, 16) 8.8305×1014 1.4051×1014 1.4051×1014

     | Show Table
    DownLoad: CSV
    Table 5.  The L errors of different variables v,p,u solved by BRI scheme.
    (N,M) E(v) E(p) E(u)
    (4, 4) 3.9671 ×102 3.7965×103 3.7965×103
    (8, 8) 1.8400×106 5.3447×107 5.3447×107
    (12, 12) 4.3919×1012 1.3874×1012 1.3874×1012
    (16, 16) 1.5009×1013 2.2794×1014 2.2794×1014

     | Show Table
    DownLoad: CSV
    Table 6.  The L errors, convergence rates of different variables v,p,u solved by FD scheme.
    (N,M) E(v) Order E(p) Order E(u) Order
    (8, 8) 4.3345 ×101 2.4250 ×102 2.4250 ×102
    (16, 16) 1.0308×101 2.0721 5.7603×103 2.8271 5.7603×103 2.8271
    (32, 32) 2.5497×102 2.0154 1.4259×103 2.0738 1.4259×103 2.0738
    (64, 64) 6.3546×103 2.0045 3.5574×104 2.0030 3.5574×104 2.0030

     | Show Table
    DownLoad: CSV

    Figures 4, 5 show the exact solutions, numerical solutions and error estimates for the state and adjoint functions when selecting Chebyshev nodes M=N=48. From these images, it is obvious that two barycentric interpolation collocation schemes are effective and stable with higher accuracy compared with the FD method.

    Figure 4.  Exact solution images, numerical solution images and error images of v, p solved by BLI scheme for Problem 2 (M=N=48).
    Figure 5.  Exact solution images, numerical solution images and error images of v, p solved by BRI scheme for Problem 2 (M=N=48).

    Figure 6 illustrates the convergence rates of different schemes and shows that both the BLI and BRI schemes can achieve exponential convergence.

    Figure 6.  Comparing results for space accuracy as solved by different schemes for Problem 2.

    In this paper, we proposed two barycentric interpolation collocation schemes for optimal control problems governed by the nonlinear convection-diffusion equation. The consistency analyses of the two collocation schemes have been derived. The numerical examples show that the efficiency of our proposed schemes which can achieve the higher-oder accurate solutions with fewer nodes compared with the FD method. Our future work will focus on other kinds of PDE-constrained optimal control problems such as parabolic equations, fluid equations and so on.

    This work is was supported in part by the NSF of China (No. 11701197) and the Natural Science Foundation of Fujian Province (No. 2022J01308).

    The authors declare that there is no conflict of interest.



    [1] A. Martínez, C. Rodríguez, M. E. Vázquez-Méndez, Theoretical and numerical analysis of an optimal control problem related to wastewater treatment, SIAM J. Control Optim, 38 (2000), 1534–1553. https://doi.org/10.1137/S0363012998345640 doi: 10.1137/S0363012998345640
    [2] J. Zhu, Q. C. Zeng, A mathematical formulation for optimal control of air pollution, Sci. China Earth Sci., 46 (2003), 994–1002. https://doi.org/10.1007/BF02959394 doi: 10.1007/BF02959394
    [3] Z. F. Weng, J. Z. Yang, X. L. Lu, A stabilized finite element method for the convection dominated diffusion optimal control problem, Appl. Anal., 95 (2016), 2807–2823. https://doi.org/10.1080/00036811.2015.1114606 doi: 10.1080/00036811.2015.1114606
    [4] R. Sandilya, S. Kumar, A discontinuous interpolated finite volume approximation of semilinear elliptic optimal control problems, Numer Methods Partial Differ Equ, 33 (2017), 2090–2113. https://doi.org/10.1002/num.22181 doi: 10.1002/num.22181
    [5] H. Yücel, M. Stoll, P. Benner, A discontinuous Galerkin method for optimal control problems governed by a system of convection-diffusion PDEs with nonlinear reaction terms, Comput. Math. Appl., 70 (2015), 2414–2431. https://doi.org/10.1016/j.camwa.2015.09.006 doi: 10.1016/j.camwa.2015.09.006
    [6] Z. L. Lu, L. Z. Cao, L. Li, Interpolation coefficients mixed finite element methods for general semilinear Dirichlet boundary elliptic optimal control problems, Appl. Anal., 97 (2018), 2496–2509. https://doi.org/10.1080/00036811.2017.1376319 doi: 10.1080/00036811.2017.1376319
    [7] X. Lin, Y. Chen, Y. Huang, Galerkin spectral methods for an elliptic optimal control problems with L2-norm state constraint, Commun. Math. Sci., 19 (2021), 1247–1267. https://doi.org/10.4310/CMS.2021.v19.n5.a4 doi: 10.4310/CMS.2021.v19.n5.a4
    [8] K. Porwal, P. Shakya, A finite element method for an elliptic optimal control problem with integral state constraints, Appl. Numer. Math., 169 (2021), 273–288. https://doi.org/10.1016/j.apnum.2021.07.002 doi: 10.1016/j.apnum.2021.07.002
    [9] F. Samadi, A. Heydari, S. Effati, A numerical method based on a bilinear pseudo-spectral method to solve the convection-diffusion optimal control problems, Int. J. Comput. Math., 98 (2021), 28–46. https://doi.org/10.1080/00207160.2020.1723563 doi: 10.1080/00207160.2020.1723563
    [10] F. Y. Wang, Z. Q. Zhang, Z. J. Zhou, A spectral Galerkin approximation of optimal control problem governed by fractional advection-diffusion-reaction equations, J. Comput. Appl. Math., 386 (2021), 113233.
    [11] P. G. Casanova, C. Gout, J. Zavaleta, Radial basis function methods for optimal control of the convection-diffusion equation: a numerical study, Eng. Anal. Bound. Elem., 108 (2019), 201–209. https://doi.org/10.1016/j.enganabound.2019.08.008 doi: 10.1016/j.enganabound.2019.08.008
    [12] F. J. Wang, C. M. Fan, C. Z. Zhang, J. Lin, A localized space-time method of fundamental solutions for diffusion and convection-diffusion problems, Adv. Appl. Math. Mech., 12 (2020), 940–958. https://doi.org/10.4208/aamm.OA-2019-0269 doi: 10.4208/aamm.OA-2019-0269
    [13] F. J. Wang, C. Wang, Z. T. Chen, Local knot method for 2D and 3D convection-diffusion-reaction equations in arbitrary domains, Appl. Math. Lett., 105 (2020), 106308. https://doi.org/10.1016/j.aml.2020.106308 doi: 10.1016/j.aml.2020.106308
    [14] F. J. Wang, Q. H. Zhao, Z. T. Chen, C. M. Fan, Localized Chebyshev collocation method for solving elliptic partial differential equations in arbitrary 2D domains, Appl. Math. Comput., 397 (2021), 125903. https://doi.org/10.1016/j.amc.2020.125903 doi: 10.1016/j.amc.2020.125903
    [15] M. Darehmiraki, A. Rezazadeh, A. Ahmadia, S. Salahshour, An interpolation method for the optimal control problem governed by the elliptic convection-diffusion equation, Numer. Meth. Part. D. E., 38 (2022), 137–159. https://doi.org/10.3917/enje.038.0137 doi: 10.3917/enje.038.0137
    [16] M. Darehmiraki, A. Rezazadeh, A new solution for optimal control of fractional convection-reaction-diffusion equation using rational barycentric interpolation, B. Iran. Math. Soc., 46 (2020), 1307–1340. https://doi.org/10.1007/s41980-019-00327-y doi: 10.1007/s41980-019-00327-y
    [17] R. Jiwari, Barycentric rational interpolation and local radial basis functions based numerical algorithms for multidimensional sine-Gordon equation, Numer Methods Partial Differ Equ, 37 (2021), 1965–1992. https://doi.org/10.1002/num.22636 doi: 10.1002/num.22636
    [18] H. Y. Liu, J. Huang, W. Zhang, Y. Y. Ma, Meshfree approach for solving multi-dimensional systems of Fredholm integral equations via barycentric Lagrange interpolation, Appl. Math. Comput., 346 (2019), 295–304. https://doi.org/10.1016/j.amc.2018.10.024 doi: 10.1016/j.amc.2018.10.024
    [19] J. Li, X. N. Su, J. Z. Qu, Linear barycentric rational collocation method for solving telegraph equation, Math. Meth. Appl. Sci., 44 (2021), 11720–11737. https://doi.org/10.1002/mma.7548 doi: 10.1002/mma.7548
    [20] O. Oruc, Two meshless methods based on local radial basis function and barycentric rational interpolation for solving 2D viscoelastic wave equation, Comput. Math. Appl., 79 (2020), 3272–3288. https://doi.org/10.1016/j.camwa.2020.01.025 doi: 10.1016/j.camwa.2020.01.025
    [21] Y. F. Deng, Z. F. Weng, Barycentric interpolation collocation method based on Crank-Nicolson scheme for the Allen-Cahn equation, AIMS Math., 6 (2021), 3857–3873. https://doi.org/10.3934/math.2021229 doi: 10.3934/math.2021229
    [22] W. H. Luo, T. Z. Huang, X. M. Gu, Y. Liu, Barycentric rational collocation methods for a class of nonlinear parabolic partial differential equations, Appl. Math. Lett., 68 (2017), 13–19. https://doi.org/10.1016/j.aml.2016.12.011 doi: 10.1016/j.aml.2016.12.011
    [23] S. C. Yi, L. Q. Yao, A steady barycentric Lagrange interpolation method for the 2D higher-order time fractional telegraph equation with nonlocal boundary condition with error analysis, Numer Methods Partial Differ Equ, 35 (2019), 1694–1716. https://doi.org/10.1002/num.22371 doi: 10.1002/num.22371
    [24] M. S. Floater, K. Hormann, Barycentric rational interpolation with no polesand high rates of approximation, Numer. Math., 107 (2007), 315–331. https://doi.org/10.1007/s00211-007-0093-y doi: 10.1007/s00211-007-0093-y
    [25] J. P. Berrut, M. S. Floater, G. Klein, Convergence rates of derivatives of a family of barycentric rational interpolants, Appl. Numer. Math., 61 (2011), 989–1000. https://doi.org/10.1016/j.apnum.2011.05.001 doi: 10.1016/j.apnum.2011.05.001
    [26] G. Klein, J. P. Berrut, Linear rational finite differences from derivatives of barycentric rational interpolants, SIAM J. Numer. Anal., 50 (2012), 643–656. https://doi.org/10.1137/110827156 doi: 10.1137/110827156
    [27] W. H. Luo, C. P. Li, T. Z. Huang, X. M. Gu, G. C. Wu, A high-order accurate numerical scheme for the Caputo derivative with applications to fractional diffusion problems, Numer Funct Anal Optim, 39 (2018), 600–622. https://doi.org/10.1055/a-0790-1999 doi: 10.1055/a-0790-1999
  • This article has been cited by:

    1. Mengli Yao, Zhifeng Weng, A Numerical Method Based on Operator Splitting Collocation Scheme for Nonlinear Schrödinger Equation, 2024, 29, 2297-8747, 6, 10.3390/mca29010006
  • Reader Comments
  • © 2023 the Author(s), licensee AIMS Press. This is an open access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0)
通讯作者: 陈斌, bchen63@163.com
  • 1. 

    沈阳化工大学材料科学与工程学院 沈阳 110142

  1. 本站搜索
  2. 百度学术搜索
  3. 万方数据库搜索
  4. CNKI搜索

Metrics

Article views(1921) PDF downloads(77) Cited by(1)

Figures and Tables

Figures(6)  /  Tables(6)

Other Articles By Authors

/

DownLoad:  Full-Size Img  PowerPoint
Return
Return

Catalog