# 3D Nonparametric Neural Identification

Authors: Rita Q. Fuentes; Isaac Chairez; Alexander Poznyak; Tatyana Poznyak
Journal: Journal of Control Science and Engineering
 (2012)

---

## Abstract

This paper presents the state identification study of 3D partial differential
equations (PDEs) using the differential neural networks (DNNs) approximation.
There are so many physical situations in applied mathematics
and engineering that can be described by PDEs; these models possess the
disadvantage of having many sources of uncertainties around their mathematical
representation. Moreover, to find the exact solutions of those uncertain
PDEs is not a trivial task especially if the PDE is described in two or
more dimensions. Given the continuous nature and the temporal evolution
of these systems, differential neural networks are an attractive option as nonparametric identifiers capable of estimating a 3D distributed model. The
adaptive laws for weights ensure the “practical stability” of the DNN trajectories
to the parabolic three-dimensional (3D) PDE states. To verify the
qualitative behavior of the suggested methodology, here a nonparametric
modeling problem for a distributed parameter plant is analyzed.


---

## Body

## 1. Introduction

### 1.1. 3D Partial Differential Equations

Partial differential equations (PDEs) are of vast importance in applied mathematics, physics, and engineering since so many real physical situations can be modelled by them. The dynamic description of natural phenomenons are usually described by a set of differential equations using mathematical modeling rules [1]. Almost every system described in PDE has already appeared in the one- and two-dimensional situations. Appending a third dimension ascends the dimensional ladder to its ultimate rung, in physical space at least. For instance, linear second-order 3D partial differential equations appear in many problems modeling equilibrium configurations of solid bodies, the three-dimensional wave equation governing vibrations of solids, liquids, gases, and electromagnetic waves, and the three-dimensional heat equation modeling basic spatial diffusion processes. These equations define a state representing rectangular coordinates on ℝ3. There are some basic underlying solution techniques to solve 3D PDEs: separation of variables and Green’s functions or fundamental solutions. Unfortunately, the most powerful of the planar tools, conformal mapping, does not carry over to higher dimensions. In this way, many numerical techniques solving such PDE, for example, the finite difference method (FDM) and the finite element method (FEM), have been developed (see [2, 3]). The principal disadvantage of these methods is that they require the complete mathematical knowledge of the system to define a mesh (domain discretization), where the functions are approximated locally. The construction of a mesh in two or more dimensions is a nontrivial task. Usually, in practice, only low-order approximations are employed resulting in a continuous approximation of the function across the mesh but not in its partial derivatives. The approximation discontinuities of the derivative can adversely affect the stability of the solution. However, all those methods are well defined if the PDE structure is perfectly known. Actually, the most of suitable numerical solutions could be achieved only if the PDE is linear. Nevertheless, there are not so many methods to solve or approximate the PDE solution when its structure (even in a linear case) is uncertain. This paper suggests a different numerical solution for uncertain systems (given by a 3D PDE) based on the Neural Network approach [4].


### 1.2. Application of Neural Networks to Model PDEs

Recent results show that neural networks techniques seem to be very effective to identify a wide class or systems when we have no complete model information, or even when the plant is considered as a gray box. It is well known that radial basis function neural networks (RBFNNs) and MultiLayer Perceptrons (MLPs) are considered as a powerful tool to approximate nonlinear uncertain functions (see [5]): any continuous function defined on a compact set can be approximated to arbitrary accuracy by such class of neural networks [6]. Since the solutions of interesting PDEs are uniformly continuous and the viable sets that arise in common problems are often compact, neural networks seem like ideal candidates for approximating viability problems (see [7, 8]). Neural networks may provide exact approximation for PDE solutions; however, numerical constraints avoid this possible exactness because it is almost impossible to simulate NN structures with infinite number of nodes (see [7, 9, 10]). The Differential Neural Network (DNN) approach avoids many problems related to global extremum search converting the learning process to an adequate feedback design (see [11, 12]). Lyapunov’s stability theory has been used within the neural networks framework (see [4, 11, 13]). The contributions given in this paper regard the development of a nonparametric identifier for 3D uncertain systems described by partial differential equations. The method produces an artificial mathematical model in three dimensions that is able to describe the PDEs dynamics. The required numerical algorithm to solve the non-parametric identifier was also developed.


## 1.1. 3D Partial Differential Equations

Partial differential equations (PDEs) are of vast importance in applied mathematics, physics, and engineering since so many real physical situations can be modelled by them. The dynamic description of natural phenomenons are usually described by a set of differential equations using mathematical modeling rules [1]. Almost every system described in PDE has already appeared in the one- and two-dimensional situations. Appending a third dimension ascends the dimensional ladder to its ultimate rung, in physical space at least. For instance, linear second-order 3D partial differential equations appear in many problems modeling equilibrium configurations of solid bodies, the three-dimensional wave equation governing vibrations of solids, liquids, gases, and electromagnetic waves, and the three-dimensional heat equation modeling basic spatial diffusion processes. These equations define a state representing rectangular coordinates on ℝ3. There are some basic underlying solution techniques to solve 3D PDEs: separation of variables and Green’s functions or fundamental solutions. Unfortunately, the most powerful of the planar tools, conformal mapping, does not carry over to higher dimensions. In this way, many numerical techniques solving such PDE, for example, the finite difference method (FDM) and the finite element method (FEM), have been developed (see [2, 3]). The principal disadvantage of these methods is that they require the complete mathematical knowledge of the system to define a mesh (domain discretization), where the functions are approximated locally. The construction of a mesh in two or more dimensions is a nontrivial task. Usually, in practice, only low-order approximations are employed resulting in a continuous approximation of the function across the mesh but not in its partial derivatives. The approximation discontinuities of the derivative can adversely affect the stability of the solution. However, all those methods are well defined if the PDE structure is perfectly known. Actually, the most of suitable numerical solutions could be achieved only if the PDE is linear. Nevertheless, there are not so many methods to solve or approximate the PDE solution when its structure (even in a linear case) is uncertain. This paper suggests a different numerical solution for uncertain systems (given by a 3D PDE) based on the Neural Network approach [4].


## 1.2. Application of Neural Networks to Model PDEs

Recent results show that neural networks techniques seem to be very effective to identify a wide class or systems when we have no complete model information, or even when the plant is considered as a gray box. It is well known that radial basis function neural networks (RBFNNs) and MultiLayer Perceptrons (MLPs) are considered as a powerful tool to approximate nonlinear uncertain functions (see [5]): any continuous function defined on a compact set can be approximated to arbitrary accuracy by such class of neural networks [6]. Since the solutions of interesting PDEs are uniformly continuous and the viable sets that arise in common problems are often compact, neural networks seem like ideal candidates for approximating viability problems (see [7, 8]). Neural networks may provide exact approximation for PDE solutions; however, numerical constraints avoid this possible exactness because it is almost impossible to simulate NN structures with infinite number of nodes (see [7, 9, 10]). The Differential Neural Network (DNN) approach avoids many problems related to global extremum search converting the learning process to an adequate feedback design (see [11, 12]). Lyapunov’s stability theory has been used within the neural networks framework (see [4, 11, 13]). The contributions given in this paper regard the development of a nonparametric identifier for 3D uncertain systems described by partial differential equations. The method produces an artificial mathematical model in three dimensions that is able to describe the PDEs dynamics. The required numerical algorithm to solve the non-parametric identifier was also developed.


