Duality Principles and Numerical Procedures for a Large Class of Non-convex Models in the Calculus of Variations

1. Introduction

In this section we establish a dual formulation for a large class of models in non-convex optimization.

The main duality principle is applied to double well models similar as those found in the phase transition theory.

Such results are based on the works of J.J. Telega and W.R. Bielski [2,3,15,16] and on a D.C. optimization approach developed in Toland [17].

About the other references, details on the Sobolev spaces involved are found in [1]. Related results on convex analysis and duality theory are addressed in [5,7,8,10,14].

Finally, in this text we adopt the standard Einstein convention of summing up repeated indices, unless otherwise indicated.

In order to clarify the notation, here we introduce the definition of topological dual space.

Definition 1.1

(Topological dual spaces). Let U be a Banach space. We shall define its dual topological space, as the set of all linear continuous functionals defined on U. We suppose such a dual space of U, may be represented by another Banach space $U^{*}$, through a bilinear form ${\langle ·,·\rangle }_U:U\times U^{*}\to \mathbb{R}$ (here we are referring to standard representations of dual spaces of Sobolev and Lebesgue spaces). Thus, given $f:U\to \mathbb{R}$ linear and continuous, we assume the existence of a unique $u^{*}\in U^{*}$ such that

$$\begin{array}{c}f\left(u\right)={\langle u,u^{*}\rangle }_U,\forall u\in U.\end{array}$$

(1)

The norm of f , denoted by ${\parallel f\parallel }_{U^{*}}$, is defined as

$$\begin{array}{c}{\parallel f\parallel }_{U^{*}}=\underset{u\in U}{sup}\left\{\right|{\langle u,u^{*}\rangle }_U{|:\parallel u\parallel }_U\le 1\}\equiv \parallel u^{*}{\parallel }_{U^{*}}.\end{array}$$

(2)

At this point we start to describe the primal and dual variational formulations.

2. A general duality principle non-convex optimization

In this section we present a duality principle applicable to a model in phase transition.

This case corresponds to the vectorial one in the calculus of variations.

Let $\Omega \subset {\mathbb{R}}^n$ be an open, bounded, connected set with a regular (Lipschitzian) boundary denoted by $\partial \Omega .$

Consider a functional $J:V\to \mathbb{R}$ where

$$J\left(u\right)=F(\nabla u_1,\cdots ,\nabla u_N)+G(u_1,\cdots ,u_N)-{\langle u_i,f_i\rangle }_{L^2},$$

and where

$$V=\{u=(u_1,\cdots ,u_N)\in W^{1,p}(\Omega ;{\mathbb{R}}^N):u=u_0\mathrm{on}\partial \Omega \},$$

$f\in L^2(\Omega ;{\mathbb{R}}^N),$ and $1<p<+\infty .$

We assume there exists $\alpha \in \mathbb{R}$ such that

$$\alpha =\underset{u\in V}{inf}J\left(u\right).$$

Moreover, suppose F and G are Fréchet differentiable but not necessarily convex. A global optimum point may not be attained for J so that the problem of finding a global minimum for J may not be a solution.

Anyway, one question remains, how the minimizing sequences behave close the infimum of J.

We intend to use duality theory to approximately solve such a global optimization problem.

Denoting $V_0=W_0^{1,p}(\Omega ;{\mathbb{R}}^N)$, $Y_1=Y_1^{*}=L^2(\Omega ;{\mathbb{R}}^{N\times n})$, $Y_2=Y_2^{*}=L^2(\Omega ;{\mathbb{R}}^{N\times n})$, $Y_3=Y_3^{*}=L^2(\Omega ;{\mathbb{R}}^N),$ at this point we define, $F_1:V\times V_0\to \mathbb{R}$, $G_1:V\to \mathbb{R}$, $G_2:V\to \mathbb{R}$, $G_3:V_0\to \mathbb{R}$ and $G_4:V\to \mathbb{R},$ by

$$\begin{array}{ccc}\hfill F_1(\nabla u,\nabla \varphi )& =& F(\nabla u_1+\nabla {\varphi }_1,\cdots ,\nabla u_N+\nabla {\varphi }_N)+\frac{K}{2}{\int }_{\Omega }\nabla u_j·\nabla u_{j}dx\hfill \\ & & +\frac{K_2}{2}{\int }_{\Omega }\nabla {\varphi }_j·\nabla {\varphi }_{j}dx\hfill \end{array}$$

(3)

and

$$G_1(u_1,\cdots ,u_n)=G(u_1,\cdots ,u_N)+\frac{K_1}{2}{\int }_{\Omega }{u}_j{u}_{j}dx-{\langle u_i,f_i\rangle }_{L^2},$$

$$G_2(\nabla u_1,\cdots ,\nabla u_N)=\frac{K_1}{2}{\int }_{\Omega }\nabla u_j·\nabla u_{j}dx,$$

$$G_3(\nabla {\varphi }_1,\cdots ,\nabla {\varphi }_N)=\frac{K_2}{2}{\int }_{\Omega }\nabla {\varphi }_j·\nabla {\varphi }_{j}dx,$$

and

$$G_4(u_1,\cdots ,u_N)=\frac{K_1}{2}{\int }_{\Omega }{u}_j{u}_{j}dx.$$

Define now $J_1:V\times V_0\to \mathbb{R}$,

$$J_1(u,\varphi )=F(\nabla u+\nabla \varphi )+G\left(u\right)-{\langle u_i,f_i\rangle }_{L^2}.$$

Observe that

$$\begin{array}{ccc}\hfill J_1(u,\varphi )& =& F_1(\nabla u,\nabla \varphi )+G_1\left(u\right)-G_2(\nabla u)-G_3(\nabla \varphi )-G_4\left(u\right)\hfill \\ & \le & F_1(\nabla u,\nabla \varphi )+G_1\left(u\right)-{\langle \nabla u,z_1^{*}\rangle }_{L^2}-{\langle \nabla \varphi ,z_2^{*}\rangle }_{L^2}-{\langle u,z_3^{*}\rangle }_{L^2}\hfill \\ & & +\underset{v_1\in Y_1}{sup}\{{\langle v_1,z_1^{*}\rangle }_{L^2}-G_2\left(v_1\right)\}\hfill \\ & & +\underset{v_2\in Y_2}{sup}\{{\langle v_2,z_2^{*}\rangle }_{L^2}-G_3\left(v_2\right)\}\hfill \\ & & +\underset{u\in V}{sup}\{{\langle u,z_3^{*}\rangle }_{L^2}-G_4\left(u\right)\}\hfill \\ & =& F_1(\nabla u,\nabla \varphi )+G_1\left(u\right)-{\langle \nabla u,z_1^{*}\rangle }_{L^2}-{\langle \nabla \varphi ,z_2^{*}\rangle }_{L^2}-{\langle u,z_3^{*}\rangle }_{L^2}\hfill \\ & & +G_2^{*}\left(z_1^{*}\right)+G_3^{*}\left(z_2^{*}\right)+G_4^{*}\left(z_3^{*}\right)\hfill \\ & =& J_1^{*}(u,\varphi ,z^{*}),\hfill \end{array}$$

(4)

$\forall u\in V,\varphi \in V_0,z^{*}=(z_1^{*},z_2^{*},z_3^{*})\in Y^{*}=Y_1^{*}\times Y_2^{*}\times Y_3^{*}$.

From the general results in [17], we may infer that

$$\begin{array}{ccc}\hfill \underset{(u,\varphi )\in V\times V_0}{inf}J(u,\varphi )& =& \underset{(u,\varphi ,z^{*})\in V\times V_0\times Y^{*}}{inf}{J}_1^{*}(u,\varphi ,z^{*}).\hfill \end{array}$$

(5)

On the other hand

$$\underset{u\in V}{inf}J\left(u\right)\ge \underset{(u,\varphi )\in V\times V_0}{inf}{J}_1(u,\varphi ).$$

From these last two results we may obtain

$$\underset{u\in V}{inf}J\left(u\right)\ge \underset{(u,\varphi ,z^{*})\in V\times V_0\times Y^{*}}{inf}{J}_1^{*}(u,\varphi ,z^{*}).$$

Moreover, from standards results on convex analysis, we may have

$$\begin{array}{ccc}\hfill \underset{u\in V}{inf}{J}_1^{*}(u,\varphi ,z^{*})& =& \underset{u\in V}{inf}\{F_1(\nabla u,\nabla \varphi )+G_1\left(u\right)\hfill \\ & & -{\langle \nabla u,z_1^{*}\rangle }_{L^2}-{\langle \nabla \varphi ,z_2^{*}\rangle }_{L^2}-{\langle u,z_3^{*}\rangle }_{L^2}\hfill \\ & & +G_2^{*}\left(z_1^{*}\right)+G_3^{*}\left(z_2^{*}\right)+G_4^{*}\left(z_3^{*}\right)\}\hfill \\ & =& \underset{(v_1^{*},v_2^{*})\in C^{*}}{sup}\{-F_1^{*}(v_1^{*}+z_1^{*},\nabla \varphi )-G_1^{*}(v_2^{*}+z_3^{*})-{\langle \nabla \varphi ,z_2^{*}\rangle }_{L^2}\hfill \\ & & +G_2^{*}\left(z_1^{*}\right)+G_3^{*}\left(z_2^{*}\right)+G_4^{*}\left(z_3^{*}\right)\},\hfill \end{array}$$

(6)

where

$$C^{*}=\{v^{*}=(v_1^{*},v_2^{*})\in Y_1^{*}\times Y_3^{*}:-\mathrm{div}{\left(v_1^{*}\right)}_i+{\left(v_2^{*}\right)}_i=\mathbf{0},\forall i\in \{1,\cdots ,N\}\},$$

$$F_1^{*}(v_1^{*}+z_1^{*},\nabla \varphi )=\underset{v_1\in Y_1}{sup}\{{\langle v_1,z_1^{*}+v_1^{*}\rangle }_{L^2}-F_1(v_1,\nabla \varphi )\},$$

and

$$G_1^{*}(v_2^{*}+z_2^{*})=\underset{u\in V}{sup}\{{\langle u,v_2^{*}+z_2^{*}\rangle }_{L^2}-G_1\left(u\right)\}.$$

Thus, defining

$$J_2^{*}(\varphi ,z^{*},v^{*})=F_1^{*}(v_1^{*}+z_1^{*},\nabla \varphi )-G_1^{*}(v_2^{*}+z_3^{*})-{\langle \nabla \varphi ,z_2^{*}\rangle }_{L^2}+G_2^{*}\left(z_1^{*}\right)+G_3^{*}\left(z_2^{*}\right)+G_4^{*}\left(z_3^{*}\right),$$

we have got

$$\begin{array}{ccc}\hfill \underset{u\in V}{inf}J\left(u\right)& \ge & \underset{(u,\varphi )\in V\times V_0}{inf}{J}_1(u,\varphi )\hfill \\ & =& \underset{(u,\varphi ,z^{*})\in V\times V_0\times Y^{*}}{inf}{J}_1^{*}(u,\varphi ,z^{*})\hfill \\ & =& \underset{z^{*}\in Y^{*}}{inf}\left\{\underset{\varphi \in V_0}{inf}\left\{\underset{v^{*}\in C^{*}}{sup}{J}_2^{*}(\varphi ,z^{*},v^{*})\right\}\right\}.\hfill \end{array}$$

(7)

Finally, observe that

$$\begin{array}{ccc}& & \underset{u\in V}{inf}J\left(u\right)\hfill \\ & \ge & \underset{z^{*}\in Y^{*}}{inf}\left\{\underset{\varphi \in V_0}{inf}\left\{\underset{v^{*}\in C^{*}}{sup}{J}_2^{*}(\varphi ,z^{*},v^{*})\right\}\right\}\hfill \\ & \ge & \underset{v^{*}\in C^{*}}{sup}\left\{\underset{(z^{*},\varphi )\in Y^{*}\times V_0}{inf}{J}_2^{*}(\varphi ,z^{*},v^{*})\right\}.\hfill \end{array}$$

(8)

This last variational formulation corresponds to a concave relaxed formulation in $v^{*}$ concerning the original primal formulation.

3. Another duality principle for a simpler related model in phase transition with a respective numerical example

In this section we present another duality principle for a related model in phase transition.

Let $\Omega =[0,1]\subset \mathbb{R}$ and consider a functional $J:V\to \mathbb{R}$ where

$$J\left(u\right)=\frac{1}{2}{\int }_{\Omega }{({\left(u^{\prime }\right)}^2-1)}^{2}dx+\frac{1}{2}{\int }_{\Omega }{u}^{2}dx-{\langle u,f\rangle }_{L^2},$$

and where

$$V=\{u\in W^{1,4}(\Omega ):u\left(0\right)=0\mathrm{and}u\left(1\right)=1/2\}$$

and $f\in L^2(\Omega ).$

A global optimum point is not attained for J so that the problem of finding a global minimum for J has no solution.

Anyway, one question remains, how the minimizing sequences behave close the infimum of J.

We intend to use duality theory to approximately solve such a global optimization problem.

Denoting $V_0=W_0^{1,4}(\Omega )$, at this point we define, $F:V\to \mathbb{R}$ and $F_1:V\times V_0\to \mathbb{R}$ by

$$F\left(u\right)=\frac{1}{2}{\int }_{\Omega }{({\left(u^{\prime }\right)}^2-1)}^{2}dx,$$

and

$$F_1(u,\varphi )=\frac{1}{2}{\int }_{\Omega }{({(u^{\prime }+{\varphi }^{\prime })}^2-1)}^{2}dx.$$

Observe that

$$F\left(u\right)\ge \underset{\varphi \in V_0}{inf}{F}_1(u,\varphi ),\forall u\in V.$$

In order to restrict the action of $\varphi$ on the region where the primal functional is non-convex, we redefine a not relabeled

$$V_0=\left\{\varphi \in W_0^{1,4}(\Omega ):{\left({\varphi }^{\prime }\right)}^2-1\le 0,\mathrm{in}\Omega \right\}$$

and define also

$$F_2:V\times V_0\to \mathbb{R},$$

$$F_3:V\times V_0\to \mathbb{R}$$

and

$$G:V\times V_0\to \mathbb{R}$$

by

$$F_2(u,\varphi )=\frac{1}{2}{\int }_{\Omega }{({(u^{\prime }+{\varphi }^{\prime })}^2-1)}^{2}dx+\frac{1}{2}{\int }_{\Omega }{u}^{2}dx-{\langle u,f\rangle }_{L^2},$$

$$\begin{array}{ccc}\hfill F_3(u,\varphi )& =& F_2(u,\varphi )+\frac{K}{2}{\int }_{\Omega }{\left(u^{\prime }\right)}^{2}dx\hfill \\ & & +\frac{K_1}{2}{\int }_{\Omega }{\left({\varphi }^{\prime }\right)}^{2}dx\hfill \end{array}$$

(9)

and

$$\begin{array}{ccc}\hfill G(u,\varphi )& =& \frac{K}{2}{\int }_{\Omega }{\left(u^{\prime }\right)}^{2}dx\hfill \\ & & +\frac{K_1}{2}{\int }_{\Omega }{\left({\varphi }^{\prime }\right)}^{2}dx\hfill \end{array}$$

(10)

Denoting $Y=Y^{*}=L^2(\Omega )$ we also define the polar functional $G^{*}:Y^{*}\times Y^{*}\to \mathbb{R}$ by

$$G^{*}(v^{*},v_0^{*})=\underset{(u,\varphi )\in V\times V_0}{sup}\{{\langle u,v^{*}\rangle }_{L^2}+{\langle \varphi ,v_0^{*}\rangle }_{L^2}-G(u,\varphi )\}.$$

Observe that

$$\underset{u\in U}{inf}J\left(u\right)\ge \underset{((u,\varphi ),(v^{*},v_0^{*}))\in V\times V_0\times {\left[Y^{*}\right]}^2}{inf}\{G^{*}(v^{*},v_0^{*})-{\langle u,v^{*}\rangle }_{L^2}-{\langle \varphi ,v_0^{*}\rangle }_{L^2}+F_3(u,\varphi )\}.$$

With such results in mind, we define a relaxed primal dual variational formulation for the primal problem, represented by $J_1^{*}:V\times V_0\times {\left[Y^{*}\right]}^2\to \mathbb{R}$, where

$$J_1^{*}(u,\varphi ,v^{*},v_0^{*})=G^{*}(v^{*},v_0^{*})-{\langle u,v^{*}\rangle }_{L^2}-{\langle \varphi ,v_0^{*}\rangle }_{L^2}+F_3(u,\varphi ).$$

Having defined such a functional, we may obtain numerical results by solving a sequence of convex auxiliary sub-problems, through the following algorithm (in order to obtain the concerning critical points, at first we have neglected the constraint ${\left({\varphi }^{\prime }\right)}^2-1\le 0\mathrm{in}\Omega$).

• Set $K\approx 0.1$ and $K_1=120.0$ and $0<\epsilon \ll 1.$
• Choose $(u_1,{\varphi }_1)\in V\times V_0,$ such that $\parallel u_1{\parallel }_{1,\infty }<1$ and $\parallel {\varphi }_1{\parallel }_{1,\infty }<1.$
• Set $n=1.$
• Calculate $(v_n^{*},{\left(v_0^{*}\right)}_n)$ solution of the system of equations:

  $$\frac{\partial J_1^{*}(u_n,{\varphi }_n,v_n^{*},{\left(v_0^{*}\right)}_n)}{\partial v^{*}}=\mathbf{0}$$

  and

  $$\frac{\partial J_1^{*}(u_n,{\varphi }_n,v_n^{*},{\left(v_0^{*}\right)}_n)}{\partial v_0^{*}}=\mathbf{0},$$

  that is

  $$\frac{\partial G^{*}(v_n^{*},{\left(v_0^{*}\right)}_n)}{\partial v^{*}}-u_n=0$$

  and

  $$\frac{\partial G^{*}(v_n^{*},{\left(v_0^{*}\right)}_n)}{\partial v_0^{*}}-{\varphi }_n=0$$

  so that

  $$v_n^{*}=\frac{\partial G(u_n,{\varphi }_n)}{\partial u}$$

  and

  $${\left(v_0^{*}\right)}_n^{*}=\frac{\partial G(u_n,{\varphi }_n)}{\partial \varphi }$$
• Calculate $(u_{n+1},{\varphi }_{n+1})$ by solving the system of equations:

  $$\frac{\partial J_1^{*}(u_{n+1},{\varphi }_{n+1},v_n^{*},{\left(v_0^{*}\right)}_n)}{\partial u}=\mathbf{0}$$

  and

  $$\frac{\partial J_1^{*}(u_{n+1},{\varphi }_{n+1},v_n^{*},{\left(v_0^{*}\right)}_n)}{\partial \varphi }=\mathbf{0}$$

  that is

  $$-v_n^{*}+\frac{\partial F_3(u_{n+1},{\varphi }_{n+1})}{\partial u}=\mathbf{0}$$

  and

  $$-{\left(v_0^{*}\right)}_n+\frac{\partial F_3(u_{n+1},{\varphi }_{n+1})}{\partial \varphi }=\mathbf{0}$$
• If $max\{\parallel u_n-u_{n+1}{\parallel }_{\infty },\parallel {\varphi }_{n+1}-{\varphi }_n{\parallel }_{\infty }\}\le \epsilon$, then stop, else set $n:=n+1$ and go to item d.

At this point, we present the corresponding software in MAT-LAB, in finite differences and based on the one-dimensional version of the generalized method of lines.

Here the software.

***********************

• clear all
  m8=300;
  d=1/m8;
  K=0.1;
  K1=120;
  for i=1:m8
  $uo(i,1)=i^2*d/2;$
  vo(i,1)=i*d/10;
  yo(i,1)=sin(i*d*pi)/2;
  end;
  k=1;
  b12=1.0;
  while $(b12>{10}^{-4.3})$ and $(k<230000)$
  k=k+1;
  for i=1:m8-1
  duo(i,1)=(uo(i+1,1)-uo(i,1))/d;
  dvo(i,1)=(vo(i+1,1)-vo(i,1))/d;
  end;
  m9=zeros(2,2);
  m9(1,1)=1;
  i=1;
  $f1=6*{(duo(i,1)+dvo(i,1))}^2-2;$
  m80(1,1,i)=-f1-K;
  m80(1,2,i)=-f1;
  m80(2,1,i)=-f1;
  m80(2,2,i)=-f1-K1;
  $y11(1,i)=K*(uo(i+1,1)-2*uo(i,1))/d^2-yo(i,1);$
  $y11(2,i)=K1*(vo(i+1,1)-2*vo(i,1))/d^2;$
  $m12=2*m80(:,:,i)-m9*d^2;$
  m50(:,:,i)=m80(:,:,i)*inv(m12);
  z(:,i)=inv(m12)*y11(:,i)*$d^2$;
  for i=2:m8-1
  $f1=6*{(duo(i,1)+dvo(i,1))}^2-2$;
  m80(1,1,i)=-f1-K;
  m80(1,2,i)=-f1;
  m80(2,1,i)=-f1;
  m80(2,2,i)=-f1-K1;
  $y11(1,i)=K*(uo(i+1,1)-2*uo(i,1)+uo(i-1,1))/d^2-yo(i,1);$
  $y11(2,i)=K1*(vo(i+1,1)-2*vo(i,1)+vo(i-1,1))/d^2;$
  $m12=2*m80(:,:,i)-m9*d^2-m80(:,:,i)*m50(:,:,i-1);$
  m50(:,:,i)=inv(m12)*m80(:,:,i);
  $z(:,i)=inv\left(m12\right)*(y11(:,i)*d^2+m80(:,:,i)*z(:,i-1));$
  end;
  U(1,m8)=1/2;
  U(2,m8)=0.0;
  for i=1:m8-1
  U(:,m8-i)=m50(:,:,m8-i)*U(:,m8-i+1)+z(:,m8-i);
  end;
  for i=1:m8
  u(i,1)=U(1,i);
  v(i,1)=U(2,i);
  end;
  b12=max(abs(u-uo))
  uo=u;
  vo=v;
  u(m8/2,1)
  end;
  for i=1:m8
  y(i)=i*d;
  end;
  plot(y,uo)
  **************************************

For the case in which $f\left(x\right)=0$, we have obtained numerical results for $K=0.1$ and $K_1=120$. For such a concerning solution $u_0$ obtained, please see Figure 1. For the case in which $f\left(x\right)=sin\left(\pi x\right)/2$, we have obtained numerical results also for $K=0.1$ and $K_1=120$. For such a concerning solution $u_0$ obtained, please see Figure 2.

Remark 3.1.

Observe that the solutions obtained are approximate critical points. They are not, in a classical sense, the global solutions for the related optimization problems. Indeed, such solutions reflect the average behavior of weak cluster points for concerning minimizing sequences.

3.1. A general proposal for relaxation

Let $\Omega \subset {\mathbb{R}}^n$ be an open, bounded and connected set with a regular (Lipschitzian) boundary denoted by $\partial \Omega .$

Consider a functional $J:V\to \mathbb{R}$ where

$$J\left(u\right)=F(\nabla u)+G\left(u\right)-{\langle u,f_1\rangle }_{L^2},$$

where

$$V=\left\{u\in W^{1,4}(\Omega ;{\mathbb{R}}^N):u=u_0\mathrm{on}\partial \Omega \right\},$$

$$u_0\in C^1(\Omega ;{\mathbb{R}}^N),$$

$f_1\in L^2(\Omega ;{\mathbb{R}}^N)$, $G:V\to \mathbb{R}$ is convex and Fréchet differentiable, and

$$F(\nabla u)={\int }_{\Omega }f(\nabla u)dx,$$

where $f:{\mathbb{R}}^{N\times n}\to \mathbb{R}$ is also Fréchet differentiable.

Assume there exists $\widehat{N}\in \mathbb{N}$ such that

$$W_h\equiv \left\{y\in {\mathbb{R}}^{N\times n}:f^{**}\left(y\right)<f\left(y\right)\right\}={\cup }_{j=1}^{\widehat{N}}{W}_j$$

where for each $j\in \{1,\cdots ,\widehat{N}\}$$W_j\subset {\mathbb{R}}^{N\times n}$ is an open connected set such that $\partial W_j$ is regular. We also suppose

$$\overline{W_j}\cap \overline{W_k}=\varnothing ,\forall j\ne k.$$

Define

$${\widehat{W}}_j=\left\{v_j\in W_0^{1,4}(\Omega ;{\mathbb{R}}^N);\nabla v_j\left(x\right)\in W_j,\mathrm{a}.\mathrm{e}.\mathrm{in}\Omega \right\}$$

and define also

$$W=\left\{v=(v_1,\cdots ,v_{\widehat{N}}):v_j\in {\widehat{W}}_j\forall j\in \{1,\cdots ,\widehat{N}\}\mathrm{and}\mathrm{supp}{v}_j\cap \mathrm{supp}{v}_k=\varnothing ,\forall j\ne k\right\}.$$

At this point we define

