Mathematical Formulation

Note

I cannot use hard mathematical terms and prove the existence of a solution on this page. Here I will focus on essential information.

Please go check the publication: Behzad Azmi, Andrea Petrocchi, Stefan Volkwein, Parameter optimization for elliptic-parabolic systems by an adaptive trust-region reduced basis method, Advances in Applied Mechanics, Elsevier, Volume 59, 2024, Pages 109-145, https://doi.org/10.1016/bs.aams.2024.07.001.

Pre-print available at https://arxiv.org/abs/2307.12723.

Important

Academic publishing companies are thieves and do not deserve our money. If you do not have access to the paper, send me an e-mail.

The geometry

Let \(T > 0\) be the finite time horizon, \(\Omega := (0, L)\) a one-dimensional interval (called in general terms domain) and the space-time domain \(Q_T := (0, T) \times (0, L)\).

The equations

We consider the following parameter-dependent parabolic-elliptic coupled system for the two state variables \(y, q: Q_T \to \mathbb{R}\) satisfying

\[ \begin{align}\begin{aligned}y_t (t,x) - \mu_1 \frac{\partial}{\partial x} \left( \kappa_1 (x) \frac{\partial}{\partial x} y (t,x) \right) - \mu_2 f \left( y(t,x), q(t,x) \right) &= 0 \qquad \text{f.a.a. } (t,x) \in Q_T,\\-\mu_3 \frac{\partial}{\partial x} \left( \kappa_2 (x) \frac{\partial}{\partial x} q (t,x) \right) + \mu_4 f \left( y(t,x), q(t,x) \right) &= 0 \qquad \text{f.a.a. } (t,x) \in Q_T\end{aligned}\end{align} \]

The functions \(\kappa_1\) and \(\kappa_2\) are called (space-dependent) diffusion coefficients and \(\boldsymbol{\mu} = [\mu_1, \mu_2, \mu_3, \mu_4]\) is a control vector.

The term “f.a.a.” means “for almost all”. This means, “for all values except those in a null set” (measure theory, fun stuff).

The boundary and initial conditions

For the state \(y\) we have homogeneous Neumann boundary conditions:

\[y_x(t,0) = y_x(t,L) = 0\quad \text{f.a.a. }t\in(0,T).\]

For \(q\) we have inhomogeneous Dirichlet-Neumann mixed boundary conditions:

\[q(t,0) = 0\text{ f.a.a. }t\in(0,T), \quad\mu_3 \kappa_2(L) q_x(t,L) = u(t)\text{ f.a.a. }t\in(0,T).\]

The time function \(u\) is called input function, it is user-defined and can change the behaviour of the two states greatly.

The initial conditions is

\[y(0,x) =y_\circ(x)\quad \text{f.a.a. }x\in\Omega.\]

The assumptions

The functions \(\kappa_1\) and \(\kappa_2\) belong to the space \(C^{0,1} (\bar\Omega)\). Namely, they are Lipschitz continuous in the closure of \(\Omega\), where we use the Hölder notation \(C^{0,\alpha}\). In simpler terms, the functions are continuous and their derivatives are dominated by the constant 1.

The initial condition \(y_\circ\) belongs to \(H^1 (\Omega)\), i.e. the Sobolev space \(W^{1,2} (\Omega)\). This means the function and its derivative are in \(L^2 (\Omega)\).

Tip

If you do not understand what I am talking about, go read a maths book.

The input function \(u\) satisfies \(u_\mathsf{a} (t) \le u (t) \le u_\mathsf{b} (t)\) for all \(t \in [0,T]\) with \(u_\mathsf{a}, u_\mathsf{b} \in L^\infty (\Omega)\). This guarantees that u is always finite and does not explode in time.

The nonlinearity is defined as \(f(\mathrm y, \mathrm q) := \sqrt{\mathrm y} \sinh(\mathrm q)\) for \(\mathrm y\in\mathbb R_{\ge}:= \{ s \in \mathbb{R} \, \vert \, s \ge 0 \}\) and \(\mathrm q \in \mathbb R\).

The weak formulation

The weak formulation is extremely helpful to define the finite element-based solution.

With \(V=H^1(\Omega)\), \(H=L^2(\Omega)\), \(V_\circ=\{ v \in V : v(0) = 0\}\), the weak formulation reads:

\[ \begin{align}\begin{aligned}\frac{\mathrm{d}}{\mathrm{d}t} \langle y (t), \varphi \rangle_H + \mu_1 a_1 (y(t), \varphi) - \mu_2 \langle g [y(t), q(t)], \varphi \rangle_{V', V} &= 0 \qquad \forall \varphi \in V,\\y(0) &= y_\circ \qquad \text{in } H,\\\mu_3 a_2 (q(t), \psi) + \mu_4 \langle g [y(t), q(t)], \psi \rangle_{V_\circ', V_\circ} - \langle b(t), \psi \rangle_{V'_\circ, V_\circ} &= 0 \qquad \forall \psi \in V_\circ,\end{aligned}\end{align} \]

where:

  • \(y(t)=y(t,\cdot) \in V\) and \(q(t)=q(t,\cdot) \in V_\circ\) f.a.a. \(t\),

  • \(\varphi \in V\) and \(\psi \in V_\circ\) are the test functions,

  • \(a_i\) for \(i=1,2\) are the diffusion operators,

  • \(g\) is the coupling operators,

  • \(\langle b(t), \psi \rangle_{V_\circ', V_\circ} = u(t) \psi(L)\) for \(\psi \in V_\circ\) and f.a.a. \(t\).

The solution space

In the paper mentioned above it is proven after 33 painful pages that there exist a solution \((y, q) \in \mathcal{Y}_T \times \mathcal{Q}_T\) to the weak formulation, where:

\[\mathcal{Y}_T := W(0, T; V, V') \cap C(\bar{Q}_T), \qquad \text{and} \qquad \mathcal{Q}_T := L^\infty (0, T; V_\circ).\]

Tip

Just go read the paper.

The time discretization

Here we just mention the Implicit Euler method, but the Crank-Nicolson method is also available in this package.

After time discretization, \(t_k = k \delta\) for \(k=0,\dots,K\), with \(\delta = T / K\), and \(K\) is the number of time steps. In the following, \(y^k\) and \(q^k\) are the approximated values of functions \(y(t_k)\) and \(q(t_k)\), respectively.

After time discretization, the solution is approximated by the following nonlinear system:

\[ \begin{align}\begin{aligned}\langle y^k, \varphi \rangle_H + \langle y^{k-1}, \varphi \rangle_H + \mu_1 \delta a_1 (y^k, \varphi) - \mu_2 \delta \langle g [y^k, q^k], \varphi \rangle_{V', V} &= 0 \qquad \forall \varphi \in V, \ k=1,\dots, K,\\y^0 &= y_\circ \qquad \text{in } H,\\\mu_3 a_2 (q^k, \psi) + \mu_4 \langle g [y^k, q^k], \psi \rangle_{V_\circ', V_\circ} - \langle b^k, \psi \rangle_{V'_\circ, V_\circ} &= 0 \qquad \forall \psi \in V_\circ, \ k=0, \dots, K.\end{aligned}\end{align} \]