## 2. 3D Finite Differences Approximation

The problem requires the proposal of a non-parametric identifier based on DNN in three dimensions. The problem here may be treated within the PDEs framework. Therefore, this section introduces the DNN approximation characteristics to reconstruct the trajectories profiles for a family of 3D PDEs.Consider the set of uncertain second-order PDEs(1)ut=f(u,ux,uxx,uy,uyy,uz,uzz,uxy,uyx,uxz,uyz)+ξ,
here uh=uh(x,y,z,t), where h represents t, x, y, z, xx, yy, zz, xy, xz, yx, yz and u=u(x,y,z,t)has n components (u(x,y,z,t)∈ℜn) defined in a domain given by [x,y,z)∈[0,1]×[0,1]×[0,1], t≥0, ξ=ξ(x,y,z,t)∈ℜn is a noise in the state dynamics. This PDE has a set of initial and boundary conditions given by(2)ux(0,y,z,t)=0∈Rn,u(x,0,z,t)=u0∈Rn,u(x,y,0,t)=u0∈Rn,u(x,y,z,0)=c∈Rn.
In (1), ut(x,y,z,t) stands for ∂u(x,y,z,t)/∂t.System (1) armed with boundary and initial conditions (2) is driven in a Hilbert space H equipped with an inner product (·,·). Let us consider a vector function g(t)∈H to be a piecewise continuous in t. By L∞(a,b;H) we denote the set of H-valued functions g such that (g(·),u) is Lebesgue measurable for all u∈H and essmaxt∈[a,b]∥g(γ,t)∥<∞, γ∈ℜn. Suppose that the nonlinear function g(γ,t) satisfies the Lipschitz condition∥g(γ,t)-g(η,t)∥≤L∥γ-η∥,(3)∀γ,η∈Br(γ0):={γ∈Rn∣‖γ-γ0‖≤r,∀t∈[t0,t1]},
where L is positive constant and ∥g∥2=(g,g) is used just to ensure that there exists some δ>0 such that the state equation γ̇=g(γ,t) with γ(t0)=γ0 has a unique solution over [t0,t0+δ] (see [14]). The norm defined above stands for the Sobolev space defined in [15] as follows.Definition 1.
LetΩ be an open set in ℜn, and let ν∈Cm(Ω). Define the norm of ν(γ) as
(4)‖ν‖m,p∶=∑0≤|α|≤m(∫Ω|Dαν(γ)|pdγ)1/p
(1≤p<∞,Dαν(γ):=(∂α/∂γα)ν(γ)). This is the Sobolev norm in which the integration is performed in the Lebesgue sense. The completion of the space of function ν(γ)∈Cm(Ω): ∥ν∥m,p<∞ with respect to ∥·∥m,p is the Sobolev space Hm,p(Ω). For p=2, the Sobolev space is a Hilbert space (see [14, 15]).Below we will use the norm (4) for the functions u(·,t) for each fixed t.


### 2.1. Numerical Approximation for Uncertain Functions

Now consider a functionh0(·) in Hm,2(Ω). By [16], h0(·) can be rewritten as(5)h0(γ,θ*)=∑i∑j∑kθijk*Ψijk(γ),θijk=∫-∞+∞h0(γ)Ψijk(γ)dγ,∀i,j,k∈Z,
where {Ψijk(γ)} are functions constituting a basis in Hm,2(Ω). Last expression is referred to as a vector function series expansion of h0(γ,θ*). Based on this series expansion, an NN takes the following mathematical structure:(6)h0(γ,θ)∶=∑i=L1L2∑j=M1M2∑k=N1N2θijkΨijk(γ)=ΘTW(γ)
that can be used to approximate a nonlinear function h0(γ,θ*)∈S with an adequate selection of integers L1, L2, M1, M2, N1, N2∈ℤ+, where(7)Θ=[θL1M1N1,…,θL1M1N2,…θL2M1N1,…θL2M2N2]⊺,W(γ)=[ΨL1M1N1,…,ΨL1M1N2,…ΨL2M1N1,…ΨL2M2N2]⊺.
Following the Stone-Weierstrass Theorem [17], if(8)ϵ(L1,L2,M1,M2,N1,N2)=h0(γ,θ*)-h0(γ,θ)
is the NN approximation error, then for any arbitrary positive constant ε there are some constants L1,L2,M1,M2,N1,N2∈ℤ such that for all x∈X⊂ℜ(9)‖ϵ(L1,L2,M1,M2,N1,N2)‖2≤ε.
The main idea behind the application of DNN [11] to approximate the 3D PDEs solution is to use a class of finite-difference methods but for uncertain nonlinear functions. So, it is necessary to construct an interior set (commonly called grid or mesh) that divides the subdomain x∈[0,1] in L equidistant sections, y∈[0,1] in M, and z∈[0,1] in N equidistant sections, each one of them (Figure 1) defined as (xi,yj,zk) in such a way that x0=y0=z0=0 and xL=yM=zN=1.Figure 1
Constructed grid for 3 dimensions.Using the mesh description, one can use the next definitions:(10)ui,j,k(t)∶=u(xi,yj,zk,t),uti,j,k(t)∶=∂u(x,y,z,t)∂t|x=xi,y=yj,z=zk,uxi,j,k(t)∶=ux(x,y,z,t)|x=xi,y=yj,z=zk,uxxi,j,k(t)∶=uxx(x,y,z,t)|x=xi,y=yj,z=zk.Analogously, we may consider the other cases (uxx, uy, uyy, uz, uzz, uxy,uxz,uyz). Using the mesh description and applying the finite-difference representation, one gets(11)uxi,j,k(t)≃ui,j,k(t)-ui-1,j,k(t)Δx,uxxi,j,k(t)≃uxi,j,k(t)-uxi-1,j,k(t)Δ2x,
and it follows for all cases such that the (Δx,Δy,Δz)-approximation of the nonlinear PDE (1) can be represented as(12)u̇i,j,k(t)=uti,j,k(t)≃Φ(ui,j,k,ui-1,j,k,ui-2,j,k,ui,j-1,k,ui,j-2,k,ui,j,k-1,ui,j,k-2,ui-1,j-1,k,ui-1,j,k-1,ui,j-1,k-1,ui-1,j-1,k-1)(i=0,…,L;j=0,…,M,k=0,…,N).


### 2.2. 3D Approximation for Uncertain PDE