$$h_5(u\left(x\right),v\left(x\right))=\left\{\begin{array}{cc}f(\nabla u\left(x\right)+\nabla v_j\left(x\right)),\hfill & \mathrm{if}\nabla u\left(x\right)\in W_j,\hfill \\ f(\nabla u(x\left)\right),\hfill & \mathrm{if}\nabla u\left(x\right)\notin W_h,\hfill \end{array}\right.$$

(11)

and

$$H\left(u\right)=\underset{v\in W_u}{inf}{\int }_{\Omega }{h}_5(u,v)dx,$$

where

$$W_u=\{v\in W:\nabla u\left(x\right)+\nabla v_j\left(x\right)\in W_j,\mathrm{if}\nabla u\left(x\right)\in W_j,\mathrm{a}.\mathrm{e}.\mathrm{in}\Omega ,\forall j\in \{1,\cdots ,\widehat{N}\}\}.$$

Moreover, we propose the relaxed functional

$$J_1\left(u\right)=H\left(u\right)+G\left(u\right)-{\langle u,f_1\rangle }_{L^2}.$$

Observe that clearly

$$\underset{u\in V}{inf}{J}_1\left(u\right)\le \underset{u\in V}{inf}J\left(u\right).$$

4. A convex dual variational formulation for a third similar model

In this section we present another duality principle for a third related model in phase transition.

Let $\Omega =[0,1]\subset \mathbb{R}$ and consider a functional $J:V\to \mathbb{R}$ where

$$J\left(u\right)=\frac{1}{2}{\int }_{\Omega }min\{{(u^{\prime }-1)}^2,{(u^{\prime }+1)}^2\}dx+\frac{1}{2}{\int }_{\Omega }{u}^{2}dx-{\langle u,f\rangle }_{L^2},$$

and where

$$V=\{u\in W^{1,2}(\Omega ):u\left(0\right)=0\mathrm{and}u\left(1\right)=1/2\}$$

and $f\in L^2(\Omega ).$

A global optimum point is not attained for J so that the problem of finding a global minimum for J has no solution.

Anyway, one question remains, how the minimizing sequences behave close to the infimum of J.

We intend to use the duality theory to solve such a global optimization problem in an appropriate sense to be specified.

At this point we define, $F:V\to \mathbb{R}$ and $G:V\to \mathbb{R}$ by

$$\begin{array}{ccc}\hfill F\left(u\right)& =& \frac{1}{2}{\int }_{\Omega }min\{{(u^{\prime }-1)}^2,{(u^{\prime }+1)}^2\}dx\hfill \\ & =& \frac{1}{2}{\int }_{\Omega }{\left(u^{\prime }\right)}^{2}dx-{\int }_{\Omega }\left|u^{\prime }\right|dx+1/2\hfill \\ & \equiv & F_1\left(u^{\prime }\right),\hfill \end{array}$$

(12)

and

$$G\left(u\right)=\frac{1}{2}{\int }_{\Omega }{u}^{2}dx-{\langle u,f\rangle }_{L^2}.$$

Denoting $Y=Y^{*}=L^2(\Omega )$ we also define the polar functional $F_1^{*}:Y^{*}\to \mathbb{R}$ and $G^{*}:Y^{*}\to \mathbb{R}$ by

$$\begin{array}{ccc}\hfill F_1^{*}\left(v^{*}\right)& =& \underset{v\in Y}{sup}\{{\langle v,v^{*}\rangle }_{L^2}-F_1\left(v\right)\}\hfill \\ & =& \frac{1}{2}{\int }_{\Omega }{\left(v^{*}\right)}^{2}dx+{\int }_{\Omega }\left|v^{*}\right|dx,\hfill \end{array}$$

(13)

and

$$\begin{array}{ccc}\hfill G^{*}\left({\left(v^{*}\right)}^{\prime }\right)& =& \underset{u\in V}{sup}\{-{\langle u^{\prime },v^{*}\rangle }_{L^2}-G\left(u\right)\}\hfill \\ & =& \frac{1}{2}{\int }_{\Omega }{({\left(v^{*}\right)}^{\prime }+f)}^{2}dx-\frac{1}{2}{v}^{*}\left(1\right).\hfill \end{array}$$

(14)

Observe this is the scalar case of the calculus of variations, so that from the standard results on convex analysis, we have

$$\underset{u\in V}{inf}J\left(u\right)=\underset{v^{*}\in Y^{*}}{max}\{-F_1^{*}\left(v^{*}\right)-G^{*}(-{\left(v^{*}\right)}^{\prime })\}.$$

Indeed, from the direct method of the calculus of variations, the maximum for the dual formulation is attained at some ${\widehat{v}}^{*}\in Y^{*}$.

Moreover, the corresponding solution $u_0\in V$ is obtained from the equation

$$u_0=\frac{\partial G\left({\left({\widehat{v}}^{*}\right)}^{\prime }\right)}{\partial {\left(v^{*}\right)}^{\prime }}={\left({\widehat{v}}^{*}\right)}^{\prime }+f.$$

Finally, the Euler-Lagrange equations for the dual problem stands for

$$\left\{\begin{array}{cc}{\left(v^{*}\right)}^{\prime \prime }+f^{\prime }-v^{*}-\mathrm{sign}\left(v^{*}\right)=0,\hfill & \mathrm{in}\Omega ,\hfill \\ {\left(v^{*}\right)}^{\prime }\left(0\right)+f\left(0\right)=0,{\left(v^{*}\right)}^{\prime }\left(1\right)+f\left(1\right)=1/2,\hfill \end{array}\right.$$

(15)

where $\mathrm{sign}\left(v^{*}\left(x\right)\right)=1$ if $v^{*}\left(x\right)>0,$$\mathrm{sign}\left(v^{*}\left(x\right)\right)=-1,$ if $v^{*}\left(x\right)<0$ and

$$-1\le \mathrm{sign}\left(v^{*}\left(x\right)\right)\le 1,$$

if $v^{*}\left(x\right)=0.$

We have computed the solutions $v^{*}$ and corresponding solutions $u_0\in V$ for the cases in which $f\left(x\right)=0$ and $f\left(x\right)=sin\left(\pi x\right)/2.$

For the solution $u_0\left(x\right)$ for the case in which $f\left(x\right)=0$, please see Figure 3.

For the solution $u_0\left(x\right)$ for the case in which $f\left(x\right)=sin\left(\pi x\right)/2$, please see Figure 4.

Remark 4.1.

Observe that such solutions $u_0$ obtained are not the global solutions for the related primal optimization problems. Indeed, such solutions reflect the average behavior of weak cluster points for concerning minimizing sequences.

4.1. The algorithm through which we have obtained the numerical results

In this subsection we present the software in MATLAB through which we have obtained the last numerical results.

This algorithm is for solving the concerning Euler-Lagrange equations for the dual problem, that is, for solving the equation

$$\left\{\begin{array}{cc}{\left(v^{*}\right)}^{\prime \prime }+f^{\prime }-v^{*}-\mathrm{sign}\left(v^{*}\right)=0,\hfill & \mathrm{in}\Omega ,\hfill \\ {\left(v^{*}\right)}^{\prime }\left(0\right)=0,{\left(v^{*}\right)}^{\prime }\left(1\right)=1/2.\hfill \end{array}\right.$$

(16)

Here the concerning software in MATLAB. We emphasize to have used the smooth approximation

$$|v^{*}|\approx \sqrt{{\left(v^{*}\right)}^2+e_1},$$

where a small value for $e_1$ is specified in the next lines.

*************************************

• clear all
• $m_8=800;$ (number of nodes)
• $d=1/m_8;$
• $e_1=0.00001;$
• $fori=1:m_8$
  $yo(i,1)=0.01;$
  $y_1(i,1)=sin(\pi *i/m_8)/2;$
  $end;$
• $fori=1:m_8-1$
  $d{y}_1(i,1)=(y_1(i+1,1)-y_1(i,1))/d;$
  $end;$
• $fork=1:3000$ (we have fixed the number of iterations)
  $i=1;$
  $h_3=1/\sqrt{vo{(i,1)}^2+e_1};$
  $m_{12}=1+d^2*h_3+d^2;$
  $m_{50}\left(i\right)=1/m_{12};$
  $z\left(i\right)=m_{50}\left(i\right)*(d{y}_1(i,1)*d^2);$
• $fori=2:m_8-1$
  $h_3=1/\sqrt{vo{(i,1)}^2+e_1};$
  $m_{12}=2+h_3*d^2+d^2-m50(i-1);$
  $m50\left(i\right)=1/m_{12};$
  $z\left(i\right)=m_{50}\left(i\right)*(z(i-1)+d{y}_1(i,1)*d^2);$
  $end;$
• $v(m_8,1)=(d/2+z(m_8-1))/(1-m_{50}(m_8-1));$
• $fori=1:m_8-1$
  $v(m_8-i,1)=m_{50}(m_8-i)*v(m_8-i+1)+z(m_8-i);$
  $end;$
• $v(m_8/2,1)$
• $vo=v;$
  $end;$
• $fori=1:m_8-1$
  $u(i,1)=(v(i+1,1)-v(i,1))/d+y_1(i,1);$
  $end;$
• $fori=1:m_8-1$
  $x\left(i\right)=i*d;$
  $end;$
  $plot(x,u(:,1\left)\right)$

********************************

5. An improvement of the convexity conditions for a non-convex related model through an approximate primal formulation

In this section we develop an approximate primal dual formulation suitable for a large class of variational models.

Here, the applications are for the Kirchhoff-Love plate model, which may be found in Ciarlet, [11].

At this point we start to describe the primal variational formulation.

Let $\Omega \subset {\mathbb{R}}^2$ be an open, bounded, connected set which represents the middle surface of a plate of thickness h. The boundary of $\Omega$, which is assumed to be regular (Lipschitzian), is denoted by $\partial \Omega$. The vectorial basis related to the cartesian system $\{x_1,x_2,x_3\}$ is denoted by $({\mathbf{a}}_{\alpha },{\mathbf{a}}_3)$, where $\alpha =1,2$ (in general Greek indices stand for 1 or 2), and where ${\mathbf{a}}_3$ is the vector normal to $\Omega$, whereas ${\mathbf{a}}_1$ and ${\mathbf{a}}_2$ are orthogonal vectors parallel to $\Omega .$ Also, $\mathbf{n}$ is the outward normal to the plate surface.

The displacements will be denoted by

$$\widehat{\mathbf{u}}=\{{\widehat{u}}_{\alpha },{\widehat{u}}_3\}={\widehat{u}}_{\alpha }{\mathbf{a}}_{\alpha }+{\widehat{u}}_3{\mathbf{a}}_3.$$

The Kirchhoff-Love relations are

$$\begin{array}{ccc}& & {\widehat{u}}_{\alpha }(x_1,x_2,x_3)=u_{\alpha }(x_1,x_2)-x_{3}w{(x_1,x_2)}_{,\alpha }\hfill \\ & & \mathrm{and}{\widehat{u}}_3(x_1,x_2,x_3)=w(x_1,x_2).\hfill \end{array}$$

(17)

Here $-h/2\le x_3\le h/2$ so that we have $u=(u_{\alpha },w)\in U$ where

$$\begin{array}{ccc}\hfill U& =& \left\{u=(u_{\alpha },w)\in W^{1,2}(\Omega ;{\mathbb{R}}^2)\times W^{2,2}(\Omega ),\right.\hfill \\ & & u_{\alpha }=w=\frac{\partial w}{\partial \mathbf{n}}=0\left.\mathrm{on}\partial \Omega \right\}\hfill \\ & =& W_0^{1,2}(\Omega ;{\mathbb{R}}^2)\times W_0^{2,2}(\Omega ).\hfill \end{array}$$

It is worth emphasizing that the boundary conditions here specified refer to a clamped plate.

We also define the operator $\Lambda :U\to Y\times Y$, where $Y=Y^{*}=L^2(\Omega ;{\mathbb{R}}^{2\times 2})$, by

$$\Lambda \left(u\right)=\left\{\gamma \right(u),\kappa (u\left)\right\},$$

$${\gamma }_{\alpha \beta }\left(u\right)=\frac{u_{\alpha ,\beta }+u_{\beta ,\alpha }}{2}+\frac{w_{,\alpha }{w}_{,\beta }}{2},$$

$${\kappa }_{\alpha \beta }\left(u\right)=-w_{,\alpha \beta }.$$

The constitutive relations are given by

$$N_{\alpha \beta }\left(u\right)=H_{\alpha \beta \lambda \mu }{\gamma }_{\lambda \mu }\left(u\right),$$

(18)

$$M_{\alpha \beta }\left(u\right)=h_{\alpha \beta \lambda \mu }{\kappa }_{\lambda \mu }\left(u\right),$$

(19)

where: $\left\{H_{\alpha \beta \lambda \mu }\right\}$ and $\left\{h_{\alpha \beta \lambda \mu }=\frac{h^2}{12}{H}_{\alpha \beta \lambda \mu }\right\}$, are symmetric positive definite fourth order tensors. From now on, we denote $\left\{{\overline{H}}_{\alpha \beta \lambda \mu }\right\}={\left\{H_{\alpha \beta \lambda \mu }\right\}}^{-1}$ and $\left\{{\overline{h}}_{\alpha \beta \lambda \mu }\right\}={\left\{h_{\alpha \beta \lambda \mu }\right\}}^{-1}$.

Furthermore $\left\{N_{\alpha \beta }\right\}$ denote the membrane force tensor and $\left\{M_{\alpha \beta }\right\}$ the moment one. The plate stored energy, represented by $(G\circ \Lambda ):U\to \mathbb{R}$ is expressed by

$$(G\circ \Lambda )\left(u\right)=\frac{1}{2}{\int }_{\Omega }{N}_{\alpha \beta }\left(u\right){\gamma }_{\alpha \beta }\left(u\right)dx+\frac{1}{2}{\int }_{\Omega }{M}_{\alpha \beta }\left(u\right){\kappa }_{\alpha \beta }\left(u\right)dx$$

(20)

and the external work, represented by $F:U\to \mathbb{R}$, is given by

$$F\left(u\right)={\langle w,P\rangle }_{L^2}+{\langle u_{\alpha },P_{\alpha }\rangle }_{L^2},$$

(21)

where $P,P_1,P_2\in L^2(\Omega )$ are external loads in the directions ${\mathbf{a}}_3$, ${\mathbf{a}}_1$ and ${\mathbf{a}}_2$ respectively. The potential energy, denoted by $J:U\to \mathbb{R}$ is expressed by:

$$J\left(u\right)=(G\circ \Lambda )\left(u\right)-F\left(u\right)$$

Define now $J_3:\tilde{U}\to \mathbb{R}$ by

$$J_3\left(u\right)=J\left(u\right)+J_5\left(w\right).$$

where

$$J_5\left(w\right)=10{\int }_{\Omega }\frac{a^{Kbw}}{ln\left(a\right)K^{3/2}}dx+10{\int }_{\Omega }\frac{a^{-K(bw-1/100)}}{ln\left(a\right)K^{3/2}}dx.$$

In such a case for $a=2.71$, $K=185$, $b=P/\left|P\right|$ in $\Omega$ and

$$\tilde{U}={\{u\in U:\parallel w\parallel }_{\infty }\le 0.01\mathrm{and}Pw\ge 0\mathrm{a}.\mathrm{e}.\mathrm{in}\Omega \},$$

we get

$$\begin{array}{ccc}\hfill \frac{\partial J_3\left(u\right)}{\partial w}& =& \frac{\partial J\left(u\right)}{\partial w}+\frac{\partial J_5\left(u\right)}{\partial w}\hfill \\ & \approx & \frac{\partial J\left(u\right)}{\partial w}+\mathcal{O}(\pm 3.0),\hfill \end{array}$$

(22)

and

$$\begin{array}{ccc}\hfill \frac{{\partial }^2{J}_3\left(u\right)}{\partial w^2}& =& \frac{{\partial }^{2}J\left(u\right)}{\partial w^2}+\frac{{\partial }^2{J}_5\left(u\right)}{\partial w^2}\hfill \\ & \approx & \frac{{\partial }^{2}J\left(u\right)}{\partial w^2}+\mathcal{O}\left(850\right).\hfill \end{array}$$

(23)

This new functional $J_3$ has a relevant improvement in the convexity conditions concerning the previous functional J.

Indeed, we have obtained a gain in positiveness for the second variation $\frac{{\partial }^{2}J\left(u\right)}{\partial w^2},$ which has increased of order $\mathcal{O}(700-1000).$

Moreover the difference between the approximate and exact equation

$$\frac{\partial J\left(u\right)}{\partial w}=\mathbf{0}$$

is of order $\mathcal{O}(\pm 3.0)$ which corresponds to a small perturbation in the original equation for a load of $P=1500N/m^2,$ for example. Summarizing, the exact equation may be approximately solved in an appropriate sense.

6. An exact convex dual variational formulation for a non-convex primal one

In this section we develop a convex dual variational formulation suitable to compute a critical point for the corresponding primal one.

Let $\Omega \subset {\mathbb{R}}^2$ be an open, bounded, connected set with a regular (Lipschitzian) boundary denoted by $\partial \Omega .$

Consider a functional $J:V\to \mathbb{R}$ where

$$J\left(u\right)=F(u_x,u_y)-{\langle u,f\rangle }_{L^2},$$

$V=W_0^{1,2}(\Omega )$ and $f\in L^2(\Omega ).$

Here we denote $Y=Y^{*}=L^2(\Omega )$ and $Y_1=Y_1^{*}=L^2(\Omega )\times L^2(\Omega ).$

Defining

$$V_1={\{u\in V:\parallel u\parallel }_{1,\infty }\le K_1\}$$

for some appropriate $K_1>0$, suppose also F is twice Fréchet differentiable and

$$det\left\{\frac{{\partial }^{2}F(u_x,u_y)}{\partial v_1\partial v_2}\right\}\ne 0,$$

$\forall u\in V_1.$

Define now $F_1:V\to \mathbb{R}$ and $F_2:V\to \mathbb{R}$ by

$$F_1(u_x,u_y)=F(u_x,u_y)+\frac{\epsilon }{2}{\int }_{\Omega }{u}_x^{2}dx+\frac{\epsilon }{2}{\int }_{\Omega }{u}_y^{2}dx,$$

and

$$F_2(u_x,u_y)=\frac{\epsilon }{2}{\int }_{\Omega }{u}_x^{2}dx+\frac{\epsilon }{2}{\int }_{\Omega }{u}_y^{2}dx,$$

where here we denote $dx=d{x}_{1}d{x}_2.$

Moreover, we define the respective Legendre transform functionals $F_1^{*}$ and $F_2^{*}$ as

$$F_1^{*}\left(v^{*}\right)={\langle v_1,v_1^{*}\rangle }_{L^2}+{\langle v_2,v_2^{*}\rangle }_{L^2}-F_1(v_1,v_2),$$

where $v_1,v_2\in Y$ are such that

$$v_1^{*}=\frac{\partial F_1(v_1,v_2)}{\partial v_1},$$

$$v_2^{*}=\frac{\partial F_1(v_1,v_2)}{\partial v_2},$$

and

$$F_2^{*}\left(v^{*}\right)={\langle v_1,v_1^{*}+f_1\rangle }_{L^2}+{\langle v_2,v_2^{*}\rangle }_{L^2}-F_2(v_1,v_2),$$

where $v_1,v_2\in Y$ are such that

$$v_1^{*}+f_1=\frac{\partial F_2(v_1,v_2)}{\partial v_1},$$

$$v_2^{*}=\frac{\partial F_2(v_1,v_2)}{\partial v_2}.$$

Here $f_1$ is any function such that

$${\left(f_1\right)}_x=f,\mathrm{in}\Omega .$$

Furthermore, we define

$$\begin{array}{ccc}\hfill J^{*}\left(v^{*}\right)& =& -F_1^{*}\left(v^{*}\right)+F_2^{*}\left(v^{*}\right)\hfill \\ & =& -F_1^{*}\left(v^{*}\right)+\frac{1}{2\epsilon }{\int }_{\Omega }{(v_1^{*}+f_1)}^{2}dx+\frac{1}{2\epsilon }{\int }_{\Omega }{\left(v_2^{*}\right)}^{2}dx.\hfill \end{array}$$

(24)

Observe that through the target conditions

$$v_1^{*}+f_1=\epsilon u_x,$$

$$v_2^{*}=\epsilon u_y,$$

we may obtain the compatibility condition

$${(v_1^{*}+f_1)}_y-{\left(v_2^{*}\right)}_x=0.$$

Define now

$$A^{*}=\{v^{*}=(v_1^{*},v_2^{*})\in B_r(0,0)\subset Y_1^{*}:{(v_1^{*}+f_1)}_y-{\left(v_2^{*}\right)}_x=0,\mathrm{in}\Omega \},$$

for some appropriate $r>0$ such that $J^{*}$ is convex in $B_r(0,0).$

Consider the problem of minimizing $J^{*}$ subject to $v^{*}\in A^{*}.$

Assuming $r>0$ is large enough so that the restriction in r is not active, at this point we define the associated Lagrangian

$$J_1^{*}(v^{*},\phi )=J^{*}\left(v^{*}\right)+{\langle \phi ,{(v_1^{*}+f)}_y-{\left(v_2^{*}\right)}_x\rangle }_{L^2},$$

where $\phi$ is an appropriate Lagrange multiplier.

Therefore

$$\begin{array}{ccc}\hfill J_1^{*}\left(v^{*}\right)& =& -F_1^{*}\left(v^{*}\right)+\frac{1}{2\epsilon }{\int }_{\Omega }{(v_1^{*}+f_1)}^{2}dx+\frac{1}{2\epsilon }{\int }_{\Omega }{\left(v_2^{*}\right)}^{2}dx\hfill \\ & & +{\langle \phi ,{(v_1^{*}+f)}_y-{\left(v_2^{*}\right)}_x\rangle }_{L^2}.\hfill \end{array}$$

(25)

The optimal point in question will be a solution of the corresponding Euler-Lagrange equations for $J_1^{*}.$

From the variation of $J_1^{*}$ in $v_1^{*}$ we obtain

$$-\frac{\partial F_1^{*}\left(v^{*}\right)}{\partial v_1^{*}}+\frac{v_1^{*}+f}{\epsilon }-\frac{\partial \phi }{\partial y}=0.$$

(26)

From the variation of $J_1^{*}$ in $v_2^{*}$ we obtain

$$-\frac{\partial F_1^{*}\left(v^{*}\right)}{\partial v_2^{*}}+\frac{v_2^{*}}{\epsilon }+\frac{\partial \phi }{\partial x}=0.$$

(27)

From the variation of $J_1^{*}$ in $\phi$ we have

$${(v_1^{*}+f)}_y-{\left(v_2^{*}\right)}_x=0.$$

From this last equation, we may obtain $u\in V$ such that

$$v_1^{*}+f=\epsilon u_x,$$

and

$$v_2^{*}=\epsilon u_y.$$

From this and the previous extremal equations indicated we have

$$-\frac{\partial F_1^{*}\left(v^{*}\right)}{\partial v_1^{*}}+u_x-\frac{\partial \phi }{\partial y}=0,$$

and

$$-\frac{\partial F_1^{*}\left(v^{*}\right)}{\partial v_2^{*}}+u_y+\frac{\partial \phi }{\partial x}=0.$$

so that

$$v_1^{*}+f=\frac{\partial F_1(u_x-{\phi }_y,u_y+{\phi }_x)}{\partial v_1},$$

and

$$v_2^{*}=\frac{\partial F_1(u_x-{\phi }_y,u_y+{\phi }_x)}{\partial v_2}.$$

From this and equation (26) and (27) we have

$$\begin{array}{ccc}& & -\epsilon {\left(\frac{\partial F_1^{*}\left(v^{*}\right)}{\partial v_1^{*}}\right)}_x-\epsilon {\left(\frac{\partial F_1^{*}\left(v^{*}\right)}{\partial v_2^{*}}\right)}_y\hfill \\ & & +{(v_1^{*}+f_1)}_x+{\left(v_2^{*}\right)}_y\hfill \\ & =& -\epsilon u_{xx}-\epsilon u_{yy}+{\left(v_1^{*}\right)}_x+{\left(v_2^{*}\right)}_y+f=0.\hfill \end{array}$$

(28)

Replacing the expressions of $v_1^{*}$ and $v_2^{*}$ into this last equation, we have

$$-\epsilon u_{xx}-\epsilon u_{yy}+{\left(\frac{\partial F_1(u_x-{\phi }_y,u_y+{\phi }_x)}{\partial v_1}\right)}_x+{\left(\frac{\partial F_1(u_x-{\phi }_y,u_y+{\phi }_x)}{\partial v_2}\right)}_y+f=0,$$

so that

$${\left(\frac{\partial F(u_x-{\phi }_y,u_y+{\phi }_x)}{\partial v_1}\right)}_x+{\left(\frac{\partial F(u_x-{\phi }_y,u_y+{\phi }_x)}{\partial v_2}\right)}_y+f=0,\mathrm{in}\Omega .$$

(29)

Observe that if

$${\nabla }^2\phi =0$$

then there exists $\widehat{u}$ such that u and $\phi$ are also such that

$$u_x-{\phi }_y={\widehat{u}}_x$$

and

$$u_y+{\phi }_x={\widehat{u}}_y.$$

The boundary conditions for $\phi$ must be such that $\widehat{u}\in W_0^{1,2}.$

From this and equation (29) we obtain

$$\delta J\left(\widehat{u}\right)=\mathbf{0}.$$

Summarizing, we may obtain a solution $\widehat{u}\in W_0^{1,2}$ of equation $\delta J\left(\widehat{u}\right)=\mathbf{0}$ by minimizing $J^{*}$ on $A^{*}$.

Finally, observe that clearly $J^{*}$ is convex in an appropriate large ball $B_r(0,0)$ for some appropriate $r>0$

7. Another primal dual formulation for a related model

Let $\Omega \subset {\mathbb{R}}^3$ be an open, bounded and connected set with a regular boundary denoted by $\partial \Omega .$

Consider the functional $J:V\to \mathbb{R}$ where

$$\begin{array}{ccc}\hfill J\left(u\right)& =& \frac{\gamma }{2}{\int }_{\Omega }\nabla u·\nabla udx+\frac{\alpha }{2}{\int }_{\Omega }{(u^2-\beta )}^{2}dx\hfill \\ & & -{\langle u,f\rangle }_{L^2},\hfill \end{array}$$

(30)

$\alpha >0,\beta >0,\gamma >0$, $V=W_0^{1,2}(\Omega )$ and $f\in L^2(\Omega ).$

Denoting $Y=Y^{*}=L^2(\Omega )$, define now $J_1^{*}:V\times Y^{*}\to \mathbb{R}$ by

$$\begin{array}{ccc}\hfill J_1^{*}(u,v_0^{*})& =& -\frac{\gamma }{2}{\int }_{\Omega }\nabla u·\nabla udx-{\langle u^2,v_0^{*}\rangle }_{L^2}\hfill \\ & & +\frac{K_1}{2}{\int }_{\Omega }{(-\gamma {\nabla }^{2}u+2v_0^{*}u-f)}^{2}dx+{\langle u,f\rangle }_{L^2}\hfill \\ & & +\frac{1}{2\alpha }{\int }_{\Omega }{\left(v_0^{*}\right)}^{2}dx+\beta {\int }_{\Omega }{v}_0^{*}dx,\hfill \end{array}$$

(31)

Define also

$$A^{+}=\{u\in V:uf\ge 0,\mathrm{a}.\mathrm{e}.\mathrm{in}\Omega \},$$

$$V_2={\{u\in V:\parallel u\parallel }_{\infty }\le K_3\},$$

and

$$V_1=V_2\cap A^{+}$$

for some appropriate $K_3>0$ to be specified.

Moreover define

$$B^{*}=\{v_0^{*}\in Y^{*}:\parallel v_0^{*}{\parallel }_{\infty }\le K\}$$

for some appropriate $K>0$ to be specified.

Observe that, denoting

$$\phi =-\gamma {\nabla }^{2}u+2v_0^{*}u-f$$

we have

$$\frac{{\partial }^2{J}_1^{*}(u,v_0^{*})}{\partial {\left(v_0^{*}\right)}^2}=\frac{1}{\alpha }+4K_1{u}^2$$

$$\frac{{\partial }^2{J}_1^{*}(u,v_0^{*})}{\partial u^2}=\gamma {\nabla }^2-2v_0^{*}+K_1{(-\gamma {\nabla }^2+2v_0^{*})}^2$$

and

$$\frac{{\partial }^2{J}_1^{*}(u,v_0^{*})}{\partial u\partial v_0^{*}}=K_1(2\phi +2(-\gamma {\nabla }^{2}u+2v_0^{*}u))-2u$$

so that

$$\begin{array}{ccc}& & det\left\{{\delta }^2{J}_1^{*}(u,v_0^{*})\right\}\hfill \\ & =& \frac{{\partial }^2{J}_1^{*}(u,v_0^{*})}{\partial {\left(v_0^{*}\right)}^2}\frac{{\partial }^2{J}_1^{*}(u,v_0^{*})}{\partial u^2}-{\left(\frac{{\partial }^2{J}_1^{*}(u,v_0^{*})}{\partial u\partial v_0^{*}}\right)}^2\hfill \\ & =& \frac{K_1{(-\gamma {\nabla }^2+2v_0^{*})}^2}{\alpha }-\frac{\gamma {\nabla }^2+2v_0^{*}+4\alpha u^2}{\alpha }\hfill \\ & & -4K_1^2{\phi }^2-8K_1\phi (-\gamma {\nabla }^2+2v_0^{*})u+8K_1\phi u\hfill \\ & & +4K_1(-\gamma {\nabla }^{2}u+2v_0^{*}u)u.\hfill \end{array}$$

(32)

Observe now that a critical point $\phi =0$ and $(-\gamma {\nabla }^{2}u+2v_0^{*}u)u=fu\ge 0$ in $\Omega$.

Therefore, for an appropriate large $K_1>0$, also at a critical point, we have

$$\begin{array}{ccc}& & det\left\{{\delta }^2{J}_1^{*}(u,v_0^{*})\right\}\hfill \\ & =& 4K_{1}fu-\frac{{\delta }^{2}J\left(u\right)}{\alpha }+K_1\frac{{(-\gamma {\nabla }^2+2v_0^{*})}^2}{\alpha }>\mathbf{0}.\hfill \end{array}$$

(33)

Remark 7.1.

From this last equation we may observe that $J_1^{*}$ has a large region of convexity about any critical point $(u_0,{\widehat{v}}_0^{*})$, that is, there exists a large $r>0$ such that $J_1^{*}$ is convex on $B_r(u_0,{\widehat{v}}_0^{*}).$

With such results in mind, we may easily prove the following theorem.

Theorem 7.2.

Assume $K_1\gg max\{1,K,K_3\}$ and suppose $(u_0,{\widehat{v}}_0^{*})\in V_1\times B^{*}$ is such that

$$\delta J_1^{*}(u_0,{\widehat{v}}_0^{*})=\mathbf{0}.$$

Under such hypotheses, there exists $r>0$ such that $J_1^{*}$ is convex in $E^{*}=B_r(u_0,{\widehat{v}}_0^{*})\cap (V_1\times B^{*})$,

$$\delta J\left(u_0\right)=\mathbf{0},$$

and

$$-J\left(u_0\right)=J_1(u_0,{\widehat{v}}_0^{*})=\underset{(u,v_0^{*})\in E^{*}}{inf}{J}_1^{*}(u,v_0^{*}).$$

8. A third primal dual formulation for a related model

Let $\Omega \subset {\mathbb{R}}^3$ be an open, bounded and connected set with a regular boundary denoted by $\partial \Omega .$

Consider the functional $J:V\to \mathbb{R}$ where

$$\begin{array}{ccc}\hfill J\left(u\right)& =& \frac{\gamma }{2}{\int }_{\Omega }\nabla u·\nabla udx+\frac{\alpha }{2}{\int }_{\Omega }{(u^2-\beta )}^{2}dx\hfill \\ & & -{\langle u,f\rangle }_{L^2},\hfill \end{array}$$

(34)

$\alpha >0,\beta >0,\gamma >0$, $V=W_0^{1,2}(\Omega )$ and $f\in L^2(\Omega ).$

Denoting $Y=Y^{*}=L^2(\Omega )$, define now $J_1^{*}:V\times Y^{*}\times Y^{*}\to \mathbb{R}$ by

$$\begin{array}{ccc}\hfill J_1^{*}(u,v_0^{*},v_1^{*})& =& \frac{\gamma }{2}{\int }_{\Omega }\nabla u·\nabla udx+\frac{1}{2}{\int }_{\Omega }K{u}^{2}dx\hfill \\ & & -{\langle u,v_1^{*}\rangle }_{L^2}+\frac{1}{2}{\int }_{\Omega }\frac{{\left(v_1^{*}\right)}^2}{(-2v_0^{*}+K)}dx\hfill \\ & & +\frac{1}{2(\alpha +\epsilon )}{\int }_{\Omega }{(v_0^{*}-\alpha (u^2-\beta ))}^{2}dx+{\langle u,f\rangle }_{L^2}\hfill \\ & & -\frac{1}{2\alpha }{\int }_{\Omega }{\left(v_0^{*}\right)}^{2}dx-\beta {\int }_{\Omega }{v}_0^{*}dx,\hfill \end{array}$$

(35)

where $\epsilon >0$ is a small real constant.

Define also

$$A^{+}=\{u\in V:uf\ge 0,\mathrm{a}.\mathrm{e}.\mathrm{in}\Omega \},$$

$$V_2={\{u\in V:\parallel u\parallel }_{\infty }\le K_3\},$$

and

$$V_1=V_2\cap A^{+}$$

for some appropriate $K_3>0$ to be specified.

Moreover define

$$B^{*}=\{v_0^{*}\in Y^{*}:\parallel v_0^{*}{\parallel }_{\infty }\le K_4\}$$

and

$$D^{*}=\{v_1^{*}\in Y^{*}:\parallel v_1^{*}\parallel \le K_5\},$$

for some appropriate real constants $K_4,K_5>0$ to be specified.

Remark 8.1.

Define now

$$H_1(u,v_0^{*})=-\gamma {\nabla }^2+2v_0^{*}+4\alpha u^2$$

For an appropriate function (or, in a more general fashion, an appropriate bounded operator) $M_1$ define

$$B_1^{*}=\{v_0^{*}\in B^{*}:2v_0^{*}+M_1\ge {\epsilon }_1\},$$

for some small parameter ${\epsilon }_1>0.$

Moreover, define

$$E^{*}=\{u\in V_1:\sqrt{4\alpha }\left|u\right|\ge \sqrt{|M_1+\gamma {\nabla }^2|}.$$

Since for $(u,v_0^{*})\in V_1\times B_1^{*}$ we have $uf\ge 0,\mathrm{in}\Omega$, so that for $u_1,u_2\in V_1$ we have

$$\mathrm{sign}\left(u_1\right)=\mathrm{sign}\left(u_2\right)\mathrm{in}\Omega ,$$

we may infer that $E^{*}$ is a convex set.

Moreover if $(u,v_0^{*})\in E^{*}\times B_1^{*}$, then

$$\sqrt{4\alpha }\left|u\right|\ge \sqrt{|M_1+\gamma {\nabla }^2|}$$

so that

$$4\alpha u^2\ge M_1+\gamma {\nabla }^2$$

and

$$2v_0^{*}+M_1\ge {\epsilon }_1$$

so that

$$H_1(u,v_0^{*})=-\gamma {\nabla }^2+2v_0^{*}+4\alpha u^2\ge {\epsilon }_1.$$

Such a result we will be used many times in the next sections.

Observe that, defining

$$\phi =v_0^{*}-\alpha (u^2-\beta )$$

we may obtain

$$\frac{{\partial }^2{J}_1^{*}(u,v_0^{*},v_1^{*})}{\partial u^2}=-\gamma {\nabla }^2+K+\frac{\alpha }{\alpha +\epsilon }4u^2-2\phi \frac{\alpha }{\alpha +\epsilon }$$

$$\frac{{\partial }^2{J}_1^{*}(u,v_0^{*},v_1^{*})}{\partial {\left(v_1^{*}\right)}^2}=\frac{1}{-2v_0^{*}+K}$$

and

$$\frac{{\partial }^2{J}_1^{*}(u,v_0^{*},v_1^{*})}{\partial u\partial v_1^{*}}=-1$$

so that

$$\begin{array}{ccc}& & det\left\{\frac{{\partial }^2{J}_1^{*}(u,v_0^{*},v_1^{*})}{\partial u\partial v_1^{*}}\right\}\hfill \\ & =& \frac{{\partial }^2{J}_1^{*}(u,v_1^{*},v_0^{*})}{\partial {\left(v_1^{*}\right)}^2}\frac{{\partial }^2{J}_1^{*}(u,v_1^{*},v_0^{*})}{\partial u^2}-{\left(\frac{{\partial }^2{J}_1^{*}(u,v_1^{*},v_0^{*})}{\partial u\partial v_1^{*}}\right)}^2\hfill \\ & =& \frac{-\gamma {\nabla }^2+2v_0^{*}+4\frac{{\alpha }^2}{\alpha +\epsilon }{u}^2-2\frac{\alpha }{\alpha +\epsilon }\phi }{-2v_0^{*}+K}\hfill \\ & \equiv & H(u,v_0^{*}).\hfill \end{array}$$

(36)

However, at a critical point, we have $\phi =\mathbf{0}$ so that, for a fixed $v_0^{*}\in B^{*}$ we define the non-active but convex restriction

$${\left(C_1\right)}_{v_0^{*}}^{*}=\{u\in V_1:{\left(\phi \right)}^2\le \epsilon \},$$

for a small parameter $\epsilon >0.$

From such results, assuming $K\gg max\{K_3,K_4,K_5\},$ and $0<\epsilon \ll {\epsilon }_1\ll 1$, we have that

$$H(u,v_0^{*})>\mathbf{0},$$

for $v_0^{*}\in B_1^{*}$ and $u\in E^{*}\cap {\left(C_1\right)}_{v_0^{*}}^{*}.$

With such results in mind, we may easily prove the following theorem.

Theorem 8.2.

Suppose $(u_0,{\widehat{v}}_1^{*},{\widehat{v}}_0^{*})\in (E^{*}\cap {\left(C_1\right)}_{{\widehat{v}}_0^{*}}^{*})\times D^{*}\times B_1^{*}$ is such that

$$\delta J_1^{*}(u_0,{\widehat{v}}_1^{*},{\widehat{v}}_0^{*})=\mathbf{0}.$$

Under such hypotheses, we have that

$$\delta J\left(u_0\right)=\mathbf{0}$$

and

$$\begin{array}{ccc}\hfill J\left(u_0\right)& =& \underset{u\in {\left(C_1\right)}_{{\widehat{v}}_0^{*}}^{*}}{inf}J\left(u\right)\hfill \\ & =& J_1^{*}(u_0,{\widehat{v}}_1^{*},{\widehat{v}}_0^{*})\hfill \\ & =& \underset{(u,v_1^{*})\in {\left(C_1\right)}_{{\widehat{v}}_0^{*}}^{*}\times D^{*}}{inf}\left\{\underset{v_0^{*}\in B^{*}}{sup}{J}_1^{*}(u,v_1^{*},v_0^{*})\right\}\hfill \\ & =& \underset{v_0^{*}\in B^{*}}{sup}\left\{\underset{(u,v_1^{*})\in {\left(C_1\right)}_{{\widehat{v}}_0^{*}}^{*}\times D^{*}}{inf}{J}_1^{*}(u,v_1^{*},v_0^{*})\right\}.\hfill \end{array}$$

(37)

Proof. 

The proof that

$$\delta J\left(u_0\right)=\mathbf{0}$$

and

$$J\left(u_0\right)=J_1^{*}(u_0,{\widehat{v}}_1^{*},{\widehat{v}}_0^{*})$$

may be easily made similarly as in the previous sections.

Moreover, observe that for $K>0$ sufficiently large, we have

$$\frac{{\partial }^2{J}_1^{*}(u_0,{\widehat{v}}_1^{*},v_0^{*})}{\partial {\left(v_0^{*}\right)}^2}<\mathbf{0},\forall v_0^{*}\in B^{*}$$

so that this and the other hypotheses, we have also

$$J_1^{*}(u_0,{\widehat{v}}_1^{*},{\widehat{v}}_0^{*})=\underset{(u,v_1^{*})\in {\left(C_1\right)}_{{\widehat{v}}_0^{*}}^{*}\times D^{*}}{inf}{J}_1^{*}(u,v_1^{*},{\widehat{v}}_0^{*})$$

and

$$J_1^{*}(u_0,{\widehat{v}}_1^{*},{\widehat{v}}_0^{*})=\underset{v_0^{*}\in B^{*}}{sup}{J}_1^{*}(u_0,{\widehat{v}}_1^{*},v_0^{*}).$$

From this, from a standard saddle point theorem and the remaining hypotheses, we may infer that

$$\begin{array}{ccc}\hfill J\left(u_0\right)& =& J_1^{*}(u_0,{\widehat{v}}_1^{*},{\widehat{v}}_0^{*})\hfill \\ & =& \underset{(u,v_1^{*})\in {\left(C_1\right)}_{{\widehat{v}}_0^{*}}^{*}\times D^{*}}{inf}\left\{\underset{v_0^{*}\in B^{*}}{sup}{J}_1^{*}(u,v_1^{*},v_0^{*})\right\}\hfill \\ & =& \underset{v_0^{*}\in B^{*}}{sup}\left\{\underset{(u,v_1^{*})\in {\left(C_1\right)}_{{\widehat{v}}_0^{*}}^{*}\times D^{*}}{inf}{J}_1^{*}(u,v_1^{*},v_0^{*})\right\}.\hfill \end{array}$$

(38)

Moreover, observe that

$$\begin{array}{ccc}\hfill J_1^{*}(u_0,{\widehat{v}}_1^{*},{\widehat{v}}_0^{*})& =& \underset{(u,v_1^{*})\in {\left(C_1\right)}_{{\widehat{v}}_0^{*}}^{*}\times D^{*}}{inf}{J}_1^{*}(u,v_1^{*},{\widehat{v}}_0^{*})\hfill \\ & \le & \frac{\gamma }{2}{\int }_{\Omega }\nabla u·\nabla udx+\frac{K}{2}{\int }_{\Omega }{u}^{2}dx\hfill \\ & & +{\langle u^2,{\widehat{v}}_0^{*}\rangle }_{L^2}-\frac{K}{2}{\int }_{\Omega }{u}^{2}dx\hfill \\ & & -\frac{1}{2\alpha }{\int }_{\Omega }{\left({\widehat{v}}_0^{*}\right)}^{2}dx-\beta {\int }_{\Omega }{\widehat{v}}_0^{*}dx\hfill \\ & & +\frac{1}{2(\alpha +\epsilon )}{\int }_{\Omega }{({\widehat{v}}_0^{*}-\alpha (u^2-\beta ))}^{2}dx-{\langle u,f\rangle }_{L^2}\hfill \\ & \le & \underset{v_0^{*}\in Y^{*}}{sup}\left\{\frac{\gamma }{2}{\int }_{\Omega }\nabla u·\nabla udx+\langle u^2,v_0^{*}\rangle \right.\hfill \\ & & -\frac{1}{2\alpha }{\int }_{\Omega }{\left(v_0^{*}\right)}^{2}dx-\beta {\int }_{\Omega }{v}_0^{*}dx\hfill \\ & & \left.+\frac{1}{2(\alpha +\epsilon )}{\int }_{\Omega }{(v_0^{*}-\alpha (u^2-\beta ))}^{2}dx-{\langle u,f\rangle }_{L^2}\right\}\hfill \\ & =& \frac{\gamma }{2}{\int }_{\Omega }\nabla u·\nabla udx+\frac{\alpha }{2}{\int }_{\Omega }{(u^2-\beta )}^{2}dx\hfill \\ & & -{\langle u,f\rangle }_{L^2},\forall u\in {\left(C_1\right)}_{{\widehat{v}}_0^{*}}^{*}.\hfill \end{array}$$

(39)

Summarizing, we have got

$$J\left(u_0\right)=J_1^{*}(u_0,{\widehat{v}}_1^{*},{\widehat{v}}_0^{*})\le \underset{u\in {\left(C_1\right)}_{{\widehat{v}}_0^{*}}^{*}}{inf}J\left(u\right).$$

From such results, we may infer that

$$\begin{array}{ccc}\hfill J\left(u_0\right)& =& \underset{u\in {\left(C_1\right)}_{{\widehat{v}}_0^{*}}^{*}}{inf}J\left(u\right)\hfill \\ & =& J_1^{*}(u_0,{\widehat{v}}_1^{*},{\widehat{v}}_0^{*})\hfill \\ & =& \underset{(u,v_1^{*})\in {\left(C_1\right)}_{{\widehat{v}}_0^{*}}^{*}\times D^{*}}{inf}\left\{\underset{v_0^{*}\in B^{*}}{sup}{J}_1^{*}(u,v_1^{*},v_0^{*})\right\}\hfill \\ & =& \underset{v_0^{*}\in B^{*}}{sup}\left\{\underset{(u,v_1^{*})\in {\left(C_1\right)}_{{\widehat{v}}_0^{*}}^{*}\times D^{*}}{inf}{J}_1^{*}(u,v_1^{*},v_0^{*})\right\}.\hfill \end{array}$$

(40)

The proof is complete. □

9. An algorithm for a related model in shape optimization

The next two subsections have been previously published by Fabio Silva Botelho and Alexandre Molter in [5], Chapter 21.

9.1. Introduction

Consider an elastic solid which the volume corresponds to an open, bounded, connected set, denoted by $\Omega \subset {\mathbb{R}}^3$ with a regular (Lipschitzian) boundary denoted by $\partial \Omega ={\Gamma }_0\cup {\Gamma }_t$ where ${\Gamma }_0\cap {\Gamma }_t=\varnothing .$ Consider also the problem of minimizing the functional $\widehat{J}:U\times B\to \mathbb{R}$ where

$$\widehat{J}(u,t)=\frac{1}{2}{\langle u_i,f_i\rangle }_{L^2(\Omega )}+\frac{1}{2}{\langle u_i,{\widehat{f}}_i\rangle }_{L^2\left({\Gamma }_t\right)},$$

subject to

$$\left\{\begin{array}{c}{\left(H_{ijkl}\left(t\right)e_{kl}\left(u\right)\right)}_{,j}+f_i=0\mathrm{in}\Omega ,\hfill \\ \\ H_{ijkl}\left(t\right)e_{kl}\left(u\right)n_j-{\widehat{f}}_i=0,\mathrm{on}{\Gamma }_t,\forall i\in \{1,2,3\}.\hfill \end{array}\right.$$

(41)

Here $\mathbf{n}=(n_1,n_2,n_3)$ denotes the outward normal to $\partial \Omega$ and

$$U=\{u=(u_1,u_2,u_3)\in W^{1,2}(\Omega ;{\mathbb{R}}^3):u=(0,0,0)=\mathbf{0}\mathrm{on}{\Gamma }_0\},$$

$$B=\left\{t:\Omega \to [0,1]\mathrm{measurable}:{\int }_{\Omega }t\left(x\right)dx=t_1|\Omega |\right\},$$

where

$$0<t_1<1$$

and $|\Omega |$ denotes the Lebesgue measure of $\Omega .$

Moreover $u=(u_1,u_2,u_3)\in W^{1,2}(\Omega ;{\mathbb{R}}^3)$ is the field of displacements relating the cartesian system $(0,x_1,x_2,x_3)$, resulting from the action of the external loads $f\in L^2(\Omega ;{\mathbb{R}}^3)$ and $\widehat{f}\in L^2({\Gamma }_t;{\mathbb{R}}^3).$

We also define the stress tensor $\left\{{\sigma }_{ij}\right\}\in Y^{*}=Y=L^2(\Omega ;{\mathbb{R}}^{3\times 3}),$ by

$${\sigma }_{ij}\left(u\right)=H_{ijkl}\left(t\right)e_{kl}\left(u\right),$$

and the strain tensor $e:U\to L^2(\Omega ;{\mathbb{R}}^{3\times 3})$ by

$$e_{ij}\left(u\right)=\frac{1}{2}(u_{i,j}+u_{j,i}),\forall i,j\in \{1,2,3\}.$$

Finally,

$$\left\{H_{ijkl}\left(t\right)\right\}=\{t{H}_{ijkl}^0+(1-t)H_{ijkl}^1\},$$

where $H^0$ corresponds to a strong material and $H^1$ to a very soft material, intending to simulate voids along the solid structure.

The variable t is the design one, which the optimal distribution values along the structure are intended to minimize its inner work with a volume restriction indicated through the set B.

The duality principle obtained is developed inspired by the works in [2,3]. Similar theoretical results have been developed in [10], however we believe the proof here presented, which is based on the min-max theorem is easier to follow (indeed we thank an anonymous referee for his suggestion about applying the min-max theorem to complete the proof). We highlight throughout this text we have used the standard Einstein sum convention of repeated indices.

Moreover, details on the Sobolev spaces addressed may be found in [1]. In addition, the primal variational development of the topology optimization problem has been described in [10].

The main contributions of this work are to present the detailed development, through duality theory, for such a kind of optimization problems. We emphasize that to avoid the check-board standard and obtain appropriate robust optimized structures without the use of filters, it is necessary to discretize more in the load direction, in which the displacements are much larger.

9.2. Mathematical formulation of the topology optimization problem

Our mathematical topology optimization problem is summarized by the following theorem.

Theorem 9.1.

Consider the statements and assumptions indicated in the last section, in particular those refereing to Ω and the functional $\widehat{J}:U\times B\to \mathbb{R}.$

Define $J_1:U\times B\to \mathbb{R}$ by

$$J_1(u,t)=-G(e\left(u\right),t)+{\langle u_i,f_i\rangle }_{L^2(\Omega )}+{\langle u_i,{\widehat{f}}_i\rangle }_{L^2\left({\Gamma }_t\right)},$$

where

$$G(e\left(u\right),t)=\frac{1}{2}{\int }_{\Omega }{H}_{ijkl}\left(t\right)e_{ij}\left(u\right)e_{kl}\left(u\right)dx,$$

and where

$$dx=d{x}_{1}d{x}_{2}d{x}_3.$$

Define also $J^{*}:U\to \mathbb{R}$ by

$$\begin{array}{ccc}\hfill J^{*}\left(u\right)& =& \underset{t\in B}{inf}\left\{J_1(u,t)\right\}\hfill \\ & =& \underset{t\in B}{inf}\{-G(e\left(u\right),t)+{\langle u_i,f_i\rangle }_{L^2(\Omega )}+{\langle u_i,{\widehat{f}}_i\rangle }_{L^2\left({\Gamma }_t\right)}\}.\hfill \end{array}$$

(42)

Assume there exists $c_0,c_1>0$ such that

$$H_{ijkl}^0{z}_{ij}{z}_{kl}>c_0{z}_{ij}{z}_{ij}$$

and

$$H_{ijkl}^1{z}_{ij}{z}_{kl}>c_1{z}_{ij}{z}_{ij},\forall z=\left\{z_{ij}\right\}\in {\mathbb{R}}^{3\times 3},\mathrm{such}\mathrm{that}z\ne \mathbf{0}.$$

Finally, define $J:U\times B\to \mathbb{R}\cup \{+\infty \}$ by

$$J(u,t)=\widehat{J}(u,t)+Ind(u,t),$$

where

$$Ind(u,t)=\left\{\begin{array}{cc}0,\hfill & \mathrm{if}(u,t)\in A^{*},\hfill \\ +\infty ,\hfill & \mathrm{otherwise},\hfill \end{array}\right.$$

(43)

where $A^{*}=A_1\cap A_2,$

$$A_1=\{(u,t)\in U\times B:{\left({\sigma }_{ij}\left(u\right)\right)}_{,j}+f_i=0,\mathrm{in}\Omega ,\forall i\in \{1,2,3\}\}$$

and

$$A_2=\{(u,t)\in U\times B:{\sigma }_{ij}\left(u\right)n_j-{\widehat{f}}_i=0,\mathrm{on}{\Gamma }_t,\forall i\in \{1,2,3\}\}.$$

Under such hypotheses, there exists $(u_0,t_0)\in U\times B$ such that

$$\begin{array}{ccc}\hfill J(u_0,t_0)& =& \underset{(u,t)\in U\times B}{inf}J(u,t)\hfill \\ & =& \underset{\widehat{u}\in U}{sup}{J}^{*}\left(\widehat{u}\right)\hfill \\ & =& J^{*}\left(u_0\right)\hfill \\ & =& \widehat{J}(u_0,t_0)\hfill \\ & =& \underset{(t,\sigma )\in B\times C^{*}}{inf}{G}^{*}(\sigma ,t)\hfill \\ & =& G^{*}(\sigma \left(u_0\right),t_0),\hfill \end{array}$$

(44)

where

$$\begin{array}{ccc}\hfill G^{*}(\sigma ,t)& =& \underset{v\in Y}{sup}\{{\langle v_{ij},{\sigma }_{ij}\rangle }_{L^2(\Omega )}-G(v,t)\}\hfill \\ & =& \frac{1}{2}{\int }_{\Omega }{\overline{H}}_{ijkl}\left(t\right){\sigma }_{ij}{\sigma }_{kl}dx,\hfill \end{array}$$

(45)

$$\left\{{\overline{H}}_{ijkl}\left(t\right)\right\}={\left\{H_{ijkl}\left(t\right)\right\}}^{-1}$$

and $C^{*}=C_1\cap C_2,$ where

$$C_1=\{\sigma \in Y^{*}:{\sigma }_{ij,j}+f_i=0,\mathrm{in}\Omega ,\forall i\in \{1,2,3\}\}$$

and

$$C_2=\{\sigma \in Y^{*}:{\sigma }_{ij}{n}_j-{\widehat{f}}_i=0,\mathrm{on}{\Gamma }_t,\forall i\in \{1,2,3\}\}.$$

Proof. 

Observe that

$$\begin{array}{ccc}\hfill \underset{(u,t)\in U\times B}{inf}J(u,t)& =& \underset{t\in B}{inf}\left\{\underset{u\in U}{inf}J(u,t)\right\}\hfill \\ & =& \underset{t\in B}{inf}\left\{\underset{\widehat{u}\in U}{sup}\left\{\underset{u\in U}{inf}\left\{\frac{1}{2}{\int }_{\Omega }{H}_{ijkl}\left(t\right)e_{ij}\left(u\right)e_{kl}\left(u\right)dx\right.\right.\right.\hfill \\ & & +{\langle {\widehat{u}}_i,{\left(H_{ijkl}\left(t\right)e_{kl}\left(u\right)\right)}_{,j}+f_i\rangle }_{L^2(\Omega )}\hfill \\ & & \left.\left.\left.-{\langle {\widehat{u}}_i,H_{ijkl}\left(t\right)e_{kl}\left(u\right)n_j-{\widehat{f}}_i\rangle }_{L^2\left({\Gamma }_t\right)}\right\}\right\}\right\}\hfill \\ & =& \underset{t\in B}{inf}\left\{\underset{\widehat{u}\in U}{sup}\left\{\underset{u\in U}{inf}\left\{\frac{1}{2}{\int }_{\Omega }{H}_{ijkl}\left(t\right)e_{ij}\left(u\right)e_{kl}\left(u\right)dx\right.\right.\right.\hfill \\ & & -{\int }_{\Omega }{H}_{ijkl}\left(t\right)e_{ij}\left(\widehat{u}\right)e_{kl}\left(u\right)dx\hfill \\ & & \left.\left.\left.+{\langle {\widehat{u}}_i,f_i\rangle }_{L^2(\Omega )}+{\langle {\widehat{u}}_i,{\widehat{f}}_i\rangle }_{L^2\left({\Gamma }_t\right)}\right\}\right\}\right\}\hfill \\ & =& \underset{t\in B}{inf}\left\{\underset{\widehat{u}\in U}{sup}\left\{-{\int }_{\Omega }{H}_{ijkl}\left(t\right)e_{ij}\left(\widehat{u}\right)e_{kl}\left(\widehat{u}\right)dx\right.\right.\hfill \\ & & \left.\left.{\langle {\widehat{u}}_i,f_i\rangle }_{L^2(\Omega )}+{\langle {\widehat{u}}_i,{\widehat{f}}_i\rangle }_{L^2\left({\Gamma }_t\right)}\right\}\right\}\hfill \\ & =& \underset{t\in B}{inf}\left\{\underset{\sigma \in C^{*}}{inf}{G}^{*}(\sigma ,t)\right\}.\hfill \end{array}$$

(46)

Also, from this and the min-max theorem, there exist $(u_0,t_0)\in U\times B$ such that

$$\begin{array}{ccc}\hfill \underset{(u,t)\in U\times B}{inf}J(u,t)& =& \underset{t\in B}{inf}\left\{\underset{\widehat{u}\in U}{sup}{J}_1(u,t)\right\}\hfill \\ & =& \underset{u\in U}{sup}\left\{\underset{t\in B}{inf}{J}_1(u,t)\right\}\hfill \\ & =& J_1(u_0,t_0)\hfill \\ & =& \underset{t\in B}{inf}{J}_1(u_0,t)\hfill \\ & =& J^{*}\left(u_0\right).\hfill \end{array}$$

(47)

Finally, from the extremal necessary condition

$$\frac{\partial J_1(u_0,t_0)}{\partial u}=\mathbf{0}$$

we obtain

$${\left(H_{ijkl}\left(t_0\right)e_{kl}\left(u_0\right)\right)}_{,j}+f_i=0\mathrm{in}\Omega ,$$

and

$$H_{ijkl}\left(t_0\right)e_{kl}\left(u_0\right)n_j-{\widehat{f}}_i=0\mathrm{on}{\Gamma }_t,\forall i\in \{1,2,3\},$$

so that

$$G\left(e\left(u_0\right)\right)=\frac{1}{2}{\langle {\left(u_0\right)}_i,f_i\rangle }_{L^2(\Omega )}+\frac{1}{2}{\langle {\left(u_0\right)}_i,{\widehat{f}}_i\rangle }_{L^2\left({\Gamma }_t\right)}.$$

Hence $(u_0,t_0)\in A^{*}$ so that $Ind(u_0,t_0)=0$ and $\sigma \left(u_0\right)\in C^{*}.$

Moreover

$$\begin{array}{ccc}\hfill J^{*}\left(u_0\right)& =& -G\left(e\left(u_0\right)\right)+{\langle {\left(u_0\right)}_i,f_i\rangle }_{L^2(\Omega )}+{\langle {\left(u_0\right)}_i,{\widehat{f}}_i\rangle }_{L^2\left({\Gamma }_t\right)}\hfill \\ & =& G\left(e\left(u_0\right)\right)\hfill \\ & =& G\left(e\left(u_0\right)\right)+Ind(u_0,t_0)\hfill \\ & =& J(u_0,t_0)\hfill \\ & =& G^{*}(\sigma \left(u_0\right),t_0).\hfill \end{array}$$

(48)

This completes the proof. □

9.3. About a concerning algorithm and related numerical method

For numerically solve this optimization problem in question, we present the following algorithm

• Set $t_1=0.5\mathrm{in}\Omega$ and $n=1$.
• Calculate $u_n\in U$ such that

  $$J_1(u_n,t_n)=\underset{u\in U}{sup}{J}_1(u,t_n).$$
• Calculate $t_{n+1}\in B$ such that

  $$J_1(u_n,t_{n+1})=\underset{t\in B}{inf}{J}_1(u_n,t).$$
• If $\parallel t_{n+1}-t_n{\parallel }_{\infty }<{10}^{-4}$ or $n>100$ then stop, else set $n:=n+1$ and go to item b.

We have developed a software in finite differences for solving such a problem.

Here the software.

**************************************

• clear all
  global P m8 d w u v Ea Eb Lo d1 z1 m9 du1 du2 dv1 dv2 c3
  m8=27;
  m9=24;
  c3=0.95;
  d=1.0/m8;
  d1=0.5/m9;
  Ea=$210*{10}^5$; (stronger material)
  Eb=1000; (softer material simulating voids)
  w=0.30;
  P=-42000000;
  z1=(m8-1)*(m9-1);
  A3=zeros(z1,z1);
  for i=1:z1
  A3(1,i)=1.0;
  end;
  b=zeros(z1,1);
  uo=0.000001*ones(z1,1);
  u1=ones(z1,1);
  b(1,1)=c3*z1;
  for i=1:m9-1
  for j=1:m8-1
  Lo(i,j)=c3;
  end; end;
  for i=1:z1
  x1(i)=c3*z1;
  end;
  for i=1:2*m8*m9
  xo(i)=0.000;
  end;
  xw=xo;
  xv=Lo;
  for k2=1:24
  c3=0.98*c3;
  b(1,1)=c3*z1;
  k2
  b14=1.0;
  k3=0;
  while $(b14>{10}^{-3.5})$ and $(k3<5)$
  k3=k3+1;
  b12=1.0;
  k=0;
  while $(b12>{10}^{-4.0})$ and $(k<120)$
  k=k+1;
  k2
  k3
  k
  X=fminunc(’funbeam’,xo);
  xo=X;
  b12=max(abs(xw-xo));
  xw=X;
  end;
  for i=1:m9-1
  for j=1:m8-1
  $E1=Lo{(i,j)}^2*(Ea-Eb);$
  ex=du1(i,j);
  ey=dv2(i,j);
  exy=1/2*(dv1(i,j)+du2(i,j));
  $Sx=E1*(ex+w*ey)/(1-w^2);$
  $Sy=E1*(w*ex+ey)/(1-w^2);$
  Sxy=E1/(2*(1+w))*exy;
  dc3(i,j)=-(Sx*ex+Sy*ey+2*Sxy*exy);
  end;
  end;
  for i=1:m9-1
  for j=1:m8-1
  f(j+(i-1)*(m8-1))=dc3(i,j);
  end;
  end;
  for k1=1:1
  k1
  X1=linprog(f,$\left[\right]$,$\left[\right]$,A3,b,uo,u1,x1);
  x1=X1;
  end;
  for i=1:m9-1
  for j=1:m8-1
  Lo(i,j)=X1(j+(m8-1)*(i-1));
  end;
  end;
  b14=max(max(abs(Lo-xv)))
  xv=Lo;
  colormap(gray); imagesc(-Lo); axis equal; axis tight; axis off;pause(1e-6)
  end;
  end;

****************************************************

Here the auxiliary Function ’funbeam’

function S=funbeam(x)

global P m8 d w u v Ea Eb Lo d1 m9 du1 du2 dv1 dv2

for i=1:m9

for j=1:m8

u(i,j)=x(j+(m8)*(i-1));

v(i,j)=x(m8*m9+(i-1)*m8+j);

end;

end;

for i=1:m9

end;

u(m9-1,1)=0;

v(m9-1,1)=0;

u(m9-1,m8-1)=0;

v(m9-1,m8-1)=0;

for i=1:m9-1

for j=1:m8-1

du1(i,j)=(u(i,j+1)-u(i,j))/d;

du2(i,j)=(u(i+1,j)-u(i,j))/d1;

dv1(i,j)=(v(i,j+1)-v(i,j))/d;

dv2(i,j)=(v(i+1,j)-v(i,j))/d1;

end;

end;

S=0;

for i=1:m9-1

for j=1:m8-1

$E1=Lo{(i,j)}^3*Ea+(1-Lo{(i,j)}^3)*Eb;$

ex=du1(i,j);

ey=dv2(i,j);

exy=1/2*(dv1(i,j)+du2(i,j));

$Sx=E1*(ex+w*ey)/(1-w^2);$

$Sy=E1*(w*ex+ey)/(1-w^2);$

Sxy=E1/(2*(1+w))*exy;

S=S+1/2*(Sx*ex+Sy*ey+2*Sxy*exy);

end;

end;

S=S*d*d1-P*v(2,(m8)/3)*d*d1;

***********************************************

For a two dimensional beam of dimensions $1m\times 0.5m$ and $t_1=0.63$ we have obtained the following results:

• Case A: For the optimal shape for a clamped beam at left (cantilever) and load $P=-4{10}^{6}N\mathbf{j}$ at $(x,y)=(1,0.25)$, please Figure 5.
• Case B :For the optimal shape for a simply supported beam at $(0,0)$ and $(1,0)$ and load $P=-4{10}^{6}N\mathbf{j}$ at $(x,y)=(1/3,0.5)$, please Figure 6.
  In the first case the mesh was $28\times 24$. In the second one the mesh was $27\times 24$

10. A duality principle for a general vectorial case in the calculus of variations

In this section we develop a duality principle for a general vectorial case in variational optimization.

Let $\Omega \subset {\mathbb{R}}^3$ be an open, bounded and connected set with a regular (Lipschitzian) boundary denoted by $\partial \Omega$. Let $J:V\to \mathbb{R}$ be a functional where

$$J\left(u\right)=G(\nabla u_1,\cdots ,\nabla u_N)-{\langle u,f\rangle }_{L^2},$$

where

$$V=W_0^{1,2}(\Omega ;{\mathbb{R}}^N)$$

and

$$f=(f_1,\cdots ,f_N)\in L^2(\Omega ;{\mathbb{R}}^N).$$

Here we have denoted $u=(u_1,\cdots ,u_N)\in V$ and

$${\langle u,f\rangle }_{L^2}={\langle u_i,f_i\rangle }_{L^2},$$

so that we may also denote

$$J\left(u\right)=G(\nabla u)-{\langle u,f\rangle }_{L^2}.$$

Assume

$$G(\nabla u)={\int }_{\Omega }g(\nabla u)dx$$

where $g:{\mathbb{R}}^{3N}\to \mathbb{R}$ is a differentiable function such that

$$g\left(y\right)\to +\infty$$

as $\left|y\right|\to \infty$. Moreover, suppose there exists $\alpha \in \mathbb{R}$ such that

$$\alpha =\underset{u\in V}{inf}J\left(u\right).$$

It is well known that

$$\begin{array}{ccc}\hfill \alpha & =& \underset{u\in V}{inf}J\left(u\right)\hfill \\ & =& \underset{u\in V}{inf}{J}^{**}\left(u\right)\hfill \\ & =& \underset{u\in V}{inf}\{{(G\circ \nabla )}^{**}\left(u\right)-{\langle u,f\rangle }_{L^2}\}.\hfill \end{array}$$

(49)

Under some mild hypotheses, from convexity, we have that

$$\begin{array}{ccc}& & \underset{u\in V}{inf}\{{(G\circ \nabla )}^{**}\left(u\right)-{\langle u,f\rangle }_{L^2}\}\hfill \\ & =& \underset{v^{*}\in A^{*}}{sup}\{-{(G\circ \nabla )}^{*}(-div{v}^{*})\}=-{(G\circ \nabla )}^{*}\left(f\right),\hfill \end{array}$$

(50)

where

$$A^{*}=\{v^{*}\in Y=Y^{*}=L^2(\Omega ;{\mathbb{R}}^{3N}):div{v}^{*}+f=0\}.$$

Now observe that the restriction $v=\nabla u$ for some $u\in V$ is equivalent to the restriction

$$\mathrm{curl}{v}_i=\mathbf{0},\mathrm{in}\Omega$$

where $v=\left\{v_i\right\}={\left\{v_{ij}\right\}}_{j=1}^3$, $\forall i\in \{1,\cdots ,N\},$ with appropriate boundary conditions, so that with an appropriate Lagrange multiplier $\varphi =\left\{{\varphi }_i\right\}$, we obtain

$$\begin{array}{ccc}\hfill {(G\circ \nabla )}^{*}(-div{v}^{*})& =& \underset{u\in V}{sup}\{{\langle u,-div{v}^{*}\rangle }_{L^2}-G(\nabla u)\}\hfill \\ & =& \underset{u\in V}{sup}\{{\langle \nabla u,v^{*}\rangle }_{L^2}-G(\nabla u)\}\hfill \\ & \le & \underset{\varphi \in Y^{*}}{inf}\left\{\underset{v\in Y}{sup}\{{\langle v,v^{*}\rangle }_{L^2}-G\left(v\right)+{\langle \varphi ,\mathrm{curl}v\rangle }_{L^2}\right\}\hfill \\ & =& \underset{\varphi \in Y^{*}}{inf}{G}^{*}(v^{*}+\mathrm{curl}\varphi ).\hfill \end{array}$$

(51)

where we have denoted

$$\mathrm{curl}v=\left\{\mathrm{curl}{v}_i\right\}$$

and

$$\mathrm{curl}\varphi =\left\{\mathrm{curl}{\varphi }_i\right\}.$$

Joining the pieces, we have got

$$\begin{array}{ccc}\hfill \underset{u\in V}{inf}J\left(u\right)& =& \underset{u\in V}{inf}\{G(\nabla u)-{\langle u,f\rangle }_{L^2}\}\hfill \\ & \ge & \underset{(v^{*},\varphi )\in A^{*}\times Y^{*}}{sup}\{-G^{*}(v^{*}+\mathrm{curl}\varphi )\},\hfill \end{array}$$

(52)

where we recall that $Y=Y^{*}=L^2(\Omega ;{\mathbb{R}}^{3N}).$

We emphasize such a dual formulation in $(v^{*},\varphi )$ is convex (in fact concave).

11. A note on the Galerkin Functional

Let $\Omega \subset {\mathbb{R}}^3$ be an open, bounded and connected set with a regular (Lipschitzian) boundary denoted by $\partial \Omega$.

Consider the functional $J:V\to \mathbb{R}$ where

$$\begin{array}{ccc}\hfill J\left(u\right)& =& \frac{\gamma }{2}{\int }_{\Omega }\nabla u·\nabla udx+\frac{\alpha }{4}{\int }_{\Omega }{u}^{4}dx\hfill \\ & & -\frac{\beta }{2}{\int }_{\Omega }{u}^{2}dx-{\langle u,f\rangle }_{L^2}\hfill \end{array}$$

(53)

Here $V=W_0^{1,2}(\Omega )$, $\gamma >0,\alpha >0,\beta >0.$

We denote also

$$Y=Y^{*}=L^2(\Omega ).$$

At this point we define

$$A^{+}=\{u\in V:uf\ge 0,\mathrm{in}\Omega \},$$

$$V_2={\{u\in V:\parallel u\parallel }_{\infty }\le K_3\},$$

for some appropriate real constant $K_3>0$ and

$$V_1=A^{+}\cap V_2.$$

Observe that

$$J^{\prime }\left(u\right)=-\gamma {\nabla }^{2}u+\alpha u^3-\beta -f,$$

so that we define the Galerkin functional $J_1:V\to \mathbb{R}$ by

$$J_1\left(u\right)=\frac{1}{2}{\parallel J^{\prime }\left(u\right)\parallel }_2^2=\frac{1}{2}{\int }_{\Omega }{(-\gamma {\nabla }^{2}u+\alpha u^3-\beta u-f)}^{2}dx.$$

From this, we get

$$\begin{array}{ccc}\hfill \frac{{\partial }^2{J}_1\left(u\right)}{\partial u^2}& =& (-\gamma \nabla u+\alpha u^3-\beta u-f)6\alpha u\hfill \\ & & +{(-\gamma {\nabla }^2+3\alpha u^2-\beta )}^2.\hfill \end{array}$$

(54)

Define now

$${\phi }_2={(-\gamma {\nabla }^{2}u+\alpha u^3-\beta u-f)}^2.$$

At this point, for an appropriate small real constant ${\epsilon }_1>0$ and bounded constant operator $M_1>{\epsilon }_1$, we set the intended non-active restriction

$$\sqrt{3\alpha }\left|u\right|\ge \sqrt{|M_1+\gamma {\nabla }^2+\beta |},$$

and define

$$B_1=\{u\in V_1:\sqrt{3\alpha }\left|u\right|\ge \sqrt{|M_1+\gamma {\nabla }^2+\beta |}\}.$$

Observe that since for $u\in V_1$ we have $uf\ge 0$ in $\Omega$ so that if $u_1,u_2\in V_1$ then

$$\mathrm{sign}\left(u_1\right)=\mathrm{sign}\left(u_2\right),\mathrm{in}\Omega ,$$

we may infer that $B_1$ is a convex set.

Furthermore, if $u\in B_1$, then

$$\sqrt{3\alpha }\left|u\right|\ge \sqrt{|M_1+\gamma {\nabla }^2+\beta |},$$

so that

$$3\alpha u^2\ge M_1+\gamma {\nabla }^2+\beta ,$$

and hence

$${\delta }^{2}J\left(u\right)=-\gamma {\nabla }^2+3\alpha u^2-\beta \ge M_1>{\epsilon }_1>\mathbf{0}.$$

For a small parameter $\epsilon >0$ we define the intended non-active restriction

$${\phi }_2\le \epsilon ,\mathrm{in}\Omega ,$$

and define

$$B_2=\{u\in V_1:{\phi }_2\le \epsilon ,\mathrm{in}\Omega \}.$$

Observe that for $\alpha >0$ and $\beta >0$ sufficiently large ${\phi }_2$ is convex in $V_1$ (positive definite Hessian) so that $B_2$ is a convex set. Assuming $0<\epsilon \ll {\epsilon }_1\ll 1$, define $B_3=B_1\cap B_2$, which is a convex set.

Summarizing, if $u\in B_3$, then

$${\delta }^2{J}_1\left(u\right)\ge \mathbf{0}.$$

With such results in mind, we define the following convex optimization problem for finding a critical point of J.

Minimize

$$J_1\left(u\right)=\frac{1}{2}{\parallel J^{\prime }\left(u\right)\parallel }_2^2=\frac{1}{2}{\int }_{\Omega }{(-\gamma {\nabla }^{2}u+\alpha u^3-\beta u-f)}^{2}dx,$$

subject to

$$u\in B_3.$$

Observe that a critical point $u_0\in B_3$ of $J_1$, from such a concerning convexity of $J_1$ on the convex set $B_1$, is also such that

$$J\left(u_0\right)=\underset{u\in B_3}{min}{J}_1\left(u\right).$$

Finally, we may also define the convex optimization problem of minimizing

$$\begin{array}{ccc}\hfill J_3\left(u\right)& =& K_1{J}_1\left(u\right)+J\left(u\right)\hfill \\ & =& \frac{K_1}{2}{\int }_{\Omega }{(-\gamma {\nabla }^{2}u+\alpha u^3-\beta u-f)}^{2}dx\hfill \\ & & +\frac{\gamma }{2}{\int }_{\Omega }\nabla u·\nabla udx+\frac{\alpha }{4}{\int }_{\Omega }{u}^{4}dx\hfill \\ & & -\frac{\beta }{2}{\int }_{\Omega }{u}^{2}dx-{\langle u,f\rangle }_{L^2},\hfill \end{array}$$

(55)

subject to

$$u\in B_3.$$

Here $K_1>0$ is a large real constant.

Such a functional $J_3$ is also convex on $B_3$ so that a critical point $u_0\in B_3$ of J is also a critical point of $J_3$, and thus

$$J_3\left(u_0\right)=\underset{u\in B_3}{min}{J}_3\left(u\right).$$

12. A note on the Legendre-Galerkin functional

Let $\Omega \subset {\mathbb{R}}^3$ be an open, bounded and connected set with a regular (Lipschitzian) boundary denoted by $\partial \Omega$.

Consider the functional $J:V\to \mathbb{R}$ where

$$\begin{array}{ccc}\hfill J\left(u\right)& =& \frac{\gamma }{2}{\int }_{\Omega }\nabla u·\nabla udx+\frac{\alpha }{4}{\int }_{\Omega }{u}^{4}dx\hfill \\ & & -\frac{\beta }{2}{\int }_{\Omega }{u}^{2}dx-{\langle u,f\rangle }_{L^2}\hfill \end{array}$$

(56)

Here $V=W_0^{1,2}(\Omega )$, $\gamma >0,\alpha >0,\beta >0.$

We denote also

$$Y=Y^{*}=L^2(\Omega )$$

and $F_1:V\to \mathbb{R}$, $F_2:V\to \mathbb{R}$ and $F_3:V\to \mathbb{R}$ by

$$F_1\left(u\right)=\frac{\gamma }{2}{\int }_{\Omega }\nabla u·\nabla udx,$$

$$F_2\left(u\right)=\frac{\alpha }{4}{\int }_{\Omega }{u}^{4}dx,$$

$$F_3\left(u\right)=\frac{\beta }{2}{\int }_{\Omega }{u}^{2}dx.$$

Moreover, we define $F_1^{*},F_2^{*},F_3^{*}:Y^{*}\to \mathbb{R}$ by

$$\begin{array}{ccc}\hfill F_1^{*}\left(v_1^{*}\right)& =& \underset{u\in V}{sup}\{{\langle u,v_1^{*}\rangle }_{L^2}-F_1\left(u\right)\}\hfill \\ & =& \frac{1}{2}{\int }_{\Omega }\frac{{\left(v_1^{*}\right)}^2}{-\gamma {\nabla }^2}dx,\hfill \end{array}$$

(57)

$$\begin{array}{ccc}\hfill F_2^{*}\left(v_2^{*}\right)& =& \underset{u\in V}{sup}\{{\langle u,v_2^{*}\rangle }_{L^2}-F_2\left(u\right)\}\hfill \\ & =& \frac{3}{4}{\int }_{\Omega }\frac{{\left(v_2^{*}\right)}^{4/3}}{{\alpha }^{1/3}}dx,\hfill \end{array}$$

(58)

$$\begin{array}{ccc}\hfill F_3^{*}\left(v_3^{*}\right)& =& \underset{u\in V}{sup}\{{\langle u,v_3^{*}\rangle }_{L^2}-F_3\left(u\right)\}\hfill \\ & =& \frac{1}{2\beta }{\int }_{\Omega }{\left(v_3^{*}\right)}^{2}dx.\hfill \end{array}$$

(59)

Observe now that these three last suprema are attained through the equations,

$$v_1^{*}=\frac{\partial F_1\left(u\right)}{\partial u}=-\gamma {\nabla }^{2}u,$$

$$v_2^{*}=\frac{\partial F_2\left(u\right)}{\partial u}=\alpha u^3$$

$$v_3^{*}=\frac{\partial F_3\left(u\right)}{\partial u}=\beta u.$$

From such results, at a critical point, we obtain the following compatibility conditions

$$u=\frac{v_1^{*}}{-\gamma {\nabla }^2}={\left(\frac{v_2^{*}}{\beta }\right)}^{1/3}=\frac{v_3^{*}}{\beta }.$$

From such relations we have

$$\frac{v_1^{*}}{-\gamma {\nabla }^2}=\frac{v_3^{*}}{\beta },$$

and

$$v_2^{*}=\alpha {\left(\frac{v_3^{*}}{\beta }\right)}^3,$$

so that

$$v_1^{*}=-\gamma {\nabla }^2\left(\frac{v_3^{*}}{\beta }\right),$$

and

$$v_2^{*}=\alpha {\left(\frac{v_3^{*}}{\beta }\right)}^3.$$

Moreover, we define the functional $F_4^{*}:Y^{*}\to \mathbb{R},$ by

$$F_4^{*}\left(v^{*}\right)=\underset{u\in V}{sup}\{{\langle u,v_1^{*}+v_2^{*}-v_3^{*}\rangle }_{L^2}-{\langle u,f\rangle }_{L^2}\}.$$

Therefore

$$F_4^{*}\left(v^{*}\right)=\left\{\begin{array}{cc}0,\hfill & \mathrm{if}{v}_1^{*}+v_2^{*}-v_3^{*}-f=0,\mathrm{in}\Omega ,\hfill \\ +\infty ,\hfill & \mathrm{otherwise}.\hfill \end{array}\right.$$

(60)

Hence, a critical point of J corresponds to the solution of the following system of equations

$$v_1^{*}=-\gamma {\nabla }^2\left(\frac{v_3^{*}}{\beta }\right),$$

$$v_2^{*}=\alpha {\left(\frac{v_3^{*}}{\beta }\right)}^3,$$

and

$$v_1^{*}+v_2^{*}-v_3^{*}-f=0,\mathrm{in}\Omega .$$

From this last equation we may obtain

$$v_1^{*}=-v_2^{*}+v_3^{*}+f,$$

so that the final equations to be solved are

$$-v_2^{*}+v_3^{*}+f+\gamma {\nabla }^2\left(\frac{v_3^{*}}{\beta }\right)=0$$

and

$$v_2^{*}-\alpha {\left(\frac{v_3^{*}}{\beta }\right)}^3=0,\mathrm{in}\Omega ,$$

with the boundary conditions

$$u=\frac{v_3^{*}}{\beta }=0,\mathrm{on}\partial \Omega .$$

With such results in mind, we define the Legendre-Galerkin functional $J^{*}:{\left[Y^{*}\right]}^2\to \mathbb{R}$, where

$$\begin{array}{ccc}\hfill J^{*}\left(v^{*}\right)& =& \frac{1}{2}{\int }_{\Omega }{\left(-v_2^{*}+v_3^{*}+f+\frac{\gamma {\nabla }^2{v}_3^{*}}{\beta }\right)}^{2}dx\hfill \\ & & +\frac{1}{2}{\int }_{\Omega }{\left(v_2^{*}-\alpha {\left(\frac{v_3^{*}}{\beta }\right)}^3\right)}^{2}dx.\hfill \end{array}$$

(61)

At this point, defining

$$\phi =v_2^{*}-\alpha {\left(\frac{v_3^{*}}{\beta }\right)}^3,$$

we obtain

$$\frac{{\partial }^2{J}^{*}\left(v^{*}\right)}{\partial {\left(v_2^{*}\right)}^2}=2;$$

$$\frac{{\partial }^2{J}^{*}\left(v^{*}\right)}{\partial {\left(v_3^{*}\right)}^2}={\left(-1-\frac{\gamma {\nabla }^2}{\beta }\right)}^2+\frac{9{\alpha }^2{\left(v_3^{*}\right)}^4}{{\beta }^6}+\mathcal{O}\left(\phi \right),$$

$$\frac{{\partial }^2{J}^{*}\left(v^{*}\right)}{\partial v_2^{*}\partial v_3^{*}}=\frac{-3\alpha {\left(v_3^{*}\right)}^2}{{\beta }^3}+\left(-1-\frac{\gamma {\nabla }^2}{\beta }\right).$$

From such results we may infer that

$$\begin{array}{ccc}\hfill det\left(\frac{{\partial }^2{J}^{*}\left(v^{*}\right)}{\partial v_2^{*}\partial v_3^{*}}\right)& =& \frac{{\partial }^2{J}^{*}\left(v^{*}\right)}{\partial {\left(v_2^{*}\right)}^2}\frac{{\partial }^2{J}^{*}\left(v^{*}\right)}{\partial {\left(v_3^{*}\right)}^2}-{\left(\frac{{\partial }^2{J}^{*}\left(v^{*}\right)}{\partial v_2^{*}\partial v_3^{*}}\right)}^2\hfill \\ & =& {\left(-1-\frac{\gamma {\nabla }^2}{\beta }+3\alpha \frac{{\left(v_3^{*}\right)}^2}{{\beta }^3}\right)}^2+\mathcal{O}\left(\phi \right)\hfill \end{array}$$

(62)

Observe that a critical point $\phi =0$ so that ${\delta }^2{J}^{*}\left(v^{*}\right)>\mathbf{0}$ at a neighborhood of any critical point.

At this point we define

$$A^{+}=\left\{v^{*}=(v_2^{*},v_3^{*})\in {\left[Y^{*}\right]}^2:\frac{v_3^{*}}{\beta }f\ge 0,\mathrm{in}\Omega \right\},$$

$$D^{*}=\{v^{*}=(v_2^{*},v_3^{*})\in {\left[Y^{*}\right]}^2:\parallel v^{*}{\parallel }_{\infty }\le K\},$$

for an appropriate real constant $K>0$.

Define now $E^{*}=A^{+}\cap D^{*}$,

$$C_1^{*}=\{v^{*}=(v_2^{*},v_3^{*})\in E^{*}:{\phi }^2\le \epsilon ,\mathrm{in}\Omega \},$$

for a small real constant $\epsilon >0,$

$$C_2^{*}=\left\{v^{*}=(v_2^{*},v_3^{*})\in E^{*}:\left(-1-\frac{\gamma {\nabla }^2}{\beta }+3\alpha \frac{{\left(v_3^{*}\right)}^2}{{\beta }^3}\right)\ge {\epsilon }_1\right\},$$

and

$$C^{*}=C_1^{*}\cap C_2^{*}.$$

Similarly as done in the previous section, we may prove that $C^{*}$ is a convex set.

Furthermore, for $0<\epsilon \ll {\epsilon }_1\ll 1,$ we have that $J^{*}$ is convex on $C^{*}$.

Summarizing, we may define the following convex optimization problem to obtain a critical point of the primal functional J,

$$\mathrm{Minimize}{J}^{*}(v_2^{*},v_3^{*})\mathrm{subject}\mathrm{to}{v}^{*}=(v_2^{*},v_3^{*})\in C^{*}.$$

We call $J^{*}$ the Legendre-Galerkin functional associated to J.

12.1. Numerical examples

We have obtained numerical solutions for two one-dimensional examples.

• For $\gamma =1.0,$$\alpha =3.0,$$\beta =30.0,$$f\equiv 10,\mathrm{in}\Omega =[0,1].$
  For the respective solution please see Figure 7.
• For $\gamma =0.01,$$\alpha =3.0,$$\beta =30.0,$$f\equiv 10,\mathrm{in}\Omega =[0,1].$
  For the respective solution please see Figure 8.

13. A general concave dual variational formulation for global optimization

Let $\Omega \subset {\mathbb{R}}^3$ be an open, bounded and connected set a regular (Lipschitzian) boundary denoted by $\partial \Omega .$

Consider a functional $J:V\to \mathbb{R}$ where

$$J\left(u\right)=G\left(u\right)-{\langle u,f\rangle }_{L^2},\forall u\in V.$$

Here $V=W_0^{1,2}(\Omega )$, $f\in L^2(\Omega )$ and we also denote $Y=Y^{*}=L^2(\Omega ).$

Assume there exists $\alpha \in \mathbb{R}$ such that

$$\alpha =\underset{u\in V}{inf}J\left(u\right).$$

Furthermore, suppose G is three times Fréchet differentiable and there exists $K>0$ such that

$$\frac{{\partial }^{2}G\left(u\right)}{\partial u^2}+K>\mathbf{0},\forall u\in V.$$

Define now $J_1:V\times Y\to \mathbb{R}$ where,

$$J_1(u,v)=G_1(u,v)+F\left(u\right),$$

where

$$G_1(u,v)=G\left(v\right)-\frac{\epsilon }{2}{\int }_{\Omega }{v}^{2}dx+\frac{K}{2}{\int }_{\Omega }{(v-u)}^{2}dx,$$

and

$$F\left(u\right)=\frac{\epsilon }{2}{\int }_{\Omega }{u}^{2}dx-{\langle u,f\rangle }_{L^2}.$$

Moreover, we define the polar functionals $G_1^{*}:Y^{*}\times V\to \mathbb{R}$ and $F^{*}:Y^{*}\to \mathbb{R}$, where

$$\begin{array}{ccc}\hfill G_1^{*}(v^{*},u)& =& \underset{v\in Y}{sup}\{{\langle v,v^{*}\rangle }_{L^2}-G_1(u,v)\}\hfill \\ & =& -G_{K_{\epsilon }}^{*}(v^{*}+Ku)+\frac{K}{2}{\int }_{\Omega }{u}^{2}dx,\hfill \end{array}$$

(63)

$$G_{K_{\epsilon }}^{*}(v^{*}+Ku)=\underset{v\in Y}{sup}\left\{{\langle v,v^{*}\rangle }_{L^2}-G\left(v\right)-\frac{K}{2}{\int }_{\Omega }{v}^{2}dx+\frac{\epsilon }{2}{\int }_{\Omega }{v}^{2}dx\right\},$$

and

$$\begin{array}{ccc}\hfill F^{*}(-v^{*})& =& \underset{u\in V}{sup}\{-{\langle u,v^{*}\rangle }_{L^2}-F\left(u\right)\}\hfill \\ & =& \frac{1}{2\epsilon }{\int }_{\Omega }{(v^{*}-f)}^{2}dx.\hfill \end{array}$$

(64)

At this point we define the functional $J_2^{*}:Y^{*}\times V\to \mathbb{R}$ by

$$J_2^{*}(v^{*},u)=-G_{K_{\epsilon }}^{*}(v^{*}+Ku)+\frac{K}{2}{\int }_{\Omega }{u}^{2}dx-F^{*}(-v^{*}).$$

With such results in mind we define

$$V_1={\{u\in V:\parallel u\parallel }_{\infty }\le K_3\},$$

and

$$D^{*}=\{v^{*}\in Y^{*}:\parallel v^{*}{\parallel }_{\infty }\le K_4\},$$

for appropriated real constants $K_3>0$ and $K_4>0.$

Moreover, we define also the penalized functional $J_3^{*}:Y^{*}\times V\to \mathbb{R}$ where

$$J_3^{*}(v^{*},u)=J_2^{*}(v^{*},u)-\frac{K_1}{2}{\int }_{\Omega }{\left(v^{*}-\frac{\partial G\left(u\right)}{\partial u}+\epsilon u\right)}^{2}dx.$$

Finally, we remark that for $\epsilon >0$ sufficiently small and $K_1>0$ sufficiently large, $J_3^{*}$ is concave in $D^{*}\times V_1$ around a concerning critical point. We recall that a critical point

$$v^{*}-\frac{\partial G\left(u\right)}{\partial u}+\epsilon u=0,\mathrm{in}\Omega .$$

14. A related restricted problem in phase transition

In this section we develop a convex (in fact concave) dual variational for a model similar to those found in phase transition problems.

Let $\Omega =[0,1]\subset \mathbb{R}.$ Consider the functional $J:V\to \mathbb{R}$ where

$$\begin{array}{ccc}\hfill J\left(u\right)& =& \frac{1}{2}{\int }_{\Omega }min\{{(u^{\prime }+1)}^2,{(u^{\prime }-1)}^2\}dx\hfill \\ & & +\frac{1}{2}{\int }_{\Omega }{u}^{2}dx-{\langle u,f\rangle }_{L^2}\hfill \\ & =& \frac{1}{2}{\int }_{\Omega }{\left(u^{\prime }\right)}^{2}dx-{\int }_{\Omega }\left|u^{\prime }\right|dx+1/2\hfill \\ & & +\frac{1}{2}{\int }_{\Omega }{u}^{2}dx-{\langle u,f\rangle }_{L^2}.\hfill \end{array}$$

(65)

Here

$$V=\{u\in W^{1,2}(\Omega ):u\left(0\right)=0\mathrm{and}u\left(1\right)=1/2\}.$$

We also denote $V_1=W_0^{1,2}(\Omega ),$ and $Y=Y^{*}=L^2(\Omega ).$

Furthermore, we define the functionals G and $F:V\times V_1\to \mathbb{R}$ by

$$G(u^{\prime },v^{\prime })=\frac{1}{2}{\int }_{\Omega }{(u^{\prime }+v^{\prime })}^{2}dx-{\int }_{\Omega }|u^{\prime }+v^{\prime }|dx+1/2,$$

and

$$F(u,v)=\frac{1}{2}{\int }_{\Omega }{u}^{2}dx-{\langle u,f\rangle }_{L^2}.$$

Moreover we define $J_1:V\times V_1\to \mathbb{R}$ by

$$J_1(u,v)=G(u^{\prime },v^{\prime })+F(u,v),$$

and consider the problem of minimizing $J_1$ on the set

$$A=\{(u,v)\in V\times V_1:{\left(v^{\prime }\right)}^2\le K_2,\mathrm{in}\Omega \}.$$

Already including the Lagrange multiplier $\varphi$ concerning such restrictions, we define

$$J_2(u,v,\varphi )=J_1(u,v)+\frac{1}{2}{\langle {\varphi }^2,{\left(v^{\prime }\right)}^2-K_2\rangle }_{L^2}.$$

Observe now that

$$\begin{array}{ccc}\hfill J_2(u,v,\varphi )& =& J_1(u,v)+\frac{1}{2}{\langle {\varphi }^2,{\left(v^{\prime }\right)}^2-K_2\rangle }_{L^2}\hfill \\ & =& G(u^{\prime },v^{\prime })+\frac{1}{2}{\langle {\varphi }^2,{\left(v^{\prime }\right)}^2-K_2\rangle }_{L^2}\hfill \\ & & +F(u,v)\hfill \\ & =& -{\langle u^{\prime },v_1^{*}\rangle }_{L^2}-{\langle v^{\prime },v_2^{*}\rangle }_{L^2}+G(u^{\prime },v^{\prime })\hfill \\ & & +\frac{1}{2}{\langle {\varphi }^2,{\left(v^{\prime }\right)}^2-K_2\rangle }_{L^2}\hfill \\ & & {\langle u^{\prime },v_1^{*}\rangle }_{L^2}+{\langle v^{\prime },v_2^{*}\rangle }_{L^2}+F(u,v)\hfill \\ & \ge & \underset{(v_1,v_2)\in Y\times Y}{inf}\left\{-{\langle v_1,v_1^{*}\rangle }_{L^2}-{\langle v_2,v_2^{*}\rangle }_{L^2}+G_1(v_1,v_2,\varphi )\right.\hfill \\ & & \left.+\frac{1}{2}{\langle {\varphi }^2,{\left(v_2\right)}^2-K_2\rangle }_{L^2}\right\}\hfill \\ & & +\underset{(u,v)\in V\times V_1}{inf}\{{\langle u^{\prime },v_1^{*}\rangle }_{L^2}+{\langle v^{\prime },v_2^{*}\rangle }_{L^2}+F(u,v)\}\hfill \\ & =& -G_1^{*}(v_1^{*},v_2^{*},\varphi )-{\tilde{F}}^{*}(v_1^{*},v_2^{*}),\forall (u,v)\in V\times V_1,(v_1^{*},v_2^{*},\varphi )\in {\left[Y^{*}\right]}^3,\hfill \end{array}$$

(66)

where

$$G_1(u^{\prime },v^{\prime },\varphi )=G(u^{\prime },v^{\prime })+\frac{1}{2}{\langle {\varphi }^2,{\left(v^{\prime }\right)}^2-K_2\rangle }_{L^2}.$$

Also,

$$\begin{array}{ccc}\hfill G_1^{*}(v_1^{*},v_2^{*},\varphi )& =& \underset{(v_1,v_2)\in Y\times Y}{sup}\{{\langle v_1,v_1^{*}\rangle }_{L^2}+{\langle v_1,v_1^{*}\rangle }_{L^2}-G_1(v_1,v_2,\varphi )\}\hfill \\ & =& \frac{1}{2}{\int }_{\Omega }{\left(v_1^{*}\right)}^{2}dx\hfill \\ & & +{\int }_{\Omega }\left|v_1^{*}\right|dx+\frac{1}{2}{\int }_{\Omega }\frac{{(v_1^{*}-v_2^{*})}^2}{{\varphi }^2}\hfill \\ & & +\frac{K_2}{2}{\int }_{\Omega }{\varphi }^{2}dx,\hfill \end{array}$$

(67)

where

$${\tilde{F}}^{*}\left(v^{*}\right)=\left\{\begin{array}{cc}\frac{1}{2}{\int }_{\Omega }{({\left(v_1^{*}\right)}^{\prime }+f)}^{2}dx-v_1^{*}\left(1\right)u\left(1\right),\hfill & \mathrm{if}{\left(v_2^{*}\right)}^{\prime }=0,\mathrm{in}\Omega ,\hfill \\ +\infty ,\hfill & \mathrm{otherwise}.\hfill \end{array}\right.$$

(68)

From this we may infer that $v_2^{*}=c,\mathrm{in}\Omega ,$ for some $c\in \mathbb{R}.$

Summarizing, denoting $v^{*}=(v_1^{*},v_2^{*})=(v_1^{*},c)$, and

$$J^{*}(v^{*},\varphi )=-G_1^{*}(v^{*},\varphi )-{\tilde{F}}^{*}\left(v^{*}\right)$$

we have got

$$\underset{(u,v)\in A}{inf}{J}_1(u,v)\ge \underset{(v^{*},\varphi )\in Y^{*}\times \mathbb{R}\times Y^{*}}{sup}{J}^{*}\left(v^{*}\varphi \right).$$

We have developed numerical results by maximizing the dual functional $J^{*}$ for two examples, namely.

• Example A: In this case, we consider $f\left(x\right)=cos\left(\pi x\right)/2$, $K_2={10}^{-4}$.
  For the optimal

  $$u_0={\left(v_1^{*}\right)}^{\prime }+f,$$

  please see Figure 9.
• Example B: In this case, we consider $f\left(x\right)=cos\left(\pi x\right)/2$, $K_2=30$.
  For the optimal

  $$u_0={\left(v_1^{*}\right)}^{\prime }+f,$$

  please see Figure 10.

15. One more dual variational formulation

In this section we develop one more dual variational formulation for a related model.

Let $\Omega =[0,1]\subset \mathbb{R}$ and consider the functional $J:V\to \mathbb{R}$ defined by

$$J\left(u\right)=\frac{1}{2}{\int }_{\Omega }{({\left(u^{\prime }\right)}^2-1)}^{2}dx+\frac{1}{2}{\int }_{\Omega }{u}^{2}dx-{\langle u,f\rangle }_{L^2},$$

where

$$V=\{u\in W^{1,4}(\Omega ):u\left(0\right)=0\mathrm{and}u\left(1\right)=1/2\}.$$

We define also the relaxed functional $J_1:V\times V_0\to \mathbb{R}$, already including a concerning restriction and corresponding non-negative Lagrange multiplier ${\Lambda }^2$, where

$$J_1(u,v,\Lambda )=\frac{1}{2}{\int }_{\Omega }{({(u^{\prime }+v^{\prime })}^2-1)}^{2}dx+\frac{1}{2}{\int }_{\Omega }{u}^{2}dx-{\langle u,f\rangle }_{L^2}+{\langle {\Lambda }^2,{\left(v^{\prime }\right)}^2-K\rangle }_{L^2}.$$

where

$$V_0=\{v\in W_0^{1,4}(\Omega ):{\left(v^{\prime }\right)}^2-K\le 0\mathrm{in}\Omega \}.$$

Observe that

$$\begin{array}{ccc}& & \frac{1}{2}{\int }_{\Omega }{({(u^{\prime }+v^{\prime })}^2-1)}^{2}dx+\frac{1}{2}{\int }_{\Omega }{u}^{2}dx-{\langle u,f\rangle }_{L^2}+{\langle {\Lambda }^2,{\left(v^{\prime }\right)}^2-K\rangle }_{L^2}\hfill \\ & =& -{\langle v_0^{*},{(u^{\prime }+v^{\prime })}^2-1\rangle }_{L^2}+\frac{1}{2}{\int }_{\Omega }{({(u^{\prime }+v^{\prime })}^2-1)}^{2}dx\hfill \\ & & +{\langle v_0^{*},{(u^{\prime }+v^{\prime })}^2-1\rangle }_{L^2}+{\langle {\Lambda }^2,{\left(v^{\prime }\right)}^2-K\rangle }_{L^2}-{\langle u^{\prime },v_1^{*}\rangle }_{L^2}-{\langle v^{\prime },v_2^{*}\rangle }_{L^2}\hfill \\ & & +{\langle u^{\prime },v_1^{*}\rangle }_{L^2}+{\langle v^{\prime },v_2^{*}\rangle }_{L^2}+\frac{1}{2}{\int }_{\Omega }{u}^{2}dx-{\langle u,f\rangle }_{L^2}\hfill \\ & \ge & \underset{w\in Y}{inf}\left\{-{\langle v_0^{*},w\rangle }_{L^2}+\frac{1}{2}{\int }_{\Omega }{\left(w\right)}^{2}dx\right\}\hfill \\ & & \underset{(v_1,v_2)\in Y\times Y}{inf}\left\{{\langle v_0^{*},{(v_1+v_2)}^2-1\rangle }_{L^2}+{\langle {\Lambda }^2,{\left(v_2\right)}^2-K\rangle }_{L^2}-{\langle v_1,v_1^{*}\rangle }_{L^2}-{\langle v_2,v_2^{*}\rangle }_{L^2}\right\}\hfill \\ & & +\underset{(u,v)\in V\times V_0}{inf}\left\{{\langle u^{\prime },v_1^{*}\rangle }_{L^2}+{\langle v^{\prime },v_2^{*}\rangle }_{L^2}+\frac{1}{2}{\int }_{\Omega }{u}^{2}dx-{\langle u,f\rangle }_{L^2}\right\}\hfill \\ & =& -\frac{1}{2}{\int }_{\Omega }{\left(v_0^{*}\right)}^{2}dx-{\int }_{\Omega }{v}_0^{*}dx\hfill \\ & & -\frac{1}{4}{\int }_{\Omega }\frac{{\left(v_1^{*}\right)}^2}{v_0^{*}}dx-\frac{1}{2}{\int }_{\Omega }\frac{{(v_1^{*}-v_2^{*})}^2}{2{\Lambda }^2}dx\hfill \\ & & -\frac{1}{2}{\int }_{\Omega }{({\left(v_1^{*}\right)}^{\prime }+f)}^{2}dx-\frac{1}{2}{\int }_{\Omega }K{\Lambda }^{2}dx+v_1^{*}\left(1\right)u\left(1\right).\hfill \end{array}$$

(69)

Here, we highlight $v_2^{*}=c\in \mathbb{R}\mathrm{in}\Omega$, for some real constant c.

Hence, denoting

$$\begin{array}{ccc}\hfill J_1^{*}(v^{*},\Lambda )& =& -\frac{1}{2}{\int }_{\Omega }{\left(v_0^{*}\right)}^{2}dx-{\int }_{\Omega }{v}_0^{*}dx\hfill \\ & & -\frac{1}{4}{\int }_{\Omega }\frac{{\left(v_1^{*}\right)}^2}{v_0^{*}}dx-\frac{1}{2}{\int }_{\Omega }\frac{{(v_1^{*}-v_2^{*})}^2}{2{\Lambda }^2}dx\hfill \\ & & -\frac{1}{2}{\int }_{\Omega }{({\left(v_1^{*}\right)}^{\prime }+f)}^{2}dx-\frac{1}{2}{\int }_{\Omega }K{\Lambda }^{2}dx+v_1^{*}\left(1\right)u\left(1\right)\hfill \end{array}$$

(70)

and

$$J_2(u,v)=\frac{1}{2}{\int }_{\Omega }{({(u^{\prime }+v^{\prime })}^2-1)}^{2}dx+\frac{1}{2}{\int }_{\Omega }{u}^{2}dx-{\langle u,f\rangle }_{L^2},$$

we have obtained

$$\underset{(u,v)\in V\times V_0}{inf}{J}_2(u,v)\}\ge \underset{(v^{*},\Lambda )\in A^{*}\times \left[Y^{*}\right]\times \mathbb{R}\times Y^{*}}{sup}{J}_1^{*}(v^{*},\Lambda ).$$

Finally, for

$$A^{*}=\{v_0^{*}\in Y^{*}:v_0^{*}\ge \epsilon \mathrm{in}\Omega \}$$

we emphasize $J_1^{*}$ is concave on $A^{*}\times \left[Y^{*}\right]\times \mathbb{R}\times Y^{*}$.

Here $\epsilon >0$ is a small regularizing real constant.

Remark 15.1.

The constraint ${\left(v^{\prime }\right)}^2-K\le 0,\mathrm{in}\Omega$ is included to restrict the action of v on the region where the primal functional is non-convex, through an appropriate constant $K>0.$

16. A duality principle for a related relaxed formulation concerning the vectorial approach in the calculus of variations

In this section we develop a duality principle for a related vectorial model in the calculus of variations.

Let $\Omega \subset {\mathbb{R}}^n$ be an open, bounded and connected set with a regular (Lipschitzian) boundary denoted by $\partial \Omega =\Gamma .$

For $1<p<+\infty$, consider a functional $J:V\to \mathbb{R}$ where

$$J\left(u\right)=G(\nabla u)+F\left(u\right)-{\langle u,f\rangle }_{L^2},$$

where

$$V=\left\{u\in W^{1,p}(\Omega ;{\mathbb{R}}^N):u=u_0\mathrm{on}\partial \Omega \right\},$$

$u_0\in C^1(\overline{\Omega };{\mathbb{R}}^N)$ and $f\in L^2(\Omega ;{\mathbb{R}}^N).$

We assume $G:Y\to \mathbb{R}$ and $F:V\to \mathbb{R}$ are Fréchet differentiable and F is also convex.

Also

$$G(\nabla u)={\int }_{\Omega }g(\nabla u)dx,$$

where $g:{\mathbb{R}}^{N\times n}\to \mathbb{R}$ it is supposed to be Fréchet differentiable. Here we have denoted $Y=L^p(\Omega ;{\mathbb{R}}^{N\times n})$.

We define also $J_1:V\times Y_1\to \mathbb{R}$ by

$$J_1(u,\varphi )=G_1(\nabla u+{\nabla }_y\varphi )+F\left(u\right)-{\langle u,f\rangle }_{L^2},$$

where

$$Y_1=W^{1,p}(\Omega \times \Omega ;{\mathbb{R}}^N)$$

and

$$G_1(\nabla u+{\nabla }_y\varphi )=\frac{1}{|\Omega |}{\int }_{\Omega }{\int }_{\Omega }g(\nabla u\left(x\right)+{\nabla }_y\varphi (x,y))dxdy.$$

Moreover, we define the relaxed functional $J_2:V\to \mathbb{R}$ by

$$J_2\left(u\right)=\underset{\varphi \in V_0}{inf}{J}_1(u,\varphi ),$$

where

$$V_0=\{\varphi \in Y_1:\varphi (x,y)=0,\mathrm{on}\Omega \times \partial \Omega \}.$$

Now observe that

$$\begin{array}{ccc}\hfill J_1(u,\varphi )& =& G_1(\nabla u+{\nabla }_y\varphi )+F\left(u\right)-{\langle u,f\rangle }_{L^2}\hfill \\ & =& -\frac{1}{|\Omega |}{\int }_{\Omega }{\int }_{\Omega }{v}^{*}(x,y)·(\nabla u+{\nabla }_y\varphi (x,y))dydx+G_1(\nabla u+{\nabla }_y\varphi )\hfill \\ & & +\frac{1}{|\Omega |}{\int }_{\Omega }{\int }_{\Omega }{v}^{*}(x,y)·(\nabla u+{\nabla }_y\varphi (x,y))dydx+F\left(u\right)-{\langle u,f\rangle }_{L^2}\hfill \\ & \ge & \underset{v\in Y_2}{inf}\left\{-\frac{1}{|\Omega |}{\int }_{\Omega }{\int }_{\Omega }{v}^{*}(x,y)·v(x,y)dydx+G_1\left(v\right)\right\}\hfill \\ & & +\underset{(v,\varphi )\in V\times V_0}{inf}\left\{\frac{1}{|\Omega |}{\int }_{\Omega }{\int }_{\Omega }{v}^{*}(x,y)·(\nabla u+{\nabla }_y\varphi (x,y))dydx+F\left(u\right)-{\langle u,f\rangle }_{L^2}\right\}\hfill \\ & =& -G_1^{*}\left(v^{*}\right)-F^{*}\left({\mathrm{div}}_x\left(\frac{1}{|\Omega |}{\int }_{\Omega }{v}^{*}(x,y)dy\right)+f\right)\hfill \\ & & +\frac{1}{|\Omega |}{\int }_{\partial \Omega }\left({\int }_{\Omega }{v}^{*}(x,y)dy\right)\otimes \mathbf{n}{u}_{0}d\Gamma ,\hfill \end{array}$$

(71)

$\forall (u,\varphi )\in V\times V_0,v^{*}\in A^{*}$, where

$$A^{*}=\{v^{*}\in Y_2^{*}:{\mathrm{div}}_y{v}^{*}(x,y)=0,\mathrm{in}\Omega \}.$$

Here we have denoted

$$G_1^{*}\left(v^{*}\right)=\underset{v\in Y_2}{sup}\left\{\frac{1}{|\Omega |}{\int }_{\Omega }{\int }_{\Omega }{v}^{*}(x,y)·v(x,y)dydx-G_1\left(v\right)\right\},$$

where $Y_2=L^p(\Omega \times \Omega ;{\mathbb{R}}^{N\times n}),$$Y_2^{*}=L^q(\Omega \times \Omega ;{\mathbb{R}}^{N\times n}),$ and where

$$\frac{1}{p}+\frac{1}{q}=1.$$

Furthermore, for $v^{*}\in A^{*}$, we have

$$\begin{array}{ccc}& & F^{*}\left({\mathrm{div}}_x\left(\frac{1}{|\Omega |}{\int }_{\Omega }{v}^{*}(x,y)dy\right)+f\right)-\frac{1}{|\Omega |}{\int }_{\partial \Omega }\left({\int }_{\Omega }{v}^{*}(x,y)dy\right)\otimes \mathbf{n}{u}_{0}d\Gamma \hfill \\ & =& \underset{(v,\varphi )\in V\times V_0}{sup}\left\{-\frac{1}{|\Omega |}{\int }_{\Omega }{\int }_{\Omega }{v}^{*}(x,y)·(\nabla u+{\nabla }_y\varphi (x,y))dydx-F\left(u\right)+{\langle u,f\rangle }_{L^2}\right\},\hfill \end{array}$$

(72)

Therefore, denoting $J_3^{*}:Y_2^{*}\to \mathbb{R}$ by

$$J_3^{*}\left(v^{*}\right)=-G_1^{*}\left(v^{*}\right)-F^{*}\left({\mathrm{div}}_x\left({\int }_{\Omega }{v}^{*}(x,y)dy\right)+f\right)+\frac{1}{|\Omega |}{\int }_{\partial \Omega }\left({\int }_{\Omega }{v}^{*}(x,y)dy\right)\otimes \mathbf{n}{u}_{0}d\Gamma ,$$

we have got

$$\underset{u\in V}{inf}{J}_2\left(u\right)\ge \underset{v^{*}\in A^{*}}{sup}{J}_3^{*}\left(v^{*}\right).$$

Finally, we highlight such a dual functional $J_3^{*}$ is convex (in fact concave).

17. A model in superconductivity through an eigenvalue approach

In this section we intend to model superconductivity through a two phase eigenvalue approach.

Let $\Omega =[0,5]\subset \mathbb{R}$ be a straight wire corresponding to a one-dimensional super-conducting sample.

Consider the functional $J:V\times V\times \mathbb{R}\to \mathbb{R}$ where

$$\begin{array}{ccc}\hfill J(u,v,E)& =& \frac{{\gamma }_1}{2}{\int }_{\Omega }\nabla u·\nabla udx+\frac{{\alpha }_1}{2}{\int }_{\Omega }{\left|u\right|}^{4}dx\hfill \\ & & -\frac{{\omega }^2}{2}{\int }_{\Omega }{\left|u\right|}^{2}dx\hfill \\ & & +\frac{{\gamma }_2}{2}{\int }_{\Omega }\nabla v·\nabla vdx+\frac{{\alpha }_2}{2}{\int }_{\Omega }{\left|v\right|}^{4}dx\hfill \\ & & -\frac{{\omega }_1^2}{2K_3^2}{\int }_{\Omega }{\left|v\right|}^{2}dx\hfill \\ & & -\frac{E}{2}\left({\int }_{\Omega }{\left(\right|u|}^2+{\left|v\right|}^2)dx-m_T\right).\hfill \end{array}$$

(73)

Here, in atomic units, $m_T$ is the total electronic charge, $V=W_0^{1,2}(\Omega )$ and we set ${\alpha }_1={10}^4$ corresponding to higher self-interacting energy which is related to a normal phase. We also set ${\alpha }_2={10}^{-1}$ corresponding to a lower self-interacting energy which is related to a super-conducting phase and respective super-currents.

Moreover, we set ${\gamma }_1={\gamma }_2=1,$ and initially $\omega =1.8$ which is gradually decreased to $\omega =1.0$.

Furthermore, we define

$$|{\varphi }_N{|}^2=\frac{{\left|u\right|}^2}{{\left|u\right|}^2+{\left|v\right|}^2}$$

and

$$|{\varphi }_S{|}^2=\frac{{\left|v\right|}^2}{{\left|u\right|}^2+{\left|v\right|}^2}$$

where ${\varphi }_N$ corresponds to a normal phase and ${\varphi }_S$ to a super-conducting one.

At this point we observe that the temperature $T=T(x,t)$ is proportional the frequency $\omega /\left(2\pi \right)$ of vibration for the normal phase.

We start the process with $\omega =1.8$ which in atomic units corresponds to a higher temperature and gradually decreases it to the value $\omega =1.0$

Between $\omega =1.2$ and $\omega =1.0$ the system changes from an almost total normal phase to an almost total super-conducting phase, as expected.

We highlight that the temperature is proportional to the vibrational kinetics energy

$$E_1\left(t\right)=\frac{1}{2}{\int }_{\Omega }{\left|u\right|}^2\frac{\partial {\mathbf{r}}_N(x,t)}{\partial t}·\frac{\partial {\mathbf{r}}_N(x,t)}{\partial t}dx$$

so that for

$${\mathbf{r}}_N(x,t)=e^{i\omega t}{\mathbf{w}}_5\left(x\right)$$

and for a suitable vectorial function ${\mathbf{w}}_5$, we have

$$T\propto E_1\propto {\omega }^2$$

so that we may model the decreasing of temperature T through the decreasing of ${\omega }^2$.

For $\omega =1.8$, for the corresponding normal phase ${\varphi }_N$ and super-conducting phase ${\varphi }_S$, please se Figure 11 and Figure 12, respectively.

For $\omega =1.0$, for the corresponding normal phase ${\varphi }_N$ and super-conducting phase ${\varphi }_S$, please se Figure 13 and Figure 14, respectively.

Finally, we have set ${\omega }_1/K_3\approx 1$ which for large ${\omega }_1$ corresponds to the super-currents.

18. A simplified qualitative many body model for the hydrogen nuclear fusion

In this section we develop a qualitative simple model for the hydrogen nuclear fusion.

Let $\Omega ={[0,L]}^3\subset {\mathbb{R}}^3$ be a box in which is confined a gas comprised by an amount of ionized deuterium and tritium isotopes of hydrogen.

Though a suitable increasing in temperature, we intend to develop the following nuclear reaction

$${\mathrm{Deuterium}}^{+}+{\mathrm{Tritium}}^{+}\to {\mathrm{Hellion}}^{++}+\mathrm{Neutron}\left(\mathrm{energetic}\right).$$

We recall that the ionized Deuterium atom comprises a proton and a neutron and the ionized Tritium atom comprises a proton and two neutrons.

Under certain conditions and at a suitable high temperature the ionized Deuterium and Tritium atoms react chemically resulting in an ionized Hellion atom, comprised by two protons and two neutrons and resulting also in one more single energetic neutron. We emphasize the higher kinetics neutron energy level has many potential practical applications, including its conversion in electric energy.

At this point we denote by $m_D,m_T,m_{H_e}\mathrm{and}{m}_N$ the masses of the ionized Deuterium, Tritium and Hellion atoms, and the single neutron, respectively.

Therefore, we have the following mass relation

$$m_D+m_T=m_{H_e}+m_N.$$

To simplify our analysis, in such a chemical reaction, denoting the total masses of ionized Deuterium, Tritium, Hellion and single Neutrons by ${\left(m_D\right)}_T,{\left(m_T\right)}_T,{\left(m_{H_e}\right)}_T$ and ${\left(m_N\right)}_T$ we assume there is a real constant $c>0$ such that

$${\left(m_D\right)}_T=c{m}_D,{\left(m_T\right)}_T=c{m}_T,{\left(m_{H_e}\right)}_T=c{m}_{H_e},{\left(m_N\right)}_T=c{m}_N.$$

With such statements and definitions in mind, we define the following functional J, where

$$J(\varphi ,\mathbf{r})=J({\varphi }_D,{\varphi }_T,{\varphi }_{H_e},{\varphi }_N,\mathbf{r})=G(\nabla \varphi )+F\left(\varphi \right)+E_c(\varphi ,\mathbf{r}),$$

where, in a simplified many body context,

$$|{\varphi }_D{(x,y)|}^2=|{\varphi }_p^D{\left(y\right)|}^2+|{\varphi }_N^D{(x,y)|}^2{\left|{\varphi }_p^D\left(y\right)\right|}^2\frac{1}{m_p},$$

$$|{\varphi }_T{(x,y)|}^2=|{\varphi }_p^T{\left(y\right)|}^2+\left(\right|{\varphi }_{N_1}^T{(x,y)|}^2+|{\varphi }_{N_2}^T{(x,y)|}^2\left)\right|{\varphi }_p^T{\left(y\right)|}^2\frac{1}{m_p},$$

$$|{\varphi }_{H_e}{(x,y)|}^2=|{\varphi }_{2P}^{H_e}{\left(y\right)|}^2+\left(\right|{\varphi }_{N_1}^{H_e}{(x,y)|}^2+|{\varphi }_{N_2}^{H_e}{(x,y)|}^2\left)\right|{\varphi }_{2P}^{H_e}{\left(y\right)|}^2\frac{1}{2m_p},$$

$${\varphi }_N={\varphi }_N\left(x\right).$$

Here $x,y\in \Omega \subset {\mathbb{R}}^3$ refers to the particle densities.

Furthermore, we assume ${\gamma }_p^D>0,{\gamma }_p^T>0,{\gamma }_N^D>0,{\gamma }_{N_1}^T>0,{\gamma }_{N_2}^T>0,$${\gamma }_{2p}^{H_e}>0,{\gamma }_{N_1}^{H_e}>0,$ ${\gamma }_{N_2}^{H_e}>0,{\gamma }_N>0,$ and ${\alpha }_D>0,{\alpha }_T>0,{\alpha }_{H_e}>0,{\alpha }_N>0$, ${\alpha }_{DT}>0,{\alpha }_{H_{e}N}>0,$ so that

$$\begin{array}{ccc}\hfill G(\nabla \varphi )& =& \frac{{\gamma }_p^D}{2}{\int }_{\Omega }(\nabla {\varphi }_p^D)·(\nabla {\varphi }_p^D)dy\hfill \\ & & +\frac{{\gamma }_N^D}{2}{\int }_{\Omega }(\nabla {\varphi }_N^D)·(\nabla {\varphi }_N^D)dxdy\hfill \\ & & \frac{{\gamma }_p^T}{2}{\int }_{\Omega }(\nabla {\varphi }_p^T)·(\nabla {\varphi }_p^T)dy\hfill \\ & & +\frac{{\gamma }_{N_1}^T}{2}{\int }_{\Omega }(\nabla {\varphi }_{N_1}^T)·(\nabla {\varphi }_{N_1}^T)dxdy\hfill \\ & & +\frac{{\gamma }_{N_2}^T}{2}{\int }_{\Omega }(\nabla {\varphi }_{N_2}^T)·(\nabla {\varphi }_{N_2}^T)dxdy\hfill \\ & & +\frac{{\gamma }_{2p}^{H_e}}{2}{\int }_{\Omega }(\nabla {\varphi }_{2p}^{H_e})·(\nabla {\varphi }_{2p}^{H_e})dy\hfill \\ & & +\frac{{\gamma }_{N_1}^{H_e}}{2}{\int }_{\Omega }(\nabla {\varphi }_{N_1}^{H_e})·(\nabla {\varphi }_{N_1}^{H_e})dxdy\hfill \\ & & +\frac{{\gamma }_{N_2}^{H_e}}{2}{\int }_{\Omega }(\nabla {\varphi }_{N_2}^{H_e})·(\nabla {\varphi }_{N_2}^{H_e})dxdy\hfill \\ & & +\frac{{\gamma }_N}{2}{\int }_{\Omega }(\nabla {\varphi }_N)·(\nabla {\varphi }_N)dx,\hfill \end{array}$$

(74)

and,

$$\begin{array}{ccc}\hfill F\left(\varphi \right)& =& \frac{{\alpha }_D}{2}{\int }_{\Omega }\frac{|{\varphi }_D(x-{\xi }_1,y-{\xi }_2){|}^2{\left|{\varphi }_D({\xi }_1,{\xi }_2)\right|}^2}{|(x,y)-({\xi }_1,{\xi }_2)|}dxdyd{\xi }_{1}d{\xi }_2\hfill \\ & & +\frac{{\alpha }_T}{2}{\int }_{\Omega }\frac{|{\varphi }_T(x-{\xi }_1,y-{\xi }_2){|}^2{\left|{\varphi }_T({\xi }_1,{\xi }_2)\right|}^2}{|(x,y)-({\xi }_1,{\xi }_2)|}dxdyd{\xi }_{1}d{\xi }_2\hfill \\ & & +\frac{{\alpha }_{DT}}{2}{\int }_{\Omega }\frac{|{\varphi }_D(x-{\xi }_1,y-{\xi }_2){|}^2{\left|{\varphi }_T({\xi }_1,{\xi }_2)\right|}^2}{|(x,y)-({\xi }_1,{\xi }_2)|}dxdyd{\xi }_{1}d{\xi }_2\hfill \\ & & +\frac{{\alpha }_{H_e}}{2}{\int }_{\Omega }\frac{|{\varphi }_{H_e}(x-{\xi }_1,y-{\xi }_2){|}^2{\left|{\varphi }_{H_e}({\xi }_1,{\xi }_2)\right|}^2}{|(x,y)-({\xi }_1,{\xi }_2)|}dxdyd{\xi }_{1}d{\xi }_2\hfill \\ & & +\frac{{\alpha }_N}{2}{\int }_0^{t_f}{\int }_{\Omega }\frac{|{\varphi }_N{(x-\xi )|}^2{\left|{\varphi }_N\left(\xi \right)\right|}^2}{|x-\xi |}dxd\xi \hfill \\ & & +\sum _{j=1}^2\frac{{\alpha }_{H_{e}N}}{2}{\int }_{\Omega }\frac{|{\varphi }_{H_e}(x_1-{\xi }_1,y-{\xi }_2){|}^2{\left|{\varphi }_N\left({\xi }_j\right)\right|}^2}{|(x,y)-({\xi }_1,{\xi }_2)|}dxdyd{\xi }_{1}d{\xi }_2\hfill \end{array}$$

(75)

and the kinetics energy is expressed by

$$\begin{array}{ccc}\hfill E_c(\varphi ,\mathbf{r})& =& \frac{1}{2}{\int }_{\Omega }{\left|{\varphi }_D\right|}^2\frac{\partial {\mathbf{r}}_D}{\partial t}·\frac{\partial {\mathbf{r}}_D}{\partial t}dxdy\hfill \\ & & +\frac{1}{2}{\int }_{\Omega }{\left|{\varphi }_T\right|}^2\frac{\partial {\mathbf{r}}_T}{\partial t}·\frac{\partial {\mathbf{r}}_T}{\partial t}dxdy\hfill \\ & & +\frac{1}{2}{\int }_{\Omega }{\left|{\varphi }_{H_e}\right|}^2\frac{\partial {\mathbf{r}}_{H_e}}{\partial t}·\frac{\partial {\mathbf{r}}_{H_e}}{\partial t}dxdy\hfill \\ & & +\frac{1}{2}{\int }_{\Omega }{\left|{\varphi }_N\right|}^2\frac{\partial {\mathbf{r}}_N}{\partial t}·\frac{\partial {\mathbf{r}}_N}{\partial t}dxdy,\hfill \end{array}$$

(76)

where we also assume

$${\mathbf{r}}_D\approx e^{i\omega t}{\mathbf{w}}_5(x,y),$$

$${\mathbf{r}}_T\approx e^{i\omega t}{\mathbf{w}}_6(x,y),$$

so that considering such a vibrational motion, the temperature T is proportional to ${\omega }^2$, that is

$$T\propto {\omega }^2.$$

Therefore, an increasing in T corresponds to a proportional increasing in ${\omega }^2.$

Summarizing, we have supposed

$$E_c(\varphi ,\mathbf{r})\approx \frac{1}{2}{\omega }^2{\int }_{\Omega }|{\varphi }_D{|}^2+|{\varphi }_T{|}^{2}dx{C}_1+\frac{1}{2}{\omega }_1^2{\int }_{\Omega }{\left|{\varphi }_N\right|}^{2}dx{C}_2,$$

so that we represent the increasing in T through an increasing in ${\omega }^2.$

Moreover, we denote by $m_N$ the mass of a single neutron and by $m_p$ the mass of a single proton.

Thus, denoting also by ${\lambda }_1,{\lambda }_2$ the proportion of non-reacted and reacted masses respectively, we have the following constraints.

• $${\int }_{\Omega }{\left|{\varphi }_N^D(x,y)\right|}^{2}dx=m_N,$$
• $${\int }_{\Omega }{\left|{\varphi }_{N_1}^T(x,y)\right|}^{2}dx=m_N,$$
• $${\int }_{\Omega }{\left|{\varphi }_{N_2}^T(x,y)\right|}^{2}dx=m_N,$$
• $${\int }_{\Omega }{\left|{\varphi }_{N_1}^{H_e}(x,y)\right|}^{2}dx=m_N,$$
• $${\int }_{\Omega }{\left|{\varphi }_{N_2}^{H_e}(x,y)\right|}^{2}dx=m_N,$$
• $${\int }_{\Omega }{\left|{\varphi }_p^D\left(y\right)\right|}^{2}dy={\lambda }_{1}c{m}_p,$$
• $${\int }_{\Omega }{\left|{\varphi }_p^T\left(y\right)\right|}^{2}dy={\lambda }_{1}c{m}_p,$$
• $${\int }_{\Omega }{\left|{\varphi }_{2P}^{H_e}\left(y\right)\right|}^{2}dy={\lambda }_2\left(2c{m}_p\right),$$

Similar constraints are valid corresponding to the charge of a single proton.

We have also the following complementing constraints,

• $${\int }_{\Omega }{\left|{\varphi }_D\right|}^{2}dxdy={\lambda }_1{\left(m_D\right)}_T,$$
• $${\int }_{\Omega }{\left|{\varphi }_T\right|}^{2}dxdy={\lambda }_1{\left(m_T\right)}_T,$$
• $${\int }_{\Omega }{\left|{\varphi }_{H_e}\right|}^{2}dxdy={\lambda }_2{\left(m_{H_e}\right)}_T,$$
• $${\int }_{\Omega }{\left|{\varphi }_N\right|}^{2}dxdy={\lambda }_2{\left(m_N\right)}_T,$$
• $${\lambda }_1+{\lambda }_2=1.$$

With such results and statements in mind and simplifying the interacting terms, we re-define the functional J now denoting it by $J_1$, here already including the Lagrange multipliers concerning the constraints, where

$$\begin{array}{ccc}\hfill J_1(\varphi ,\omega ,E,\lambda )& =& \frac{{\gamma }_p^D}{2}{\int }_{\Omega }(\nabla {\varphi }_p^D)·(\nabla {\varphi }_p^D)dy\hfill \\ & & +\frac{{\gamma }_N^D}{2}{\int }_{\Omega }(\nabla {\varphi }_N^D)·(\nabla {\varphi }_N^D)dxdy\hfill \\ & & \frac{{\gamma }_p^T}{2}{\int }_{\Omega }(\nabla {\varphi }_p^T)·(\nabla {\varphi }_p^T)dy\hfill \\ & & +\frac{{\gamma }_{N_1}^T}{2}{\int }_{\Omega }(\nabla {\varphi }_{N_1}^T)·(\nabla {\varphi }_{N_1}^T)dxdy\hfill \\ & & +\frac{{\gamma }_{N_2}^T}{2}{\int }_{\Omega }(\nabla {\varphi }_{N_2}^T)·(\nabla {\varphi }_{N_2}^T)dxdy\hfill \\ & & +\frac{{\gamma }_{2p}^{H_e}}{2}{\int }_{\Omega }(\nabla {\varphi }_{2p}^{H_e})·(\nabla {\varphi }_{2p}^{H_e})dy\hfill \\ & & +\frac{{\gamma }_{N_1}^{H_e}}{2}{\int }_{\Omega }(\nabla {\varphi }_{N_1}^{H_e})·(\nabla {\varphi }_{N_1}^{H_e})dxdy\hfill \\ & & +\frac{{\gamma }_{N_2}^{H_e}}{2}{\int }_{\Omega }(\nabla {\varphi }_{N_2}^{H_e})·(\nabla {\varphi }_{N_2}^{H_e})dxdy\hfill \\ & & +\frac{{\gamma }_N}{2}{\int }_{\Omega }(\nabla {\varphi }_N)·(\nabla {\varphi }_N)dx\hfill \\ & & +\frac{{\alpha }_D}{2}{\int }_{\Omega }|{\varphi }_D{|}^{4}dx+\frac{{\alpha }_T}{2}{\int }_{\Omega }{\left|{\varphi }_T\right|}^{4}dx\hfill \\ & & +\frac{{\alpha }_{H_e}}{2}{\int }_{\Omega }|{\varphi }_{H_e}{|}^{4}dx+\frac{{\alpha }_N}{2}{\int }_{\Omega }{\left|{\varphi }_N\right|}^{4}dx\hfill \\ & & -{\omega }^2{\int }_{\Omega }\left(\right|{\varphi }_D{|}^2+|{\varphi }_T{|}^2)dx\hfill \\ & & -{\omega }_1^2{\int }_{\Omega }{\left|{\varphi }_N\right|}^{2}dx+J_{Aux},\hfill \end{array}$$

(77)

where the functional $J_{Aux}$ stands for

$$\begin{array}{ccc}\hfill J_{Aux}& =& -{\int }_{\Omega }{\left(E_N^D\right)}_5\left(y\right)\left({\int }_{\Omega }{\left|{\varphi }_N^D(x,y)\right|}^{2}dx-m_N\right)dy\hfill \\ & & -{\int }_{\Omega }{\left(E_{N_1}^T\right)}_6\left(y\right)\left({\int }_{\Omega }{\left|{\varphi }_{N_1}^T(x,y)\right|}^{2}dx-m_N\right)dy\hfill \\ & & -{\int }_{\Omega }{\left(E_{N_2}^T\right)}_7\left(y\right)\left({\int }_{\Omega }{\left|{\varphi }_{N_2}^T(x,y)\right|}^{2}dx-m_N\right)dy\hfill \\ & & -{\int }_{\Omega }{\left(E_{N_1}^{H_e}\right)}_8\left(y\right)\left({\int }_{\Omega }{\left|{\varphi }_{N_1}^{H_e}(x,y)\right|}^{2}dx-m_N\right)dy\hfill \\ & & -{\int }_{\Omega }{\left(E_{N_2}^{H_e}\right)}_9\left(y\right)\left({\int }_{\Omega }{\left|{\varphi }_{N_2}^{H_e}(x,y)\right|}^{2}dx-m_N\right)dy\hfill \\ & & -{\left(E_D\right)}_2\left({\int }_{\Omega }{\left|{\varphi }_p^D\left(y\right)\right|}^{2}dy-{\lambda }_{1}c{m}_p\right)\hfill \\ & & -{\left(E_T\right)}_3\left({\int }_{\Omega }{\left|{\varphi }_p^T\left(y\right)\right|}^{2}dy-{\lambda }_{1}c{m}_p\right)\hfill \\ & & -{\left(E_{H_e}\right)}_3\left({\int }_{\Omega }{\left|{\varphi }_{2P}^{H_e}(x,y)\right|}^{2}dy-{\lambda }_{2}2c{m}_p\right)\hfill \\ & & -E_5\left({\int }_{\Omega }{\left|{\varphi }_D\right|}^{2}dxdy-{\lambda }_1{\left(m_D\right)}_T\right)\hfill \\ & & -E_6\left({\int }_{\Omega }{\left|{\varphi }_T\right|}^{2}dxdy-{\lambda }_1{\left(m_T\right)}_T\right)\hfill \\ & & -E_7\left({\int }_{\Omega }{\left|{\varphi }_{H_e}\right|}^{2}dxdy-{\lambda }_2{\left(m_{H_e}\right)}_T\right)\hfill \\ & & -E_8\left({\int }_{\Omega }{\left|{\varphi }_N\right|}^{2}dxdy-{\lambda }_2{\left(m_N\right)}_T\right)\hfill \\ & & -E_9({\lambda }_1+{\lambda }_2-1).\hfill \end{array}$$

(78)

Remark 18.1.

In order to obtain consistent results it is necessary to set

$$({\alpha }_N,{\alpha }_{H_e})\gg ({\alpha }_D,{\alpha }_T).$$

In such a case, a higher temperature corresponding to a large ${\omega }^2$, though such a nuclear reaction, will result in a small ${\lambda }_1$ and a higher kinetics energy for the neutron field, corresponding to a large ${\omega }_1^2$ and ${\lambda }_2$ closer to 1.

19. A more detailed mathematical description of the hydrogen nuclear fusion

In this section we develop in more details another model for the hydrogen nuclear fusion.

Let $\Omega \subset {\mathbb{R}}^3$ be an open, bounded and connected set with a regular (Lipschitzian) boundary denoted by $\partial \Omega .$

Here such a set $\Omega$ stands for a control volume in which an ionized gas (plasma) flows. Such a gas comprises ionized Deuterium and Tritium atoms intended, through a suitable higher temperature, to chemically react resulting in atoms of Hellion and a field of single energetic Neutrons.

Symbolically such a reaction stands for

$${\mathrm{Deuterium}}^{+}+{\mathrm{Tritium}}^{+}\to {\mathrm{Hellion}}^{++}+\mathrm{Neutron}\left(\mathrm{energetic}\right).$$

We recall that the ionized Deuterium atom is comprised by a proton and a neutron and the ionized Tritium atom is comprised by a proton and two neutrons.

Moreover, the ionized Hellion atom is comprised by two protons and two neutrons.

As previously mentioned, resulting from such a chemical reaction up surges also an energetic neutron which the higher kinetics energy has a great variety of applications, including its conversion in electric energy.

We highlight the model here presented includes electric and magnetic fields and the corresponding potential ones.

Denoting by t the time on the interval $[0,t_f],$ at this point we define the following density functions:

• For the Deuterium field

  $$|{\varphi }_D{(x,y,t)|}^2=|{\varphi }_p^D{(y,t)|}^2+|{\varphi }_N^D{(x,y,t)|}^2{\left|{\varphi }_p^D(y,t)\right|}^2\frac{1}{m_p},$$
• For the Tritium field

  $$|{\varphi }_T{(x,y,t)|}^2=|{\varphi }_p^T{(y,t)|}^2+\left(\right|{\varphi }_{N_1}^T{(x,y,t)|}^2+|{\varphi }_{N_2}^T{(x,y,t)|}^2\left)\right|{\varphi }_p^T{(y,t)|}^2\frac{1}{m_p},$$
• For the Hellion field

  $$|{\varphi }_{H_e}{(x,y,t)|}^2=|{\varphi }_{2p}^{H_e}{(y,t)|}^2+\left(\right|{\varphi }_{N_1}^{H_e}{(x,y,t)|}^2+|{\varphi }_{N_2}^{H_e}{(x,y,t)|}^2\left)\right|{\varphi }_{2p}^{H_e}{(y,t)|}^2\frac{1}{2m_p},$$
• For the Neutron field

  $${\varphi }_N={\varphi }_N(x,t),$$
• For the electronic field resulting from the ionization

  $${\varphi }_e={\varphi }_e(x,y,t).$$

Furthermore, we define also the related densities

• $${\rho }_D(y,t)={\int }_{\Omega }{\left|{\varphi }_D(x,y,t)\right|}^{2}dx,$$
• $${\rho }_T(y,t)={\int }_{\Omega }{\left|{\varphi }_T(x,y,t)\right|}^{2}dx,$$

  $${\rho }_{H_e}(y,t)={\int }_{\Omega }{\left|{\varphi }_{H_e}(x,y,t)\right|}^{2}dx,$$

  $${\rho }_N(x,t)={\left|{\varphi }_N(x,t)\right|}^2,$$

  $${\rho }_e(y,t)={\int }_{\Omega }{\left|{\varphi }_e(x,y,t)\right|}^{2}dx.$$

For the chemical reaction in question we consider that one unit of mass of fractional proportion ${\alpha }_D$ of ionized Deuterium and ${\alpha }_T$ of ionized Tritium results in one unit of mass of fractional proportion ${\alpha }_{H_e}$ of ionized Hellion and ${\alpha }_N$ of neutrons.

Symbolic, this stands for

$$1={\alpha }_D+{\alpha }_T={\alpha }_{H_e}+{\alpha }_N.$$

Concerning the control volume $\Omega$ in question and related surface control $\partial \Omega ,$ we assume such a volume has an initial (fot $t=0$) amount of ionized Deuterium of ${\left(m_D\right)}_0$ and an initial amount of ionized Tritium of ${\left(m_T\right)}_0.$ The initial amount of ionized Hellion and single neutrons are supposed to be zero.

On the other hand, about the surface control $\partial \Omega$, we assume there is a part ${\Omega }_1\subset \partial \Omega$ for which is allowed the entrance and exit of Deuterium and Tritium ionized atoms.

We assume also there is another part $\partial {\Omega }_2\subset \partial \Omega$ such that $\partial {\Omega }_1\cap \partial {\Omega }_2=\varnothing$ for which is allowed only the exit of ionized Hellion atoms and neutrons, but not their entrance.

In $\partial {\Omega }_2$ is allowed the exit only (not the entrance) of ionized Deuterium and Tritium atoms.

Indeed, we assume the following relations for the masses:

• $${\left(m_{H_e,N}\right)}_T\left(t\right)=m_{H_e,N}\left(t\right)+{\int }_0^t{\int }_{\partial {\Omega }_2}({\rho }_{H_e}(x,\tau )+{\rho }_N(x,\tau ))\mathbf{u}·\mathbf{n}dSd\tau ,$$
• $$m_{H_e,N}\left(t\right)=m_{H_e}\left(t\right)+m_N\left(t\right),$$
• $$m_{H_e}\left(t\right)={\int }_{\Omega }{\rho }_{H_e}(x,t)dx,$$
• $$m_N\left(t\right)={\int }_{\Omega }{\rho }_N(x,t)dx,$$
• $$\left(m_D\right)\left(t\right)={\left(m_D\right)}_0-{\int }_0^t{\int }_{\partial {\Omega }_1\cup \partial {\Omega }_2}\left({\rho }_D(x,\tau )\right)\mathbf{u}·\mathbf{n}dSd\tau -{\alpha }_D{\left(m_{H_e,N}\right)}_T\left(t\right),$$
• $$\left(m_T\right)\left(t\right)={\left(m_T\right)}_0-{\int }_0^t{\int }_{\partial {\Omega }_1\cup \partial {\Omega }_2}\left({\rho }_T(x,\tau )\right)\mathbf{u}·\mathbf{n}dSd\tau -{\alpha }_T{\left(m_{H_e,N}\right)}_T\left(t\right),$$
• $${\left(m_e\right)}_T\left(t\right)=m_e\left(t\right)+{\int }_0^t{\int }_{\partial {\Omega }_2}\left({\rho }_e(x,\tau )\right)\mathbf{u}·\mathbf{n}dSd\tau ,$$
• $$m_e\left(t\right)={\int }_{\Omega }{\rho }_e(x,t)dx.$$
• $$m_e\left(t\right)={\int }_{\Omega }|{\varphi }_p^D{(x,t)|}^{2}dx\frac{m_e}{m_p}+{\int }_{\Omega }|{\varphi }_p^T{(x,t)|}^{2}dx\frac{m_e}{m_p}+{\int }_{\Omega }{\left|{\varphi }_{2p}^{H_e}(x,t)\right|}^{2}dx\frac{m_e}{m_p}.$$

Here $\mathbf{n}$ denotes the outward normal vectorial fields to the concerning surfaces.

Having clarified such masses relations, we define the functional

$$J(\varphi ,\rho ,\mathbf{r},\mathbf{u},\mathbf{E},\mathbf{A},\mathbf{B})$$

where

$$J=G(\nabla u)+F\left(\varphi \right)+E_c(\varphi ,\mathbf{r})+F_1+F_2+F_3,$$

and where we assume ${\gamma }_p^D>0,{\gamma }_p^T>0,{\gamma }_N^D>0,{\gamma }_{N_1}^T>0,{\gamma }_{N_2}^T>0,$${\gamma }_{2p}^{H_e}>0,{\gamma }_{N_1}^{H_e}>0,$ ${\gamma }_{N_2}^{H_e}>0,{\gamma }_N>0,{\gamma }_e>0$ and ${\alpha }_D>0,{\alpha }_T>0,{\alpha }_{H_e}>0,{\alpha }_N>0$, ${\alpha }_{DT}>0,{\alpha }_{H_{e}N}>0,{\alpha }_{e,e}>0,{\alpha }_{H_e,e}<0$ so that

$$\begin{array}{ccc}\hfill G(\nabla \varphi )& =& \frac{{\gamma }_p^D}{2}{\int }_0^{t_f}{\int }_{\Omega }(\nabla {\varphi }_p^D)·(\nabla {\varphi }_p^D)dydt\hfill \\ & & +\frac{{\gamma }_N^D}{2}{\int }_0^{t_f}{\int }_{\Omega }(\nabla {\varphi }_N^D)·(\nabla {\varphi }_N^D)dxdydt\hfill \\ & & \frac{{\gamma }_p^T}{2}{\int }_0^{t_f}{\int }_{\Omega }(\nabla {\varphi }_p^T)·(\nabla {\varphi }_p^T)dydt\hfill \\ & & +\frac{{\gamma }_{N_1}^T}{2}{\int }_0^{t_f}{\int }_{\Omega }(\nabla {\varphi }_{N_1}^T)·(\nabla {\varphi }_{N_1}^T)dxdydt\hfill \\ & & +\frac{{\gamma }_{N_2}^T}{2}{\int }_0^{t_f}{\int }_{\Omega }(\nabla {\varphi }_{N_2}^T)·(\nabla {\varphi }_{N_2}^T)dxdydt\hfill \\ & & +\frac{{\gamma }_{2p}^{H_e}}{2}{\int }_0^{t_f}{\int }_{\Omega }(\nabla {\varphi }_{2p}^{H_e})·(\nabla {\varphi }_{2p}^{H_e})dydt\hfill \\ & & +\frac{{\gamma }_{N_1}^{H_e}}{2}{\int }_0^{t_f}{\int }_{\Omega }(\nabla {\varphi }_{N_1}^{H_e})·(\nabla {\varphi }_{N_1}^{H_e})dxdydt\hfill \\ & & +\frac{{\gamma }_{N_2}^{H_e}}{2}{\int }_0^{t_f}{\int }_{\Omega }(\nabla {\varphi }_{N_2}^{H_e})·(\nabla {\varphi }_{N_2}^{H_e})dxdydt\hfill \\ & & +\frac{{\gamma }_N}{2}{\int }_0^{t_f}{\int }_{\Omega }(\nabla {\varphi }_N)·(\nabla {\varphi }_N)dxdt\hfill \\ & & +\frac{{\gamma }_e}{2}{\int }_0^{t_f}{\int }_{\Omega }(\nabla {\varphi }_e)·(\nabla {\varphi }_e)dxdydt,\hfill \end{array}$$

(79)

and

$$\begin{array}{ccc}\hfill F\left(\varphi \right)& =& \frac{{\alpha }_D}{2}{\int }_0^{t_f}{\int }_{\Omega }\frac{|{\varphi }_D(x-{\xi }_1,y-{\xi }_2,t){|}^2{\left|{\varphi }_D({\xi }_1,{\xi }_2,t)\right|}^2}{|(x,y)-({\xi }_1,{\xi }_2)|}dxdyd{\xi }_{1}d{\xi }_{2}dt\hfill \\ & & +\frac{{\alpha }_T}{2}{\int }_0^{t_f}{\int }_{\Omega }\frac{|{\varphi }_T(x-{\xi }_1,y-{\xi }_2,t){|}^2{\left|{\varphi }_T({\xi }_1,{\xi }_2,t)\right|}^2}{|(x,y)-({\xi }_1,{\xi }_2)|}dxdyd{\xi }_{1}d{\xi }_{2}dt\hfill \\ & & +\frac{{\alpha }_{DT}}{2}{\int }_0^{t_f}{\int }_{\Omega }\frac{|{\varphi }_D(x-{\xi }_1,y-{\xi }_2,t){|}^2{\left|{\varphi }_T({\xi }_1,{\xi }_2,t)\right|}^2}{|(x,y)-({\xi }_1,{\xi }_2)|}dxdyd{\xi }_{1}d{\xi }_{2}dt\hfill \\ & & +\frac{{\alpha }_{H_e}}{2}{\int }_0^{t_f}{\int }_{\Omega }\frac{|{\varphi }_{H_e}(x-{\xi }_1,y-{\xi }_2,t){|}^2{\left|{\varphi }_{H_e}({\xi }_1,{\xi }_2,t)\right|}^2}{|(x,y)-({\xi }_1,{\xi }_2)|}dxdyd{\xi }_{1}d{\xi }_{2}dt\hfill \\ & & +\frac{{\alpha }_N}{2}{\int }_0^{t_f}{\int }_{\Omega }\frac{|{\varphi }_N{(x-\xi ,t)|}^2{\left|{\varphi }_N\left(\xi \right)\right|}^2}{|x-\xi ,t|}dxd\xi dt\hfill \\ & & +\sum _{j=1}^2\frac{{\alpha }_{H_{e}N}}{2}{\int }_0^{t_f}{\int }_{\Omega }\frac{|{\varphi }_{H_e}(x_1-{\xi }_1,y-{\xi }_2,t){|}^2{\left|{\varphi }_N({\xi }_j,t)\right|}^2}{|(x,y)-({\xi }_1,{\xi }_2)|}dxdyd{\xi }_{1}d{\xi }_{2}dt\hfill \\ & & +\frac{{\alpha }_{H_e,e}}{2}{\int }_0^{t_f}{\int }_{\Omega }\frac{|{\varphi }_{H_e}(x-{\xi }_1,y-{\xi }_2,t){|}^2{\left|{\varphi }_e({\xi }_1,{\xi }_2,t)\right|}^2}{|(x,y)-({\xi }_1,{\xi }_2)|}dxdyd{\xi }_{1}d{\xi }_{2}dt\hfill \\ & & +\frac{{\alpha }_{e,e}}{2}{\int }_0^{t_f}{\int }_{\Omega }\frac{|{\varphi }_e(x-{\xi }_1,y-{\xi }_2,t){|}^2{\left|{\varphi }_e({\xi }_1,{\xi }_2,t)\right|}^2}{|(x,y)-({\xi }_1,{\xi }_2)|}dxdyd{\xi }_{1}d{\xi }_{2}dt\hfill \end{array}$$

(80)

and the internal kinetics energy is expressed by

$$\begin{array}{ccc}\hfill E_c(\varphi ,\mathbf{r})& =& \frac{1}{2}{\int }_0^{t_f}{\int }_{\Omega }{\left|{\varphi }_D\right|}^2\frac{\partial {\mathbf{r}}_D}{\partial t}·\frac{\partial {\mathbf{r}}_D}{\partial t}dxdydt\hfill \\ & & +\frac{1}{2}{\int }_0^{t_f}{\int }_{\Omega }{\left|{\varphi }_T\right|}^2\frac{\partial {\mathbf{r}}_T}{\partial t}·\frac{\partial {\mathbf{r}}_T}{\partial t}dxdydt\hfill \\ & & +\frac{1}{2}{\int }_0^{t_f}{\int }_{\Omega }{\left|{\varphi }_{H_e}\right|}^2\frac{\partial {\mathbf{r}}_{H_e}}{\partial t}·\frac{\partial {\mathbf{r}}_{H_e}}{\partial t}dxdydt\hfill \\ & & +\frac{1}{2}{\int }_0^{t_f}{\int }_{\Omega }{\left|{\varphi }_N\right|}^2\frac{\partial {\mathbf{r}}_N}{\partial t}·\frac{\partial {\mathbf{r}}_N}{\partial t}dxdydt\hfill \\ & & +\frac{1}{2}{\int }_0^{t_f}{\int }_{\Omega }{\left|{\varphi }_e\right|}^2\frac{\partial {\mathbf{r}}_e}{\partial t}·\frac{\partial {\mathbf{r}}_e}{\partial t}dxdydt,\hfill \end{array}$$

(81)

Here it is worth highlighting we have approximated the initially discrete set of indices s of particles as a continuous positive real variable s.

Moreover,

$$F_1=\frac{1}{4\pi }{\int }_0^{t_f}{\parallel \mathrm{curl}\mathbf{A}-{\mathbf{B}}_0\parallel }_{2}dt,$$

$$\begin{array}{ccc}\hfill F_2& =& {\int }_0^{t_f}{\int }_{\Omega }{\mathbf{E}}_{ind}·K_p{\left|{\varphi }_p^D\right|}^2\left(\mathbf{u}+\frac{\partial {\mathbf{r}}_D}{\partial t}\right)dxdydt\hfill \\ & & +{\int }_0^{t_f}{\int }_{\Omega }{\mathbf{E}}_{ind}·K_p{\left|{\varphi }_p^T\right|}^2\left(\mathbf{u}+\frac{\partial {\mathbf{r}}_T}{\partial t}\right)dxdydt\hfill \\ & & +{\int }_0^{t_f}{\int }_{\Omega }{\mathbf{E}}_{ind}·K_p{\left|{\varphi }_{2p}^{H_e}\right|}^2\left(\mathbf{u}+\frac{\partial {\mathbf{r}}_{H_e}}{\partial t}\right)dxdydt\hfill \\ & & +{\int }_0^{t_f}{\int }_{\Omega }{\mathbf{E}}_{ind}·K_e{\left|{\varphi }_e\right|}^2\left(\mathbf{u}+\frac{\partial {\mathbf{r}}_e}{\partial t}\right)dxdydt,\hfill \end{array}$$

(82)

where $K_p$ and $K_e$ are appropriate real constants related to the respective charges.

Here $\mathbf{u}=(u_1,u_2,u_3)$ is the fluid velocity field and

$${\mathbf{r}}_D,{\mathbf{r}}_T,{\mathbf{r}}_{H_e},{\mathbf{r}}_N,{\mathbf{r}}_e$$

are fields of displacements for the corresponding atom fields.

Also $\mathbf{A}$ denotes the magnetic potential, ${\mathbf{B}}_0$ an external magnetic field and $\mathbf{B}$ is the total magnetic field.

Moreover, ${\mathbf{E}}_{ind}$ is an induced electric field.

Finally,

$$\begin{array}{ccc}\hfill F_3& =& \frac{C_D}{2}{\int }_0^{t_f}{\int }_{\Omega }{\nabla }_{(x,y)}{\mathbf{r}}_D·{\nabla }_{(x,y)}{\mathbf{r}}_{D}dxdydt+\frac{C_T}{2}{\int }_0^{t_f}{\int }_{\Omega }{\nabla }_{(x,y)}{\mathbf{r}}_T·{\nabla }_{(x,y)}{\mathbf{r}}_{T}dxdydt\hfill \\ & & \frac{C_{H_e}}{2}{\int }_0^{t_f}{\int }_{\Omega }{\nabla }_{(x,y)}{\mathbf{r}}_{H_e}·{\nabla }_{(x,y)}{\mathbf{r}}_{H_e}dxdydt+\frac{C_N}{2}{\int }_0^{t_f}{\int }_{\Omega }{\nabla }_{(x,y)}{\mathbf{r}}_N·{\nabla }_{(x,y)}{\mathbf{r}}_{N}dxdydt\hfill \\ & & \frac{C_e}{2}{\int }_0^{t_f}{\int }_{\Omega }{\nabla }_{(x,y)}{\mathbf{r}}_e·{\nabla }_{(x,y)}{\mathbf{r}}_{e}dxdydt,\hfill \end{array}$$

(83)

for appropriate real positive constants $C_D{C}_T,C_{H_e},C_N,C_e.$

Such a functional J is subject to the following constraints:

• The momentum conservation equation for the fluid motion

  $$\rho \left(\frac{\partial u_k}{\partial t}+u_j\frac{\partial u_k}{\partial x_j}\right)=\rho f_k-\frac{\partial P}{\partial x_k}+{\tau }_{kj,j}+{\left(F_E\right)}_k+{\left(F_M\right)}_k,$$
  $\forall k\in \{1,2,3\}.$
  Here $\rho ={\rho }_D+{\rho }_T+{\rho }_{H_e}+{\rho }_N+{\rho }_e$ is the total density and P is the fluid pressure field.
  Furthermore,

  $${\tau }_{ij}=\mu \left(\frac{\partial u_i}{\partial x_j}+\frac{\partial u_j}{\partial x_i}-\frac{2}{3}{\delta }_{ij}\sum _{k=1}^3\frac{\partial u_k}{\partial x_k}\right),$$
  $\forall i,j\in \{1,2,3\},$

  $${\mathbf{F}}_E=\left\{{\left(F_E\right)}_k\right\}=\left(K_p\left(\right|{\varphi }_p^D{|}^2+|{\varphi }_p^T{|}^2+|{\varphi }_{2p}^{H_e}{|}^2)+K_e{\int }_{\Omega }{\left|{\varphi }_e\right|}^{2}dx\right)\mathbf{E},$$

  and

  $$\begin{array}{ccc}\hfill {\mathbf{F}}_M& =& \left\{{\left(F_M\right)}_k\right\}\hfill \\ & =& \left(K_p\left(|{\varphi }_p^D{|}^2\left(\mathbf{u}+\frac{\partial {\mathbf{r}}_D}{\partial t}\right)\right.\right.\hfill \\ & & |{\varphi }_p^T{|}^2\left(\mathbf{u}+\frac{\partial {\mathbf{r}}_T}{\partial t}\right)\hfill \\ & & \left.+|{\varphi }_{2p}^{H_e}{|}^2\left(\mathbf{u}+\frac{\partial {\mathbf{r}}_{H_e}}{\partial t}\right)\right)\hfill \\ & & \left.+K_e{\left|{\varphi }_e\right|}^2\left(\mathbf{u}+\frac{\partial {\mathbf{r}}_e}{\partial t}\right)\right)\times \mathbf{B}.\hfill \end{array}$$

  (84)
• Mass conservation equation:

  $$\frac{\partial \rho }{\partial t}+\mathrm{div}\left(\rho \mathbf{u}\right)=0.$$
• Energy equation

  $$\rho \frac{De}{Dt}+P\left(\mathrm{div}\mathbf{u}\right)=\frac{\partial Q}{\partial t}-\mathrm{div}\mathbf{q},$$

  where we assume the Fourier law

  $$\mathbf{q}=-K\nabla T,$$

  where $T=T(x,t)$ is the scalar field of temperature.
  Also,

  $$\begin{array}{ccc}\hfill e& =& \frac{\rho }{2}\mathbf{u}·\mathbf{u}+\frac{{\rho }_D}{2}\frac{\partial {\mathbf{r}}_D}{\partial t}·\frac{\partial {\mathbf{r}}_D}{\partial t}\hfill \\ & & +\frac{{\rho }_T}{2}\frac{\partial {\mathbf{r}}_T}{\partial t}·\frac{\partial {\mathbf{r}}_T}{\partial t}\hfill \\ & & +\frac{{\rho }_{H_e}}{2}\frac{\partial {\mathbf{r}}_{H_e}}{\partial t}·\frac{\partial {\mathbf{r}}_{H_e}}{\partial t}\hfill \\ & & +\frac{{\rho }_N}{2}\frac{\partial {\mathbf{r}}_N}{\partial t}·\frac{\partial {\mathbf{r}}_N}{\partial t}\hfill \\ & & +\frac{{\rho }_e}{2}\frac{\partial {\mathbf{r}}_e}{\partial t}·\frac{\partial {\mathbf{r}}_e}{\partial t}\hfill \end{array}$$

  (85)

  and

  $$\frac{De}{Dt}=\frac{\partial e}{\partial t}+u_j\frac{\partial e}{\partial x_j}.$$
• $$P=F_7(\rho ,T),$$

  for an appropriate scalar function $F_7$.
• Mass relations

  (a)

  $$m_D\left(t\right)={\int }_{\Omega }{\rho }_D(x,t)dx,$$

  (b)

  $$m_T\left(t\right)={\int }_{\Omega }{\rho }_T(x,t)dx,$$

  (c)

  $$m_{He}\left(t\right)={\int }_{\Omega }{\rho }_{H_e}(x,t)dx,$$

  (d)

  $$m_N\left(t\right)={\int }_{\Omega }{\rho }_N(x,t)dx,$$

  (e)

  $$m_e\left(t\right)={\int }_{\Omega }{\rho }_e(x,t)dx,$$

  where,

  (a)

  $${\left(m_{H_e,N}\right)}_T\left(t\right)=m_{H_e,N}\left(t\right)+{\int }_0^t{\int }_{\partial {\Omega }_2}\left({\rho }_{H_e}(x,\tau )\right)\mathbf{u}·\mathbf{n}dSd\tau ,$$

  (b)

  $$m_{H_e,N}\left(t\right)=m_{H_e}\left(t\right)+m_N\left(t\right),$$

  (c)

  $$\left(m_D\right)\left(t\right)={\left(m_D\right)}_0-{\int }_0^t{\int }_{\partial {\Omega }_1\cup \partial {\Omega }_2}\left({\rho }_D(x,\tau )\right)\mathbf{u}·\mathbf{n}dSd\tau -{\alpha }_D{\left(m_{H_e,N}\right)}_T\left(t\right),$$

  (d)

  $$\left(m_T\right)\left(t\right)={\left(m_T\right)}_0-{\int }_0^t{\int }_{\partial {\Omega }_1\cup \partial {\Omega }_2}\left({\rho }_T(x,\tau )\right)\mathbf{u}·\mathbf{n}dSd\tau -{\alpha }_T{\left(m_{H_e,N}\right)}_T\left(t\right),$$

  (e)

  $${\left(m_e\right)}_T\left(t\right)=m_e\left(t\right)+{\int }_0^t{\int }_{\partial {\Omega }_2}\left({\rho }_T(x,\tau )\right)\mathbf{u}·\mathbf{n}dSd\tau .$$

  (f)

  $$m_e\left(t\right)={\int }_{\Omega }|{\varphi }_p^D{(x,t)|}^{2}dx\frac{m_e}{m_p}+{\int }_{\Omega }|{\varphi }_p^T{(x,t)|}^{2}dx\frac{m_e}{m_p}+{\int }_{\Omega }{\left|{\varphi }_{2p}^{H_e}(x,t)\right|}^{2}dx\frac{m_e}{m_p}.$$

• Other mass constraints

  (a)

  $${\int }_{\Omega }{\left|{\varphi }_N^D(x,y,t)\right|}^{2}dx=m_N,$$

  (b)

  $${\int }_{\Omega }{\left|{\varphi }_{N_1}^T(x,y,t)\right|}^{2}dx=m_N,$$

  (c)

  $${\int }_{\Omega }{\left|{\varphi }_{N_2}^T(x,y,t)\right|}^{2}dx=m_N,$$

  (d)

  $${\int }_{\Omega }{\left|{\varphi }_{N_1}^{H_e}(x,y,t)\right|}^{2}dx=m_N,$$

  (e)

  $${\int }_{\Omega }{\left|{\varphi }_{N_2}^{H_e}(x,y,t)\right|}^{2}dx=m_N.$$

• For the induced electric field, we must have

  $$\begin{array}{ccc}& & \mathrm{curl}{\mathbf{E}}_{ind}+\frac{1}{c}\mathrm{curl}\left({\widehat{K}}_p{\left|{\varphi }_p^D\right|}^2\left(\mathbf{u}+\frac{\partial {\mathbf{r}}_D}{\partial t}\right)\right.\hfill \\ & & +{\widehat{K}}_p{\left|{\varphi }_p^T\right|}^2\left(\mathbf{u}+\frac{\partial {\mathbf{r}}_T}{\partial t}\right)\hfill \\ & & +{\widehat{K}}_p{\left|{\varphi }_{2p}^{H_e}\right|}^2\left(\mathbf{u}+\frac{\partial {\mathbf{r}}_{H_e}}{\partial t}\right)\hfill \\ & & +\left.{\widehat{K}}_e{\int }_{\Omega }{\left|{\varphi }_e(x,y,t)\right|}^2\left(\mathbf{u}(y,t)+\frac{\partial {\mathbf{r}}_e(x,y,t)}{\partial t}dx\right)\right)\hfill \\ & & \times \left(\mathrm{curl}\mathbf{A}-{\mathbf{B}}_0\right)-\frac{1}{c}\frac{\partial }{\partial t}\left(\mathrm{curl}\mathbf{A}-{\mathbf{B}}_0\right)=\mathbf{0},\hfill \end{array}$$

  (86)

  where ${\widehat{K}}_p$ and ${\widehat{K}}_e$ are appropriate real constants related to the respective charges.
• A Maxwell equation:

  $$\mathrm{div}\mathbf{B}=0,$$

  where

  $$\mathbf{B}={\mathbf{B}}_0-\mathrm{curl}\mathbf{A}.$$
• Another Maxwell equation:

  $$\mathrm{div}\mathbf{E}=4\pi \left(K_p\left(\right|{\varphi }_p^D{|}^2+|{\varphi }_p^T{|}^2+|{\varphi }_{2p}^{H_e}{|}^2)+K_e{\int }_{\Omega }{\left|{\varphi }_e(x,y,t)\right|}^{2}dx\right),$$

  where the total electric field $\mathbf{E}$ stands for

  $$\mathbf{E}={\mathbf{E}}_{ind}+{\mathbf{E}}_{\rho },$$

  and where generically denoting

  $$F\left(\varphi \right)={\int }_{\Omega }{f}_5(\varphi ,x,\xi )dxd\xi ,$$

  we have also

  $${\mathbf{E}}_{\rho }=\left\{{\int }_{\Omega }\frac{\partial f_5(\varphi ,x,\xi )}{\partial x_k}d\xi \right\}.$$

At this point we generically denote

$${\langle h_1,h_2\rangle }_{L^2}={\int }_0^{t_f}{\int }_{\Omega }{h}_1{h}_{2}dxdydt.$$

Thus, already including the Lagrange multipliers concerning the restrictions indicated, the extended functional $J_3$ stands for

$$\begin{array}{ccc}& & J_3=J_3(\varphi ,\mathbf{u},\mathbf{r},P,\mathbf{A},\mathbf{B},\mathbf{E},\Lambda ,E)\hfill \\ & =& G(\nabla \varphi )+F\left(\varphi \right)+E_c(\varphi ,\mathbf{r})+F_1+F_2+F_3\hfill \\ & & +{〈{\Lambda }_k,\rho \left(\frac{\partial u_k}{\partial t}+u_j\frac{\partial u_k}{\partial x_j}\right)-\rho f_k+\frac{\partial P}{\partial x_k}-{\tau }_{kj,j}-{\left(F_E\right)}_k-{\left(F_M\right)}_k〉}_{L^2}\hfill \\ & & +{〈{\Lambda }_4,\frac{\partial \rho }{\partial t}+\mathrm{div}\left(\rho \mathbf{u}\right)〉}_{L^2}+J_{Au{x}_1}+J_{Au{x}_2}+J_{Au{x}_3}+J_{Au{x}_4},\hfill \end{array}$$

(87)

where,

$$\begin{array}{ccc}\hfill J_{Au{x}_1}& =& {〈{\Lambda }_5,\rho \frac{De}{Dt}+P\left(\mathrm{div}\mathbf{u}\right)-\frac{\partial Q}{\partial t}+\mathrm{div}\mathbf{q}〉}_{L^2}\hfill \\ & & +{〈{\Lambda }_6,P-F_7(\rho ,T)〉}_{L^2},\hfill \end{array}$$

(88)

$$\begin{array}{ccc}\hfill J_{Au{x}_2}& =& {〈{\Lambda }_7,m_D\left(t\right)-{\int }_{\Omega }{\rho }_D(x,t)dx〉}_{L^2}\hfill \\ & & +{〈{\Lambda }_8,m_T\left(t\right)-{\int }_{\Omega }{\rho }_T(x,t)dx〉}_{L^2}\hfill \\ & & {〈{\Lambda }_9,m_{H_e}\left(t\right)-{\int }_{\Omega }{\rho }_{H_e}(x,t)dx〉}_{L^2}\hfill \\ & & {〈{\Lambda }_{10},m_N\left(t\right)-{\int }_{\Omega }{\rho }_N(x,t)dx〉}_{L^2}\hfill \\ & & {〈{\Lambda }_{11},m_e\left(t\right)-{\int }_{\Omega }{\rho }_e(x,t)dx〉}_{L^2},\hfill \end{array}$$

(89)

$$\begin{array}{ccc}\hfill J_{Au{x}_3}& =& -{\int }_0^{t_f}{\int }_{\Omega }{\left(E_N^D\right)}_5(y,t)\left({\int }_{\Omega }{\left|{\varphi }_N^D(x,y,t)\right|}^{2}dx-m_N\right)dydt\hfill \\ & & -{\int }_0^{t_f}{\int }_{\Omega }{\left(E_{N_1}^T\right)}_6(y,t)\left({\int }_{\Omega }{\left|{\varphi }_{N_1}^T(x,y,t)\right|}^{2}dx-m_N\right)dydt\hfill \\ & & -{\int }_0^{t_f}{\int }_{\Omega }{\left(E_{N_2}^T\right)}_7(y,t)\left({\int }_{\Omega }{\left|{\varphi }_{N_2}^T(x,y,t)\right|}^{2}dx-m_N\right)dydt\hfill \\ & & -{\int }_0^{t_f}{\int }_{\Omega }{\left(E_{N_1}^{H_e}\right)}_8(y,t)\left({\int }_{\Omega }{\left|{\varphi }_{N_1}^{H_e}(x,y,t)\right|}^{2}dx-m_N\right)dydt\hfill \\ & & -{\int }_0^{t_f}{\int }_{\Omega }{\left(E_{N_2}^{H_e}\right)}_9(y,t)\left({\int }_{\Omega }{\left|{\varphi }_{N_2}^{H_e}(x,y,t)\right|}^{2}dx-m_N\right)dydt,\hfill \end{array}$$

(90)

$$\begin{array}{ccc}\hfill J_{Au{x}_4}& =& 〈{\Lambda }_{12},\mathrm{curl}{\mathbf{E}}_{ind}\hfill \\ & & +\frac{1}{c}\mathrm{curl}\left({\widehat{K}}_p{\left|{\varphi }_p^D\right|}^2\left(\mathbf{u}+\frac{\partial {\mathbf{r}}_D}{\partial t}\right)\right.\hfill \\ & & +{\widehat{K}}_p{\left|{\varphi }_p^T\right|}^2\left(\mathbf{u}+\frac{\partial {\mathbf{r}}_T}{\partial t}\right)\hfill \\ & & +{\widehat{K}}_p{\left|{\varphi }_{2p}^{H_e}\right|}^2\left(\mathbf{u}+\frac{\partial {\mathbf{r}}_{H_e}}{\partial t}\right)\hfill \\ & & +\left.{\widehat{K}}_e{\int }_{\Omega }{\left|{\varphi }_e(x,y,t)\right|}^2\left(\mathbf{u}(y,t)+\frac{\partial {\mathbf{r}}_e(x,y,t)}{\partial t}dx\right)\right)\hfill \\ & & {\times \left(\mathrm{curl}\mathbf{A}-{\mathbf{B}}_0\right)-\frac{1}{c}\frac{\partial }{\partial t}\left(\mathrm{curl}\mathbf{A}-{\mathbf{B}}_0\right)〉}_{L^2}\hfill \\ & & +{〈{\Lambda }_{13},\mathrm{div}\mathbf{B}〉}_{L^2}\hfill \\ & & +{〈{\Lambda }_{14},\mathrm{div}\mathbf{E}-4\pi \left(K_p\left(\right|{\varphi }_p^D{|}^2+|{\varphi }_p^T{|}^2+|{\varphi }_{2p}^{H_e}{|}^2)+K_e{\int }_{\Omega }{\left|{\varphi }_e(x,y,t)\right|}^{2}dx\right)〉}_{L^2}.\hfill \end{array}$$

(91)

Here we recall the following definitions and relations:

• For the Deuterium field

  $$|{\varphi }_D{(x,y,t)|}^2=|{\varphi }_p^D{(y,t)|}^2+|{\varphi }_N^D{(x,y,t)|}^2{\left|{\varphi }_p^D(y,t)\right|}^2\frac{1}{m_p},$$
• For the Tritium field

  $$|{\varphi }_D{(x,y,t)|}^2=|{\varphi }_p^D{(y,t)|}^2+\left(\right|{\varphi }_{N_1}^D{(x,y,t)|}^2+|{\varphi }_{N_2}^D{(x,y,t)|}^2\left)\right|{\varphi }_p^D{(y,t)|}^2\frac{1}{m_p},$$
• For the Hellion field

  $$|{\varphi }_{H_e}{(x,y,t)|}^2=|{\varphi }_{2p}^{H_e}{(y,t)|}^2+\left(\right|{\varphi }_{N_1}^{H_e}{(x,y,t)|}^2+|{\varphi }_{N_2}^{H_e}{(x,y,t)|}^2\left)\right|{\varphi }_{2p}^{H_e}{(y,t)|}^2\frac{1}{2m_p},$$
• For the Neutron field

  $${\varphi }_N={\varphi }_N(x,t),$$
• For the electronic field resulting from the ionization

  $${\varphi }_e={\varphi }_e(x,y,t).$$

• $${\rho }_D(y,t)={\int }_{\Omega }{\left|{\varphi }_D(x,y,t)\right|}^{2}dx,$$
• $${\rho }_T(y,t)={\int }_{\Omega }{\left|{\varphi }_T(x,y,t)\right|}^{2}dx,$$

  $${\rho }_{H_e}(y,t)={\int }_{\Omega }{\left|{\varphi }_{H_e}(x,y,t)\right|}^{2}dx,$$

  $${\rho }_N(x,t)={\left|{\varphi }_N(x,t)\right|}^2,$$

  $${\rho }_e(y,t)={\int }_{\Omega }{\left|{\varphi }_e(x,y,t)\right|}^{2}dx.$$

Also,

$$\rho ={\rho }_D+{\rho }_T+{\rho }_{H_e}+{\rho }_N+{\rho }_e,$$

• $${\left(m_{H_e,N}\right)}_T\left(t\right)=m_{H_e,N}\left(t\right)+{\int }_0^t{\int }_{\partial {\Omega }_2}({\rho }_{H_e}(x,\tau )+{\rho }_N(x,\tau ))\mathbf{u}·\mathbf{n}dSd\tau ,$$
• $$m_{H_e,N}\left(t\right)=m_{H_e}\left(t\right)+m_N\left(t\right),$$
• $$m_{H_e}\left(t\right)={\int }_{\Omega }{\rho }_{H_e}(x,t)dx,$$
• $$m_N\left(t\right)={\int }_{\Omega }{\rho }_N(x,t)dx,$$
• $$\left(m_D\right)\left(t\right)={\left(m_D\right)}_0-{\int }_0^t{\int }_{\partial {\Omega }_1\cup \partial {\Omega }_2}\left({\rho }_D(x,\tau )\right)\mathbf{u}·\mathbf{n}dSd\tau -{\alpha }_D{\left(m_{H_e,N}\right)}_T\left(t\right),$$
• $$\left(m_T\right)\left(t\right)={\left(m_T\right)}_0-{\int }_0^t{\int }_{\partial {\Omega }_1\cup \partial {\Omega }_2}\left({\rho }_T(x,\tau )\right)\mathbf{u}·\mathbf{n}dSd\tau -{\alpha }_T{\left(m_{H_e,N}\right)}_T\left(t\right),$$
• $${\left(m_e\right)}_T\left(t\right)=m_e\left(t\right)-{\int }_0^t{\int }_{\partial {\Omega }_2}\left({\rho }_e(x,\tau )\right)\mathbf{u}·\mathbf{n}dSd\tau ,$$
• $$m_e\left(t\right)={\int }_{\Omega }{\rho }_e(x,t)dx.$$
• $$m_e\left(t\right)={\int }_{\Omega }|{\varphi }_p^D{(x,t)|}^{2}dx\frac{m_e}{m_p}+{\int }_{\Omega }|{\varphi }_p^T{(x,t)|}^{2}dx\frac{m_e}{m_p}+{\int }_{\Omega }{\left|{\varphi }_{2p}^{H_e}(x,t)\right|}^{2}dx\frac{m_e}{m_p}.$$

Finally,

$$\mathbf{E}={\mathbf{E}}_{ind}+{\mathbf{E}}_{\rho },$$

and where generically denoting

$$F\left(\varphi \right)={\int }_{\Omega }{f}_5(\varphi ,x,\xi )dxd\xi ,$$

we have also

$${\mathbf{E}}_{\rho }=\left\{{\int }_{\Omega }\frac{\partial f_5(\varphi ,x,\xi )}{\partial x_k}d\xi \right\}.$$

and,

$$\mathbf{B}={\mathbf{B}}_0-\mathrm{curl}\mathbf{A}.$$

20. A final mathematical description of the hydrogen nuclear fusion

In this section we develop in even more details another model for the hydrogen nuclear fusion.

Let $\Omega \subset {\mathbb{R}}^3$ be an open, bounded and connected set with a regular (Lipschitzian) boundary denoted by $\partial \Omega .$

Here such a set $\Omega$ stands for a control volume in which an ionized gas (plasma) flows. Such a gas comprises ionized Deuterium and Tritium atoms intended, through a suitable higher temperature, to chemically react resulting in atoms of Hellion and a field of single energetic Neutrons.

Symbolically such a reaction stands for

$${\mathrm{Deuterium}}^{+}+{\mathrm{Tritium}}^{+}\to {\mathrm{Hellion}}^{++}+\mathrm{Neutron}\left(\mathrm{energetic}\right).$$

We recall that the ionized Deuterium atom is comprised by a proton and a neutron and the ionized Tritium atom is comprised by a proton and two neutrons.

Moreover, the ionized Hellion atom is comprised by two protons and two neutrons.

As previously mentioned, resulting from such a chemical reaction up surges also an energetic neutron which the higher kinetics energy has a great variety of applications, including its conversion in electric energy.

We highlight the model here presented includes electric and magnetic fields and the corresponding potential ones.

Denoting by t the time on the interval $[0,t_f],$ at this point we define the following density functions:

• For a single Deuterium atom indexed by s:

  $$|{\varphi }_D{(x,y,t,s)|}^2=|{\varphi }_p^D{(y,t,s)|}^2+|{\varphi }_N^D{(x,y,t,s)|}^2{\left|{\varphi }_p^D(y,t,s)\right|}^2\frac{1}{m_p},$$
• For a single Tritium atom indexed by s:

  $$|{\varphi }_T{(x,y,t,s)|}^2=|{\varphi }_p^T{(y,t,s)|}^2+\left(\right|{\varphi }_{N_1}^T{(x,y,t,s)|}^2+|{\varphi }_{N_2}^T{(x,y,t,s)|}^2\left)\right|{\varphi }_p^T{(y,t,s)|}^2\frac{1}{m_p},$$
• For a single Hellion atom indexed by s:

  $$|{\varphi }_{H_e}{(x,y,t,s)|}^2=|{\varphi }_{2p}^{H_e}{(y,t,s)|}^2+\left(\right|{\varphi }_{N_1}^{H_e}{(x,y,t,s)|}^2+|{\varphi }_{N_2}^{H_e}{(x,y,t,s)|}^2\left)\right|{\varphi }_{2p}^{H_e}{(y,t,s)|}^2\frac{1}{2m_p},$$
• For the Neutron field:

  $${\varphi }_N={\varphi }_N(x,t,s),$$
• For the electronic field resulting from the ionization

  $${\varphi }_e={\varphi }_e(x,y,t,s).$$

Furthermore, we define also the related densities

• $${\rho }_D(y,t)={\int }_0^{N_D\left(t\right)}{\int }_{\Omega }{\left|{\varphi }_D(x,y,t,s)\right|}^{2}dxds,$$
• $${\rho }_T(y,t)={\int }_0^{N_T\left(t\right)}{\int }_{\Omega }{\left|{\varphi }_T(x,y,t,s)\right|}^{2}dxds,$$

  $${\rho }_{H_e}(y,t)={\int }_0^{N_{H_e}\left(t\right)}{\int }_{\Omega }{\left|{\varphi }_{H_e}(x,y,t,s)\right|}^{2}dxds,$$

  $${\rho }_N(x,t)={\int }_0^{N_N\left(t\right)}{\left|{\varphi }_N(x,t,s)\right|}^{2}ds,$$

  $${\rho }_e(y,t)={\int }_0^{N_e\left(t\right)}{\int }_{\Omega }{\left|{\varphi }_e(x,y,t,s)\right|}^{2}dxds.$$

For the chemical reaction in question we consider that one unit of mass of fractional proportion ${\alpha }_D$ of ionized Deuterium and ${\alpha }_T$ of ionized Tritium results in one unit of mass of fractional proportion ${\alpha }_{H_e}$ of ionized Hellion and ${\alpha }_N$ of neutrons.

Symbolically, this stands for

$$1={\alpha }_D+{\alpha }_T={\alpha }_{H_e}+{\alpha }_N.$$

Concerning the control volume $\Omega$ in question and related surface control $\partial \Omega ,$ we assume such a volume has an initial (fot $t=0$) amount of ionized Deuterium of ${\left(m_D\right)}_0$ and an initial amount of ionized Tritium of ${\left(m_T\right)}_0.$ The initial amount of ionized Hellion and single neutrons are supposed to be zero.

On the other hand, about the surface control $\partial \Omega$, we assume there is a part ${\Omega }_1\subset \partial \Omega$ for which is allowed the entrance and exit of Deuterium and Tritium ionized atoms.

We assume also there is another part $\partial {\Omega }_2\subset \partial \Omega$ such that $\partial {\Omega }_1\cap \partial {\Omega }_2=\varnothing$ for which is allowed only the exit of ionized Hellion atoms and neutrons, but not their entrance.

In $\partial {\Omega }_2$ is allowed the exit only (not the entrance) of ionized Deuterium and Tritium atoms.

Indeed, we assume the following relations for the masses:

• $${\left(m_{H_e,N}\right)}_T\left(t\right)=m_{H_e,N}\left(t\right)+{\int }_0^t{\int }_{\partial {\Omega }_2}({\rho }_{H_e}(x,\tau )+{\rho }_N(x,\tau ))\mathbf{u}·\mathbf{n}dSd\tau ,$$
• $$m_{H_e,N}\left(t\right)=m_{H_e}\left(t\right)+m_N\left(t\right),$$
• $$m_{H_e}\left(t\right)={\int }_{\Omega }{\rho }_{H_e}(x,t)dx,$$
• $$m_N\left(t\right)={\int }_{\Omega }{\rho }_N(x,t)dx,$$
• $$\left(m_D\right)\left(t\right)={\left(m_D\right)}_0-{\int }_0^t{\int }_{\partial {\Omega }_1\cup \partial {\Omega }_2}\left({\rho }_D(x,\tau )\right)\mathbf{u}·\mathbf{n}dSd\tau -{\alpha }_D{\left(m_{H_e,N}\right)}_T\left(t\right),$$
• $$\left(m_T\right)\left(t\right)={\left(m_T\right)}_0-{\int }_0^t{\int }_{\partial {\Omega }_1\cup \partial {\Omega }_2}\left({\rho }_T(x,\tau )\right)\mathbf{u}·\mathbf{n}dSd\tau -{\alpha }_T{\left(m_{H_e,N}\right)}_T\left(t\right),$$
• $${\left(m_e\right)}_T\left(t\right)=m_e\left(t\right)+{\int }_0^t{\int }_{\partial {\Omega }_2}\left({\rho }_e(x,\tau )\right)\mathbf{u}·\mathbf{n}dSd\tau ,$$
• $$m_e\left(t\right)={\int }_{\Omega }{\rho }_e(x,t)dx.$$
• $$\begin{array}{ccc}\hfill m_e\left(t\right)& =& {\int }_0^{N_D\left(t\right)}{\int }_{\Omega }|{\varphi }_p^D{(y,t,s)|}^{2}dyds\frac{m_e}{m_p}+{\int }_0^{N_T\left(t\right)}{\int }_{\Omega }{\left|{\varphi }_p^T(y,t,s)\right|}^{2}dyds\frac{m_e}{m_p}\hfill \\ & & +{\int }_0^{N_p\left(t\right)}{\int }_{\Omega }{\left|{\varphi }_{2p}^{H_e}(y,t,s)\right|}^{2}dyds\frac{m_e}{m_p}.\hfill \end{array}$$

  (92)

Here $\mathbf{n}$ denotes the outward normal vectorial fields to the concerning surfaces.

Having clarified such masses relations, denoting by $N_D\left(t\right)N_T\left(t\right),N_{H_e}\left(t\right),N_N\left(t\right),N_e\left(t\right)$ the respective indexed number of particles at time t, we define the functional

$$J(\varphi ,\rho ,\mathbf{r},\mathbf{u},\mathbf{E},\mathbf{A},\mathbf{B},\{N_D,N_T,N_{H_e},N_N,N_e\})$$

where

$$J=G(\nabla u)+F\left(\varphi \right)+E_c(\varphi ,\mathbf{r})+F_1+F_2+F_3+F_4,$$

and where we assume ${\gamma }_p^D>0,{\gamma }_p^T>0,{\gamma }_N^D>0,{\gamma }_{N_1}^T>0,{\gamma }_{N_2}^T>0,$${\gamma }_{2p}^{H_e}>0,{\gamma }_{N_1}^{H_e}>0,$ ${\gamma }_{N_2}^{H_e}>0,{\gamma }_N>0,{\gamma }_e>0$ and ${\alpha }_D>0,{\alpha }_T>0,{\alpha }_{H_e}>0,{\alpha }_N>0$, ${\alpha }_{DT}>0,{\alpha }_{H_{e}N}>0,{\alpha }_{e,e}>0,{\alpha }_{H_e,e}<0$ so that

$$\begin{array}{ccc}\hfill G(\nabla \varphi )& =& \frac{{\gamma }_p^D}{2}{\int }_0^{t_f}{\int }_0^{N_D\left(t\right)}{\int }_{\Omega }(\nabla {\varphi }_p^D)·(\nabla {\varphi }_p^D)dydsdt\hfill \\ & & +\frac{{\gamma }_N^D}{2}{\int }_0^{t_f}{\int }_0^{N_D\left(t\right)}{\int }_{\Omega }(\nabla {\varphi }_N^D)·(\nabla {\varphi }_N^D)dxdydsdt\hfill \\ & & \frac{{\gamma }_p^T}{2}{\int }_0^{t_f}{\int }_0^{N_T\left(t\right)}{\int }_{\Omega }(\nabla {\varphi }_p^T)·(\nabla {\varphi }_p^T)dydsdt\hfill \\ & & +\frac{{\gamma }_{N_1}^T}{2}{\int }_0^{t_f}{\int }_0^{N_T\left(t\right)}{\int }_{\Omega }(\nabla {\varphi }_{N_1}^T)·(\nabla {\varphi }_{N_1}^T)dxdydsdt\hfill \\ & & +\frac{{\gamma }_{N_2}^T}{2}{\int }_0^{t_f}{\int }_0^{N_T\left(t\right)}{\int }_{\Omega }(\nabla {\varphi }_{N_2}^T)·(\nabla {\varphi }_{N_2}^T)dxdydsdt\hfill \\ & & +\frac{{\gamma }_{2p}^{H_e}}{2}{\int }_0^{t_f}{\int }_0^{N_{H_e}\left(t\right)}{\int }_{\Omega }(\nabla {\varphi }_{2p}^{H_e})·(\nabla {\varphi }_{2p}^{H_e})dydsdt\hfill \\ & & +\frac{{\gamma }_{N_1}^{H_e}}{2}{\int }_0^{t_f}{\int }_0^{N_{H_e}\left(t\right)}{\int }_{\Omega }(\nabla {\varphi }_{N_1}^{H_e})·(\nabla {\varphi }_{N_1}^{H_e})dxdydsdt\hfill \\ & & +\frac{{\gamma }_{N_2}^{H_e}}{2}{\int }_0^{t_f}{\int }_0^{N_{H_e}\left(t\right)}{\int }_{\Omega }(\nabla {\varphi }_{N_2}^{H_e})·(\nabla {\varphi }_{N_2}^{H_e})dxdydsdt\hfill \\ & & +\frac{{\gamma }_N}{2}{\int }_0^{t_f}{\int }_0^{N_N\left(t\right)}{\int }_{\Omega }(\nabla {\varphi }_N)·(\nabla {\varphi }_N)dxdsdt\hfill \\ & & +\frac{{\gamma }_e}{2}{\int }_0^{t_f}{\int }_0^{N_e\left(t\right)}{\int }_{\Omega }(\nabla {\varphi }_e)·(\nabla {\varphi }_e)dxdydsdt,\hfill \end{array}$$

(93)

and

$$\begin{array}{ccc}& & F\left(\varphi \right)=\hfill \\ & & \frac{{\alpha }_D}{2}{\int }_0^{t_f}{\int }_0^{N_D\left(t\right)}{\int }_0^{N_D\left(t\right)}{\int }_{\Omega }\frac{|{\varphi }_D(x-{\xi }_1,y-{\xi }_2,t,s-s_1){|}^2{\left|{\varphi }_D({\xi }_1,{\xi }_2,t,s_1)\right|}^2}{|(x,y)-({\xi }_1,{\xi }_2)|}dxdyd{\xi }_{1}d{\xi }_{2}dsd{s}_{1}dt\hfill \\ & & +\frac{{\alpha }_T}{2}{\int }_0^{t_f}{\int }_0^{N_T\left(t\right)}{\int }_0^{N_T\left(t\right)}{\int }_{\Omega }\frac{|{\varphi }_T(x-{\xi }_1,y-{\xi }_2,t,s-s_1){|}^2{\left|{\varphi }_T({\xi }_1,{\xi }_2,t,s_1)\right|}^2}{|(x,y)-({\xi }_1,{\xi }_2)|}dxdyd{\xi }_{1}d{\xi }_{2}dsd{s}_{1}dt\hfill \\ & & +\frac{{\alpha }_{DT}}{2}{\int }_0^{t_f}{\int }_0^{N_D\left(t\right)}{\int }_0^{N_T\left(t\right)}{\int }_{\Omega }\frac{|{\varphi }_D(x-{\xi }_1,y-{\xi }_2,t,s-s_1){|}^2{\left|{\varphi }_T({\xi }_1,{\xi }_2,t,s_1)\right|}^2}{|(x,y)-({\xi }_1,{\xi }_2)|}dxdyd{\xi }_{1}d{\xi }_{2}dt\hfill \\ & & +\frac{{\alpha }_{H_e}}{2}{\int }_0^{t_f}{\int }_0^{N_{H_e}\left(t\right)}{\int }_0^{N_{H_e}\left(t\right)}{\int }_{\Omega }\frac{|{\varphi }_{H_e}(x-{\xi }_1,y-{\xi }_2,t,s-s_1){|}^2{\left|{\varphi }_{H_e}({\xi }_1,{\xi }_2,t)\right|}^2}{|(x,y)-({\xi }_1,{\xi }_2,s_1)|}dxdyd{\xi }_{1}d{\xi }_{2}dsd{s}_{1}dt\hfill \\ & & +\frac{{\alpha }_N}{2}{\int }_0^{t_f}{\int }_0^{N_N\left(t\right)}{\int }_0^{N_N\left(t\right)}{\int }_{\Omega }\frac{|{\varphi }_N(x-\xi ,t,s-s_1){|}^2{\left|{\varphi }_N(\xi ,t,s_1)\right|}^2}{|x-\xi |}dxd\xi dsd{s}_{1}dt\hfill \\ & & +\sum _{j=1}^2\frac{{\alpha }_{H_{e}N}}{2}{\int }_0^{t_f}{\int }_0^{N_{H_e}\left(t\right)}{\int }_0^{N_D\left(t\right)}{\int }_{\Omega }\frac{|{\varphi }_{H_e}(x_1-{\xi }_1,y-{\xi }_2,t){|}^2{\left|{\varphi }_N({\xi }_j,t)\right|}^2}{|(x,y)-({\xi }_1,{\xi }_2)|}dxdyd{\xi }_{1}d{\xi }_{2}dsd{s}_{1}dt\hfill \\ & & +\frac{{\alpha }_{H_e,e}}{2}{\int }_0^{t_f}{\int }_0^{N_{H_e}\left(t\right)}{\int }_0^{N_e\left(t\right)}{\int }_{\Omega }\frac{|{\varphi }_{H_e}(x-{\xi }_1,y-{\xi }_2,t,s-s_1){|}^2{\left|{\varphi }_e({\xi }_1,{\xi }_2,t,s_1)\right|}^2}{|(x,y)-({\xi }_1,{\xi }_2)|}dxdyd{\xi }_{1}d{\xi }_{2}dsd{s}_{1}dt\hfill \\ & & +\frac{{\alpha }_{e,e}}{2}{\int }_0^{t_f}{\int }_0^{N_e\left(t\right)}{\int }_0^{N_e\left(t\right)}{\int }_{\Omega }\frac{|{\varphi }_e(x-{\xi }_1,y-{\xi }_2,t,s-s_1){|}^2{\left|{\varphi }_e({\xi }_1,{\xi }_2,t,s_1)\right|}^2}{|(x,y)-({\xi }_1,{\xi }_2)|}dxdyd{\xi }_{1}d{\xi }_{2}dsd{s}_{1}dt\hfill \end{array}$$

(94)

and the internal kinetics energy is expressed by

$$\begin{array}{ccc}\hfill E_c(\varphi ,\mathbf{r})& =& \frac{1}{2}{\int }_0^{t_f}{\int }_0^{N_D\left(t\right)}{\int }_{\Omega }{\left|{\varphi }_D\right|}^2\frac{\partial {\mathbf{r}}_D}{\partial t}·\frac{\partial {\mathbf{r}}_D}{\partial t}dxdydsdt\hfill \\ & & +\frac{1}{2}{\int }_0^{t_f}{\int }_0^{N_T\left(t\right)}{\int }_{\Omega }{\left|{\varphi }_T\right|}^2\frac{\partial {\mathbf{r}}_T}{\partial t}·\frac{\partial {\mathbf{r}}_T}{\partial t}dxdydsdt\hfill \\ & & +\frac{1}{2}{\int }_0^{t_f}{\int }_0^{N_{H_e}\left(t\right)}{\int }_{\Omega }{\left|{\varphi }_{H_e}\right|}^2\frac{\partial {\mathbf{r}}_{H_e}}{\partial t}·\frac{\partial {\mathbf{r}}_{H_e}}{\partial t}dxdydsdt\hfill \\ & & +\frac{1}{2}{\int }_0^{t_f}{\int }_0^{N_N\left(t\right)}{\int }_{\Omega }{\left|{\varphi }_N\right|}^2\frac{\partial {\mathbf{r}}_N}{\partial t}·\frac{\partial {\mathbf{r}}_N}{\partial t}dxdydsdt\hfill \\ & & +\frac{1}{2}{\int }_0^{t_f}{\int }_0^{N_e\left(t\right)}{\int }_{\Omega }{\left|{\varphi }_e\right|}^2\frac{\partial {\mathbf{r}}_e}{\partial t}·\frac{\partial {\mathbf{r}}_e}{\partial t}dxdydsdt,\hfill \end{array}$$

(95)

Moreover,

$$F_1=\frac{1}{4\pi }{\int }_0^{t_f}{\parallel \mathrm{curl}\mathbf{A}-{\mathbf{B}}_0\parallel }_{2}dt,$$

$$\begin{array}{ccc}\hfill F_2& =& {\int }_0^{t_f}{\int }_0^{N_D\left(t\right)}{\int }_{\Omega }{\mathbf{E}}_{ind}·K_p{\left|{\varphi }_p^D\right|}^2\left(\mathbf{u}+\frac{\partial {\mathbf{r}}_D}{\partial t}\right)dxdydsdt\hfill \\ & & +{\int }_0^{t_f}{\int }_0^{N_T\left(t\right)}{\int }_{\Omega }{\mathbf{E}}_{ind}·K_p{\left|{\varphi }_p^T\right|}^2\left(\mathbf{u}+\frac{\partial {\mathbf{r}}_T}{\partial t}\right)dxdydsdt\hfill \\ & & +{\int }_0^{t_f}{\int }_0^{N_{H_e}\left(t\right)}{\int }_{\Omega }{\mathbf{E}}_{ind}·K_p{\left|{\varphi }_{2p}^{H_e}\right|}^2\left(\mathbf{u}+\frac{\partial {\mathbf{r}}_{H_e}}{\partial t}\right)dxdydsdt\hfill \\ & & +{\int }_0^{t_f}{\int }_0^{N_e\left(t\right)}{\int }_{\Omega }{\mathbf{E}}_{ind}·K_e{\left|{\varphi }_e\right|}^2\left(\mathbf{u}+\frac{\partial {\mathbf{r}}_e}{\partial t}\right)dxdydsdt,\hfill \end{array}$$

(96)

where $K_p$ and $K_e$ are appropriate real constants related to the respective charges.

Here $\mathbf{u}=(u_1,u_2,u_3)$ is the fluid velocity field and

$${\mathbf{r}}_D,{\mathbf{r}}_T,{\mathbf{r}}_{H_e},{\mathbf{r}}_N,{\mathbf{r}}_e$$

are fields of displacements for the corresponding particle fields.

Also $\mathbf{A}$ denotes the magnetic potential, ${\mathbf{B}}_0$ an external magnetic field and $\mathbf{B}$ is the total magnetic field.

Moreover, ${\mathbf{E}}_{ind}$ is an induced electric field.

Also,

$$\begin{array}{ccc}\hfill F_3& =& \frac{C_D}{2}{\int }_0^{t_f}{\int }_0^{N_D\left(t\right)}{\int }_{\Omega }{\nabla }_{(x,y)}{\mathbf{r}}_D·{\nabla }_{(x,y)}{\mathbf{r}}_{D}dxdydsdt\hfill \\ & & +\frac{C_T}{2}{\int }_0^{t_f}{\int }_0^{N_T\left(t\right)}{\int }_{\Omega }{\nabla }_{(x,y)}{\mathbf{r}}_T·{\nabla }_{(x,y)}{\mathbf{r}}_{T}dxdydsdt\hfill \\ & & +\frac{C_{H_e}}{2}{\int }_0^{t_f}{\int }_0^{N_{H_e}\left(t\right)}{\int }_{\Omega }{\nabla }_{(x,y)}{\mathbf{r}}_{H_e}·{\nabla }_{(x,y)}{\mathbf{r}}_{H_e}dxdydsdt\hfill \\ & & +\frac{C_N}{2}{\int }_0^{t_f}{\int }_0^{N_N\left(t\right)}{\int }_{\Omega }{\nabla }_{(x,y)}{\mathbf{r}}_N·{\nabla }_{(x,y)}{\mathbf{r}}_{N}dxdydsdt\hfill \\ & & \frac{C_e}{2}{\int }_0^{t_f}{\int }_0^{N_e\left(t\right)}{\int }_{\Omega }{\nabla }_{(x,y)}{\mathbf{r}}_e·{\nabla }_{(x,y)}{\mathbf{r}}_{e}dxdydsdt,\hfill \end{array}$$

(97)

for appropriate real positive constants $C_D{C}_T,C_{H_e},C_N,C_e.$

Finally,

$$\begin{array}{ccc}\hfill F_4& =& \frac{{\epsilon }_D}{2}{\int }_0^{t_f}{\left(\frac{\partial N_D\left(t\right)}{\partial t}\right)}^{2}dt+\frac{{\epsilon }_T}{2}{\int }_0^{t_f}{\left(\frac{\partial N_D\left(t\right)}{\partial t}\right)}^{2}dt\hfill \\ & & +\frac{{\epsilon }_N}{2}{\int }_0^{t_f}{\left(\frac{\partial N_N\left(t\right)}{\partial t}\right)}^{2}dt+\frac{{\epsilon }_{H_e}}{2}{\int }_0^{t_f}{\left(\frac{\partial N_{H_e}\left(t\right)}{\partial t}\right)}^{2}dt\hfill \\ & & +\frac{{\epsilon }_e}{2}{\int }_0^{t_f}{\left(\frac{\partial N_e\left(t\right)}{\partial t}\right)}^{2}dt,\hfill \end{array}$$

(98)

where ${\epsilon }_D,{\epsilon }_T,{\epsilon }_N,{\epsilon }_{H_e},{\epsilon }_e$ are small real positive constants.

Such a functional J is subject to the following constraints:

• The momentum conservation equation for the fluid motion

  $$\rho \left(\frac{\partial u_k}{\partial t}+u_j\frac{\partial u_k}{\partial x_j}\right)=\rho f_k-\frac{\partial P}{\partial x_k}+{\tau }_{kj,j}+{\left(F_E\right)}_k+{\left(F_M\right)}_k,$$
  $\forall k\in \{1,2,3\}.$
  Here $\rho ={\rho }_D+{\rho }_T+{\rho }_{H_e}+{\rho }_N+{\rho }_e$ is the total density and P is the fluid pressure field.
  Furthermore,

  $${\tau }_{ij}=\mu \left(\frac{\partial u_i}{\partial x_j}+\frac{\partial u_j}{\partial x_i}-\frac{2}{3}{\delta }_{ij}\sum _{k=1}^3\frac{\partial u_k}{\partial x_k}\right),$$
  $\forall i,j\in \{1,2,3\},$

  $$\begin{array}{ccc}& & {\mathbf{F}}_E=\left\{{\left(F_E\right)}_k\right\}=\hfill \\ & & \left(K_p\left({\int }_0^{N_D\left(t\right)}|{\varphi }_p^D{|}^{2}ds+{\int }_0^{N_T\left(t\right)}|{\varphi }_p^T{|}^{2}ds+{\int }_0^{N_{H_e}\left(t\right)}{\left|{\varphi }_{2p}^{H_e}\right|}^{2}ds\right)+K_e{\int }_0^{N_e\left(t\right)}{\left|{\varphi }_e\right|}^{2}ds\right)\mathbf{E},\hfill \end{array}$$

  and

  $$\begin{array}{ccc}\hfill {\mathbf{F}}_M& =& \left\{{\left(F_M\right)}_k\right\}\hfill \\ & =& \left(K_p\left({\int }_0^{N_D\left(t\right)}{\left|{\varphi }_p^D\right|}^2\left(\mathbf{u}+\frac{\partial {\mathbf{r}}_D}{\partial t}\right)ds\right.\right.\hfill \\ & & {\int }_0^{N_T\left(t\right)}{\left|{\varphi }_p^T\right|}^2\left(\mathbf{u}+\frac{\partial {\mathbf{r}}_T}{\partial t}\right)ds\hfill \\ & & \left.+{\int }_0^{N_{H_e}\left(t\right)}{\left|{\varphi }_{2p}^{H_e}\right|}^2\left(\mathbf{u}+\frac{\partial {\mathbf{r}}_{H_e}}{\partial t}\right)ds\right)\hfill \\ & & \left.+K_e{\int }_0^{N_e\left(t\right)}{\left|{\varphi }_e\right|}^2\left(\mathbf{u}+\frac{\partial {\mathbf{r}}_e}{\partial t}\right)ds\right)\times \mathbf{B}.\hfill \end{array}$$

  (99)
• Mass conservation equation:

  $$\frac{\partial \rho }{\partial t}+\mathrm{div}\left(\rho \mathbf{u}\right)=0.$$
• Energy equation

  $$\rho \frac{De}{Dt}+P\left(\mathrm{div}\mathbf{u}\right)=\frac{\partial Q}{\partial t}-\mathrm{div}\mathbf{q},$$

  where we assume the Fourier law

  $$\mathbf{q}=-K\nabla T,$$

  where $T=T(x,t)$ is the scalar field of temperature.
  Also,

  $$\begin{array}{ccc}\hfill e& =& \frac{\rho }{2}\mathbf{u}·\mathbf{u}+\frac{{\rho }_D}{2}\frac{\partial {\mathbf{r}}_D}{\partial t}·\frac{\partial {\mathbf{r}}_D}{\partial t}\hfill \\ & & +\frac{{\rho }_T}{2}\frac{\partial {\mathbf{r}}_T}{\partial t}·\frac{\partial {\mathbf{r}}_T}{\partial t}\hfill \\ & & +\frac{{\rho }_{H_e}}{2}\frac{\partial {\mathbf{r}}_{H_e}}{\partial t}·\frac{\partial {\mathbf{r}}_{H_e}}{\partial t}\hfill \\ & & +\frac{{\rho }_N}{2}\frac{\partial {\mathbf{r}}_N}{\partial t}·\frac{\partial {\mathbf{r}}_N}{\partial t}\hfill \\ & & +\frac{{\rho }_e}{2}\frac{\partial {\mathbf{r}}_e}{\partial t}·\frac{\partial {\mathbf{r}}_e}{\partial t}\hfill \end{array}$$

  (100)

  and

  $$\frac{De}{Dt}=\frac{\partial e}{\partial t}+u_j\frac{\partial e}{\partial x_j}.$$
• $$P=F_7(\rho ,T),$$

  for an appropriate scalar function $F_7$.
• Mass relations

  (a)

  $$m_D\left(t\right)={\int }_{\Omega }{\rho }_D(x,t)dx,$$

  (b)

  $$m_T\left(t\right)={\int }_{\Omega }{\rho }_T(x,t)dx,$$

  (c)

  $$m_{He}\left(t\right)={\int }_{\Omega }{\rho }_{H_e}(x,t)dx,$$

  (d)

  $$m_N\left(t\right)={\int }_{\Omega }{\rho }_N(x,t)dx,$$

  (e)

  $$m_e\left(t\right)={\int }_{\Omega }{\rho }_e(x,t)dx,$$

  where,

  (a)

  $${\left(m_{H_e,N}\right)}_T\left(t\right)=m_{H_e,N}\left(t\right)+{\int }_0^t{\int }_{\partial {\Omega }_2}\left({\rho }_{H_e}(x,\tau )\right)\mathbf{u}·\mathbf{n}dSd\tau ,$$

  (b)

  $$m_{H_e,N}\left(t\right)=m_{H_e}\left(t\right)+m_N\left(t\right),$$

  (c)

  $$\left(m_D\right)\left(t\right)={\left(m_D\right)}_0-{\int }_0^t{\int }_{\partial {\Omega }_1\cup \partial {\Omega }_2}\left({\rho }_D(x,\tau )\right)\mathbf{u}·\mathbf{n}dSd\tau -{\alpha }_D{\left(m_{H_e,N}\right)}_T\left(t\right),$$

  (d)

  $$\left(m_T\right)\left(t\right)={\left(m_T\right)}_0-{\int }_0^t{\int }_{\partial {\Omega }_1\cup \partial {\Omega }_2}\left({\rho }_T(x,\tau )\right)\mathbf{u}·\mathbf{n}dSd\tau -{\alpha }_T{\left(m_{H_e,N}\right)}_T\left(t\right),$$

  (e)

  $${\left(m_e\right)}_T\left(t\right)=m_e\left(t\right)+{\int }_0^t{\int }_{\partial {\Omega }_2}\left({\rho }_T(x,\tau )\right)\mathbf{u}·\mathbf{n}dSd\tau .$$

  (f)

  $$\begin{array}{ccc}\hfill m_e\left(t\right)& =& {\int }_0^{N_D\left(t\right)}{\int }_{\Omega }|{\varphi }_p^D{(y,t,s)|}^{2}dydyds\frac{m_e}{m_p}+{\int }_0^{N_T\left(t\right)}{\int }_{\Omega }{\left|{\varphi }_p^T(y,t,s)\right|}^{2}dyds\frac{m_e}{m_p}\hfill \\ & & +{\int }_0^{N_p\left(t\right)}{\int }_{\Omega }{\left|{\varphi }_{2p}^{H_e}(y,t,s)\right|}^{2}dyds\frac{m_e}{m_p}.\hfill \end{array}$$

  (101)

• Other mass constraints

  (a)

  $${\int }_{\Omega }{\left|{\varphi }_N^D(x,y,t,s)\right|}^{2}dx=m_N,$$

  (b)

  $${\int }_{\Omega }{\left|{\varphi }_{N_1}^T(x,y,t,s)\right|}^{2}dx=m_N,$$

  (c)

  $${\int }_{\Omega }{\left|{\varphi }_{N_2}^T(x,y,t,s)\right|}^{2}dx=m_N,$$

  (d)

  $${\int }_{\Omega }{\left|{\varphi }_{N_1}^{H_e}(x,y,t,s)\right|}^{2}dx=m_N,$$

  (e)

  $${\int }_{\Omega }{\left|{\varphi }_{N_2}^{H_e}(x,y,t,s)\right|}^{2}dx=m_N,$$

  (f)

  $${\int }_{\Omega }{\left|{\varphi }_p^D(x,t,s)\right|}^{2}dx=m_p,$$

  (g)

  $${\int }_{\Omega }{\left|{\varphi }_p^T(x,t,s)\right|}^{2}dx=m_p,$$

  (h)

  $${\int }_{\Omega }{\left|{\varphi }_{2p}^{H_e}(x,t,s)\right|}^{2}dx=2m_p,$$

• $$m_D\left(t\right)=m_p{N}_D\left(t\right)+m_N{N}_D\left(t\right)$$

  $$m_T\left(t\right)=m_p{N}_T\left(t\right)+m_N{N}_T\left(t\right),$$

  $$m_{H_e}\left(t\right)=2m_p{N}_{H_e}\left(t\right)+2m_N{N}_{H_e}\left(t\right),$$

  $$m_e\left(t\right)=m_e{N}_D\left(t\right)+m_e{N}_T\left(t\right)+2m_e{N}_{H_e}\left(t\right).$$
• For the induced electric field, we must have

  $$\begin{array}{ccc}& & \mathrm{curl}{\mathbf{E}}_{ind}+\frac{1}{c}\mathrm{curl}\left({\widehat{K}}_p{\int }_0^{N_D\left(t\right)}{\left|{\varphi }_p^D\right|}^2\left(\mathbf{u}+\frac{\partial {\mathbf{r}}_D}{\partial t}\right)ds\right.\hfill \\ & & +{\widehat{K}}_p{\int }_0^{N_T\left(t\right)}{\left|{\varphi }_p^T\right|}^2\left(\mathbf{u}+\frac{\partial {\mathbf{r}}_T}{\partial t}\right)ds\hfill \\ & & +{\widehat{K}}_p{\int }_0^{N_{H_e}\left(t\right)}{\left|{\varphi }_{2p}^{H_e}\right|}^2\left(\mathbf{u}+\frac{\partial {\mathbf{r}}_{H_e}}{\partial t}\right)ds\hfill \\ & & +\left.{\widehat{K}}_e{\int }_0^{N_e\left(t\right)}{\int }_{\Omega }{\left|{\varphi }_e(x,y,t,s)\right|}^2\left(\mathbf{u}(y,t)+\frac{\partial {\mathbf{r}}_e(x,y,t)}{\partial t}dx\right)ds\right)\hfill \\ & & \times \left(\mathrm{curl}\mathbf{A}-{\mathbf{B}}_0\right)-\frac{1}{c}\frac{\partial }{\partial t}\left(\mathrm{curl}\mathbf{A}-{\mathbf{B}}_0\right)=\mathbf{0},\hfill \end{array}$$

  (102)

  where ${\widehat{K}}_p$ and ${\widehat{K}}_e$ are appropriate real constants related to the respective charges.
• A Maxwell equation:

  $$\mathrm{div}\mathbf{B}=0,$$

  where

  $$\mathbf{B}={\mathbf{B}}_0-\mathrm{curl}\mathbf{A}.$$
• Another Maxwell equation:

  $$\begin{array}{ccc}\hfill \mathrm{div}\mathbf{E}& =& 4\pi \left(K_p\left({\int }_0^{N_D\left(t\right)}|{\varphi }_p^D{|}^{2}ds+{\int }_0^{N_T\left(t\right)}|{\varphi }_p^T{|}^{2}ds+{\int }_0^{N_{H_e}\left(t\right)}{\left|{\varphi }_{2p}^{H_e}\right|}^{2}ds\right)\right.\hfill \\ & & \left.+K_e{\int }_0^{N_e\left(t\right)}{\int }_{\Omega }{\left|{\varphi }_e(x,y,t,s)\right|}^{2}dxds\right),\hfill \end{array}$$

  (103)

  where the total electric field $\mathbf{E}$ stands for

  $$\mathbf{E}={\mathbf{E}}_{ind}+{\mathbf{E}}_{\rho },$$

  and where generically denoting

  $$F\left(\varphi \right)={\int }_{\Omega }{f}_5(\varphi ,x,t\xi ,s)dxd\xi ds,$$

  we have also

  $${\mathbf{E}}_{\rho }=\left\{{\int }_{\Omega }\frac{\partial f_5(\varphi ,x,t,\xi ,s)}{\partial x_k}d\xi ds\right\}.$$

At this point we generically denote

$${\langle h_1,h_2\rangle }_{L^2}={\int }_0^{t_f}{\int }_{\Omega }{h}_1{h}_{2}dxdydt.$$

Thus, already including the Lagrange multipliers concerning the restrictions indicated, the extended functional $J_3$ stands for

$$\begin{array}{ccc}& & J_3=J_3(\varphi ,\mathbf{u},\mathbf{r},P,\mathbf{A},\mathbf{B},\mathbf{E},\Lambda ,E,\{N_D,N_T,N_{H_e},N_N,N_e\})\hfill \\ & =& G(\nabla \varphi )+F\left(\varphi \right)+E_c(\varphi ,\mathbf{r})+F_1+F_2+F_3+F_4\hfill \\ & & +{〈{\Lambda }_k,\rho \left(\frac{\partial u_k}{\partial t}+u_j\frac{\partial u_k}{\partial x_j}\right)-\rho f_k+\frac{\partial P}{\partial x_k}-{\tau }_{kj,j}-{\left(F_E\right)}_k-{\left(F_M\right)}_k〉}_{L^2}\hfill \\ & & +{〈{\Lambda }_4,\frac{\partial \rho }{\partial t}+\mathrm{div}\left(\rho \mathbf{u}\right)〉}_{L^2}+J_{Au{x}_1}+J_{Au{x}_2}+J_{Au{x}_3}+J_{Au{x}_4}+J_{Au{x}_5},\hfill \end{array}$$

(104)

where,

$$\begin{array}{ccc}\hfill J_{Au{x}_1}& =& {〈{\Lambda }_5,\rho \frac{De}{Dt}+P\left(\mathrm{div}\mathbf{u}\right)-\frac{\partial Q}{\partial t}+\mathrm{div}\mathbf{q}〉}_{L^2}\hfill \\ & & +{〈{\Lambda }_6,P-F_7(\rho ,T)〉}_{L^2},\hfill \end{array}$$

(105)

$$\begin{array}{ccc}\hfill J_{Au{x}_2}& =& {〈{\Lambda }_7,m_D\left(t\right)-{\int }_{\Omega }{\rho }_D(x,t)dx〉}_{L^2}\hfill \\ & & +{〈{\Lambda }_8,m_T\left(t\right)-{\int }_{\Omega }{\rho }_T(x,t)dx〉}_{L^2}\hfill \\ & & {〈{\Lambda }_9,m_{H_e}\left(t\right)-{\int }_{\Omega }{\rho }_{H_e}(x,t)dx〉}_{L^2}\hfill \\ & & {〈{\Lambda }_{10},m_N\left(t\right)-{\int }_{\Omega }{\rho }_N(x,t)dx〉}_{L^2}\hfill \\ & & {〈{\Lambda }_{11},m_e\left(t\right)-{\int }_{\Omega }{\rho }_e(x,t)dx〉}_{L^2},\hfill \end{array}$$

(106)

$$\begin{array}{ccc}\hfill J_{Au{x}_3}& =& -{\int }_0^{t_f}{\int }_{\Omega }{\left(E_N^D\right)}_5(y,t,s)\left({\int }_{\Omega }{\left|{\varphi }_N^D(x,y,t,s)\right|}^{2}dx-m_N\right)dydt\hfill \\ & & -{\int }_0^{t_f}{\int }_{\Omega }{\left(E_{N_1}^T\right)}_6(y,t,s)\left({\int }_{\Omega }{\left|{\varphi }_{N_1}^T(x,y,t,s)\right|}^{2}dx-m_N\right)dydt\hfill \\ & & -{\int }_0^{t_f}{\int }_{\Omega }{\left(E_{N_2}^T\right)}_7(y,t,s)\left({\int }_{\Omega }{\left|{\varphi }_{N_2}^T(x,y,t,s)\right|}^{2}dx-m_N\right)dydt\hfill \\ & & -{\int }_0^{t_f}{\int }_{\Omega }{\left(E_{N_1}^{H_e}\right)}_8(y,t,s)\left({\int }_{\Omega }{\left|{\varphi }_{N_1}^{H_e}(x,y,t,s)\right|}^{2}dx-m_N\right)dydt\hfill \\ & & -{\int }_0^{t_f}{\int }_{\Omega }{\left(E_{N_2}^{H_e}\right)}_9(y,t,s)\left({\int }_{\Omega }{\left|{\varphi }_{N_2}^{H_e}(x,y,t,s)\right|}^{2}dx-m_N\right)dydt,\hfill \\ & & -{\int }_0^{t_f}{\int }_{\Omega }\left(E_p^D\right)(t,s)\left({\int }_{\Omega }{\left|{\varphi }_p^D(y,t,s)\right|}^{2}dy-m_p\right)dsdt,\hfill \\ & & -{\int }_0^{t_f}{\int }_{\Omega }\left(E_p^T\right)(t,s)\left({\int }_{\Omega }{\left|{\varphi }_p^T(y,t,s)\right|}^{2}dy-m_p\right)dsdt,\hfill \\ & & -{\int }_0^{t_f}{\int }_{\Omega }\left(E_{2p}^{H_e}\right)(t,s)\left({\int }_{\Omega }{\left|{\varphi }_{2p}^{H_e}(y,t,s)\right|}^{2}dy-2m_p\right)dsdt,\hfill \end{array}$$

(107)

$$\begin{array}{ccc}\hfill J_{Au{x}_4}& =& 〈{\Lambda }_{12},\mathrm{curl}{\mathbf{E}}_{ind}\hfill \\ & & +\frac{1}{c}\mathrm{curl}\left({\widehat{K}}_p{\int }_0^{N_D\left(t\right)}{\left|{\varphi }_p^D\right|}^2\left(\mathbf{u}+\frac{\partial {\mathbf{r}}_D}{\partial t}\right)ds\right.\hfill \\ & & +{\widehat{K}}_p{\int }_0^{N_T\left(t\right)}{\left|{\varphi }_p^T\right|}^2\left(\mathbf{u}+\frac{\partial {\mathbf{r}}_T}{\partial t}\right)ds\hfill \\ & & +{\widehat{K}}_p{\int }_0^{N_{H_e}\left(t\right)}{\left|{\varphi }_{2p}^{H_e}\right|}^2\left(\mathbf{u}+\frac{\partial {\mathbf{r}}_{H_e}}{\partial t}\right)ds\hfill \\ & & +\left.{\widehat{K}}_e{\int }_0^{N_e\left(t\right)}{\int }_{\Omega }{\left|{\varphi }_e(x,y,t,s)\right|}^2\left(\mathbf{u}(y,t)+\frac{\partial {\mathbf{r}}_e(x,y,t,s)}{\partial t}dx\right)ds\right)\hfill \\ & & {\times \left(\mathrm{curl}\mathbf{A}-{\mathbf{B}}_0\right)-\frac{1}{c}\frac{\partial }{\partial t}\left(\mathrm{curl}\mathbf{A}-{\mathbf{B}}_0\right)〉}_{L^2}\hfill \\ & & +{〈{\Lambda }_{13},\mathrm{div}\mathbf{B}〉}_{L^2}\hfill \\ & & +〈{\Lambda }_{14},\mathrm{div}\mathbf{E}-4\pi \left(K_p\left({\int }_0^{N_D\left(t\right)}|{\varphi }_p^D{|}^{2}ds+{\int }_0^{N_T\left(t\right)}|{\varphi }_p^T{|}^{2}ds+{\int }_0^{N_{H_e}\left(t\right)}{\left|{\varphi }_{2p}^{H_e}\right|}^{2}ds\right)\right.\hfill \\ & & {\left.+K_e{\int }_{\Omega }{\left|{\varphi }_e\right|}^{2}dxds\right)〉}_{L^2}.\hfill \end{array}$$

(108)

$$\begin{array}{ccc}\hfill J_{Au{x}_5}& =& {\langle {\Lambda }_{15},m_D\left(t\right)-(m_p{N}_D\left(t\right)+m_N{N}_D\left(t\right))\rangle }_{L^2}\hfill \\ & & +{\langle {\Lambda }_{16},m_T\left(t\right)-(m_p{N}_T\left(t\right)+m_N{N}_T\left(t\right))\rangle }_{L^2}\hfill \\ & & +{\langle {\Lambda }_{17},m_{H_e}\left(t\right)-(2m_p{N}_{H_e}\left(t\right)+2m_N{N}_{H_e}\left(t\right))\rangle }_{L^2}\hfill \\ & & +{\langle {\Lambda }_{18},m_e\left(t\right)-(m_e{N}_D\left(t\right)+m_e{N}_T\left(t\right)+2m_e{N}_{H_e}\left(t\right))\rangle }_{L^2}.\hfill \end{array}$$

(109)

Here we recall the following definitions and relations:

• For the Deuterium field

  $$|{\varphi }_D{(x,y,t,s)|}^2=|{\varphi }_p^D{(y,t,s)|}^2+|{\varphi }_N^D{(x,y,t,s)|}^2{\left|{\varphi }_p^D(y,t,s)\right|}^2\frac{1}{m_p},$$
• For the Tritium field

  $$|{\varphi }_T{(x,y,t,s)|}^2=|{\varphi }_p^T{(y,t,s)|}^2+\left(\right|{\varphi }_{N_1}^T{(x,y,t,s)|}^2+|{\varphi }_{N_2}^T{(x,y,t,s)|}^2\left)\right|{\varphi }_p^D{(y,t,s)|}^2\frac{1}{m_p},$$
• For the Hellion field

  $$|{\varphi }_{H_e}{(x,y,t,s)|}^2=|{\varphi }_{2p}^{H_e}{(y,t,s)|}^2+\left(\right|{\varphi }_{N_1}^{H_e}{(x,y,t,s)|}^2+|{\varphi }_{N_2}^{H_e}{(x,y,t,s)|}^2\left)\right|{\varphi }_{2p}^{H_e}{(y,t,s)|}^2\frac{1}{2m_p},$$
• For the Neutron field

  $${\varphi }_N={\varphi }_N(x,t,s),$$
• For the electronic field resulting from the ionization

  $${\varphi }_e={\varphi }_e(x,y,t,s).$$

• $${\rho }_D(y,t)={\int }_0^{N_D\left(t\right)}{\int }_{\Omega }{\left|{\varphi }_D(x,y,t,s)\right|}^{2}dxds,$$
• $${\rho }_T(y,t)={\int }_0^{N_T\left(t\right)}{\int }_{\Omega }{\left|{\varphi }_T(x,y,t,s)\right|}^{2}dxds,$$

  $${\rho }_{H_e}(y,t)={\int }_0^{N_{H_e}\left(t\right)}{\int }_{\Omega }{\left|{\varphi }_{H_e}(x,y,t,s)\right|}^{2}dxds,$$

  $${\rho }_N(x,t)={\int }_0^{N_N\left(t\right)}{\left|{\varphi }_N(x,t,s)\right|}^{2}ds,$$

  $${\rho }_e(y,t)={\int }_0^{N_e\left(t\right)}{\int }_{\Omega }{\left|{\varphi }_e(x,y,t,s)\right|}^{2}dxds.$$

Also,

$$\rho ={\rho }_D+{\rho }_T+{\rho }_{H_e}+{\rho }_N+{\rho }_e,$$

• $${\left(m_{H_e,N}\right)}_T\left(t\right)=m_{H_e,N}\left(t\right)+{\int }_0^t{\int }_{\partial {\Omega }_2}({\rho }_{H_e}(x,\tau )+{\rho }_N(x,\tau ))\mathbf{u}·\mathbf{n}dSd\tau ,$$
• $$m_{H_e,N}\left(t\right)=m_{H_e}\left(t\right)+m_N\left(t\right),$$
• $$m_{H_e}\left(t\right)={\int }_{\Omega }{\rho }_{H_e}(x,t)dx,$$
• $$m_N\left(t\right)={\int }_{\Omega }{\rho }_N(x,t)dx,$$
• $$\left(m_D\right)\left(t\right)={\left(m_D\right)}_0-{\int }_0^t{\int }_{\partial {\Omega }_1\cup \partial {\Omega }_2}\left({\rho }_D(x,\tau )\right)\mathbf{u}·\mathbf{n}dSd\tau -{\alpha }_D{\left(m_{H_e,N}\right)}_T\left(t\right),$$
• $$\left(m_T\right)\left(t\right)={\left(m_T\right)}_0-{\int }_0^t{\int }_{\partial {\Omega }_1\cup \partial {\Omega }_2}\left({\rho }_T(x,\tau )\right)\mathbf{u}·\mathbf{n}dSd\tau -{\alpha }_T{\left(m_{H_e,N}\right)}_T\left(t\right),$$
• $${\left(m_e\right)}_T\left(t\right)=m_e\left(t\right)-{\int }_0^t{\int }_{\partial {\Omega }_2}\left({\rho }_e(x,\tau )\right)\mathbf{u}·\mathbf{n}dSd\tau ,$$
• $$m_e\left(t\right)={\int }_{\Omega }{\rho }_e(x,t)dx.$$
• $$\begin{array}{ccc}\hfill m_e\left(t\right)& =& {\int }_0^{N_D\left(t\right)}{\int }_{\Omega }|{\varphi }_p^D{(y,t,s)|}^{2}dydyds\frac{m_e}{m_p}+{\int }_0^{N_T\left(t\right)}{\int }_{\Omega }{\left|{\varphi }_p^T(y,t,s)\right|}^{2}dyds\frac{m_e}{m_p}\hfill \\ & & +{\int }_0^{N_p\left(t\right)}{\int }_{\Omega }{\left|{\varphi }_{2p}^{H_e}(y,t,s)\right|}^{2}dyds\frac{m_e}{m_p}.\hfill \end{array}$$

  (110)

Finally,

$$\mathbf{E}={\mathbf{E}}_{ind}+{\mathbf{E}}_{\rho },$$

and where generically denoting

$$F\left(\varphi \right)={\int }_{\Omega }{f}_5(\varphi ,x,t,\xi ,s)dxd\xi ds,$$

we have also

$${\mathbf{E}}_{\rho }=\left\{{\int }_{\Omega }\frac{\partial f_5(\varphi ,x,t,\xi ,s)}{\partial x_k}d\xi ds\right\}.$$

and,

$$\mathbf{B}={\mathbf{B}}_0-\mathrm{curl}\mathbf{A}.$$

21. A qualitative modeling for a general phase transition process

In this section we develop a general qualitative modeling for a phase transition process.

Let $\Omega \subset {\mathbb{R}}^3$ be an open, bounded and connected set with a regular (Lipschitzian) boundary denoted by $\partial \Omega .$

Such a set $\Omega$ is supposed to a be a fixed volume in which an amount of mass of a substance A with a density function u will develop phase a transition for another phase with corresponding density function $v.$ The total mass $m_T$ is suppose to be kept constant throughout such a process.

We model such transition in phase through a functional $J:V\times V\to \mathbb{R}$ where

$$\begin{array}{ccc}\hfill J(u,v)& =& \frac{{\gamma }_1}{2}{\int }_{\Omega }\nabla u·\nabla udx+\frac{{\alpha }_1}{2}{\int }_{\Omega }{u}^{4}dx\hfill \\ & & \frac{{\gamma }_2}{2}{\int }_{\Omega }\nabla v·\nabla vdx+\frac{{\alpha }_2}{2}{\int }_{\Omega }{v}^{4}dx\hfill \\ & & -\frac{1}{2}{\int }_{\Omega }{\omega }^2(u^2+v^2)dx-\frac{E}{2}\left({\int }_{\Omega }(u^2+v^2)dx-m_T\right).\hfill \end{array}$$

(111)

Here ${\gamma }_1>0,{\gamma }_2>0,{\alpha }_1>0,{\alpha }_2>0$ and $V=W^{1,2}(\Omega ).$

The phases corresponding to u and v are connected through a Lagrange multiplier E, which represents the chemical potential of the chemical process in question.

We assume the temperature is directly proportional to the internal kinetics $E_C$ energy where

$$E_C=\frac{1}{2}{\int }_{\Omega }{u}^2\frac{\partial {\mathbf{r}}_u}{\partial t}·\frac{\partial {\mathbf{r}}_u}{\partial t}dx.$$

For a internal vibrational motion, we assume approximately

$${\mathbf{r}}_u\approx e^{i\omega t}{\mathbf{w}}_5\left(x\right),$$

for an appropriate frequency $\omega$ and vectorial function ${\mathbf{w}}_5.$

Thus, the temperature $T=T(x,t)$ is indeed proportional to ${\omega }^2$, that is, symbolically, we may write

$$T\propto E_1\propto {\omega }^2.$$

Therefore, we start with the system with a phase corresponding to $u\approx 1$ and $v\approx 0$ at $\omega =1$. Gradually increasing the temperature to a corresponding $\omega =15$, we obtain a transition to a phase corresponding to $u\approx 0$ and $v\approx 1$.

At this point, we also define the index normalized corresponding densities

$${\varphi }_u=\frac{u^2}{u^2+v^2}$$

and

$${\varphi }_v=\frac{v^2}{u^2+v^2}.$$

Finally, we have obtained some numerical results for the following parameters:

$\Omega =[0,1]\subset \mathbb{R},$${\gamma }_1={\gamma }_2=1$, $\alpha =0.1$, ${\alpha }_2={10}^3.$

• We start with $\omega =1$ corresponding to ${\varphi }_u\approx 1$ and ${\varphi }_v\approx 0$ in $\Omega$.
  For the corresponding solutions ${\varphi }_u$ and ${\varphi }_v$, please see Figure 15 and Figure 16, respectively.
• We end the process with $\omega =15$ corresponding to ${\varphi }_u\approx 0$ and ${\varphi }_v\approx 1$ in $\Omega$.
  For the corresponding solutions ${\varphi }_u$ and ${\varphi }_v$, please see Figure 17 and Figure 18, respectively.