By simple adding and subtracting the corresponding terms, one can describe (1) as(13)ut=Au+V°1σu+V°2φ1ux+V°3γ1uxxV°4φ2uy+V°5γ2uyy+V°6φ3uz+V°7γ3uzz+V°8ψ1uxy+V°9ψ2uyzV°10ψ3uxz+V°11σ2uxyz+f̃,
where ut=ut(x,y,z,t), u=u(x,y,z,t), σ=σ(x,y,z), ux=ux(x,y,z,t), uxx=uxx(x,y,z,t), the same for σ2, φi, γi, ψi, uy, uyy, uz, uzz, uxy, uyz, uxz, and uxyz (i=1,3¯), A∈ℜn×n, V°1∈ℜn×s1, V°2∈ℜn×s2, V°3∈ℜn×s3, f̃=f̃(x,y,z,t), V°4∈ℜn×s4, V°5∈ℜn×s5, V°6∈ℜn×s6, V°7∈ℜn×s7, V°8∈ℜn×s8, V°9∈ℜn×s9, V°10∈ℜn×s10.Heref̃(x,t)∈ℜn represents the modelling error term, A and V°k (k=1,6¯) any constant matrices and the set of sigmoidal functions have the corresponding size (σ(x,y,z)∈ℜs1, φ1(x,y,z)∈ℜs2, γ1(x,y,z)∈ℜs3, φ2(x,y,z)∈ℜs4, γ2(x,y,z)∈ℜs5, ψ1(x,y,z)∈ℜs6, φ3(x,y,z)∈ℜs7, γ3(x,y,z)∈ℜs8, ψ2(x,y,z)∈ℜs9, ψ3(x,y,z)∈ℜs10, and σ2(x,y,z)∈ℜs11) and are known as the neural network set of activation functions. These functions obey the following sector conditions:(14)‖σ1(x,y,z)-σ1(x′,y′,z′)‖2≤Lσ1(‖x-x′‖2+‖y-y′‖2+‖z-z′‖2),‖σ2(x,y,z)-σ2(x′,y′,z′)‖2≤Lσ2(‖x-x′‖2+‖y-y′‖2+‖z-z′‖2),‖ϕs(x,y,z)-ϕs(x′,y′,z′)‖2≤Lϕs(‖x-x′‖2+‖y-y′‖2+‖z-z′‖2)
which are bounded in x, y, and z, that is,(15)‖σ(⋅)‖2≤σ(l-1)+,‖φl(⋅)‖2≤φl+,‖γl(⋅)‖2≤γl+,‖ψl(⋅)‖≤ψl+,l=1,3.
Following the methodology of DNN [11] and applying the same representation to (12), we get for each i∈(1,…,L), j∈(1,…,M), k∈(1,…,N) the following robust adaptive non-parametric identifier:(16)uti,j,k(t)=Ai,j,kui,j,k(t)+∑p=111W°pi,j,kϕsU(t)+f̃i,j,k(t).
In the sum we have that s=1,3¯, ϕ represents functions σ, φ, γ, ψ, and U can be taken as the corresponding ui,j,k, ui-1,j,k, ui-2,j,k, ui,j-1,k, ui,j-2,k, ui,j,k-1, ui,j,k-2, ui-1,j-1,k, ui,j-1,k-1, ui-1,j,k-1, ui-1,j-1,k-1. In this equation the term f̃i,j,k(t), which is usually recognized as the modeling error, satisfies the following identify, and here, it has been omitted the dependence to xi,yj,zk of each sigmoidal function:(17)f̃i,j,k(t)∶=Φ(ui,j,k,ui-1,j,k,ui-2,j,k,ui,j-1,k,ui,j-2,k,ui,j,k-1,ui,j,k-2,ui-1,j-1,k,ui-1,j,k-1,ui,j-1,k-1,ui-1,j-1,k-1)-Ai,j,kui,j,k(t)-∑p=111W°pi,j,kϕsU(t),
where W°pi∈ℜn×sp,p=1,11¯,ϕs, ϕ represents functions σ,φ,γ,ψ and s=1,3¯,U(t) represents the corresponding (ui,j,k, ui-1,j,k, ui-2,j,k, ui,j-1,k, ui,j-2,k, ui,j,k-1, ui,j,k-2, ui-1,j-1,k, ui-1,j,k-1, ui,j-1,k-1, ui-1,j-1,k-1).We will assume that the modeling error terms satisfy the following.Assumption 2.
The modeling error is absolutely bounded inΩ:
(18)‖f̃i,j,k‖2≤f1i,j,k.Assumption 3.
The error modeling gradient(19)∇sf̃(x,y,z,t)|s=si:=∇sf̃i,j,k
is bounded as∥∇sf̃i,j,k∥2≤fri,j,k, wheres=x,y,z andfri,j,k(r=1,3¯) are constants.


## 2.1. Numerical Approximation for Uncertain Functions

Now consider a functionh0(·) in Hm,2(Ω). By [16], h0(·) can be rewritten as(5)h0(γ,θ*)=∑i∑j∑kθijk*Ψijk(γ),θijk=∫-∞+∞h0(γ)Ψijk(γ)dγ,∀i,j,k∈Z,
where {Ψijk(γ)} are functions constituting a basis in Hm,2(Ω). Last expression is referred to as a vector function series expansion of h0(γ,θ*). Based on this series expansion, an NN takes the following mathematical structure:(6)h0(γ,θ)∶=∑i=L1L2∑j=M1M2∑k=N1N2θijkΨijk(γ)=ΘTW(γ)
that can be used to approximate a nonlinear function h0(γ,θ*)∈S with an adequate selection of integers L1, L2, M1, M2, N1, N2∈ℤ+, where(7)Θ=[θL1M1N1,…,θL1M1N2,…θL2M1N1,…θL2M2N2]⊺,W(γ)=[ΨL1M1N1,…,ΨL1M1N2,…ΨL2M1N1,…ΨL2M2N2]⊺.
Following the Stone-Weierstrass Theorem [17], if(8)ϵ(L1,L2,M1,M2,N1,N2)=h0(γ,θ*)-h0(γ,θ)
is the NN approximation error, then for any arbitrary positive constant ε there are some constants L1,L2,M1,M2,N1,N2∈ℤ such that for all x∈X⊂ℜ(9)‖ϵ(L1,L2,M1,M2,N1,N2)‖2≤ε.
The main idea behind the application of DNN [11] to approximate the 3D PDEs solution is to use a class of finite-difference methods but for uncertain nonlinear functions. So, it is necessary to construct an interior set (commonly called grid or mesh) that divides the subdomain x∈[0,1] in L equidistant sections, y∈[0,1] in M, and z∈[0,1] in N equidistant sections, each one of them (Figure 1) defined as (xi,yj,zk) in such a way that x0=y0=z0=0 and xL=yM=zN=1.Figure 1
Constructed grid for 3 dimensions.Using the mesh description, one can use the next definitions:(10)ui,j,k(t)∶=u(xi,yj,zk,t),uti,j,k(t)∶=∂u(x,y,z,t)∂t|x=xi,y=yj,z=zk,uxi,j,k(t)∶=ux(x,y,z,t)|x=xi,y=yj,z=zk,uxxi,j,k(t)∶=uxx(x,y,z,t)|x=xi,y=yj,z=zk.Analogously, we may consider the other cases (uxx, uy, uyy, uz, uzz, uxy,uxz,uyz). Using the mesh description and applying the finite-difference representation, one gets(11)uxi,j,k(t)≃ui,j,k(t)-ui-1,j,k(t)Δx,uxxi,j,k(t)≃uxi,j,k(t)-uxi-1,j,k(t)Δ2x,
and it follows for all cases such that the (Δx,Δy,Δz)-approximation of the nonlinear PDE (1) can be represented as(12)u̇i,j,k(t)=uti,j,k(t)≃Φ(ui,j,k,ui-1,j,k,ui-2,j,k,ui,j-1,k,ui,j-2,k,ui,j,k-1,ui,j,k-2,ui-1,j-1,k,ui-1,j,k-1,ui,j-1,k-1,ui-1,j-1,k-1)(i=0,…,L;j=0,…,M,k=0,…,N).


## 2.2. 3D Approximation for Uncertain PDE

By simple adding and subtracting the corresponding terms, one can describe (1) as(13)ut=Au+V°1σu+V°2φ1ux+V°3γ1uxxV°4φ2uy+V°5γ2uyy+V°6φ3uz+V°7γ3uzz+V°8ψ1uxy+V°9ψ2uyzV°10ψ3uxz+V°11σ2uxyz+f̃,
where ut=ut(x,y,z,t), u=u(x,y,z,t), σ=σ(x,y,z), ux=ux(x,y,z,t), uxx=uxx(x,y,z,t), the same for σ2, φi, γi, ψi, uy, uyy, uz, uzz, uxy, uyz, uxz, and uxyz (i=1,3¯), A∈ℜn×n, V°1∈ℜn×s1, V°2∈ℜn×s2, V°3∈ℜn×s3, f̃=f̃(x,y,z,t), V°4∈ℜn×s4, V°5∈ℜn×s5, V°6∈ℜn×s6, V°7∈ℜn×s7, V°8∈ℜn×s8, V°9∈ℜn×s9, V°10∈ℜn×s10.Heref̃(x,t)∈ℜn represents the modelling error term, A and V°k (k=1,6¯) any constant matrices and the set of sigmoidal functions have the corresponding size (σ(x,y,z)∈ℜs1, φ1(x,y,z)∈ℜs2, γ1(x,y,z)∈ℜs3, φ2(x,y,z)∈ℜs4, γ2(x,y,z)∈ℜs5, ψ1(x,y,z)∈ℜs6, φ3(x,y,z)∈ℜs7, γ3(x,y,z)∈ℜs8, ψ2(x,y,z)∈ℜs9, ψ3(x,y,z)∈ℜs10, and σ2(x,y,z)∈ℜs11) and are known as the neural network set of activation functions. These functions obey the following sector conditions:(14)‖σ1(x,y,z)-σ1(x′,y′,z′)‖2≤Lσ1(‖x-x′‖2+‖y-y′‖2+‖z-z′‖2),‖σ2(x,y,z)-σ2(x′,y′,z′)‖2≤Lσ2(‖x-x′‖2+‖y-y′‖2+‖z-z′‖2),‖ϕs(x,y,z)-ϕs(x′,y′,z′)‖2≤Lϕs(‖x-x′‖2+‖y-y′‖2+‖z-z′‖2)
which are bounded in x, y, and z, that is,(15)‖σ(⋅)‖2≤σ(l-1)+,‖φl(⋅)‖2≤φl+,‖γl(⋅)‖2≤γl+,‖ψl(⋅)‖≤ψl+,l=1,3.
Following the methodology of DNN [11] and applying the same representation to (12), we get for each i∈(1,…,L), j∈(1,…,M), k∈(1,…,N) the following robust adaptive non-parametric identifier:(16)uti,j,k(t)=Ai,j,kui,j,k(t)+∑p=111W°pi,j,kϕsU(t)+f̃i,j,k(t).
In the sum we have that s=1,3¯, ϕ represents functions σ, φ, γ, ψ, and U can be taken as the corresponding ui,j,k, ui-1,j,k, ui-2,j,k, ui,j-1,k, ui,j-2,k, ui,j,k-1, ui,j,k-2, ui-1,j-1,k, ui,j-1,k-1, ui-1,j,k-1, ui-1,j-1,k-1. In this equation the term f̃i,j,k(t), which is usually recognized as the modeling error, satisfies the following identify, and here, it has been omitted the dependence to xi,yj,zk of each sigmoidal function:(17)f̃i,j,k(t)∶=Φ(ui,j,k,ui-1,j,k,ui-2,j,k,ui,j-1,k,ui,j-2,k,ui,j,k-1,ui,j,k-2,ui-1,j-1,k,ui-1,j,k-1,ui,j-1,k-1,ui-1,j-1,k-1)-Ai,j,kui,j,k(t)-∑p=111W°pi,j,kϕsU(t),
where W°pi∈ℜn×sp,p=1,11¯,ϕs, ϕ represents functions σ,φ,γ,ψ and s=1,3¯,U(t) represents the corresponding (ui,j,k, ui-1,j,k, ui-2,j,k, ui,j-1,k, ui,j-2,k, ui,j,k-1, ui,j,k-2, ui-1,j-1,k, ui-1,j,k-1, ui,j-1,k-1, ui-1,j-1,k-1).We will assume that the modeling error terms satisfy the following.Assumption 2.
The modeling error is absolutely bounded inΩ:
(18)‖f̃i,j,k‖2≤f1i,j,k.Assumption 3.
The error modeling gradient(19)∇sf̃(x,y,z,t)|s=si:=∇sf̃i,j,k
is bounded as∥∇sf̃i,j,k∥2≤fri,j,k, wheres=x,y,z andfri,j,k(r=1,3¯) are constants.


## 3. DNN Identification for Distributed Parameters Systems

### 3.1. DNN Identifier Structure

Based on the DNN methodology [11], consider the DNN identifier(20)ddtûi,j,k=Ai,j,kûi,j,k+∑p=111Wpi,j,k(t)ϕsÛ
for all i=0,…,L; û-1(t)=û-2(t)=0, where ϕ represents activation functions σ, φ, γ, and ψ, s=1,3¯, Û is each one of the states ûi,j,k, ûi-1,j,k, ûi-2,j,k, ûi,j-1,k, ûi,j-2,k, ûi,j,k-1, ûi,j,k-1, ûi,j,k-2, ûi-1,j-1,k, ûi-1,j,k-1, ûi-1,j-1,k-1, and Ai,j,k∈ℜn×n is a constant matrix to be selected, ûi,j,k(t) is the estimate of ui,j,k(t). Obviously that proposed methodology implies the designing of individual DNN identifier for each point xi, yj, zk. The collection of such identifiers will constitute a DNN net containing N×M connected DNN identifiers working in parallel. Here σ1(xi,yj,zk), φ1(xi,yj,zk), φ2(xi,yj,zk), φ3(xi,yj,zk), γ1(xi,yj,zk), γ2(xi,yj,zk), γ3(xi,yj,zk), ψ1(xi,yj,zk), ψ2(xi,yj,zk), ψ3(xi,yj,zk), and σ2(xi,yj,zk) are the NN activation vectors. This means that the applied DNN-approximation significantly simplifies the specification of σ1(·,·), φ1(·,·), φ2(·,·), φ3(·,·), γ1(·,·), γ2(·,·), γ3(·,·) and ψ1(·,·), ψ2(·,·), ψ3(·,·), σ2(·,·) which now are constant for any xi, yj, zk fixed.


### 3.2. Learning Laws for Identifier’s Weights

For eachi=0,…,L,j=0,…,M,k=0,…,N, define the vector-functions defining the error between the trajectories produced by the model and the DNN-identifier as well as their derivatives with respect to x, y, and z for each i,j,k:(21)ũi,j,k(t)∶=ûi,j,k(t)-ui,j,k(t),ũsi,j,k(t)∶=ûsi,j,k(t)-usi,j,k(t),s=x,y,z.
Let Wri,j,k(t)∈ℜn, r=1,11¯ be time-variant matrices. These matrices satisfy the following nonlinear matrix differential equations: (22)Ẇri,j,k(t)∶=ddtW̅ri,j,k(t)=-αW̃ri,j,k(t)-Kr-1Pi,j,kûi,j,k(Û(i,j,k),r)⊺(Ωr(xi,yj,zj))⊺-Kr-1S1i,j,kũxi,j,k(Û(i,j,k),r)⊺(Ωxr(xi,yj,zj))⊺-Kr-1S2i,j,kũyi,j,k(Û(i,j,k),r)⊺(Ωyr(xi,yj,zj))⊺-Kr-1S3i,j,kũzi,j,k(Û(i,j,k),r)⊺(Ωzr(xi,yj,zj))⊺,
where Ωh(xi,yj,zj)=Sl(xi,yj,zj) (h=1,10¯, l=1,3¯) represents the corresponding sigmoidal functions σl(xi,yj,zj), φl(xi,yj,zj), γl(xi,yj,zj), and ψl(xi,yj,zj). Here(23)Û(i,j,k),1(t)=ûi,j,k(t),Ũ(i,j,k),2(t)=ûi-1,j,k(t),Ũ(i,j,k),3(t)=ûi-2,j,k(t),Ũ(i,j,k),4(t)=ûi,j-1,k(t),Ũ(i,j,k),5(t)=ûi,j-2,k(t),Ũ(i,j,k),6(t)=ûi-1,j-1,k(t),Ũ(i,j,k),7(t)=ûi-1,j-1,k(t),Ũ(i,j,k),8(t)=ûi-1,j-1,k(t),Ũ(i,j,k),9(t)=ûi-1,j-1,k(t),Ũ(i,j,k),10(t)=ûi-1,j-1,k(t),Ũ(i,j,k),11(t)=ûi-1,j-1,k-1(t)
with positive matrices Kr(r=1,11¯) and Pi,j,k, S1i,j,k, and S2i,j,k (i=0,N¯;j=0,M¯) which are positive definite solutions (Pi,j>0, S1i,j>0, S2i,j>0) and S3i,j,k (i=0,N¯;j=0,M¯, k=0,L¯) of the algebraic Riccati equations defined as follows:(24)Ric(Pi,j,k)∶=Pi,j,kAi,j,k+[Ai,j,k]⊺Pi,j,k+Pi,j,kRPi,j,kPi,j,k+QPi,j,k=0,Ric(S1i,j,k)∶=S1i,j,kAi,j,k+[Ai,j]⊺S1i,j,k+S1i,j,kRS1i,j,kS1i,j,k+QS1i,j,k=0,Ric(S2i,j,k)∶=S2i,j,kAi,j,k+[Ai,j,k]⊺S2i,j,k+S2i,j,kRS2i,j,kS2i,j,k+QS2i,j,k=0,Ric(S3i,j,k)∶=S3i,j,kAi,j,k+[Ai,j,k]⊺S3i,j,k+S3i,j,kRS3i,j,kS3i,j,k+QS3i,j,k=0,
where each RBi,j,k has the form(25)RBi,j,k∶=∑r=111W°ri,j,kΛr+b(W°ri,j,k)⊺+Λb,
where B can be P, S1, S2, and S3 and b=(7,14,21,28).MatricesQBi,j,k have the form(26)QBi,j,k∶=‖Ωm1(xi,yj,zk)‖Λ1-12+‖Ωm2(xi-1,yj,zk)‖Λ2-12+‖Ωm3(xi-2,yj,zk)‖Λ3-12+‖Ωm4(xi,yj-1,zk)‖Λ4-12+‖Ωm5(xi,yj-2,zk)‖Λ5-12+‖Ωm6(xi,yj,zk-1)‖Λ6-12+‖Ωm7(xi,yj,zk-2)‖Λ7-12+‖Ωm8(xi-1,yj-1,zk)‖Λ8-12+‖Ωm9(xi,yj-1,zk-1)‖Λ9-12+‖Ωm10(xi-1,yj,zk-1)‖Λ10-12+‖Ωm11(xi-1,yj-1,zk-1)‖Λ11-12+QBi,j,k,
where B can be P, S1, S2, or S3, representing the partial derivative; for S1 it is with respect to x, for S2 with respect to y, and for S3 with respect to z, and Λl-1 (l=1,46¯), where(27)Ωxr(xi,yj,zj)∶=ddxΩr(x,y,z)|x=xi,y=yj,z=zk,Ωyr(xi,yj,zj)∶=ddyΩr(x,y,z)|x=xi,y=yj,z=zk,Ωzr(xi,yj,zj)∶=ddzΩr(x,y,z)|x=xi,y=yj,z=zk.
Here, W̃ki,j,k(t):=Wki,j,k(t)-W°ki,j,k, k=1,11¯.Thespecial class of Riccati equation(28)PA+ATP+PRP+Q=0
has a unique positive solution P if and only if the four conditions given in [11] (page 65, chapter 2 Nonlinear System Identification: Differential Learning) are fulfilled (see [11]): (1) matrix A is stable, (2) pair (A,R1/2) is controllable, (3) pair (Q1/2,A) is observable, and (4) matrices (A,Q,R) should be selected in such a way to satisfy the following inequality:(29)14(A⊺R-1-R-1A)R(A⊺R-1-R-1A)⊺+Q≤A⊺R-1A
which restricts the largest eigenvalue of R guarantying the existence of a unique positive solution. The main result obtained in this part is in the practical stability framework.


## 3.1. DNN Identifier Structure

Based on the DNN methodology [11], consider the DNN identifier(20)ddtûi,j,k=Ai,j,kûi,j,k+∑p=111Wpi,j,k(t)ϕsÛ
for all i=0,…,L; û-1(t)=û-2(t)=0, where ϕ represents activation functions σ, φ, γ, and ψ, s=1,3¯, Û is each one of the states ûi,j,k, ûi-1,j,k, ûi-2,j,k, ûi,j-1,k, ûi,j-2,k, ûi,j,k-1, ûi,j,k-1, ûi,j,k-2, ûi-1,j-1,k, ûi-1,j,k-1, ûi-1,j-1,k-1, and Ai,j,k∈ℜn×n is a constant matrix to be selected, ûi,j,k(t) is the estimate of ui,j,k(t). Obviously that proposed methodology implies the designing of individual DNN identifier for each point xi, yj, zk. The collection of such identifiers will constitute a DNN net containing N×M connected DNN identifiers working in parallel. Here σ1(xi,yj,zk), φ1(xi,yj,zk), φ2(xi,yj,zk), φ3(xi,yj,zk), γ1(xi,yj,zk), γ2(xi,yj,zk), γ3(xi,yj,zk), ψ1(xi,yj,zk), ψ2(xi,yj,zk), ψ3(xi,yj,zk), and σ2(xi,yj,zk) are the NN activation vectors. This means that the applied DNN-approximation significantly simplifies the specification of σ1(·,·), φ1(·,·), φ2(·,·), φ3(·,·), γ1(·,·), γ2(·,·), γ3(·,·) and ψ1(·,·), ψ2(·,·), ψ3(·,·), σ2(·,·) which now are constant for any xi, yj, zk fixed.


## 3.2. Learning Laws for Identifier’s Weights

For eachi=0,…,L,j=0,…,M,k=0,…,N, define the vector-functions defining the error between the trajectories produced by the model and the DNN-identifier as well as their derivatives with respect to x, y, and z for each i,j,k:(21)ũi,j,k(t)∶=ûi,j,k(t)-ui,j,k(t),ũsi,j,k(t)∶=ûsi,j,k(t)-usi,j,k(t),s=x,y,z.
Let Wri,j,k(t)∈ℜn, r=1,11¯ be time-variant matrices. These matrices satisfy the following nonlinear matrix differential equations: (22)Ẇri,j,k(t)∶=ddtW̅ri,j,k(t)=-αW̃ri,j,k(t)-Kr-1Pi,j,kûi,j,k(Û(i,j,k),r)⊺(Ωr(xi,yj,zj))⊺-Kr-1S1i,j,kũxi,j,k(Û(i,j,k),r)⊺(Ωxr(xi,yj,zj))⊺-Kr-1S2i,j,kũyi,j,k(Û(i,j,k),r)⊺(Ωyr(xi,yj,zj))⊺-Kr-1S3i,j,kũzi,j,k(Û(i,j,k),r)⊺(Ωzr(xi,yj,zj))⊺,
where Ωh(xi,yj,zj)=Sl(xi,yj,zj) (h=1,10¯, l=1,3¯) represents the corresponding sigmoidal functions σl(xi,yj,zj), φl(xi,yj,zj), γl(xi,yj,zj), and ψl(xi,yj,zj). Here(23)Û(i,j,k),1(t)=ûi,j,k(t),Ũ(i,j,k),2(t)=ûi-1,j,k(t),Ũ(i,j,k),3(t)=ûi-2,j,k(t),Ũ(i,j,k),4(t)=ûi,j-1,k(t),Ũ(i,j,k),5(t)=ûi,j-2,k(t),Ũ(i,j,k),6(t)=ûi-1,j-1,k(t),Ũ(i,j,k),7(t)=ûi-1,j-1,k(t),Ũ(i,j,k),8(t)=ûi-1,j-1,k(t),Ũ(i,j,k),9(t)=ûi-1,j-1,k(t),Ũ(i,j,k),10(t)=ûi-1,j-1,k(t),Ũ(i,j,k),11(t)=ûi-1,j-1,k-1(t)
with positive matrices Kr(r=1,11¯) and Pi,j,k, S1i,j,k, and S2i,j,k (i=0,N¯;j=0,M¯) which are positive definite solutions (Pi,j>0, S1i,j>0, S2i,j>0) and S3i,j,k (i=0,N¯;j=0,M¯, k=0,L¯) of the algebraic Riccati equations defined as follows:(24)Ric(Pi,j,k)∶=Pi,j,kAi,j,k+[Ai,j,k]⊺Pi,j,k+Pi,j,kRPi,j,kPi,j,k+QPi,j,k=0,Ric(S1i,j,k)∶=S1i,j,kAi,j,k+[Ai,j]⊺S1i,j,k+S1i,j,kRS1i,j,kS1i,j,k+QS1i,j,k=0,Ric(S2i,j,k)∶=S2i,j,kAi,j,k+[Ai,j,k]⊺S2i,j,k+S2i,j,kRS2i,j,kS2i,j,k+QS2i,j,k=0,Ric(S3i,j,k)∶=S3i,j,kAi,j,k+[Ai,j,k]⊺S3i,j,k+S3i,j,kRS3i,j,kS3i,j,k+QS3i,j,k=0,
where each RBi,j,k has the form(25)RBi,j,k∶=∑r=111W°ri,j,kΛr+b(W°ri,j,k)⊺+Λb,
where B can be P, S1, S2, and S3 and b=(7,14,21,28).MatricesQBi,j,k have the form(26)QBi,j,k∶=‖Ωm1(xi,yj,zk)‖Λ1-12+‖Ωm2(xi-1,yj,zk)‖Λ2-12+‖Ωm3(xi-2,yj,zk)‖Λ3-12+‖Ωm4(xi,yj-1,zk)‖Λ4-12+‖Ωm5(xi,yj-2,zk)‖Λ5-12+‖Ωm6(xi,yj,zk-1)‖Λ6-12+‖Ωm7(xi,yj,zk-2)‖Λ7-12+‖Ωm8(xi-1,yj-1,zk)‖Λ8-12+‖Ωm9(xi,yj-1,zk-1)‖Λ9-12+‖Ωm10(xi-1,yj,zk-1)‖Λ10-12+‖Ωm11(xi-1,yj-1,zk-1)‖Λ11-12+QBi,j,k,
where B can be P, S1, S2, or S3, representing the partial derivative; for S1 it is with respect to x, for S2 with respect to y, and for S3 with respect to z, and Λl-1 (l=1,46¯), where(27)Ωxr(xi,yj,zj)∶=ddxΩr(x,y,z)|x=xi,y=yj,z=zk,Ωyr(xi,yj,zj)∶=ddyΩr(x,y,z)|x=xi,y=yj,z=zk,Ωzr(xi,yj,zj)∶=ddzΩr(x,y,z)|x=xi,y=yj,z=zk.
Here, W̃ki,j,k(t):=Wki,j,k(t)-W°ki,j,k, k=1,11¯.Thespecial class of Riccati equation(28)PA+ATP+PRP+Q=0
has a unique positive solution P if and only if the four conditions given in [11] (page 65, chapter 2 Nonlinear System Identification: Differential Learning) are fulfilled (see [11]): (1) matrix A is stable, (2) pair (A,R1/2) is controllable, (3) pair (Q1/2,A) is observable, and (4) matrices (A,Q,R) should be selected in such a way to satisfy the following inequality:(29)14(A⊺R-1-R-1A)R(A⊺R-1-R-1A)⊺+Q≤A⊺R-1A
which restricts the largest eigenvalue of R guarantying the existence of a unique positive solution. The main result obtained in this part is in the practical stability framework.


## 4. Practical Stability and Stabilization

The following definition and proposition are needed for the main results of the paper. Consider the following ODE nonlinear system:(30)żt=g(zt,vt)+ϖt
with zt∈ℜn, and vt∈ℜm, and ϖt an external perturbation or uncertainty such that ∥ϖt∥2≤ϖ+.Definition 4 (Practical Stability).
Assume that a time intervalT and a fixed function vt*∈ℜm over T are given. Given ε>0, the nonlinear system (30) is said to be ε-practically stable over T under the presence of ϖt if there exists a δ>0 (δ depends on ε and the interval T) such that zt∈B[0,ε], for all 0≤t, whenever zt0∈B[0,δ].Similarly to the Lyapunov stability theory for nonlinear systems, it was applied the aforementioned direct method for theε-practical stability of nonlinear systems using-practical Lyapunov-like functions under the presence of external perturbations and model uncertainties. Note that these functions have properties differing significantly from the usual Lyapunov functions in classic stability theory.The subsequent proposition requires the following lemma.Lemma 5.
Let a nonnegative functionV(t) satisfying the following differential inequality:
(31)V̇(t)≤-αV(t)+β,
where α>0 and β≥0. Then
(32)[1-μV(t)]+⟶0,
with μ=β/α and the function [·]+ defined as
(33)[z]+:={zifz≥0,0ifz<0.Proof.
The proof of this lemma can be obtained directly by the application of the Gronwall-Bellman Lemma.Proposition 6.
Given a time intervalT and a function v(·) over a continuously differentiable real-valued function V(z,t) satisfying V(0,t)=0, for all tϵT, is said to be ε-practical Lyapunov-like function over T under v if there exists a constant α>0 such that
(34)V̇(z,t)≤-αV(z,t)+H(ϖ+),
with H a bounded nonnegative nonlinear function with upper bound H+. Moreover, the trajectories of zt belong to the zone ε:=H+/α when t→∞. In this proposition V̇(zt,t) denotes the derivative of V(z,t) along zt, that is, V̇(z,t)=Vz(z,t)·(g(zt,vt)+ϖt)+Vt(x,t).Proof.
The proof follows directly from Lemma5.Definition 7.
Given a time intervalT and a function v(·) over T, nonlinear system (30) is ε-practically stable, T under v if there exists an ε-practical Lyapunov-like function V(x,t) over T under v.


## 5. Identification Problem Formulation

The state identification problem for nonlinear systems (13) analyzed in this study, could be now stated as follows.Problem.
For the nonlinear system, given by the vector PDE (20), to study the quality of the DNN identifier supplied with the adjustment (learning) laws (22), estimate the upper bound of the identification error δ given by
(35)δ∶=lim⁡t→∞¯∑i=0L∑j=0M∑k=0N‖ûi,j,k(t)-ui,j,k(t)‖Pi,j,k2
(with Pi,j,k from (24)) and, if it is possible, to reduce to its lowest possible value, selecting free parameters participating into the DNN identifier (A,Kr,r=1,11¯).This implicates that the reduction of the identification errorδ means that the differential neural network has converged to the solution of the 3D PDE; this can be observed in the matching of the DNN to the PDE state.


## 6. Main Result

The main result of this paper is presented in the following theorem.Theorem 8.
Consider the nonlinear model (1) (i=0,…,L; j=0,…,M, k=0,…,N), given by the system of PDEs with uncertainties (perturbations) in the states, under the border conditions (2). Let us also suppose that the DNN-identifier is given by (20) which parameters are adjusted by the learning laws (22) with parameters αmi,j,k(i=0,…,L;j=0,…,M;k=0,…,N). If positive definite matrices Q1i,j,k, Q2i,j,k, and Q3i,j,k provide the existence of positive solutions Pi,j, S1i,j,k, S2i,j,k, and S3i,j,k(i=0,…,L) to the Riccati equations (24), then for all i=0,…,L; j=0,…,M; k=0,…,N the following ρ-upper bound:
(36)lim⁡t→∞¯∑i=0L∑j=0M×∑k=0N[‖ûi,j,k(t)-ui,j,k(t)‖Pi,j,k2+‖ûxi,j,k(t)-uxi,j,k(t)‖S1i,j,k+‖ûyi,j,k(t)-uyi,j,k(t)‖S2i,j,k+‖ûzi,j,k(t)-uzi,j,k(t)‖S3i,j,k]≤μ
is ensured with μ0:=LMNα-1/2 and
(37)μ∶=μ0max⁡i,j,kΨ,Ψ∶=ϖ∑s=13fsi,j,k≥(f̃i,j,k(t))TΛ7-1f̃i,j,k(t)+(f̃xi,j,k(t))TΛ14-1f̃xi,j,k(t)+(f̃yi,j,k(t))TΛ21-1f̃yi,j,k(t)+(f̃zi,j,k(t))TΛ28-1f̃zi,j,k(t),ϖ∶=max⁡{λmax⁡(Λ7-1),λmax⁡(Λ14-1),λmax⁡(Λ21-1),λmax⁡(Λ28-1)}.
Moreover, the weights Wr,tr=1,11¯ remain bounded being proportional to μ, that is, ∥Wr,t∥≤Krμ, r=1,11¯.Proof.
The proof is given in the appendix.


## 7. Simulation Results

### 7.1. Numerical Example

Below, the numerical smulation shows the qualitative illustration for a benchmark system. Therefore, consider the following three-dimensional PDE described as follows:(38)ut(x,y,z,t)=-c1uxx(x,y,z,t)-c2uyy(x,y,z,t)-c3uzz(x,y,z,t)+ξ(x,y,z,t),
where c1=c2=c3=0.15. It is assumed that there is access to discrete measurements of the state u(x,y,z,t) along the whole domain, which is feasible in practice. ξ(x,y,z,t) is a noise in the state dynamics. This model will be used just to generate the data to test the 3D identifier based on DNN. Boundary conditions and initial conditions were selected as follows:(39)u(0,0,0,t)=rand(1),ux(0,0,0,t)=0,uy(0,0,0,t)=0,uz(0,0,0,t)=0.
The trajectories of the model can be seen in Figure 3 as well as the estimated state, produced by the DNN identifier. The efficiency of the identification process provided by the suggested DNN algorithm shown in Figure 4.In Figures5 and 6 there are shown the trajectories of the PDE and the DNN for the coordinate z and the index error of the PDE and the DNN at 10 seconds for coordinates x, z, respectively (Figure 2).Figure 2
Numerical trajectory produced by the mathematical model described by 3D partial differential equation at time 10 s along the whole domain.Figure 3
Numerical trajectory produced by the neural network described reproducing the model at time 10 s along the whole domain.Figure 4
Error of the approximation. Difference between 3D PDE numerical solution and 3D DNN approximation at time 10 s along the whole domain.Figure 5
Comparison between 3D PDE trajectories and 3D DNN approximation for coordinatez at time 10 s.Figure 6
Index quadratic error for the whole time with coordinatesx, y, and z fixed.


### 7.2. Tumour Growth Example

The mathematical model of the brain tumour growth is presented in this section based on the results of [18]. Here the diffusion coefficient is considered as a constant. Let consider the following three-dimensional parabolic equation of the tumour growth described as follows:(40)ut(x,y,z,t)=-Pux(x,y,z,t)-Ruy(x,y,z,t)-Suz(x,y,z,t)+Quxx(x,y,z,t)+Quyy(x,y,z,t)+Quzz(x,y,z,t)+Γ-Lu(x,y,z,t),
where u(x,y,z,t) is the growth rate of a brain tumour, Q is the diffusion coefficient, W=(P,R,S) is the drift velocity field, Γ=Γ(u) is the proliferation coefficient, and L=L(u) is the decay coefficient of cells. It is assumed that there is access to discrete measurements of the state u(x,y,z,t) along the whole domain, which is feasible in practice by PET-CT (Positron emission tomography-computed tomography) technology. The domain 0≤x≤1, 0≤y≤1, and 0≤z≤1 This model will be used just to generate the data to test the 3D identifier based on DNN. Boundary conditions and initial conditions were selected as follows:(41)u(0,0,0,t)=200±20μ,ux(0,0,0,t)=0,uy(0,0,0,t)=0,uz(0,0,0,t)=0.
The trajectories of the model and the estimate state produced by the DNN identifier can be seen in Figure 7. The dissimilarity between both trajectories depends on the learning period required for adjusting the DNN identifier. The error between trajectories produced by the model and the proposed identifier is close to zero almost for all x, y, z and all t that shows the efficiency of the identification process provided by the suggested DNN algorithm is shown in Figure 8.Numerical trajectory produced by the mathematical model described by 3D partial differential equation at time 46 days along the whole domain.
(a)(b)Figure 8
Error of the approximation. Difference between 3D PDE numerical solution and 3D DNN approximation at time 46 days along the whole domain.In Figures9 and 10 there are shown the trajectories of the PDE and the DNN for the coordinate z and the index error of the PDE and the DNN at 46 days for coordinates x, z, respectively.Figure 9
Comparison between 3D PDE trajectories and 3D DNN approximation for coordinatez at time 46 days.Figure 10
Index quadratic error for the whole time with coordinatesx, y, and z fixed.


## 7.1. Numerical Example

Below, the numerical smulation shows the qualitative illustration for a benchmark system. Therefore, consider the following three-dimensional PDE described as follows:(38)ut(x,y,z,t)=-c1uxx(x,y,z,t)-c2uyy(x,y,z,t)-c3uzz(x,y,z,t)+ξ(x,y,z,t),
where c1=c2=c3=0.15. It is assumed that there is access to discrete measurements of the state u(x,y,z,t) along the whole domain, which is feasible in practice. ξ(x,y,z,t) is a noise in the state dynamics. This model will be used just to generate the data to test the 3D identifier based on DNN. Boundary conditions and initial conditions were selected as follows:(39)u(0,0,0,t)=rand(1),ux(0,0,0,t)=0,uy(0,0,0,t)=0,uz(0,0,0,t)=0.
The trajectories of the model can be seen in Figure 3 as well as the estimated state, produced by the DNN identifier. The efficiency of the identification process provided by the suggested DNN algorithm shown in Figure 4.In Figures5 and 6 there are shown the trajectories of the PDE and the DNN for the coordinate z and the index error of the PDE and the DNN at 10 seconds for coordinates x, z, respectively (Figure 2).Figure 2
Numerical trajectory produced by the mathematical model described by 3D partial differential equation at time 10 s along the whole domain.Figure 3
Numerical trajectory produced by the neural network described reproducing the model at time 10 s along the whole domain.Figure 4
Error of the approximation. Difference between 3D PDE numerical solution and 3D DNN approximation at time 10 s along the whole domain.Figure 5
Comparison between 3D PDE trajectories and 3D DNN approximation for coordinatez at time 10 s.Figure 6
Index quadratic error for the whole time with coordinatesx, y, and z fixed.


## 7.2. Tumour Growth Example

The mathematical model of the brain tumour growth is presented in this section based on the results of [18]. Here the diffusion coefficient is considered as a constant. Let consider the following three-dimensional parabolic equation of the tumour growth described as follows:(40)ut(x,y,z,t)=-Pux(x,y,z,t)-Ruy(x,y,z,t)-Suz(x,y,z,t)+Quxx(x,y,z,t)+Quyy(x,y,z,t)+Quzz(x,y,z,t)+Γ-Lu(x,y,z,t),
where u(x,y,z,t) is the growth rate of a brain tumour, Q is the diffusion coefficient, W=(P,R,S) is the drift velocity field, Γ=Γ(u) is the proliferation coefficient, and L=L(u) is the decay coefficient of cells. It is assumed that there is access to discrete measurements of the state u(x,y,z,t) along the whole domain, which is feasible in practice by PET-CT (Positron emission tomography-computed tomography) technology. The domain 0≤x≤1, 0≤y≤1, and 0≤z≤1 This model will be used just to generate the data to test the 3D identifier based on DNN. Boundary conditions and initial conditions were selected as follows:(41)u(0,0,0,t)=200±20μ,ux(0,0,0,t)=0,uy(0,0,0,t)=0,uz(0,0,0,t)=0.
The trajectories of the model and the estimate state produced by the DNN identifier can be seen in Figure 7. The dissimilarity between both trajectories depends on the learning period required for adjusting the DNN identifier. The error between trajectories produced by the model and the proposed identifier is close to zero almost for all x, y, z and all t that shows the efficiency of the identification process provided by the suggested DNN algorithm is shown in Figure 8.Numerical trajectory produced by the mathematical model described by 3D partial differential equation at time 46 days along the whole domain.
(a)(b)Figure 8
Error of the approximation. Difference between 3D PDE numerical solution and 3D DNN approximation at time 46 days along the whole domain.In Figures9 and 10 there are shown the trajectories of the PDE and the DNN for the coordinate z and the index error of the PDE and the DNN at 46 days for coordinates x, z, respectively.Figure 9
Comparison between 3D PDE trajectories and 3D DNN approximation for coordinatez at time 46 days.Figure 10
Index quadratic error for the whole time with coordinatesx, y, and z fixed.


## 8. Conclusions

The adaptive DNN method proposed here solves the problem of non-parametric identification of nonlinear systems (with uncertainties) given by 3D uncertain PDE. Practical stability for the identification process is demonstrated based on the Lyapunov-like analysis. The upper bound of the identification error is explicitly constructed. Numerical examples demonstrate the estimation efficiency of the suggested methodology.



---

*Source: 618403-2012-01-24.xml*