text:
\begin{document}

\maketitle

\begin{abstract}

In this paper a new observer is introduced to estimate the Crystal Size Distribution (CSD) only from the measurements of the solute concentration, temperature
\startlucas
and a model of the growth rate.
No model of the nucleation rate is needed.
This approach is based on the use of a Kazantzis-Kravaris/Luenberger observer which exponentially estimates functionals of the CSD.
Then, the full state is estimated by means of a Tikhonov regularization procedure.
Numerical simulations are provided.
Our approach relies on an infinite-dimensional observer, contrarily to the usual moment based observers.
\stoplucas

\end{abstract}

{\small
Keywords:
Observers, Crystallization,
\startlucas
Crystal Size Distribution,
\stoplucas
Tikhonov regularization, Solute concentration
}

\setlength{\baselineskip}{1.6\baselineskip}
\section{Introduction}

Crystallization is one of the oldest and major processes used in industry (chemical, pharmaceutical, food, etc.) to produce, purify or separate solid compounds or products \cite{kleinboistelle}.
This unit operation aims to produce solid crystals with well defined specifications including (among others) the Crytal Size Distribution (CSD) which is of critical importance.
At the industrial scale, the CSD is neither well controlled nor monitored during the crystallization process and a grinding step is usually performed before delivering the final product.
Hence,
a key point is the real-time control and the ``\emph{online}" monitoring of the CSD during the crystallization operation in order to avoid the grinding step.
Unfortunately,
the CSD is nowadays not directly measurable in real-time by existing sensors.
Nevertheless,
the Process Analytical Technologies (PATs) allow us to get access to real time information such as the solute concentration based on the Attenuated Total Reflectance Fourier Transform InfraRed spectroscopy (ATR-FTIR) and the Chord Length Distribution based on the Focused Beam Reflectance Measurement (\fbrmp).
\startlucas
 \fbrm is commonly employed in pharmaceutical industries to detect online some process deviations.
 Note however that its use for online estimation of the CSD remains challenging \cite{Poblete,scheler2013ray}.

Obtaining
an online
CSD estimation from the solute concentration, temperature and a growth rate model only is an interesting problem which is the purpose of the present paper.
To obtain online state estimation for dynamical systems from measured outputs, control engineers usually employ asymptotic state  observers.
Designing state observer for complex dynamical systems is an active research area. Some studies are devoted to nonlinear dynamics (see \cite{Bernard_Book_2019} and  references therein). Others consider infinite dimensional systems (see for instance \cite{TW}).

Designing observers for a batch crystallization process has been addressed in recent years by several researchers (see for instance \cite{mesbah2011real, mesbah2011comparison,motz2008state, nagy2003robust, porru2017monitoring, Uccheddu, vissers2012model}).

These studies are based
on spatial discretization of the PDE (see \cite{porru2017monitoring} or \cite{vissers2012model})
or more frequently on moments analysis (see survey \cite{mesbah2011comparison}).

However, the moment based approaches suffer from several drawbacks: (i) the numerical moment values exhibit a very large difference in their order of magnitude \cite{Thevenin}. Consequently, a small numerical/experimental error has a significant impact on the moments estimation. The error will then be propagated during the computation of the moment transport equations \cite{Jakobsen};
(ii) the recovering of the CSD from a finite number of its moments is still an open area of research in mathematics and highly dependant of the moments quality \cite{Lebaz, Rohani}.
\stoplucas
\startlucas
Furthermore, the
\stoplucas
moment based observers depend on the knowledge of the nucleation rate which may be tricky to model.
\startlucas
All these
\stoplucas
reasons may explain the difficulty to develop
an efficient algorithm following this route.

In this paper,
another approach is adopted to describe the CSD without using its moments.
Moreover,
a key point of the proposed approach is that no information is needed on the nucleation rate.

In the first part of the paper, we recall the \startlucas Kazantzis-Kravaris/Luenberger \stoplucas methodology for observer design which has been introduced in its original form in \cite{luenberger1964observing} for linear systems and adapted for nonlinear dynamics in \startlucas
\cite{andrieu2006existence}  and \stoplucas \cite{kazantzis1998nonlinear}.
\startlucas
Section~\ref{secappli}
\stoplucas

is devoted to the modeling of the batch crystallization process.
Finally some simulation results are presented
\startlucas
in Section~\ref{secsimu}
\stoplucas
which highlight the practical interest of the suggested methodology.

\section{Kazantzis-Kravaris/Luenberger observer methodology}
\label{secluen}

\startlucas
In control engineering, the algorithm which is employed to reconstruct online missing data on a partially measured dynamical process is named an observer.
Observer design for nonlinear dynamical systems is a very active research area and has been the subject of numerous studies in the past 40 years (see \cite{Bernard_Book_2019} for a recent survey on this topic).
\stoplucas
One very efficient way to design observers for finite-dimensional linear as well as nonlinear dynamical systems
\startlucas is the Kazantzis-Karavis/Lunberger (KKL) observer methodology. It
\stoplucas
is an approach which follows Luenberger original idea in \cite{luenberger1964observing} where asymptotic observers for linear systems were introduced for the first time.
The approach of \cite{luenberger1964observing} is somehow different from the way \textit{Luenberger observers} are introduced nowadays (which follows Luenberger second paper on observers \cite{luenberger1971introduction}).
It has been recently employed by Kazantzis and Kravaris in \cite{kazantzis1998nonlinear} to design local observers for finite-dimensional nonlinear dynamics based on Luenberger original idea in a local version. We hope it may also be fruitful in the context of distributed parameter systems.
In this section we adapt to the context of linear abstract Cauchy problem the \startlucas KKL \stoplucas observer methodology which
is a two-steps design procedure.
In the first step we consider the problem of estimating a function of the state.
The reconstruction of the overall state of the system is obtained in the second step by inverting this function.

\subsection{Step 1: reconstruction of a function of the state}

As shown in \cite{andrieu2006existence} and \cite{bernard2019luenberger}, it is always possible to exponentially estimate a function of the state, even for nonlinear finite-dimensional systems,
\startlucas
that will carry enough information about the state to estimate it in Step 2.
\stoplucas

In order to do so,
it is sufficient to introduce an auxiliary dynamical system fed by the measured output such that its solutions provide an estimation of this function of the state.

\startlucas
Let $\RR$ (resp. $\RR_+$, $\RR_-$ and $\RR^*$) denotes the set of real numbers (resp. non-negative, non-positive and non-zero real numbers).
Let $\XR$ be a Hilbert space. We denote $\mathcal L (\XR)$ the set of all endomorphisms of $\XR$ and $\mathcal L (\XR, \RR)$ the space of linear forms from $\XR$ to $\RR$. Consider the abstract Cauchy problem on $\XR$
\stoplucas

\begin{equation}\label{eq_AbsCauchyProb}
\dot \psi = \FR \psi, \quad\psi(0) = \psi_0,
\end{equation}
where
$\FR: D(\FR)\subset\XR\rightarrow \XR$ is a linear operator which is the generator of a strongly continuous semi-group denoted $(T_t)_{t\in\RR_+}$ in $\mathcal L (\XR)$
\startlucas
and $\psi_0\in D(\FR)$.
\stoplucas
Let $\rho(\FR) = \{\lambda \in \CC \mid (\FR-\lambda \Id)^{-1}\in \mathcal L (\XR)\}$ denote  the resolvent set of $\FR$.
Moreover, consider a bounded output operator
\begin{equation}\label{eq_OutputAbstract}
y = \HR \psi    ,
\end{equation}

\startlucas
where $\HR\in\mathcal L (\XR, \RR)$ is bounded.
\stoplucas
Following the \startlucas KKL
\stoplucas
methodology, we obtain the following proposition.
\begin{proposition}\label{propobs}
For all $\lambda$ in $\rho(\FR)\cap \startlucas\RR_-^*\stoplucas$, let $\TR_\lambda$ in $\mathcal L(\XR,\RR)$ be the operator defined as
$$
\TR_{\lambda}: \XR\ni \psi \mapsto \HR(\FR-\lambda \Id)^{-1}\psi\in \RR\ .
$$
Then, the dynamical system
\begin{equation}\label{eq_Obs}
\dot z_\lambda = \lambda z_\lambda + y,
\end{equation}
is an
\startlucas
exponential
\stoplucas
observer for $\TR_{\lambda}\psi$. More precisely, for all $(\psi_0,z_0)$ in $D(\FR)\times\RR$, it yields for all $t\geq0$

\startlucas
\begin{align}
\TR_{\lambda}(T_t \psi_0)- z_\lambda(t) = \exp(\lambda t)\left(\TR_{\lambda}(\psi_0)- z_0\right).
\label{eq_vitconv}
\end{align}
\stoplucas
where $z_\lambda:\RR_+\rightarrow\RR$ is the solution of system (\ref{eq_Obs}) when $y$ is given by (\ref{eq_OutputAbstract}) and initiated from $z_0$.
\end{proposition}
\begin{proof}
Let $\psi_0$ be in $D(\FR)$.
Equations \eqref{eq_AbsCauchyProb}--\eqref{eq_Obs} yield

\begin{align*}
\frac{d}{dt}\big(\TR_{\lambda}(T_t \psi_0)- z_\lambda(t)\big)
&=  \TR_{\lambda}(\FR T_t \psi_0) - \lambda z_\lambda(t) - \HR T_t \psi_0
\\

&=  \TR_{\lambda}(\FR - \lambda \Id) T_t \psi_0  +  \lambda \big(\TR_{\lambda}(T_t \psi_0)  -  z_\lambda(t) \big) - \HR T_t \psi_0
\\

&= \lambda (\TR_{\lambda}(T_t\psi_0) -  z_\lambda(t)),
\end{align*}

where the last equality follows since $\TR_{\lambda}(\FR - \lambda \Id) = \HR$.
\startlucas
Hence, \eqref{eq_vitconv} follows by integrating in time the former equation.
\stoplucas

\startlucas
Keeping in mind that $\lambda$ is negative in Proposition~\ref{propobs}, \eqref{eq_vitconv} implies
\begin{equation}\label{eq_Conv}
\lim_{t\rightarrow +\infty} \left|\TR_{\lambda}(T_t \psi_0)- z_\lambda(t)\right| = 0 .
\end{equation}
This ends the proof.
\stoplucas
\end{proof}
\startlucas
\begin{remark}
The operator $\TR_{\lambda}$ is solution to the Sylvester equation~:
\begin{equation}
\FR\TR_{\lambda} = \lambda \TR_\lambda + \HR .
\end{equation}
We recognize here the algebraic equation which was already given in Luenberger seminal paper \cite{luenberger1964observing} and which becomes a nonlinear partial differential equation in \cite{andrieu2006existence}.
\end{remark}
\stoplucas
\subsection{Step 2: reconstruction of the entire state of the system}
According to step 1, we can easily estimate $\TR_{\lambda}\psi$ for all $\lambda$ in $\rho(\FR)\cap \startlucas\RR_-^*\stoplucas$ via the observer system (\ref{eq_Obs}).
The idea of the \startlucas KKL
\stoplucas
observer methodology is to consider the mapping
\(
\TR:\XR\mapsto \RR^p
\)
given by $\psi \mapsto (\TR_{\lambda_1} \psi, \dots, \TR_{\lambda_p} \psi)$

which will be

\startlucas
exponentially
\stoplucas
estimated along the trajectory of \eqref{eq_AbsCauchyProb} via a bench of observers of the form (\ref{eq_Obs}).
To solve the estimation problem, the question is to solve the inverse problem
\begin{equation}\label{E:equality-constraint}
\TR \hat \psi = z
\end{equation}
with the unknown  $\hat \psi$ in $\XR$.
\startlucas
Let $\Im\TR=\{\TR x\mid x\in\XR\}$ be the image of $\TR$.
Then \eqref{E:equality-constraint} admits a solution only if $z$ is in $\Im\TR$.
\stoplucas

The former condition is in general too restrictive due to the fact that $z$ is only an estimation of $\TR \psi$ and consequently may not be in $\Im\TR$.

A solution to overcome this problem is to replace the equality constraint \eqref{E:equality-constraint} by the minimization problem
\begin{equation}\label{pb:min_ill}
\textrm{minimize~}\|\TR \psi - z\|^2 \textrm{~subject to~} \psi \in \XR.
\end{equation}
The set of solutions of \eqref{pb:min_ill} is denoted by
\(
\argmin_{\psi \in \XR}\|\TR \psi - z\|^2\).

The following results can be found in \cite[Chapter 4]{kern2016methodes}:
\begin{itemize}
    \item If $\TR$ is injective, then \eqref{pb:min_ill} has at most one solution.
    \item If $z\in \Im\TR \oplus \left(\Im\TR\right)^\perp$, then the set $\argmin_{\psi \in \XR}\|\TR \psi - z\|^2$ is closed, convex and non-empty (in particular \eqref{pb:min_ill} admits at least one solution).
    \item If $\TR$ is bijective and admits a left inverse denoted $\TR^{-1}$ then the unique solution of \eqref{pb:min_ill} is $\hat \psi(t)=\TR^{-1}z(t)$.
\end{itemize}

For finite dimensional systems, the injectivity of $\TR$ is directly linked with an observability property of the dynamical system \eqref{eq_AbsCauchyProb}.
Indeed, following \cite{luenberger1964observing}, it can be shown that for finite dimensional systems,
if the pair $(\FR,\HR)$ is observable (in the sense that the Kalman observability matrix is full rank), then picking $p$ pairwise distinct $\lambda_i$'s with $p = \dim \XR$, then the obtained mapping $\TR$ is invertible.
In that case, a \startlucas KKL
\stoplucas
 observer can simply be obtained as
$$
\dot z =\begin{bmatrix}\lambda_1 & & \\
&\ddots&\\
&&\lambda_p\end{bmatrix}  z +\begin{bmatrix}1\\\vdots\\1\end{bmatrix} y\ ,\quad \hat \psi(t) = \TR^{-1} z(t).
$$

In the infinite dimensional case, there is no hope that a finite number of $\lambda_i$'s can give all the information allowing to reconstruct the entire system's state. It is an open question to know if observability properties of the infinite dimensional system allows to obtain that a certain countable set of $\lambda_i$'s gives enough information.
In that case, the minimization problem \eqref{pb:min_ill} becomes ill-posed.
\startlucas
In other words, it may have no solution, or numerous solutions, and its solutions may depend on $z$ in a non-continuous manner.
In our case, $z$ is an estimation of $\TR\psi$. Hence, the minimizer $\hat\psi$ of \eqref{pb:min_ill} may be very different of the real state $\psi$.
A typical approach to overcome this problem is to consider a Tikhonov regularization method in which the optimization problem \eqref{pb:min_ill} is slightly modified.
This new minimization problem will have solutions close to the former one, and be well-posed.
We recall the next proposition on which this method is based for the convenience of the reader.
\stoplucas
\begin{proposition}[\cite{kern2016methodes}, Proposition \startlucas 6\stoplucas.1]
Given $\TR$ in $\mathcal L(\XR,\RR^p)$ and $\delta>0$,
the minimization problem
\begin{equation}\label{pb:min_well}
\textrm{minimize~} \|\TR \psi - z\|^2+ \delta \|\psi \|^2 \textrm{~subject to~} \psi \in \XR
\end{equation}
admits a unique solution and is a continuous function of $z$ in $\RR^p$.
\end{proposition}

The choice of the \emph{regularization parameter} $\delta$ is of the uttermost importance. Indeed, as $\delta\to+\infty$, the solution of \eqref{pb:min_well} goes to zero. Conversely, as $\delta\to0$, the problem \eqref{pb:min_well} goes closer to the ill-posed problem \eqref{pb:min_ill}. A compromise shall be made. In practice, $\delta$ is chosen experimentally, and is linked to the measure confidence: the more uncertain is the output $y$, the bigger is $\delta$.

We also have the following theoretical result, that describes what happens when $\delta$ goes to zero.

\begin{theorem}[\protect{\cite[Theorem \startlucas 6\stoplucas.1]{kern2016methodes}}]
Let $\TR$ in $\mathcal L(\XR,\RR^p)$ and $z \in \Im \TR$. Let $\psi_0\in \XR$ and $\psi$ the solution of \eqref{pb:min_ill} closest to $\psi_0$. Let $(z_n)_{n\in \NN}$ be a sequence in $\RR^p$ converging to $z$.
Let $\varepsilon_n = |z_n-z|$. Let $(\delta_n)_{n\in \NN}$ be a sequence of regularization parameters converging to zero. For any $n\in \NN$, let $\psi_n$ be the solution of the problem \eqref{pb:min_well} associated to $z_n$ and $\delta_n$. Then,
\begin{itemize}
\item $\displaystyle|\TR \psi_n-z_n|\underset{n\rightarrow+\infty}{\longrightarrow}0$;
\item
if
\(\displaystyle
\frac{\varepsilon_n}{\delta_n}\underset{n\rightarrow+\infty}{\longrightarrow}0,
\)
then
\(|\TR\psi_n-z_n| = \OO(\varepsilon_n)\) and \(
\psi_n\underset{n\rightarrow+\infty}{\longrightarrow}\psi;
\)

\item
if
\(\displaystyle
\frac{\varepsilon_n}{\delta_n^2}\underset{n\rightarrow+\infty}{\longrightarrow}0
\)
and
\(
\psi\in\left(\Im\TR\right)^*,
\)
then
\(
|\TR\psi_n-z_n| = \OO(\varepsilon_n^2)
\)
and
\(
|\psi_n-\psi| = \OO(\varepsilon_n).
\)

\end{itemize}
\end{theorem}
\startlucas
To summarize, a possible observer design of a given abstract Cauchy problem in the form \eqref{eq_AbsCauchyProb} is then given by
\begin{equation}\label{eq_Obs_KKL_Tikhonov}
\left\{
\begin{aligned}
&\dot z =\begin{bmatrix}\lambda_1 & & \\
&\ddots&\\
&&\lambda_p\end{bmatrix}  z +\begin{bmatrix}1\\\vdots\\1\end{bmatrix} y
\\
&\hat \psi(t)
= \argmin_{\psi \in \XR}\left\{\|\TR \psi - z(t)\|^2
+ \delta \|\psi\|^2\right\}, \quad \delta >0
\\
&\TR
= \left(\TR_{\lambda_1}, \dots, \TR_{\lambda_p}\right), \quad \TR_{\lambda_i}  = \HR(\FR - \lambda_i \Id)^{-1}, \quad
i=1,\dots, p,
\end{aligned}\right.
\end{equation}
where the $\lambda_i$ are  pairwise distinct elements of $\rho(\FR)\cap \startlucas\RR_-^*\stoplucas$.

Note that despite the fact that this dynamical system is well-defined for all abstract Cauchy problem in the form \eqref{eq_AbsCauchyProb}, its convergence to the real state is not guaranteed \textit{a priori} and may be linked to observability properties.

\stoplucas

\section{Application to the crystallization process}\label{secappli}

\subsection{Modeling the batch crystallization process}

In this subsection, a dynamical model representing the batch crystallization process and the concentration measurement is given.

\subsubsection{Population balance equation}\label{secbilan}

In a first step a batch crystallization process is modeled.
\startlucas
We assume that the size of the crystals is described by a scalar parameter $x$ (in $m$).
For example, if the crystals are spherical, then $x$ may represent their diameter.
\stoplucas

Let us denote $\etat(t, \cdot)$
the Number Density Function in terms of the crystal size (NDF) at time $t$ \startlucas
(in $\mathrm{m}^{-1}.\mathrm{m}^{-3}$),
\stoplucas
so that
$\int_{x_1}^{x_2}\etat(t, x)dx$ is the total amount of crystals in the reactor at time $t$ with size between $x_1$ and $x_2$
\startlucas
(in $\mathrm{m}^{-3}$).
\stoplucas
Let  $[t_0, t_1]$ be the time window in which the crystallization process occurs.

We assume that the crystals never reach a specific maximal allowable size
denoted  $\xmax$ during the experiment.
Typically this size can simply be related to the size of the reactor in which the crystallization occurs:
\begin{equation}\label{eq_EtatMax}
\forall t \in [t_0, t_1],\quad \etat(t, \xmax) = 0.
\end{equation}
We assume that all the crystals appear at the same positive size $\xmin$,  and we denote $\cont(t)$ the appearance rate of new crystals at size $\xmin$ at time $t$.
\startlucas
The function $\cont$ quantifies the nucleation rate in the reactor.
In \cite{nagy2013recent}, for instance, an expression for the function $\cont$ is given as
\(
\cont(t) = {R_n(t)}/{\vit(t)} ,
\)
where $R_n$ is the overall rate of nucleation expressed in $\mathrm{s^{-1}.m^{-3}}$ and $G$ is the growth rate of the crystals in $\mathrm{m.s^{-1}}$.
Note however that in our approach for NDF observation,
we don't need to know precisely this expression. We do not use any model of $u$, and assume this quantity to be unknown.
\stoplucas
We have
\startlucas
\begin{align}
\forall t \in [t_0, t_1],\quad \etat(t, \xmin) = \cont(t).
\end{align}
\stoplucas

Finally, let  $\vit(t, x) \geq 0$ be the growth rate of the crystals in $\mathrm{m.s^{-1}}$, in other words the rate at which a crystal of size $x$ grows at time $t$.
The population balance
\startlucas
leads to
\stoplucas

\begin{equation}
\forall t \in (t_0, t_1),\quad \forall x\in (\xmin, \xmax),\quad \partial_t \etat(t, x) + \partial_x (\vit \etat)(t, x) = 0.
\label{trans1}
\end{equation}
By considering the McCabe assumption, we assume that $\vit$ does not depend on $x$.
In that case equation \eqref{trans1} becomes
\begin{equation}
\forall t \in (t_0, t_1),\quad \forall x\in (\xmin, \xmax),\quad \partial_t\etat(t, x) + \vit(t)\partial_x\etat(t, x) = 0.
\label{trans}
\end{equation}
Equation (\ref{trans}) is a time-varying one-dimensional transport equation.

We also assume that at the beginning of the experiment, some seed particles are in the reactor.
This yields
\begin{align*}
\forall x \in [\xmin, \xmax], \quad \etat(t_0, x) = \etat_0(x).
\end{align*}

To summarize, the model of the NDF in a batch crystallization process is:
\begin{equation}
\begin{aligned}
\begin{cases}
\partial_t \etat(t, x) = -\vit(t)\partial_x \etat(t, x) & \forall t \in (t_0, t_1), \ \forall x \in (\xmin, \xmax)\\
\etat(t_0, x) = \etat_0(x) & \forall x \in [\xmin, \xmax]\\
\etat(t, \xmin) = \cont(t) &\forall t \in [t_0, t_1]
\end{cases}
\end{aligned}
\label{systbilan}
\end{equation}
to which is added the other boundary condition in (\ref{eq_EtatMax}) which in our case is seen as a knowledge on the particular solution we wish to estimate.

The following theorem states that system \eqref{systbilan} of this model admits weak solution in $L^2$ space and strong solution in $H^1$.
\begin{theorem}\label{th:exist}
Let $t_1 > t_0 \geq 0,\ \xmax > \xmin\geq 0,\ \etat_0 \in H^1(\xmin, \xmax),\ \cont \in H^1(t_0, t_1),\ \vit \in C^0([t_0, t_1]; \RR^*_+)$. Assume that $\cont(t_0) = \etat_0(\xmin)$. Then system \eqref{systbilan} admits a unique solution
\begin{align*}
\etat \in C^0([t_0, t_1]; H^1(\xmin, \xmax)) \cap C^1([t_0, t_1]; L^2(\xmin, \xmax)).
\end{align*}
Moreover, for all $(t, x) \in [t_0, t_1] \times [\xmin, \xmax]$,
\begin{align}

\etat(t, x)=
\begin{cases}
\etat_0 (x - \intg(t)) &\text{ if } x - \xmin \geq \intg(t) \\
\cont \circ \intg^{-1} (\intg(t) - x + \xmin) &\text{ else.}
\end{cases}
\label{eqf}
\end{align}
where $\intg: [t_0, t_1] \ni t \mapsto \int_{t_0}^t\vit(\tau)\dd \tau$.

\end{theorem}
The proof of this theorem can be found in \cite[Theorem 2.4]{Coron} in the case $\vit = 1$, and can be easily adapted by means of a time reparametrization.

It is worth noticing that this theorem does not take into account hypothesis \eqref{eq_EtatMax}.
However, the following proposition holds.
\begin{proposition}
Assume that the hypothesis of Theorem \ref{th:exist} are satisfied.
Assume that there exists $\bar{x}\in[\xmin, \xmax)$ such that $\etat_0(x) = 0$ for all $x\in[\bar{x}, \xmax]$.
If
\begin{align}
    \bar{x} + \intg(t_1) < \xmax,
\end{align}
then $\psi(t, x) = 0$ for all $t\in[t_0, t_1]$ and all $x\in\left[\bar{x} + \intg(t_1), \xmax \right]$.
\end{proposition}
\begin{proof}
Let $t\in[t_0, t_1]$ and $x\in\left[\bar{x} + \intg(t_1), \xmax \right]$. Then
\begin{align*}
    x - \intg(t_1) \geq \bar{x} \geq \xmin.
\end{align*}
Consequently, according to \eqref{eqf},
\(
    \etat(t, x) = \etat_0(x - \intg(t)) = 0.
\)
\end{proof}

Hence, one must choose $t_1$ small enough so that the particles did not reach the size $\xmax$. Roughly speaking, this means that the observer that we are going to design must estimate the state in small time, \emph{i.e.} before the particles reach the size $\xmax$.

In the following, we always assume that \eqref{eq_EtatMax} is satisfied.

\subsection{Concentration sensor modeling}

In the considered batch crystallization process, the measured outputs are the temperature and the solute concentration denoted $\cste(t)$.
These two measurements allow to obtain online estimation of the growth rate (\emph{i.e.} $\vit$) and the third moment of the NDF.

\subsubsection{\texorpdfstring{Estimation of $\vit$}{Estimation of G}}
\label{estG}

The knowledge of the temperature and the solute concentration allows to obtain some approximation of the growth rate $\vit$ .
Indeed, following \cite{Uccheddu}, a model of $\vit$ can be given for all time $t\in [t_0, t_1]$ by
\begin{equation}
\vit(t) = k_g \frac{\cste(t) - C^*(t)}{C^*(t)}
\label{eqG}
\end{equation}
where
\begin{itemize}

\item $k_g$ is a known growth rate parameter (in $\mathrm{m.s^{-1}}$),
\item $C^*(t)$ is the solubility at time $t$ (in $\mathrm{kg}$ of solute per $\mathrm{kg}$ of solvent),
\item $\cste(t)$ is the solute concentration at time $t$ (in $\mathrm{kg}$ of solute per $\mathrm{kg}$ of solvent).
\end{itemize}
Since $C^*(t)$ depends on the temperature at time $t$ ,

the growth rate $\vit$ of the crystals can be estimated online with the available sensors.\\
\startlucas
Other model expressions of $G$ are available in the literature, for more details one may refer to \cite{Mersmann,Mullin}.
\stoplucas

\subsubsection{Estimation of the third moment of the NDF}
It is possible to link the solute concentration with the NDF.
Indeed, for each $t\in \RR_+$, let $\cs(t)$ (in $\mathrm{kg}$ of solid per $\mathrm{kg}$ of solvent) be the solid concentration in the reactor at time $t$,
in other words, the ratio between the total crystals mass in the reactor at time $t$ and the solvent mass.
Let $\rho_s$ (in $\mathrm{kg}$.$\mathrm{m^{-3}}$) be the density of the solute in solid phase and $M_e$ the solvent mass (in $\mathrm{kg}$).
It yields:
\begin{align*}
\cs(t) = \frac{\rho_s}{M_e}V_s(t)
\end{align*}
where $V_s(t)$ is the volume (in $\mathrm{m^3}$) occupied by the crystals at time $t$.
\startlucas
Then the volume of a crystal with size $x$ (in $\mathrm{m}$) is simply $V = k_v x^3$
where $k_v$ is a volumetric shape factor (see \emph{e.g.} \cite{Hulburtkatz,Randolph}). For example, $k_v = \pi/6$ for spherical crystals.
The total volume of the crystals is then
\stoplucas

\begin{align*}
V_s(t) = k_v\int_{\xmin}^{\xmax} \etat(t, x) x^3 \dd x
\end{align*}
Hence, the solid concentration in the reactor can be expressed as follows.
\begin{align}
\forall t \in [t_0, t_1],\quad \cs(t) = \frac{\rho_s k_v}{M_e}\int_{\xmin}^{\xmax} \etat(t, x) x^3 \dd x.
\end{align}

Assume moreover that $\rho_s$ is a known parameter.
This implies that we can associate  to system \eqref{systbilan}  the measurement $\mes$ defined as
\begin{equation}
\forall t \in [t_0, t_1],\quad \mes(t) = \int_{\xmin}^{\xmax} \etat(t, x) x^3 \dd x.
\label{eqy}
\end{equation}

From there, the observation problem we intend to solve is the following.
\textit{From the knowledge of the output function $y(t)$ and the growth rate, give an online estimation of the NDF.}

The purpose of the  Section \ref{Sec_Obs_Cryst} is to propose a novel algorithm to solve this problem based on the
\startlucas KKL
\stoplucas
 methodology which was described in Section \ref{secluen}.
Note however that before implementing the observer an observability analysis can be carried out.

\subsection{About the observability of the crystallization model}

In this section, we study how the third moment may help us to estimate the NDF.
First, we have the following result.
\begin{proposition}\label{prop_pos}
Let $\tau\in(t_0, t_1]$.
\startlucas
Assume that there exists $\mu>0$ such that $\vit(t)\geq\mu$ for all $t\in[t_0, t_1]$.
\stoplucas

Then for all $\mes\in C^0(t_0, \tau)$, there exists at most one function $\cont\in H^4(t_0, \tau)$ such that the solution $\etat$ of \eqref{systbilan} given by $\cont$ and $\etat_0 = 0$ satisfies
\begin{align*}
    \mes(t) = \int_{\xmin}^{\xmax} \etat(t, x)x^3\dd x,\qquad \forall t\in[t_0, \tau].
\end{align*}
\end{proposition}
In other words,
Proposition \ref{prop_pos} states that the map $\cont\mapsto \mes$ is injective, where $\mes$ denotes the third moment of the solution of \eqref{systbilan} with null initial condition.
Its proof can be found in appendix.

Hence, one can hope that our method may reconstruct $\etat$ from $\mes$, at least when the initial condition is zero (\emph{i.e.} there is no crystals at the beginning of the process).

However,
one can wonder what happens if the initial condition is not zero. Can we still reconstruct the state from the measurement of its third moment and the knowledge of its dynamics? In other words, is the map $\etat\mapsto \mes$ injective?
If yes,
then one can hope that our algorithm is robust, so that the estimation of the state converges to the actual NDF.

\startlucas
Unfortunately,
the answer is no. Indeed, we have the following proposition which is a slight modification of \cite[Theorem 3.2.3]{vissers2012model} that we state in our own context only. For the convenience of the reader, its proof is given in appendix.
\stoplucas

\begin{proposition}
\label{prop_neg}
Let $\tau\in(t_0, t_1]$.
\startlucas
Assume that there exists $\mu>0$ such that $\vit(t)\geq\mu$ for all $t\in[t_0, t_1]$.
\stoplucas

There exist infinitely many solutions of \eqref{systbilan} with different initial conditions and boundary conditions that have the same third moment $\mes\in C^0(t_0, \tau)$.
\end{proposition}
We shall say that system \eqref{systbilan} with measurement \eqref{eqy} is not observable.
\startlucas
Thus, we cannot guarantee that our estimation of the NDF converges to the actual NDF.
Despite this fact, our methodology should be able to reconstruct partially the actual NDF.
Indeed, the linear function that maps the NDF to its third moment has rank 4.
The image of this state-output mapping is sometimes called the \emph{observable part} of the system (see \emph{e.g.} \cite{haine2014recovering}), due to the fact that an observer shall estimate at least the projection of the actual state on this subspace.{}
See the proof in appendix for more details.
\stoplucas

\subsection{A dynamical observer from the concentration measurement}
\label{Sec_Obs_Cryst}

Eventhough it was shown in the former section that the system is not observable, in this subsection, we show how the
\startlucas KKL
\stoplucas
observer approach can be employed on the considered model.
Following the procedure given in Section \ref{secluen}, we
consider $\lambda$ a negative real number and the dynamical system
\begin{equation}\label{eq_ObsCrist}
\dot z = \lambda z + y\ .
\end{equation}
We must find a mapping $\TR_\lambda$ which is estimated by this dynamical equation.
We have the following proposition.
\begin{proposition}
Let $\TR_\lambda : C^1([t_0,t_1]; L^2(\xmin, \xmax))\mapsto C^1([t_0,t_1])$ be the functional defined as
\begin{equation}\label{defT}
    \TR_\lambda(\etat) : t \mapsto \int_{\xmin}^{\xmax}a(t,x)\etat(t,x)\dd x
\end{equation}
where
\begin{equation}

\begin{cases}
\partial_t a(t, x) + \vit(t) \partial_x a(t, x) = \lambda a(t, x) + x^3 & \forall t \in (t_0, t_1), \ \forall x \in (\xmin, \xmax)\\
a(t_0, x) = 0 &\forall x \in [\xmin, \xmax]\\
a(t, \xmin) = 0 &\forall t \in [t_0, t_1],
\end{cases}
\label{edpa}
\end{equation}
then, along the solution of \eqref{eq_ObsCrist} which satisfies \eqref{eq_EtatMax} where $y$ is given in \eqref{eqy}, it yields for all $t$ in $[t_0,t_1]$
\begin{equation}\label{eq_ObTCrist}
\startlucas
\TR_\lambda(\etat)(t) - z(t) = \exp(\lambda t) (\TR_\lambda(\etat)(0) - z_0).
\stoplucas

\end{equation}
\end{proposition}
\begin{proof}
Using \eqref{systbilan} and an integration by parts yields
\begin{align*}
\frac{d}{dt}(\TR_\lambda(\etat)(t) - z(t))
&= \int_{\xmin}^{\xmax}\partial_t a(t, x) \etat(t,x) dx
    - \int_{\xmin}^{\xmax}\vit(t)  a(t,x) \partial_x \etat(t,x)  dx
\\
& \qquad -\lambda z(t) -\int_{\xmin}^{\xmax}x^3 \etat(t,x) dx
\\
&= \int_{\xmin}^{\xmax}\partial_t a(t, x) \etat(t,x) dx
    + \int_{\xmin}^{\xmax}\vit(t)  \partial_x a(t,x)  \etat(t,x)  dx \\
& \qquad - G(t)\big[a(t,x)\etat(t,x) \big]_{\xmin}^{\xmax}
    -\lambda z(t) -  \int_{\xmin}^{\xmax}x^3 \etat(t,x) dx .
\end{align*}
Hence, with \eqref{edpa} and also the boundary condition in  \eqref{systbilan} and \eqref{eq_EtatMax}, this implies
$$
\frac{d}{dt}(\TR_\lambda(\etat)(t) - z(t)) =  \lambda (\TR_\lambda(\etat)(t) - z(t)).
$$
By integrating in time the former equation, we obtain \eqref{eq_ObTCrist}.
\end{proof}

Consequently, for each negative $\lambda$ we exponentially estimate the functional $\TR_\lambda \etat(t)$.
It is interesting to remark that no information on the nucleation rate is needed to obtain this estimation.
The state observer is given as \eqref{eq_Obs_KKL_Tikhonov}, after a choice of the regularization parameter $\delta$.

\section{Numerical simulations}\label{secsimu}
In this subsection

\startlucas
numerical simulations are
\stoplucas
carried out.

\startlucas
Let $(\mathbb{x}_j)_{1\leq j\leq N_x}$ be a uniform discretization of the space interval $(x_{\min}, x_{\max})$ with space step $\Delta x$ and $(\mathbb{t}_k)_{1\leq k\leq N_t}$ be a uniform discretization of the time interval $(t_0, t_1)$ with time step $\Delta t$.
We fix $N_x = N_t = 100$.
Let $(\lambda_j)_{1\leq i\leq p}$ be the considered negative values of $\lambda$.
An approximation of $(\TR_{\lambda_i}\psi)(\mathbb{t}_k)$ is given by $\Delta x \sum_{j=1}^{N_x} a_{i, j, k} \psi_{j, k}$ where $a_{i, j, k}$ is an approximation of $a_{\lambda_i}(\mathbb{t}_k, \mathbb{x}_j)$ (solution of \eqref{edpa}) and $\psi_{j, k}$ an approximation of $\psi(\mathbb{t}_k, \mathbb{x}_j)$.
\stoplucas
The transport equation which describes the crystallization process is simulated via the method of characteristics.

\startlucas
We consider the system \eqref{systbilan} with $\vit$ as in \eqref{eqG} with a null initial condition $\etat_0 = 0$ and a boundary condition similar to a truncated normal distribution reaching its maximum at $t=3$ $\mathrm{s}$ and with a compact support $[0, 6]$
\startlucas
(see Fig.~\ref{fig:u}).
\stoplucas
The unique solution of this system is drawn in Fig.~\ref{fig:sol} (solid line),
\startlucas
and the corresponding growth rate is drawn in Fig.~\ref{fig:G}.
\stoplucas

\begin{figure}
    \begin{subfigure}{.5\linewidth}
    \centering
    \includegraphics[width = \linewidth]{u.pdf}
    \caption{Boundary condition $u$ due to the nucleation}
    \label{fig:u}
    \end{subfigure}
    \begin{subfigure}{.5\linewidth}
    \centering
    \includegraphics[width = \linewidth]{G.pdf}
    \caption{Growth rate $\vit$, given by \eqref{eqG}}
    \label{fig:G}
    \end{subfigure}\par\medskip
    \centering
    \begin{subfigure}{.5\linewidth}
    \centering
    \includegraphics[width = \linewidth]{sol.pdf}
    \caption{Simulation of the NDF $\psi$}
    \label{fig:sol}
    \end{subfigure}
    \caption{Numerical simulation of the process with $(\xmin, \xmax) = (0, 10)$, $(t_0, t_1) = (0, 10)$ and $N_x = N_t = 100$.}
    \label{fig:all_sim}
\end{figure}

\subsection{Step 1: reconstruction of a function of the state}

Following the methodology developed in Section \ref{secluen}, we first try to estimate the function $\TR_\lambda(\etat)$ of the state via the dynamical system \eqref{eq_ObsCrist} for some fixed negative values of $\lambda$.
\startlucas
All along the simulation of \eqref{systbilan}, we compute $\mes$ and estimate the solution of \eqref{edpa} via the method of characteristics.

\stoplucas

We integrate the solution of \eqref{eq_ObsCrist} with the first order Euler's method.
\stoplucas
Then we plot the evolution of the relative error between $z$ and $\TR_\lambda(\etat)$ in Fig.~\ref{fig:func_state} for some values of $\lambda$. One can check that the error goes to zero as $t\to+\infty$. Moreover, the bigger is $|\lambda|$, the faster is the convergence. This is due to the exponential convergence of $z-\TR_\lambda(\etat)$ to zero given by \eqref{eq_vitconv}.
Hence, we are able to approximate any function $\TR_\lambda(\etat)$ of the state.

Now, we can move to the second part of the methodology of Section \ref{secluen}.

\begin{figure}[htbp]
\centering

\begin{tabular}{r@{}l}
\includegraphics[width = 0.5\textwidth]{lambda01v2.pdf}&
\includegraphics[width = 0.5\textwidth]{lambda1v2.pdf}\\
\includegraphics[width = 0.5\textwidth]{lambda10v2.pdf}&
\includegraphics[width = 0.5\textwidth]{lambda100v2.pdf}
\end{tabular}

\caption{\small Convergence of $\TR_\lambda(\etat) - z$ to zero for different values of $\lambda$.
We choose $z(0) \neq 0$ arbitrarily. The bigger is $|\lambda|$, the faster is the convergence. By means of a linear regression, one can estimate the convergence rate of the relative error to zero:
$
\mathcal{O}(e^{-7.4t})\text{ if } \lambda = -0,1,\quad
\mathcal{O}(e^{-8.2t})\text{ if } \lambda = -1,\quad
\mathcal{O}(e^{-14.2t})\text{ if } \lambda = -10,\quad
\mathcal{O}(e^{-32.9t})\text{ if } \lambda = -100.
$

}

\label{fig:func_state}
\end{figure}

\subsection{Step 2: reconstruction of the entire state of the system}

Following Step 1, we estimate simultaneously numerous functions $\TR_{\lambda_i}(\etat)$ which correspond to different values $\lambda_i <0$.
\startlucas
These estimations are denoted $z_{\lambda_i}$.
The aim of this section is to estimate the state $\psi$ from the knowledge of $(z_{\lambda_i})_{1\leq i\leq p}$.

Then, we choose a regularization parameter $\delta>0$ and solve the discrete version of the quadratic minimization problem \eqref{pb:min_well} at each time step, that is for each time $\mathbb{t_k}$, find $\big(\psi_{j, k}\big)_{1\leq j\leq N_x}$ that
\begin{equation}
    \textrm{minimize~} \left\|\Delta x \big(a_{i, j, k}\big)_{1\leq i\leq p, 1\leq j\leq N_x}\ \cdot\ \big(\psi_{j, k}\big)_{1\leq j\leq N_x} - \big(z_{\lambda_i}(\mathbb{t}_k)\big)_{1\leq i\leq p}\right\|^2+ \delta \left\| \big(\psi_{j, k}\big)_{1\leq j\leq N_x}\right\|^2.
\end{equation}

This is a quadratic minimization problem, which we solve via an interior-point method (see \emph{e.g.} \cite[Chapter III.11]{BV}.
We need to fix an initial condition to apply this algorithm.
Following a continuation method, we choose as an initial condition at time $\mathbb{t}_k$ the minimum value obtained at time $\mathbb{t}_{k-1}$, transported during a time $\Delta t$ at speed $G(\mathbb{t}_{k-1})$.

The choice of parameters $p$, $\lambda_1,\dots,\lambda_p$ and $\delta$ and their influence are investigated in the paragraphs below.
\begin{itemize}

    \item \textbf{Choice of the $p$ and $(\lambda_i)$.}
Note that the matrix $\big(a_{i, j, k}\big)_{1\leq i\leq p, 1\leq j\leq N_x}$ may be injective only if $p\geq N_x$, that is if the discretization in $\lambda$ is thinner than in $x$. Therefore, we fix $p= 2N_x = 200$.
Moreover, even if the matrix $(a_{i, j, k})_{i, j}$ is injective, a regularization method is needed to left-inverse it. Indeed, for all $t\in(t_0, t_1)$, the operator
$$L^2(\xmin, \xmax)\ni \psi\mapsto \left(\lambda\mapsto \int_{\xmin}^{\xmax} a_{\lambda}(t, x) \psi(x)\dd x \right)\in L^2(\lambda_{\min}, \lambda_{\max})$$
is compact (as an integral operator).
Hence, even if it is injective, its inverse is not continuous.
The matrix $(a_{i, j, k})_{i, j}$ is a discretization of this operator. Then, the more the discretization is thinner, the more it is ill-conditioned. This emphasizes the necessity of using a regularization method.

In Fig.~\ref{fig:lambda}, we plot the estimation of the NDF for different values of $(\lambda_i)$.
For large values of $|\lambda|$, $z$ converges quickly to $\TR\psi$. However, it appears that functions $a_\lambda$ carry less information for large values of $|\lambda|$, so that the map $\TR$ is more difficult to inverse. This explains Fig.~\ref{fig:llarge}, on which the estimation $\hat\psi$ is worst than on Fig.~\ref{fig:ltot}.

On the contrary, for small values of $|\lambda|$, it seems that functions $a_\lambda$ carry more information, since the estimation $\hat\psi$ is similar on Fig.~\ref{fig:lsmall} and Fig.~\ref{fig:ltot} at $t=10$ s. However, we also see a peaking phenomenon (for $t\leq5$ s on Fig.~\ref{fig:lsmall}), due to the fact that $z$ is slower to converge to $\TR\psi$ than for large values of $|\lambda|$.

Thus, one must find a compromise for the choice of $(\lambda_i)$: take large values for fast convergence and avoiding peaking, and small values for efficient estimation.

\begin{figure}[htbp]

\begin{subfigure}{.5\linewidth}
\includegraphics[width = \linewidth]{delta01lambda110bruit0.pdf}
\caption{$(\lambda_i)\subset[-10, -1]$}
\label{fig:lsmall}
\end{subfigure}
\begin{subfigure}{.5\linewidth}
\includegraphics[width = \linewidth]{delta01lambda10100bruit0.pdf}
\caption{$(\lambda_i)\subset[-100, -10]$}
\label{fig:llarge}
\end{subfigure}\par\medskip
\centering
\begin{subfigure}{.5\linewidth}
\includegraphics[width = \linewidth]{delta01lambda1100bruit0.pdf}
\caption{$(\lambda_i)\subset[-100, -1]$}
\label{fig:ltot}
\end{subfigure}
\caption{Influence of $(\lambda_i)_{1\leq i\leq p}$ on the reconstruction of the CSD.}
\label{fig:lambda}
\end{figure}

    \item \textbf{Choice of the regularization parameter $\delta$.}

\end{itemize}

The regularization parameter $\delta$ must be chosen numerically, in order to find a compromise between the minimization of the norm of the state, and the minimization of the gap $\TR\psi - z$.
This compromise can be interpreted as a measurement reliability. Indeed, if the measurement has a small uncertainty, then we choose a small $\delta$. On the contrary, if the measurement is highly uncertain, then we fix a large value of $\delta$ in order to regularize the solution.

In Fig.~\ref{fig:obs}, we plot the actual NDF $\etat$ and its estimation $\hat\etat$ at different times, for different values of $\delta$, and with or without measurement noise.
Measurement noise is fixed at $2\%$ of the maximal value of the output on the time interval.
For small values of $\delta$ and/or with measurement noise, we see that a peaking phenomenon appear: this is due to a lack of regularization of the solution.
On the contrary, if $\delta$ is too large, then the minimization of the norm of the state takes too much importance in the minimization problem, and $\hat\psi$ is too attenuated.

\stoplucas

\begin{figure}[htbp]
\centering

\begin{tabular}{r@{}l}
\subcaptionbox{$\delta=0.05,\ \alpha = 0$}{
\includegraphics[width = 0.5\textwidth]{delta005lambda1100bruit0.pdf}}&
\subcaptionbox{$\delta=0.05,\ \alpha = 2\%$}{
\includegraphics[width = 0.5\textwidth]{delta005lambda1100bruit002.pdf}}\\
\subcaptionbox{$\delta=0.1,\ \alpha = 0$}{
\includegraphics[width = 0.5\textwidth]{delta01lambda1100bruit0.pdf}}&
\subcaptionbox{$\delta=0.1,\ \alpha = 2\%$}{
\includegraphics[width = 0.5\textwidth]{delta01lambda1100bruit002.pdf}}\\
\subcaptionbox{$\delta=0.2,\ \alpha = 0$}{
\includegraphics[width = 0.5\textwidth]{delta02lambda1100bruit0.pdf}}&
\subcaptionbox{$\delta=0.2,\ \alpha = 2\%$}{
\includegraphics[width = 0.5\textwidth]{delta02lambda1100bruit002.pdf}}
\end{tabular}
\caption{Influence of the regularization parameter $\delta$ and measurement noise $\alpha$ on the reconstruction of the CSD}
\label{fig:obs}
\end{figure}

\section{Conclusion}

In this paper a new observer has been introduced to estimate the Crystal Size Distribution
\startlucas
from the measurements of the solute concentration, temperature and a model of the growth rate.
No model of the nucleation rate is needed.
\stoplucas

This approach is based on the use of Kazantzis-Kravaris/Luenberger observer and a Tikhonov regularization procedure.
The numerical results obtained are promising.
Even though, the knowledge of the solute concentration alone does not allow to accurately reconstruct the full CSD as shown by our observability analysis, we believe that this tool could be used in addition with supplementary information given by other sensors (for instance \fbrmp).
This approach will be evaluated on experiments in a future research project.

\section*{Appendix: proofs of Propositions \ref{prop_pos} and \ref{prop_neg}.}

\startlucas
Assume that $\vit\geq\mu>0$.
Let $\tau\in(t_0, t_1]$ and $u\in H^4(t_0, \tau)$.
Let $\etat$ be the solution of \eqref{systbilan} with initial condition $\etat_0$ and boundary condition $\cont$.
We introduce a time reparametrization $\Tilde{t} = \int_{t_0}^tG(s)\dd s$, which is well defined since $\vit\geq\mu$.
Let $\Tilde{\psi}$, $\Tilde{\cont}$ and $\Tilde{\mes}$ be such that
$\Tilde{\psi}(\Tilde{t}) = \psi(t)$, $\Tilde{\cont}(\Tilde{t}) = \cont(t)$ and $\Tilde{\mes}(\Tilde{t}) = \mes(t)$ for all $t\in[t_0, t_1]$. Then
\begin{equation}
\begin{aligned}
\begin{cases}
\partial_{\Tilde{t}} \Tilde{\etat}(\Tilde{t}, x) = -\partial_x \Tilde{\etat}(\Tilde{t}, x)\\
\Tilde{\etat}(t_0, x) = \etat_0(x)\\
\Tilde{\etat}(\Tilde{t}, \xmin) = \Tilde{\cont}(\Tilde{t})
\end{cases}
\end{aligned}
\label{syst_vit1}
\end{equation}
and $\Tilde{\mes}(\Tilde{t}) = \int_{\xmin}^{\xmax}\Tilde{\psi}(\Tilde{t}, x)x^3\dd x$.
Since the observability properties are not affected by the time reparametrization, one can investigate observability properties of the system \eqref{syst_vit1} instead of \eqref{systbilan}.
Therefore, one can assume without loss of generality that $\vit=1$ in the rest of the proof.
Since $\cont\in H^4(t_0, \tau)$, we have $\mes\in C^4(t_0, \tau)$.
Equation \eqref{eq_EtatMax} and system \eqref{systbilan} yield
\begin{align}
\mes'
&= 3 \int_{x_{\min}}^{x_{\max}}x^2\etat(\cdot, x)\dd x -  \left[x^3\etat(\cdot, x)\right]_{x_{\min}}^{x_{\max}}
\nonumber\\
&= 3 \int_{x_{\min}}^{x_{\max}}x^2\etat(\cdot, x)\dd x +  x_{\min}^3\cont,
\label{E:y(1)}
\\
\mes^{(2)}
&= 6 \int_{x_{\min}}^{x_{\max}}x\etat(\cdot, x)\dd x + 3 \left[x^2\etat(\cdot, x)\right]_{x_{\min}}^{x_{\max}}
\startlucas
+
\stoplucas

 x_{\min}^3\cont' \nonumber\\
&= 6 \int_{x_{\min}}^{x_{\max}}x\etat(\cdot, x)\dd x + 3 x_{\min}^2\cont +  x_{\min}^3\cont',
\label{E:y(2)}
\\
\mes^{(3)}
&= 6\int_{x_{\min}}^{x_{\max}}\etat(\cdot, x)\dd x - 6 \left[x\etat(\cdot, x)\right]_{x_{\min}}^{x_{\max}} + 3\ x_{\min}^2\cont' + x_{\min}^3\cont^{(2)} \nonumber\\
&= 6\int_{x_{\min}}^{x_{\max}}\etat(\cdot, x)\dd x + 6x_{\min}\cont + 3 x_{\min}^2\cont' +  x_{\min}^3\cont^{(2)},
\label{E:y(3)}
\\
\mes^{(4)}
&= - 6 \left[\etat(\cdot, x)\right]_{x_{\min}}^{x_{\max}} + 6x_{\min}\cont' + 3 x_{\min}^2\cont^{(2)}(t) + x_{\min}^3\cont^{(3)} \nonumber\\
&= 6\cont + 6x_{\min}\cont' + 3 x_{\min}^2\cont^{(2)} +  x_{\min}^3\cont^{(3)}.
\label{E:y(4)}
\end{align}

\medskip
\begin{proof}[End of the proof of Proposition \ref{prop_pos}]
By hypothesis, $\etat_0=0$.
Consequently, Equations \eqref{E:y(1)}--\eqref{E:y(3)} yield
\begin{align*}
\begin{cases}
\mes'(t_0) = x_{\min}^3\cont(t_0)\\
\mes^{(2)}(t_0) = 3 x_{\min}^2\cont(t_0) + x_{\min}^3\cont'(t_0)\\
\mes^{(3)}(t_0) = 6 x_{\min}\cont(t_0) + 3 x_{\min}^2\cont'(t_0) +x_{\min}^3\cont^{(2)}(t_0),
\end{cases}
\end{align*}
which is a triangular system with non vanishing diagonal since $\xmin>0$. Hence $\cont(t_0)$, $\cont'(t_0)$ and $\cont^{(2)}(t_0)$ are determined by $\mes$.
Moreover, on $[t_0, \tau]$, $\cont$ satisfies Equation \eqref{E:y(4)} which is a 3rd order ordinary differential equation.
Hence, according to the Cauchy-Lipschitz theorem, there exits a unique solution $\cont$ to this problem.
Thus $\mes$ determines $\cont$ uniquely, that $\cont\mapsto\mes$ is injective.
\end{proof}

\begin{proof}[End of the proof of Proposition \ref{prop_neg}]
Substituting the boundary condition in equation \eqref{E:y(4)} with  $u=0$ yields $\mes^{(4)} = 0$ identically on $[t_0, \tau]$.
Hence $\mes$ is a polynomial function of degree less or equal than 3.
Thus the linear function that maps any solution of \eqref{systbilan} with null boundary condition to its third moment has rank 4. Since $\etat$ lies in an infinite dimensional vector space, we get by the rank-nullity theorem that its kernel is non-trivial, \emph{i.e.} the state-output map $\etat\mapsto\mes$ is not injective, and the system has a 4-dimensional observable part.
\end{proof}

Note that Proposition \ref{prop_neg} relies deeply on Hypothesis \eqref{eq_EtatMax}.
Hence the non-injectivity of the measurement is due to the fact that the system is observed on a too small time interval.
If the system was observed on $[t_0, +\infty)$, then one could show with similar argues an injectivity result.
\stoplucas

\bibliographystyle{abbrv}

\bibliography{Biblio}

\end{document}

config: arxiv


TITLE: What is your idea on "3-d trig ratios" for surface area?
QUESTION [1 upvotes]: Let's say we take a triple integral and use it to find "3d angles", which use the idea of of taking a portion of a 3-d sphere, compared to the entire sphere.
If 2-d trig ratios can be used to find side-lengths of 2-d triangles, is it possible we can use "3-d trig ratios" to find "surface area sections" of 3-d right pyramids. A "right" 3-d angle is one-fourth of a sphere with a volume of 1 cubic centimeter.
So for finding surface areas, 3-d trigonometry could be used to solve the surface areas of right pyramids by only knowing one side of a surface area.
3-d Pythagorean theorem can be used to be applied to surface areas, or at least that's what I think?
So are there any questions or comments? Is this "possible", could this require help from several people?

REPLY [2 votes]: So indeed there do exist concepts called Solid angles and there is a version of the pythagorean theorem in 3d space (not only for line segment lengths BUT ALSO for entire areas of SIDES) which is extremely elegant.
That being the issue with attempting to create a trigonmetry for angles in 3d space is that that there are infinitely many such angles! Since a solid angle can have an arbitrary number of line segments connecting to it and therefore to gain the some control over 3d geometry as 2d geometry we need to essentially define
Sin for the 4 sided pyramid
Sin for the 5 sided pyramid
Sin for the 6 sided pyramid
Etc...
As far as I know no one had done this but it is definitely something interesting to look at.
And I would imagine you would traditionally consider various ratios of side lengths as well as face areas to generate your functions (why not try to prove some funky identities?)
Now if you generalized the concept of an 'Angle' as a collection of various pieces of information (how many lines are connecting to it, what are the angles between these lines, etc...) then you can define your new ocean of trig functions over which angles are valid arguments and which aren't and see if you can find some identities.
Good luck!


TITLE: Random rotations in SO(3) and free group
QUESTION [10 upvotes]: Is it true that two random (w.r.t. Haar measure) rotations in $SO(3)$ generate a free group?

REPLY [7 votes]: See also:  Epstein - Almost all subgroups of a Lie group are free (MR).


TITLE: Simple Algebra: How has this equation been rearranged?
QUESTION [1 upvotes]: Solving for ΔY
So I'm doing economics and it's been a good 10 years since I've done any algebra. I'm having difficulty rearranging equations. In the picture it shows the equation moving from:
ΔY = MPC x (ΔY-ΔT)
to
(1-MPC) x ΔY = -MPC x ΔT
to the "Final result" in the picture. I'm unsure of the steps involved to get there. I've watched many videos and googled a lot but I still can't wrap my head around it. Could someone please take me through the steps on how we arrive at the final result?
NOTE: MPC, ΔY etc are single values, not separate characters ie M, P and C

REPLY [2 votes]: Your original equation is $$\Delta Y=MPC(\Delta Y-\Delta T)$$ Using the distributive property, you can expand the right-hand side, leading to the following result:$$\Delta Y=MPC\times\Delta Y-MPC\times\Delta T$$ Grouping all the terms that contain $\Delta Y$ gives $$\Delta Y-MPC\times\Delta Y=-MPC\times\Delta T$$ Since $\Delta Y$ is a common factor on the left-hand side, you can factorize it out as follows:
$$\Delta Y(1-MPC)=-MPC\times\Delta T$$ Since you want to make $\Delta Y$ the subject of the formula, you can divide through by $(1-MPC)$ to obtain $$\Delta Y=\frac{-MPC}{1-MPC}\times\Delta T$$ and that is your answer. I hope it helps.


TITLE: How to think about the presentation of $\Bbb Z \times\Bbb Z$
QUESTION [0 upvotes]: I keep trying to understand this basic fact that $\mathbb{Z} \times \mathbb{Z}$ has a presentation $G = (x, y \vert xy = yx)$. I am stuck trying to understand the map between tuples $(x, y)$ in $\mathbb{Z} \times \mathbb{Z}$ and the scalars that result when I generate words from $G$. Clearly the tuples $(x, y)$ and $(y, x)$ are different in $\mathbb{Z} \times \mathbb{Z}$, so where does the commutative relation come from in the presentation?
I looked at this question but the explanation "since the exponents are just adding" really did not clarify it for me.

REPLY [3 votes]: The components of this presentation notation are
$$\langle  \underbrace{ {x, y}}_{\text{names for generators of the group}} \mid
\underbrace{xy = yx}_\text{relations between the generators} \rangle$$
It has nothing to do with tuples.
Let $G$ be the group given by this presentation. To obtain an isomorphism  $f : G \to \mathbb{Z} \times \mathbb{Z}$ you might:

Choose two elements $u, v \in \mathbb{Z} \times \mathbb{Z}$ to be the image of $x$ and $y$ respectively; that is, to set $f(x) = u$ and $f(y) = v$.
Verify that all of the relations hold after substituting $u$ for $x$ and $v$ for $y$; that is, check $f(xy) = f(yx)$, or equivalently, that $u + v = v + u$.
Find a homomorphism $\mathbb{Z} \times \mathbb{Z} \to G$ that is inverse to $f$

The first bullet point is enough to define a homomorphism from the free group on two elements to $\mathbb{Z} \times \mathbb{Z}$. The second bullet point shows this homomorphism is well-defined on the group $G$. The third bullet point shows that it is an isomorphism.
Ultimately, it's not hard to guess at the choice of $u$ and $v$; the very first thing that springs to mind in response to the question "Give me a pair of elements of $\mathbb{Z} \times \mathbb{Z}$ that generate the group!" are probably going to work.


TITLE: Is the function $f(z)=e^{i\bar z^2}$ analytic?
QUESTION [0 upvotes]: Is the function $f(z)=e^{i\bar {z}^2}$ analytic?
My attempt: Let $z=x+iy$. Then $$f(z)=f(x+iy)=e^{i(x-iy)^2}=e^{2xy}\cdot e^{i(x^2-y^2)}=e^{2xy}[\cos(x^2-y^2)+i\sin(x^2-y^2)].$$ Therefore, $u(x,y)=e^{2xy}\cos(x^2-y^2)$ and $v(x,y)=e^{2xy}\sin(x^2-y^2)$. Now we want to calculate the values of $x$ and $y$ such that $u_x(x,y)=v_y(x,y)$ and $u_y(x,y)=-v_x(x,y)$. I found that those values of $x$ and $y$ have to satisfy these equations:
\begin{gather*}
x\sin(x^2-y^2)=y\cos(x^2-y^2)\\
x\cos(x^2-y^2)=-y\sin(x^2-y^2)
\end{gather*}
Clearly $(x,y)=(0,0)$ satisfies both equations, but I don't know if it is the only point of $\mathbb{C}$ that satisfies both equations.

REPLY [2 votes]: Let's solve these equations. First, note that, if $\cos(x^2 - y^2) = 0$, then $\sin(x^2 - y^2) = \pm 1$. This means the equations simplify to $\pm x = 0$ and $\mp y = 0$, and in any case, $x = y = 0$. This contradicts $\cos(x^2 - y^2) = 0$.
By rearranging the equations, we get
\begin{align*}
x \tan(x^2 - y^2) &= y \\
y \tan(x^2 - y^2) &= -x.
\end{align*}
Multiplying them respectively by a factor of $y$ and $x$, we get
$$y^2 = yx\tan(x^2 - y^2) = xy\tan(x^2 - y^2) = -x^2,$$
and hence $x^2 + y^2 = 0$. Since $x$ and $y$ are real numbers, this implies $x = y = 0$. That is, the one and only solution is $x = y = 0$, proving that $f$ is analytic nowhere.


TITLE: Expectation of a Normal distribution with dependent conditions
QUESTION [1 upvotes]: $X$  and $Z$ are two independent Normal random variables and $k$ is a constant. I am interested in the following conditional expected value
$$E(X\mid X<k,X+Z>k).$$

REPLY [2 votes]: The general approach is:
$$  E[X| X<k, X+Z>k] = \frac{E[X1_{X<k, X+Z>k}]}{P(X<k, X+Z>k)} $$
$$ E[X1_{X<k, X+Z>k}] =  \iint_D x f_{X,Z} (x,z) dx dz$$
$$ P(X<k, X+Z>k) = \iint_D f_{X,Z} (x,z) dx dz $$
where $D =\{(x,z)\in \mathbf{R}^2 \mid   x<k, z >k-x \}  $.


TITLE: Non-degenerate pairing on Néron–Severi Group
QUESTION [3 upvotes]: $\DeclareMathOperator\NS{NS}\DeclareMathOperator\Pic{Pic}$Let $X$ be a smooth projective variety over an algebraically closed field $k$. Let $\NS(X)$ denote the group $[\Pic(X)/\Pic^{\text a}(X)]\otimes_{\mathbb{Z}}\mathbb{Q}$, where $\Pic^{\text a}(X)$ denotes the divisors which are algebraically equivalent to 0. If $X$ is a surface, then the pairing on $\NS(X)$ given by $(D,D')\mapsto D\cdot D'$ is non-degenerate.
For higher dimensions, let $\Theta$ denote an ample divisor on $X$. Let $d$ be the dimension of $X$. Is it true that the pairing on $\NS(X)$ given by $(D,D')\mapsto D\cdot D'\cdot \Theta^{d-2}$ is non-degenerate?
Is there a good reference for the above result?

REPLY [1 votes]: S. Kleiman, Les théorèmes de finitude pour le foncteur de Picard, Exp. XIII of Séminaire de géométrie algébrique du Bois-Marie - 1966--1967 - Théorie des intersections et théorème de Riemann-Roch (SGA6) (P. Berthelot, A. Grothendieck, L. Illusie, eds), Lect. Notes in Math. 225, Springer, 1971, 616--666,
See Cor. 7.4 (i) p. 665.


text:
:: Properties of Trigonometric Function
::  by Takashi Mitsuishi and Yuguang Yang

environ

 vocabularies NUMBERS, XREAL_0, SUBSET_1, XXREAL_0, CARD_1, ARYTM_3, RELAT_1,
      SQUARE_1, ARYTM_1, SIN_COS, XXREAL_1, ORDINAL2, REAL_1, FDIFF_1, FUNCT_1,
      VALUED_0, NAT_1, NEWTON, BINOP_2, PARTFUN1, FUNCOP_1, SEQ_1, VALUED_1,
      FUNCT_2, SEQ_2, TARSKI, RCOMP_1, XBOOLE_0, FCONT_1, SIN_COS2, FUNCT_7,
      ORDINAL1;
 notations TARSKI, XBOOLE_0, SUBSET_1, ORDINAL1, NUMBERS, XCMPLX_0, XREAL_0,
      VALUED_0, FUNCT_1, FUNCT_2, RELSET_1, FCONT_1, FDIFF_1, NAT_1, FUNCOP_1,
      SQUARE_1, NEWTON, PARTFUN1, SEQ_2, RCOMP_1, RFUNCT_1, SIN_COS, VALUED_1,
      SEQ_1, BINOP_2, XXREAL_0, REAL_1;
 constructors PARTFUN1, FUNCOP_1, REAL_1, SQUARE_1, NAT_1, BINOP_2, SEQ_2,
      RCOMP_1, RFUNCT_1, RFUNCT_2, FCONT_1, FDIFF_1, SIN_COS, VALUED_1,
      RELSET_1, NEWTON, COMSEQ_2, NUMBERS;
 registrations XBOOLE_0, RELSET_1, NUMBERS, XXREAL_0, XREAL_0, SQUARE_1,
      MEMBERED, RCOMP_1, FDIFF_1, SIN_COS, VALUED_0, VALUED_1, FUNCT_2,
      FCONT_3, NEWTON, SEQ_2, ORDINAL1;
 requirements REAL, NUMERALS, SUBSET, BOOLE, ARITHM;
 equalities SQUARE_1, SUBSET_1;
 theorems SEQ_1, SEQ_4, SQUARE_1, RCOMP_1, FUNCT_1, FDIFF_1, FDIFF_2, FUNCT_2,
      RFUNCT_1, PARTFUN1, ROLLE, SIN_COS, NEWTON, XREAL_0, XCMPLX_0, XCMPLX_1,
      BINOP_2, FUNCOP_1, XREAL_1, VALUED_1, XXREAL_1, VALUED_0, ORDINAL1;
 schemes NAT_1, FUNCT_2;

begin :: Monotone increasing and monotone decreasing of sine and cosine.

reserve p,q,r,th,th1 for Real;
reserve n for Nat;

theorem Th1:
  p>=0 & r>=0 implies p+r>=2*sqrt(p*r)
proof
  assume that
A1: p>=0 and
A2: r>=0;
A3: (sqrt p - sqrt r)^2>=0 by XREAL_1:63;
  (sqrt p - sqrt r)^2 =(sqrt p)^2 - 2*sqrt p*sqrt r + (sqrt r)^2
    .=p-2*sqrt p*sqrt r + (sqrt r)^2 by A1,SQUARE_1:def 2
    .=p -2*sqrt p*sqrt r +r by A2,SQUARE_1:def 2
    .=p+r-2*sqrt p*sqrt r;
  then 0+2*sqrt p*sqrt r <= p+r by A3,XREAL_1:19;
  then 2*(sqrt p*sqrt r) <= p+r;
  hence thesis by A1,A2,SQUARE_1:29;
end;

theorem
  sin|].0,PI/2.[ is increasing
proof
  for th st th in ].0,PI/2.[ holds 0 < diff(sin,th)
  proof
    let th;
    assume th in ].0,PI/2.[;
    then cos.th > 0 by SIN_COS:80;
    hence thesis by SIN_COS:68;
  end;
  hence thesis by FDIFF_1:26,ROLLE:9,SIN_COS:24,68;
end;

Lm1: for th st th in ].0,PI/2.[ holds 0 < sin.th
proof
  let th;
  assume
A1: th in ].0,PI/2.[;
  then th < PI/2 by XXREAL_1:4;
  then -th>-PI/2 by XREAL_1:24;
  then
A2: -th+PI/2 > -PI/2+ PI/2 by XREAL_1:6;
  0 < th by A1,XXREAL_1:4;
  then 0-th<0;
  then -th+PI/2 < 0+PI/2 by XREAL_1:6;
  then PI/2-th in ].0,PI/2.[ by A2,XXREAL_1:4;
  then cos.(PI/2-th) > 0 by SIN_COS:80;
  hence thesis by SIN_COS:78;
end;

theorem
  sin|].PI/2,PI.[ is decreasing
proof
  for th1 st th1 in ].PI/2,PI.[ holds 0 > diff(sin,th1)
  proof
    let th1;
    assume
A1: th1 in ].PI/2,PI.[;
    then th1 < PI by XXREAL_1:4;
    then
A2: th1 - PI/2 < PI - PI/2 by XREAL_1:9;
    PI/2 < th1 by A1,XXREAL_1:4;
    then PI/2 - PI/2 < th1 - PI/2 by XREAL_1:9;
    then th1-PI/2 in ].0,PI/2.[ by A2,XXREAL_1:4;
    then sin.(th1-PI/2) > 0 by Lm1;
    then
A3: 0-sin.(th1-PI/2) < 0;
    diff(sin,(th1))=cos.(PI/2+(th1-PI/2)) by SIN_COS:68
      .=-sin.(th1-PI/2) by SIN_COS:78;
    hence thesis by A3;
  end;
  hence thesis by FDIFF_1:26,ROLLE:10,SIN_COS:24,68;
end;

theorem
  cos|].0,PI/2.[ is decreasing
proof
  for th st th in ].0,PI/2.[ holds diff(cos,th) < 0
  proof
    let th;
    assume th in ].0,PI/2.[;
    then 0 < sin.(th) by Lm1;
    then diff(cos,th) = -sin.(th) & 0-sin.(th) < 0 by SIN_COS:67;
    hence thesis;
  end;
  hence thesis by FDIFF_1:26,ROLLE:10,SIN_COS:24,67;
end;

theorem
  cos|].PI/2,PI.[ is decreasing
proof
  for th st th in ].PI/2,PI.[ holds diff(cos,th) < 0
  proof
    let th;
    assume
A1: th in ].PI/2,PI.[;
    then th < PI by XXREAL_1:4;
    then
A2: th - PI/2 < PI - PI/2 by XREAL_1:9;
    PI/2 < th by A1,XXREAL_1:4;
    then PI/2 - PI/2 < th - PI/2 by XREAL_1:9;
    then th-PI/2 in ].0,PI/2.[ by A2,XXREAL_1:4;
    then cos.(th-PI/2) > 0 by SIN_COS:80;
    then
A3: 0-cos.(th-PI/2) < 0;
    diff(cos,(th))=-sin.(PI/2+(th-PI/2)) by SIN_COS:67
      .=-cos.(th-PI/2) by SIN_COS:78;
    hence thesis by A3;
  end;
  hence thesis by FDIFF_1:26,ROLLE:10,SIN_COS:24,67;
end;

theorem
  sin|].PI,3/2*PI.[ is decreasing
proof
  for th st th in ].PI,3/2*PI.[ holds diff(sin,(th))<0
  proof
    let th such that
A1: th in ].PI,3/2*PI.[;
    th < 3/2*PI by A1,XXREAL_1:4;
    then
A2: th-PI < 3/2*PI-PI by XREAL_1:9;
    PI < th by A1,XXREAL_1:4;
    then PI-PI < th-PI by XREAL_1:9;
    then th-PI in ].0,PI/2.[ by A2,XXREAL_1:4;
    then cos.(th-PI) > 0 by SIN_COS:80;
    then
A3: 0-cos.(th-PI) < 0;
    diff(sin,(th))= cos.(PI+(th-PI)) by SIN_COS:68
      .= -cos.(th-PI) by SIN_COS:78;
    hence thesis by A3;
  end;
  hence thesis by FDIFF_1:26,ROLLE:10,SIN_COS:24,68;
end;

theorem
  sin|].3/2*PI,2*PI.[ is increasing
proof
  for th st th in ].3/2*PI,2*PI.[ holds diff(sin,(th))>0
  proof
    let th such that
A1: th in ].3/2*PI,2*PI.[;
    th < 2*PI by A1,XXREAL_1:4;
    then
A2: th-3/2*PI < 2*PI-3/2*PI by XREAL_1:9;
A3: diff(sin,(th)) = cos.(PI + (PI/2+(th-3/2*PI))) by SIN_COS:68
      .= -cos.(PI/2+(th-3/2*PI)) by SIN_COS:78
      .= -(-sin.(th-3/2*PI)) by SIN_COS:78
      .= sin.(th-3/2*PI);
    3/2*PI < th by A1,XXREAL_1:4;
    then 3/2*PI-3/2*PI < th-3/2*PI by XREAL_1:9;
    then th-3/2*PI in ].0,PI/2.[ by A2,XXREAL_1:4;
    hence thesis by A3,Lm1;
  end;
  hence thesis by FDIFF_1:26,ROLLE:9,SIN_COS:24,68;
end;

theorem
  cos|].PI,3/2*PI.[ is increasing
proof
  for th st th in ].PI,3/2*PI.[ holds diff(cos,(th))>0
  proof
    let th such that
A1: th in ].PI,3/2*PI.[;
    th < 3/2*PI by A1,XXREAL_1:4;
    then
A2: th-PI < 3/2*PI-PI by XREAL_1:9;
A3: diff(cos,(th)) = -sin.(PI+(th-PI)) by SIN_COS:67
      .= -(-sin.(th-PI)) by SIN_COS:78
      .= sin.(th-PI);
    PI < th by A1,XXREAL_1:4;
    then PI-PI < th-PI by XREAL_1:9;
    then th-PI in ].0,PI/2.[ by A2,XXREAL_1:4;
    hence thesis by A3,Lm1;
  end;
  hence thesis by FDIFF_1:26,ROLLE:9,SIN_COS:24,67;
end;

theorem
  cos|].3/2*PI,2*PI.[ is increasing
proof
  for th st th in ].3/2*PI,2*PI.[ holds diff(cos,(th))>0
  proof
    let th such that
A1: th in ].3/2*PI,2*PI.[;
    th < 2*PI by A1,XXREAL_1:4;
    then
A2: th-3/2*PI < 2*PI-3/2*PI by XREAL_1:9;
A3: diff(cos,(th)) = -sin.(PI + (PI/2+(th-3/2*PI))) by SIN_COS:67
      .= -(-sin.(PI/2+(th-3/2*PI))) by SIN_COS:78
      .= cos.(th-3/2*PI) by SIN_COS:78;
    3/2*PI < th by A1,XXREAL_1:4;
    then 3/2*PI-3/2*PI < th-3/2*PI by XREAL_1:9;
    then th-3/2*PI in ].0,PI/2.[ by A2,XXREAL_1:4;
    hence thesis by A3,SIN_COS:80;
  end;
  hence thesis by FDIFF_1:26,ROLLE:9,SIN_COS:24,67;
end;

theorem
  for n being Nat holds sin.th = sin.(2*PI*n + th)
proof
  defpred X[Nat] means for th holds sin.th = sin.(2*PI*$1 + th);
  let n be Nat;
A1: for n being Nat st X[n] holds X[n+1]
  proof
    let n be Nat such that
A2: for th holds sin.th = sin.(2*PI*n + th);
    for th holds sin.th = sin.(2*PI*(n+1) + th)
    proof
      let th;
      sin.(2*PI*(n+1) + th) = sin.((2*PI*n+th) + 2*PI)
        .= sin.(2*PI*n + th) by SIN_COS:78
        .= sin.th by A2;
      hence thesis;
    end;
    hence thesis;
  end;
A3: X[0];
  for n being Nat holds X[n] from NAT_1:sch 2(A3,A1);
  hence thesis;
end;

theorem
  for n being Nat holds cos.th = cos.(2*PI*n + th)
proof
  defpred X[Nat] means for th holds cos.th = cos.(2*PI*$1 + th);
  let n be Nat;
A1: for n being Nat st X[n] holds X[n+1]
  proof
    let n be Nat such that
A2: for th holds cos.th = cos.(2*PI*n + th);
    for th holds cos.th = cos.(2*PI*(n+1) + th)
    proof
      let th;
      cos.(2*PI*(n+1) + th) = cos.((2*PI*n+th) + 2*PI)
        .= cos.(2*PI*n + th) by SIN_COS:78
        .= cos.th by A2;
      hence thesis;
    end;
    hence thesis;
  end;
A3: X[0];
  for n being Nat holds X[n] from NAT_1:sch 2(A3,A1);
  hence thesis;
end;

begin ::Hyperbolic sine, hyperbolic cosine and hyperbolic tangent.

definition
  func sinh -> Function of REAL, REAL means
  :Def1:
  for d being Real holds it.d=(exp_R.d - exp_R.(-d))/2;
  existence
  proof
    deffunc U(Real) = In((exp_R.($1) - exp_R.(-$1))/2,REAL);
    consider f being Function of REAL, REAL such that
A1: for d be Element of REAL holds f.d = U(d) from FUNCT_2:sch 4;
    take f;
    let d be Real;
    d in REAL by XREAL_0:def 1;
    then f.d = U(d) by A1;
    hence thesis;
  end;
  uniqueness
  proof
    let f1,f2 be Function of REAL, REAL;
    assume
A2: for d being Real holds f1.d=(exp_R.d - exp_R.(-d))/2;
    assume
A3: for d being Real holds f2.d=(exp_R.d - exp_R.(-d))/2;
    for d being Element of REAL holds f1.d = f2.d
    proof
      let d be Element of REAL;
      thus f1.d=(exp_R.d - exp_R.(-d))/2 by A2
        .=f2.d by A3;
    end;
    hence f1=f2 by FUNCT_2:63;
  end;
end;

definition
  let d be object;
  func sinh(d) -> number equals
  sinh.d;
  coherence;
end;

registration
  let d be object;
  cluster sinh(d) -> real;
  coherence;
end;

definition
  func cosh -> Function of REAL, REAL means
  :Def3:
  for d being Real holds it.d=(exp_R.d + exp_R.(-d))/2;
  existence
  proof
    deffunc U(Real) = In((exp_R.$1 + exp_R.(-$1))/2,REAL);
    consider f being Function of REAL, REAL such that
A1: for d be Element of REAL holds f.d = U(d) from FUNCT_2:sch 4;
    take f;
    let d be Real;
    d in REAL by XREAL_0:def 1;
    then f.d = U(d) by A1;
    hence thesis;
  end;
  uniqueness
  proof
    let f1,f2 be Function of REAL, REAL;
    assume
A2: for d being Real holds f1.d=(exp_R.d + exp_R.(-d))/2;
    assume
A3: for d being Real holds f2.d=(exp_R.d + exp_R.(-d))/2;
    for d being Element of REAL holds f1.d = f2.d
    proof
      let d be Element of REAL;
      thus f1.d=(exp_R.d + exp_R.(-d))/2 by A2
        .=f2.d by A3;
    end;
    hence f1=f2 by FUNCT_2:63;
  end;
end;

definition
  let d be object;
  func cosh(d) -> number equals
  cosh.d;
  coherence;
end;

registration
  let d be object;
  cluster cosh(d) -> real;
  coherence;
end;

definition
  func tanh -> Function of REAL, REAL means
  :Def5:
  for d being Real
  holds it.d=(exp_R.d - exp_R.(-d))/(exp_R.d + exp_R.(-d));
  existence
  proof
    deffunc U(Real) =
      In((exp_R.$1 - exp_R.(-$1))/(exp_R.$1 + exp_R.(-$1)),REAL);
    consider f being Function of REAL, REAL such that
A1: for d be Element of REAL holds f.d = U(d) from FUNCT_2:sch 4;
    take f;
    let d be Real;
    d in REAL by XREAL_0:def 1;
    then f.d = U(d) by A1;
    hence thesis;
  end;
  uniqueness
  proof
    let f1,f2 be Function of REAL, REAL;
    assume
A2: for d being Real holds f1.d=(exp_R.(d) - exp_R.(-d))/(exp_R
    .d + exp_R.(-d));
    assume
A3: for d being Real holds f2.d=(exp_R.(d) - exp_R.(-d))/(exp_R
    .(d) + exp_R.(-d));
    for d being Element of REAL holds f1.d = f2.d
    proof
      let d be Element of REAL;
      thus f1.d=(exp_R.(d) - exp_R.(-d))/(exp_R.(d) + exp_R.(-d)) by A2
        .=f2.d by A3;
    end;
    hence f1=f2 by FUNCT_2:63;
  end;
end;

definition
  let d be object;
  func tanh(d) -> number equals
  tanh.d;
  coherence;
end;

registration
  let d be object;
  cluster tanh(d) -> real;
  coherence;
end;

theorem Th12:
  exp_R.(p+q) = exp_R.(p) * exp_R.(q)
proof
  exp_R.(p+q) = exp_R(p+q) by SIN_COS:def 23
    .= exp_R(p) * exp_R(q) by SIN_COS:50
    .= exp_R.(p) * exp_R(q) by SIN_COS:def 23
    .= exp_R.(p) * exp_R.(q) by SIN_COS:def 23;
  hence thesis;
end;

theorem Th13:
  exp_R.0 = 1 by SIN_COS:51,def 23;

theorem Th14:
  (cosh.p)^2-(sinh.p)^2=1 & (cosh.p)*(cosh.p)-(sinh.p)*(sinh.p)=1
proof
A1: (sinh.p)*(sinh.p) = ((exp_R.(p) - exp_R.(-p))/2)*(sinh.p) by Def1
    .= ((exp_R.(p) - exp_R.(-p))/2)*((exp_R.(p) - exp_R.(-p))/2) by Def1
    .= ((exp_R.(p))*(exp_R.(p))-(exp_R.(p))*(exp_R.(-p)) -(exp_R.(-p))*(
  exp_R.(p))+(exp_R.(-p))*(exp_R.(-p)))/(2*2)
    .= ((exp_R.(p+p))-(exp_R.(p))*(exp_R.(-p)) -(exp_R.(-p))*(exp_R.(p))+(
  exp_R.(-p))*(exp_R.(-p)))/4 by Th12
    .= ((exp_R.(p+p))-(exp_R.(p+-p)) -(exp_R.(-p))*(exp_R.(p))+(exp_R.(-p))*
  (exp_R.(-p)))/4 by Th12
    .= ((exp_R.(p+p))-(exp_R.(p+-p)) -(exp_R.(-p))*(exp_R.(p))+(exp_R.(-p+-p
  )))/4 by Th12
    .= ((exp_R.(p+p)) - 1 - 1 + (exp_R.(-p+-p)))/4 by Th12,Th13;
  (cosh.p)*(cosh.p) = ((exp_R.(p) + exp_R.(-p))/2)*(cosh.p) by Def3
    .= ((exp_R.(p) + exp_R.(-p))/2)*((exp_R.(p) + exp_R.(-p))/2) by Def3
    .= ((exp_R.(p))*(exp_R.(p))+(exp_R.(p))*(exp_R.(-p)) +(exp_R.(-p))*(
  exp_R.(p))+(exp_R.(-p))*(exp_R.(-p)))/(2*2)
    .= ((exp_R.(p+p))+(exp_R.(p))*(exp_R.(-p)) +(exp_R.(-p))*(exp_R.(p))+(
  exp_R.(-p))*(exp_R.(-p)))/4 by Th12
    .= ((exp_R.(p+p))+(exp_R.(p+-p)) +(exp_R.(-p))*(exp_R.(p))+(exp_R.(-p))*
  (exp_R.(-p)))/4 by Th12
    .= ((exp_R.(p+p))+(exp_R.(p+-p)) +(exp_R.(-p))*(exp_R.(p))+(exp_R.(-p+-p
  )))/4 by Th12
    .= ((exp_R.(p+p)) + 1 + 1 + (exp_R.(-p+-p)))/4 by Th12,Th13
    .= ((exp_R.(p+p)) + 2 + (exp_R.(-p+-p)))/4;
  hence thesis by A1;
end;

Lm2: cosh.(p+r)=(cosh.p)*(cosh.r) + (sinh.p)*(sinh.r)
proof
  (cosh.p)*(cosh.r) + (sinh.p)*(sinh.r) =((exp_R.(p) + exp_R.(-p))/2)*(
  cosh.r) +(sinh.p)*(sinh.r) by Def3
    .=((exp_R.(p) + exp_R.(-p))/2)*((exp_R.(r) + exp_R.(-r))/2) +(sinh.p)*(
  sinh.r) by Def3
    .=((exp_R.(p) + exp_R.(-p))/2)*((exp_R.(r) + exp_R.(-r))/2) +((exp_R.(p)
  - exp_R.(-p))/2)*(sinh.r) by Def1
    .=((exp_R.(p) + exp_R.(-p))/2)*((exp_R.(r) + exp_R.(-r))/2) +((exp_R.(p)
  - exp_R.(-p))/2)*((exp_R.(r) - exp_R.(-r))/2) by Def1
    .=( 2*((exp_R.(p))*(exp_R.(r))) + 2*((exp_R.(-p))*(exp_R.(-r))) )/4
    .=( 2*(exp_R.(p+r)) + 2*((exp_R.(-p))*(exp_R.(-r))) )/4 by Th12
    .=( 2*(exp_R.(p+r)) + 2*(exp_R.(-p+-r)) )/4 by Th12
    .=(exp_R.(p+r)+exp_R.(-(p+r)))/2
    .=cosh.(p+r) by Def3;
  hence thesis;
end;

Lm3: sinh.(p+r)=(sinh.p)*(cosh.r) + (cosh.p)*(sinh.r)
proof
  (sinh.p)*(cosh.r) + (cosh.p)*(sinh.r) =((exp_R.(p) - exp_R.(-p))/2)*(
  cosh.r) + (cosh.p)*(sinh.r) by Def1
    .=((exp_R.(p) - exp_R.(-p))/2)*((exp_R.(r) + exp_R.(-r))/2) + (cosh.p)*(
  sinh.r) by Def3
    .=((exp_R.(p) - exp_R.(-p))/2)*((exp_R.(r) + exp_R.(-r))/2) + ((exp_R.(p
  ) + exp_R.(-p))/2)*(sinh.r) by Def3
    .=((exp_R.(p) - exp_R.(-p))/2)*((exp_R.(r) + exp_R.(-r))/2) + ((exp_R.(p
  ) + exp_R.(-p))/2)*((exp_R.(r) - exp_R.(-r))/2) by Def1
    .=( 2*( (exp_R.(p))*(exp_R.(r)) )+2*( -((exp_R.(-p))*(exp_R.(-r))) ) )/4
    .=( 2*(exp_R.(p+r))+2*( -(exp_R.(-p))*(exp_R.(-r)) ) )/4 by Th12
    .=( 2*(exp_R.(p+r))+2*( -(exp_R.(-p+-r)) ) )/4 by Th12
    .=(exp_R.(p+r)-exp_R.(-(p+r)))/2
    .=sinh.(p+r) by Def1;
  hence thesis;
end;

theorem Th15:
  cosh.p <> 0 & cosh.p > 0 & cosh.0 = 1
proof
  exp_R.p > 0 & exp_R.(-p) > 0 by SIN_COS:54;
  then
A1: (exp_R.p+exp_R.(-p))/2 > 0 by XREAL_1:139;
  cosh.0 = (exp_R.0+exp_R.(-0))/2 by Def3
    .= 1 by SIN_COS:51,def 23;
  hence thesis by A1,Def3;
end;

theorem Th16:
  sinh.0 = 0
proof
  sinh.0 = (exp_R.0-exp_R.(-0))/2 by Def1
    .= 0;
  hence thesis;
end;

theorem Th17:
  tanh.p = (sinh.p)/(cosh.p)
proof
  (sinh.p)/(cosh.p) =((exp_R.(p) - exp_R.(-p))/2)/(cosh.p) by Def1
    .=((exp_R.(p) - exp_R.(-p))/2)/((exp_R.(p) + exp_R.(-p))/2) by Def3
    .=(exp_R.(p) - exp_R.(-p))/(exp_R.(p) + exp_R.(-p)) by XCMPLX_1:55
    .=tanh.p by Def5;
  hence thesis;
end;

Lm4: for a1 be Real holds r<>0 & q<>0 implies (p*q + r*a1)/(r*q + p*a1)
= (p/r + a1/q)/(1 + (p/r)*(a1/q))
proof
  let a1 be Real;
  assume that
A1: r<>0 and
A2: q<>0;
  (p*q + r*a1)/(r*q + p*a1) =( (p*q + r*a1)/(r*q) )/( (r*q + p*a1)/(r*q) )
  by A1,A2,XCMPLX_1:6,55
    .=( (p*q)/(r*q) + (r*a1)/(r*q) )/( (r*q + p*a1)/(r*q) ) by XCMPLX_1:62
    .=( (p*q)/(r*q) + (r*a1)/(r*q) )/ ( (r*q)/(r*q) + (p*a1)/(r*q) ) by
XCMPLX_1:62
    .=( p/r + (a1*r)/(q*r) )/( (r*q)/(r*q) + (p*a1)/(r*q) ) by A2,XCMPLX_1:91
    .=( p/r + a1/q )/( (r*q)/(r*q) + (p*a1)/(r*q) ) by A1,XCMPLX_1:91
    .=( p/r + a1/q )/( 1 + (p*a1)/(r*q) ) by A1,A2,XCMPLX_1:6,60
    .=( p/r + a1/q )/( 1 + (p/r)*(a1/q) ) by XCMPLX_1:76;
  hence thesis;
end;

Lm5: tanh.(p+r) = (tanh.p + tanh.r)/(1+ (tanh.p)*(tanh.r))
proof
A1: cosh.r <> 0 & cosh.p <> 0 by Th15;
  tanh.(p+r) =(sinh.(p+r))/(cosh.(p+r)) by Th17
    .=( (sinh.p)*(cosh.r) + (cosh.p)*(sinh.r) )/(cosh.(p+r)) by Lm3
    .=( (sinh.p)*(cosh.r) + (cosh.p)*(sinh.r) )/ ( (cosh.p)*(cosh.r) + (sinh
  .p)*(sinh.r) ) by Lm2
    .=( (sinh.p)/(cosh.p) + (sinh.r)/(cosh.r) )/ ( 1 + ( (sinh.p)/(cosh.p) )
  *( (sinh.r)/(cosh.r) ) ) by A1,Lm4
    .=( tanh.p + (sinh.r)/(cosh.r) )/ ( 1 + ( (sinh.p)/(cosh.p) )*( (sinh.r)
  /(cosh.r) ) ) by Th17
    .=( tanh.p + tanh.r )/ ( 1 + ( (sinh.p)/(cosh.p) )*( (sinh.r)/(cosh.r) )
  ) by Th17
    .=( tanh.p + tanh.r )/ ( 1 + ( tanh.p )*( (sinh.r)/(cosh.r) ) ) by Th17
    .=( tanh.p + tanh.r )/( 1 + ( tanh.p )*( tanh.r ) ) by Th17;
  hence thesis;
end;

theorem Th18:
  (sinh.p)^2 = 1/2*(cosh.(2*p) - 1) & (cosh.p)^2 = 1/2*(cosh.(2*p) + 1)
proof
A1: (cosh.p)^2 =( (exp_R.(p) + exp_R.(-p))/2 )*(cosh.p) by Def3
    .=( (exp_R.(p) + exp_R.(-p))/2 )*( (exp_R.(p) + exp_R.(-p))/2 ) by Def3
    .=( (exp_R.(p))*(exp_R.(p))+(exp_R.(p))*(exp_R.(-p)) +(exp_R.(-p))*(
  exp_R.(p))+(exp_R.(-p))*(exp_R.(-p)) )/4
    .=( (exp_R.(p+p))+(exp_R.(p))*(exp_R.(-p)) +(exp_R.(p))*(exp_R.(-p))+(
  exp_R.(-p))*(exp_R.(-p)) )/4 by Th12
    .=( (exp_R.(p+p))+(exp_R.(p))*(exp_R.(-p)) +(exp_R.(p))*(exp_R.(-p))+(
  exp_R.(-p+ -p)) )/4 by Th12
    .=( (exp_R.(p+p))+(exp_R.(p+ -p)) +(exp_R.(p))*(exp_R.(-p))+(exp_R.(-p+
  -p)) )/4 by Th12
    .=( (exp_R.(p+p))+(exp_R.(p+ -p)) +(exp_R.(p+ -p))+(exp_R.(-p+ -p)) )/4
  by Th12
    .=( exp_R.(2*p) + 2*(exp_R.(0)) +exp_R.(-p+-p) )/4
    .=( exp_R.(2*p) + 2*1 +exp_R.(-(2*p)) )/4 by SIN_COS:51,def 23
    .=1/2*( ( exp_R.(2*p) +exp_R.(-(2*p)) )/2 ) + (1*2)/(2*2)
    .=1/2*(cosh.(2*p))+ 1/2*((2*1/2)) by Def3
    .=1/2*( cosh.(2*p)+ 1 );
  (sinh.p)^2 =( (exp_R.(p) - exp_R.(-p))/2 )*(sinh.p) by Def1
    .=( (exp_R.(p) - exp_R.(-p))/2 )*( (exp_R.(p) - exp_R.(-p))/2 ) by Def1
    .=( (exp_R.(p))*(exp_R.(p))-(exp_R.(p))*(exp_R.(-p)) -(exp_R.(-p))*(
  exp_R.(p))+(exp_R.(-p))*(exp_R.(-p)) )/4
    .=( (exp_R.(p+p))-(exp_R.(p))*(exp_R.(-p)) -(exp_R.(p))*(exp_R.(-p))+(
  exp_R.(-p))*(exp_R.(-p)) )/4 by Th12
    .=( (exp_R.(p+p))-(exp_R.(p))*(exp_R.(-p)) -(exp_R.(p))*(exp_R.(-p))+(
  exp_R.(-p+ -p)) )/4 by Th12
    .=( (exp_R.(p+p))-(exp_R.(p+ -p)) -(exp_R.(p))*(exp_R.(-p))+(exp_R.(-p+
  -p)) )/4 by Th12
    .=( ((exp_R.(p+p))+ -(exp_R.(p+ -p)) -(exp_R.(p+ -p)) )+(exp_R.(-p+ -p))
  )/4 by Th12
    .=( exp_R.(2*p) + 2*(-(exp_R.(0))) +exp_R.(-p+-p) )/4
    .=( exp_R.(2*p) + 2*(-1) +exp_R.(-(2*p)) )/4 by SIN_COS:51,def 23
    .=1/2*( ( exp_R.(2*p) +exp_R.(-(2*p)) )/2 ) + ((-1)*2)/(2*2)
    .=1/2*(cosh.(2*p))+ 1/2*((2*(-1))/2) by Def3
    .=1/2*( cosh.(2*p) - 1 );
  hence thesis by A1;
end;

Lm6: sinh.(2*p) = 2*(sinh.p)*(cosh.p) & cosh.(2*p) = 2*(cosh.p)^2 - 1
proof
A1: 2*(cosh.p)^2 =2*( 1/2*(cosh.(2*p) + 1) ) by Th18
    .=cosh.(2*p) + 1;
  2*(sinh.p)*(cosh.p) =2*( (exp_R.(p) - exp_R.(-p))/2 )*(cosh.p) by Def1
    .=2*( (exp_R.(p) - exp_R.(-p))/2 )*( (exp_R.(p) + exp_R.(-p))/2 ) by Def3
    .=2/4*((exp_R.p)^2 - (exp_R.(-p))^2)
    .=2/4*( exp_R.(p+p)- (exp_R.(-p))*(exp_R.(-p)) ) by Th12
    .=2/4*(exp_R.(2*p)- exp_R.(-p+ -p)) by Th12
    .=( ( exp_R.(2*p)- exp_R.(-(2*p)) )/2 )*1
    .=sinh.(2*p) by Def1;
  hence thesis by A1;
end;

theorem Th19:
  cosh.(-p) = cosh.p & sinh.(-p) = -sinh.p & tanh.(-p) = -tanh.p
proof
A1: cosh.(-p) = (exp_R.(-p) + exp_R.(-(-p)))/2 by Def3
    .= cosh.p by Def3;
A2: sinh.(-p) = (exp_R.(-p) - exp_R.(-(-p)))/2 by Def1
    .= -(exp_R.(p) - exp_R.(-p))/2
    .= -sinh.p by Def1;
  then tanh.(-p) = (-sinh.p)/(cosh.(-p)) by Th17
    .= -(sinh.p)/(cosh.p) by A1,XCMPLX_1:187
    .= -tanh.p by Th17;
  hence thesis by A1,A2;
end;

Lm7: cosh.(p-r)=(cosh.p)*(cosh.r) - (sinh.p)*(sinh.r)
proof
  cosh.(p-r)=cosh.(p+ -r) .=(cosh.p)*(cosh.-r) + (sinh.p)*(sinh.-r) by Lm2
    .=(cosh.p)*(cosh.r) + (sinh.p)*(sinh.-r) by Th19
    .=(cosh.p)*(cosh.r) + (sinh.p)*(-sinh.r) by Th19
    .=(cosh.p)*(cosh.r) - ((sinh.p)*(sinh.r));
  hence thesis;
end;

theorem
  cosh.(p+r)=(cosh.p)*(cosh.r) + (sinh.p)*(sinh.r) & cosh.(p-r)=(cosh.p)
  *(cosh.r) - (sinh.p)*(sinh.r) by Lm2,Lm7;

Lm8: sinh.(p-r)=(sinh.p)*(cosh.r) - (cosh.p)*(sinh.r)
proof
  sinh.(p-r)=sinh.(p+-r) .=(sinh.p)*(cosh.-r) +(cosh.p)*(sinh.-r) by Lm3
    .=(sinh.p)*(cosh.r) + (cosh.p)*(sinh.-r) by Th19
    .=(sinh.p)*(cosh.r) + (cosh.p)*(-sinh.r) by Th19
    .=(sinh.p)*(cosh.r) -((cosh.p)*(sinh.r));
  hence thesis;
end;

theorem
  sinh.(p+r)=(sinh.p)*(cosh.r) + (cosh.p)*(sinh.r) & sinh.(p-r)=(sinh.p)
  *(cosh.r) - (cosh.p)*(sinh.r) by Lm3,Lm8;

Lm9: tanh.(p-r) = (tanh.p - tanh.r)/(1- (tanh.p)*(tanh.r))
proof
  tanh.(p-r)=tanh.(p+-r) .=(tanh.p + tanh.-r)/(1+ (tanh.p)*(tanh.-r)) by Lm5
    .=(tanh.p + -tanh.r)/(1+ (tanh.p)*(tanh.-r)) by Th19
    .=(tanh.p - tanh.r)/(1+ (tanh.p)*(-tanh.r)) by Th19
    .=(tanh.p - tanh.r)/(1 -(tanh.p*tanh.r));
  hence thesis;
end;

theorem
  tanh.(p+r) = (tanh.p + tanh.r)/(1+ (tanh.p)*(tanh.r)) & tanh.(p-r) = (
  tanh.p - tanh.r)/(1- (tanh.p)*(tanh.r)) by Lm5,Lm9;

theorem
  sinh.(2*p) = 2*(sinh.p)*(cosh.p) & cosh.(2*p) = 2*(cosh.p)^2 - 1 &
  tanh.(2*p) = (2*tanh.p)/(1+(tanh.p)^2)
proof
  tanh.(2*p)=tanh.(p+p) .=(tanh.p + tanh.p)/(1+ (tanh.p)*(tanh.p)) by Lm5
    .=(2*tanh.p)/(1+ (tanh.p)^2);
  hence thesis by Lm6;
end;

theorem Th24:
  (sinh.p)^2 - (sinh.q)^2 = (sinh.(p+q))*(sinh.(p-q)) & (sinh.(p+q
))*(sinh.(p-q)) = (cosh.p)^2 - (cosh.q)^2 & (sinh.p)^2 - (sinh.q)^2 = (cosh.p)
  ^2 - (cosh.q)^2
proof
A1: (sinh.(p+q))*(sinh.(p-q)) =((sinh.p)*(cosh.q) + (cosh.p)*(sinh.q))*(sinh
  .(p+(-q))) by Lm3
    .=((sinh.p)*(cosh.q) + (cosh.p)*(sinh.q)) *((sinh.p)*(cosh.(-q)) + (cosh
  .p)*(sinh.(-q))) by Lm3
    .=((sinh.p)*(cosh.q) + (cosh.p)*(sinh.q)) *((sinh.p)*(cosh.q) + (cosh.p)
  *(sinh.(-q))) by Th19
    .=((sinh.p)*(cosh.q) + (cosh.p)*(sinh.q)) *((sinh.p)*(cosh.q) + (cosh.p)
  *(-sinh.q)) by Th19
    .= (sinh.p)^2*(cosh.q)^2-(sinh.q)^2*(cosh.p)^2;
  then
A2: (sinh.(p+q))*(sinh.(p-q)) =(cosh.q)^2*( -((cosh.p)^2-(sinh.p)^2) ) +( ((
  cosh.p)^2*(cosh.q)^2)+ -(sinh.q)^2*(cosh.p)^2 )
    .=(cosh.q)^2*( -1 ) +( (cosh.p)^2*( (cosh.q)^2-(sinh.q)^2 ) ) by Th14
    .=(cosh.q)^2*( -1 )+ ((cosh.p)^2*1) by Th14
    .=(cosh.p)^2 -(cosh.q)^2;
  (sinh.(p+q))*(sinh.(p-q)) =(sinh.p)^2*( (cosh.q)^2 -(sinh.q)^2 ) + ((
  sinh.q)^2*(sinh.p)^2)-(sinh.q)^2*(cosh.p)^2 by A1
    .=(sinh.p)^2*1 + (sinh.p)^2*(sinh.q)^2-(sinh.q)^2*(cosh.p)^2 by Th14
    .=(sinh.p)^2 + (sinh.q)^2*(-((cosh.p)^2-(sinh.p)^2))
    .=(sinh.p)^2 + (sinh.q)^2*(-1) by Th14
    .=(sinh.p)^2 - (sinh.q)^2;
  hence thesis by A2;
end;

theorem Th25:
  (sinh.p)^2 + (cosh.q)^2 = (cosh.(p+q))*(cosh.(p-q)) & (cosh.(p+q
))*(cosh.(p-q)) = (cosh.p)^2 + (sinh.q)^2 & (sinh.p)^2 + (cosh.q)^2 = (cosh.p)
  ^2 + (sinh.q)^2
proof
A1: (cosh.(p+q))*(cosh.(p-q)) =( (cosh.p)*(cosh.q) + (sinh.p)*(sinh.q) )* (
  cosh.(p + -q)) by Lm2
    .=( ( (cosh.p)*(cosh.q) ) + ( (sinh.p)*(sinh.q) ) )* ( ( (cosh.p)*(cosh.
  -q) ) + ( (sinh.p)*(sinh.-q) ) ) by Lm2
    .=( (cosh.p)*(cosh.q) )*( (cosh.p)*(cosh.-q) ) +( (cosh.p)*(cosh.q) )*(
(sinh.p)*(sinh.-q) ) +( (sinh.p)*(sinh.q) )*( (cosh.p)*(cosh.-q) ) +( (sinh.p)*
  (sinh.q) )*( (sinh.p)*(sinh.-q) )
    .=( (cosh.p)*(cosh.q) )*( (cosh.p)*(cosh.q) ) +( (cosh.p)*(cosh.q) )*( (
sinh.-q)*(sinh.p) ) +( (sinh.p)*(sinh.q) )*( (cosh.p)*(cosh.-q) ) +( (sinh.p)*(
  sinh.q) )*( (sinh.-q)*(sinh.p) ) by Th19
    .=( (cosh.p)*(cosh.q) )^2 +( (sinh.-q)*(sinh.p) )*( (cosh.p)*(cosh.q) )
+( (sinh.p)*(sinh.q) )*( (cosh.p)*(cosh.q) ) +( (sinh.-q)*(sinh.p) )*( (sinh.p)
  *(sinh.q) ) by Th19
    .=( (cosh.p)*(cosh.q) )^2 +( (-sinh.q)*(sinh.p) )*( (cosh.p)*(cosh.q) )
+( (sinh.p)*(sinh.q) )*( (cosh.p)*(cosh.q) ) +( (sinh.-q)*(sinh.p) )*( (sinh.p)
  *(sinh.q) ) by Th19
    .=( (cosh.p)*(cosh.q) )^2 + 0 +( (-1)*(sinh.q)*(sinh.p) )*( (sinh.p)*(
  sinh.q) ) by Th19
    .=(cosh.p)^2*(cosh.q)^2 - (sinh.q)^2*(sinh.p)^2;
  then
A2: (cosh.(p+q))*(cosh.(p-q)) =(cosh.p)^2*( (cosh.q)^2-(sinh.q)^2 ) + ((cosh
  .p)^2*(sinh.q)^2) + -((sinh.q)^2*(sinh.p)^2)
    .=(cosh.p)^2*1 + ((cosh.p)^2*(sinh.q)^2) + -((sinh.q)^2*(sinh.p)^2) by Th14
    .=(cosh.p)^2+(sinh.q)^2*((cosh.p)^2-(sinh.p)^2)
    .=(cosh.p)^2+(sinh.q)^2*1 by Th14
    .=(cosh.p)^2+(sinh.q)^2;
  (cosh.(p+q))*(cosh.(p-q)) =(cosh.q)^2*((cosh.p)^2-(sinh.p)^2) +( (cosh.q
  )^2*(sinh.p)^2 + -(sinh.p)^2*(sinh.q)^2 ) by A1
    .=(cosh.q)^2*1 +( (cosh.q)^2*(sinh.p)^2-(sinh.p)^2*(sinh.q)^2 ) by Th14
    .=(cosh.q)^2+(sinh.p)^2*( (cosh.q)^2-(sinh.q)^2 )
    .=(cosh.q)^2+(sinh.p)^2*1 by Th14
    .=(cosh.q)^2+(sinh.p)^2;
  hence thesis by A2;
end;

theorem
  sinh.p + sinh.r = 2*(sinh.(p/2+r/2))*(cosh.(p/2-r/2)) & sinh.p - sinh.
  r = 2*(sinh.(p/2-r/2))*(cosh.(p/2+r/2))
proof
A1: 2*(sinh.(p/2-r/2))*(cosh.(p/2+r/2)) =2*( (sinh.(p/2))*(cosh.(r/2))-(cosh
  .(p/2))*(sinh.(r/2)) ) *(cosh.(p/2+r/2)) by Lm8
    .=2*((sinh.(p/2))*(cosh.(r/2))-(cosh.(p/2))*(sinh.(r/2))) *((cosh.(p/2))
  *(cosh.(r/2))+(sinh.(p/2))*(sinh.(r/2))) by Lm2
    .=2*( ((cosh.(r/2))*(sinh.(r/2)))*(-((cosh.(p/2))^2-(sinh.(p/2))^2)) +(
sinh.(p/2))*( (cosh.(p/2))*(cosh.(r/2))^2 ) -(cosh.(p/2))*( (sinh.(p/2))*(sinh.
  (r/2))^2 ) )
    .=2*( ((cosh.(r/2))*(sinh.(r/2)))*(-1) +(sinh.(p/2))*( (cosh.(p/2))*(
  cosh.(r/2))^2 ) -(cosh.(p/2))*( (sinh.(p/2))*(sinh.(r/2))^2 ) ) by Th14
    .=2*( ((sinh.(p/2))*(cosh.(p/2)))*((cosh.(r/2))^2-(sinh.(r/2))^2) +(-1)*
  ((cosh.(r/2))*(sinh.(r/2))) )
    .=2*( ((sinh.(p/2))*(cosh.(p/2)))*1 +(-1)*((cosh.(r/2))*(sinh.(r/2))) )
  by Th14
    .= 2*(sinh.(p/2))*(cosh.(p/2))-2*((sinh.(r/2))*(cosh.(r/2)))
    .=sinh.(2*(p/2)) - 2*(sinh.(r/2))*(cosh.(r/2)) by Lm6
    .=sinh.p - sinh.(2*(r/2)) by Lm6
    .=sinh.p - sinh.r;
  2*(sinh.(p/2+r/2))*(cosh.(p/2-r/2)) =2*((sinh.(p/2))*(cosh.(r/2))+(cosh.
  (p/2))*(sinh.(r/2))) *(cosh.(p/2-r/2)) by Lm3
    .=2*((sinh.(p/2))*(cosh.(r/2))+(cosh.(p/2))*(sinh.(r/2))) *((cosh.(p/2))
  *(cosh.(r/2))-(sinh.(p/2))*(sinh.(r/2))) by Lm7
    .=2*( ((sinh.(r/2))*(cosh.(r/2)))*((cosh.(p/2))^2-(sinh.(p/2))^2) +(sinh
.(p/2))*( (cosh.(p/2))*(cosh.(r/2))^2 ) -(cosh.(p/2))*( (sinh.(p/2))*(sinh.(r/2
  ))^2 ) )
    .=2*( ((sinh.(r/2))*(cosh.(r/2)))*1 +(sinh.(p/2))*( (cosh.(p/2))*(cosh.(
  r/2))^2 ) -(cosh.(p/2))*( (sinh.(p/2))*(sinh.(r/2))^2 ) ) by Th14
    .=2*( ((sinh.(p/2))*(cosh.(p/2)))*((cosh.(r/2))^2 -(sinh.(r/2))^2 ) +(
  sinh.(r/2))*(cosh.(r/2)) )
    .=2*( ((sinh.(p/2))*(cosh.(p/2))*1)+(sinh.(r/2))*(cosh.(r/2)) ) by Th14
    .=2*(sinh.(p/2))*(cosh.(p/2)) + 2*((sinh.(r/2))*(cosh.(r/2)))
    .=sinh.(2*(p/2)) + 2*(sinh.(r/2))*(cosh.(r/2)) by Lm6
    .=sinh.p + sinh.(2*(r/2)) by Lm6
    .=sinh.p + sinh.r;
  hence thesis by A1;
end;

theorem
  cosh.p + cosh.r = 2*(cosh.(p/2+r/2))*(cosh.(p/2-r/2)) & cosh.p - cosh.
  r = 2*(sinh.(p/2-r/2))*(sinh.(p/2+r/2))
proof
A1: 2*(sinh.(p/2-r/2))*(sinh.(p/2+r/2)) =2*((sinh.(p/2-r/2))*(sinh.(p/2+r/2) ))
    .=2*( (cosh.(p/2))^2 - (cosh.(r/2))^2 ) by Th24
    .=2*( 1/2*(cosh.(2*(p/2)) + 1)-(cosh.(r/2))^2 ) by Th18
    .=2*( 1/2*(cosh.p + 1)-1/2*(cosh.(2*(r/2)) + 1) ) by Th18
    .=cosh.p-cosh.r;
  2*(cosh.(p/2+r/2))*(cosh.(p/2-r/2)) =2*((cosh.(p/2+r/2))*(cosh.(p/2-r/2) ))
    .=2*( (sinh.(p/2))^2+(cosh.(r/2))^2 ) by Th25
    .=2*( 1/2*(cosh.(2*(p/2)) - 1)+(cosh.(r/2))^2 ) by Th18
    .=2*( 1/2*(cosh.p - 1)+1/2*(cosh.(2*(r/2)) + 1) ) by Th18
    .=cosh.r+cosh.p;
  hence thesis by A1;
end;

theorem
  tanh.p + tanh.r = (sinh.(p+r))/((cosh.p)*(cosh.r)) & tanh.p - tanh.r =
  (sinh.(p-r))/((cosh.p)*(cosh.r))
proof
A1: (sinh.(p-r))/((cosh.p)*(cosh.r)) =( (sinh.p)*(cosh.r)-(cosh.p)*(sinh.r)
  )/((cosh.p)*(cosh.r)) by Lm8
    .= ((sinh.p)*(cosh.r))/((cosh.p)*(cosh.r)) -((cosh.p)*(sinh.r))/((cosh.p
  )*(cosh.r)) by XCMPLX_1:120
    .=(sinh.p)/(cosh.p)-((cosh.p)*(sinh.r))/((cosh.p)*(cosh.r)) by Th15,
XCMPLX_1:91
    .=(sinh.p)/(cosh.p)-(sinh.r)/(cosh.r) by Th15,XCMPLX_1:91
    .=tanh.p-(sinh.r)/(cosh.r) by Th17
    .=tanh.p-tanh.r by Th17;
  (sinh.(p+r))/((cosh.p)*(cosh.r)) =( (sinh.p)*(cosh.r)+(cosh.p)*(sinh.r)
  )/((cosh.p)*(cosh.r)) by Lm3
    .= ((sinh.p)*(cosh.r))/((cosh.p)*(cosh.r)) +((cosh.p)*(sinh.r))/((cosh.p
  )*(cosh.r)) by XCMPLX_1:62
    .=(sinh.p)/(cosh.p)+((cosh.p)*(sinh.r))/((cosh.p)*(cosh.r)) by Th15,
XCMPLX_1:91
    .=(sinh.p)/(cosh.p)+(sinh.r)/(cosh.r) by Th15,XCMPLX_1:91
    .=tanh.p+(sinh.r)/(cosh.r) by Th17
    .=tanh.p+tanh.r by Th17;
  hence thesis by A1;
end;

theorem
  (cosh.p + sinh.p) |^ n = cosh.(n*p) + sinh.(n*p)
proof
  defpred X[Nat] means for p holds (cosh.p + sinh.p) |^ $1 = cosh.(
  $1*p) + sinh.($1*p);
A1: for n st X[n] holds X[n+1]
  proof
    let n such that
A2: for p holds (cosh.p + sinh.p) |^ n = cosh.(n*p) + sinh.(n*p);
    for p holds (cosh.p + sinh.p) |^ (n+1) = cosh.((n+1)*p) + sinh.((n+1)* p)
    proof
      let p;
A3:   (cosh.(n*p))*(cosh.p)+(sinh.(n*p))*(sinh.p) =((exp_R.(n*p) + exp_R.(
      -n*p))/2)*(cosh.p) +(sinh.(n*p))*(sinh.p) by Def3
        .=((exp_R.(n*p) + exp_R.(-n*p))/2)*((exp_R.p + exp_R.-p)/2) +(sinh.(
      n*p))*(sinh.p) by Def3
        .=((exp_R.(n*p) + exp_R.(-n*p))/2)*((exp_R.p + exp_R.-p)/2) +((exp_R
      .(n*p) - exp_R.(-n*p))/2)*(sinh.p) by Def1
        .=((exp_R.(n*p))/2 + (exp_R.(-n*p))/2)*((exp_R.p)/2 + (exp_R.-p)/2)
      +((exp_R.(n*p))/2 - (exp_R.(-n*p))/2)*((exp_R.p - exp_R.-p)/2) by Def1
        .=2*( ( (exp_R.(n*p))*(exp_R.p) )/(2*2) ) +2*(((exp_R.(-n*p))/2)*((
      exp_R.-p)/2))
        .=2*( exp_R.(p*n+p*1)/(2*2) ) +2*( ( (exp_R.(-n*p))*(exp_R.-p) )/(2*
      2) ) by Th12
        .=2*( exp_R.(p*(n+1))/(2*2) ) +2*( exp_R.(-n*p+-p)/(2*2) ) by Th12
        .=(exp_R.(p*(n+1))+exp_R.(-(p*(n+1))))/2
        .=cosh.(p*(n+1)) by Def3;
A4:   (cosh.(n*p))*(sinh.p)+(sinh.(n*p))*(cosh.p) =((exp_R.(n*p) + exp_R.(
      -n*p))/2)*(sinh.p) +(sinh.(n*p))*(cosh.p) by Def3
        .=((exp_R.(n*p) + exp_R.(-n*p))/2)*((exp_R.p - exp_R.-p)/2) +(sinh.(
      n*p))*(cosh.p) by Def1
        .=((exp_R.(n*p) + exp_R.(-n*p))/2)*((exp_R.p - exp_R.-p)/2) +((exp_R
      .(n*p) - exp_R.(-n*p))/2)*(cosh.p) by Def1
        .=((exp_R.(n*p))/2)*((exp_R.p)/2) - ((exp_R.(n*p))/2)*((exp_R.-p)/2)
+((exp_R.(-n*p))/2)*((exp_R.p)/2) - ((exp_R.(-n*p))/2)*((exp_R.-p)/2) +((exp_R.
      (n*p) - exp_R.(-n*p))/2)*((exp_R.p + exp_R.-p)/2) by Def3
        .=2*(( (exp_R.(n*p))*(exp_R.p) )/4 ) +2*( -( ((exp_R.(-n*p))*(exp_R.
      -p))/(2*2) ))
        .=2*( exp_R.(n*p+p)/4 ) +2*( -( ((exp_R.(-n*p))*(exp_R.-p))/4 )) by
Th12
        .=2*( exp_R.(n*p+p)/4 ) +2*( -( (exp_R.(-n*p+-p))/4 ) ) by Th12
        .=( exp_R.(p*(n+1))- exp_R.(-(p*(n+1))) )/2
        .=sinh.(p*(n+1)) by Def1;
      (cosh.p + sinh.p) |^ (n+1) =(cosh.p + sinh.p) |^ n * (cosh.p + sinh.
      p) by NEWTON:6
        .=(cosh.(n*p) + sinh.(n*p)) * (cosh.p + sinh.p) by A2
        .=sinh.((n+1)*p)+cosh.((n+1)*p) by A3,A4;
      hence thesis;
    end;
    hence thesis;
  end;
A5: X[0] by Th15,Th16,NEWTON:4;
  for n holds X[n] from NAT_1:sch 2(A5,A1);
  hence thesis;
end;

theorem Th30:
  dom sinh=REAL & dom cosh=REAL & dom tanh=REAL by FUNCT_2:def 1;

Lm10: for d being Real holds compreal.d=(-1)*d
proof
  let d be Real;
  thus compreal.d= -d by BINOP_2:def 7
    .= (-1)*d;
end;

Lm11: dom(compreal)=REAL by FUNCT_2:def 1;

Lm12: for f being PartFunc of REAL,REAL st f = compreal holds for p holds f
is_differentiable_in p & diff(f,p) = -1
proof
  reconsider f2 = REAL --> In(0,REAL) as Function of REAL,REAL;
  deffunc U(Real) = In(-$1,REAL);
  let f be PartFunc of REAL,REAL such that
A1: f = compreal;
  let p;
  consider f1 being Function of REAL,REAL such that
A2: for p be Element of REAL holds f1.p=U(p) from FUNCT_2:sch 4;
A3: dom f2 = REAL by FUNCOP_1:13;
  for hy1 be 0-convergent non-zero Real_Sequence holds (hy1")(#)(f2/*hy1) is
  convergent & lim ((hy1")(#)(f2/*hy1)) = 0
  proof
    let hy1 be 0-convergent non-zero Real_Sequence;
A4: for n being Nat holds ((hy1")(#)(f2/*hy1)).n=In(0,REAL)
    proof
      let n be Nat;
A5:   rng hy1 c= dom f2 by A3;
A6:   n in NAT by ORDINAL1:def 12;
      ((hy1")(#)(f2/*hy1)).n = (hy1".n)*((f2/*hy1).n) by SEQ_1:8
        .= (hy1.n)"*((f2/*hy1).n) by VALUED_1:10;
      then ((hy1")(#)(f2/*hy1)).n =(hy1.n)"*(f2.(hy1.n)) by A6,A5,FUNCT_2:108
        .=(hy1.n)"*(0) by FUNCOP_1:7
        .= 0;
      hence thesis;
    end;
    then
A7: (hy1")(#)(f2/*hy1) is constant by VALUED_0:def 18;
    for n holds lim ((hy1")(#)(f2/*hy1)) = 0
    proof
      let n;
      lim ((hy1")(#)(f2/*hy1)) =((hy1")(#)(f2/*hy1)).n by A7,SEQ_4:26
        .=0 by A4;
      hence thesis;
    end;
    hence thesis by A7;
  end;
  then
A8: f2 is RestFunc by FDIFF_1:def 2;
A9: ex N being Neighbourhood of p st N c= dom compreal & for r st r
  in N holds compreal.r - compreal.p = f1.(r-p) + f2.(r-p)
  proof
    take ].p-1,p+1 .[;
    for r st r in ].p-1,p+1 .[ holds compreal.r - compreal.p = f1.(r-p) +
    f2.(r-p)
    proof
      let r;
A10:  for d be Real holds f2.d = 0
      by XREAL_0:def 1,FUNCOP_1:7;
A11:  for p be Real holds f1.p = - p
      proof
        let p be Real;
        p in REAL by XREAL_0:def 1;
        then f1.p=U(p) by A2;
        hence thesis;
      end;
      compreal.r-compreal.p =(-1)*r - compreal.p by Lm10
        .=(-1)*r -(-1)*p + 0 by Lm10
        .=-(r-p) + f2.(r-p) by A10
        .=f1.(r-p) + f2.(r-p) by A11;
      hence thesis;
    end;
    hence thesis by Lm11,RCOMP_1:def 6;
  end;
     reconsider pp = p as Real;
A12: 1 in REAL by XREAL_0:def 1;
  for p holds f1.p=(-1)*p
  proof
    let p;
     reconsider pp = p as Element of REAL by XREAL_0:def 1;
    f1.p = U(pp) by A2;
    hence thesis;
  end;
  then
A13: f1 is LinearFunc by FDIFF_1:def 3;
  then f is_differentiable_in pp by A1,A8,A9,FDIFF_1:def 4;
  then diff(f,pp) = f1.1 by A1,A13,A8,A9,FDIFF_1:def 5
    .= U(1) by A2,A12;
  hence thesis by A1,A13,A8,A9,FDIFF_1:def 4;
end;

Lm13: for f being PartFunc of REAL,REAL st f = compreal holds exp_R*f
is_differentiable_in p & diff(exp_R*f,p) = (-1)*exp_R.(f.p)
proof
  let f be PartFunc of REAL,REAL such that
A1: f = compreal;
A2: p is Real & exp_R is_differentiable_in f.p by SIN_COS:65;
A3: f is_differentiable_in p by A1,Lm12;
  then diff(exp_R*f,p) =diff(exp_R,f.p)*diff(f,p) by A2,FDIFF_2:13
    .=diff(exp_R,f.p)*(-1) by A1,Lm12
    .=exp_R.(f.p)*(-1) by SIN_COS:65;
  hence thesis by A2,A3,FDIFF_2:13;
end;

Lm14: for f being PartFunc of REAL,REAL st f = compreal holds exp_R.((-1)*p)=
(exp_R*f).p
proof
  let f be PartFunc of REAL,REAL;
A1: p in REAL by XREAL_0:def 1;
  assume
A2: f = compreal;
  then exp_R.((-1)*p) =exp_R.(f.p) by Lm10
    .=(exp_R*f).p by A2,A1,FUNCT_2:15;
  hence thesis;
end;

Lm15: for f being PartFunc of REAL,REAL st f = compreal holds (exp_R-(exp_R*f)
) is_differentiable_in p & (exp_R+(exp_R*f)) is_differentiable_in p & diff(
exp_R-(exp_R*f),p)=exp_R.p+exp_R.(-p) & diff(exp_R+(exp_R*f),p)=exp_R.p-exp_R.(
-p)
proof
  let f be PartFunc of REAL,REAL;
A1: p is Real & exp_R is_differentiable_in p by SIN_COS:65;
  assume
A2: f = compreal;
  then
A3: exp_R*f is_differentiable_in p by Lm13;
  then
A4: diff(exp_R+(exp_R*f),p) =diff(exp_R,p)+diff((exp_R*f),p) by A1,FDIFF_1:13
    .=exp_R.p+diff((exp_R*f),p) by SIN_COS:65
    .=exp_R.p+ (-1)*exp_R.(f.p) by A2,Lm13
    .=exp_R.p+ (-1)*exp_R.((-1)*p) by A2,Lm10
    .=exp_R.p-exp_R.(-p);
  diff(exp_R-(exp_R*f),p) =diff(exp_R,p)-diff((exp_R*f),p) by A1,A3,FDIFF_1:14
    .=exp_R.p-diff((exp_R*f),p) by SIN_COS:65
    .=exp_R.p- (-1)*exp_R.(f.p) by A2,Lm13
    .=exp_R.p- (-1)*exp_R.((-1)*p) by A2,Lm10
    .=exp_R.p+exp_R.(-p);
  hence thesis by A1,A3,A4,FDIFF_1:13,14;
end;

Lm16: for f being PartFunc of REAL,REAL st f = compreal holds (1/2)(#)(exp_R-
exp_R*f) is_differentiable_in p & diff(((1/2)(#)(exp_R-exp_R*f)),p) = (1/2)*
diff((exp_R-(exp_R*f)),p)
proof
  let f be PartFunc of REAL,REAL;
  assume f = compreal;
  then p is Real & (exp_R-exp_R*f) is_differentiable_in p by Lm15;
  hence thesis by FDIFF_1:15;
end;

Lm17: for ff being PartFunc of REAL,REAL st ff = compreal holds sinh.p=((1/2)
(#)(exp_R-exp_R*ff)).p
proof
A1: p in REAL by XREAL_0:def 1;
  let ff be PartFunc of REAL,REAL such that
A2: ff = compreal;
A3: dom (exp_R - exp_R*ff) = dom exp_R /\ dom (exp_R*ff) & dom (exp_R*ff) =
  REAL by A2,FUNCT_2:def 1,VALUED_1:12;
A4: dom((1/2)(#)(exp_R-exp_R*ff))=REAL by A2,FUNCT_2:def 1;
  sinh.p = (exp_R.p - exp_R.(-p))/2 by Def1
    .=(1/2)*(exp_R.p - exp_R.((-1)*p))
    .=(1/2)*(exp_R.p - (exp_R*ff).p) by A2,Lm14
    .=(1/2)*((exp_R - exp_R*ff).p) by A1,A3,SIN_COS:47,VALUED_1:13
    .=((1/2)(#)(exp_R-exp_R*ff)).p by A1,A4,VALUED_1:def 5;
  hence thesis;
end;

Lm18: for ff being PartFunc of REAL,REAL st ff = compreal holds sinh = (1/2)
(#)(exp_R-exp_R*ff)
proof
  let ff be PartFunc of REAL,REAL;
  assume ff = compreal;
  then
A1: REAL=dom((1/2)(#)(exp_R-exp_R*ff)) & for p being Element of REAL st p in
  REAL holds sinh.p=((1/2)(#)(exp_R-exp_R*ff)).p by Lm17,FUNCT_2:def 1;
  REAL=dom(sinh) by FUNCT_2:def 1;
  hence thesis by A1,PARTFUN1:5;
end;

theorem Th31:
  sinh is_differentiable_in p & diff(sinh,p)=cosh.p
proof
  set ff = compreal;
A1: sinh = (1/2)(#)(exp_R-exp_R*ff) by Lm18;
  diff(sinh,p) =diff( ((1/2)(#)(exp_R-exp_R*ff)),p ) by Lm18
    .=(1/2)*diff((exp_R-(exp_R*ff)),p) by Lm16
    .=(1/2)*( exp_R.p+exp_R.(-p) ) by Lm15
    .=( exp_R.p+exp_R.(-p) )/2
    .=cosh.p by Def3;
  hence thesis by A1,Lm16;
end;

Lm19: for ff being PartFunc of REAL,REAL st ff = compreal holds (1/2)(#)(exp_R
+exp_R*ff) is_differentiable_in p & diff(((1/2)(#)(exp_R+exp_R*ff)),p) = (1/2)*
diff((exp_R+(exp_R*ff)),p)
proof
  let ff be PartFunc of REAL,REAL;
  assume ff = compreal;
  then p is Real & (exp_R+exp_R*ff) is_differentiable_in p by Lm15;
  hence thesis by FDIFF_1:15;
end;

Lm20: for ff being PartFunc of REAL,REAL st ff = compreal holds cosh.p=((1/2)
(#)(exp_R+exp_R*ff)).p
proof
A1: p in REAL by XREAL_0:def 1;
  let ff be PartFunc of REAL,REAL such that
A2: ff = compreal;
A3: dom (exp_R + exp_R*ff) = dom exp_R /\ dom (exp_R*ff) & dom (exp_R*ff) =
  REAL by A2,FUNCT_2:def 1,VALUED_1:def 1;
A4: dom((1/2)(#)(exp_R+exp_R*ff))=REAL by A2,FUNCT_2:def 1;
  cosh.p =(exp_R.p + exp_R.(-p))/2 by Def3
    .=(1/2)*(exp_R.p + exp_R.((-1)*p))
    .=(1/2)*(exp_R.p + (exp_R*ff).p) by A2,Lm14
    .=(1/2)*((exp_R + exp_R*ff).p) by A1,A3,SIN_COS:47,VALUED_1:def 1
    .=((1/2)(#)(exp_R+exp_R*ff)).p by A1,A4,VALUED_1:def 5;
  hence thesis;
end;

Lm21: for ff being PartFunc of REAL,REAL st ff = compreal holds cosh = (1/2)
(#)(exp_R+exp_R*ff)
proof
  let ff be PartFunc of REAL,REAL;
  assume ff = compreal;
  then
A1: REAL=dom((1/2)(#)(exp_R+exp_R*ff)) & for p being Element of REAL st p in
  REAL holds cosh.p=((1/2)(#)(exp_R+exp_R*ff)).p by Lm20,FUNCT_2:def 1;
  REAL=dom(cosh) by FUNCT_2:def 1;
  hence thesis by A1,PARTFUN1:5;
end;

theorem Th32:
  cosh is_differentiable_in p & diff(cosh,p)=sinh.p
proof
  reconsider ff = compreal as PartFunc of REAL,REAL;
A1: cosh = (1/2)(#)(exp_R+exp_R*ff) by Lm21;
  diff(cosh,p) =diff( ((1/2)(#)(exp_R+exp_R*ff)),p ) by Lm21
    .=(1/2)*diff((exp_R+(exp_R*ff)),p) by Lm19
    .=(1/2)*( exp_R.p-exp_R.(-p) ) by Lm15
    .=( exp_R.p-exp_R.(-p) )/2
    .=sinh.p by Def1;
  hence thesis by A1,Lm19;
end;

Lm22: sinh/cosh is_differentiable_in p & diff(sinh/cosh,p) = 1/(cosh.p)^2
proof
A1: p is Real & cosh.p<>0 by Th15;
A2: sinh is_differentiable_in p & cosh is_differentiable_in p by Th31,Th32;
  then
  diff(sinh/cosh,p) =(diff(sinh,p) * cosh.p - diff(cosh,p)*sinh.p)/(cosh.p
  )^2 by A1,FDIFF_2:14
    .=((cosh.p)*(cosh.p) - diff(cosh,p) * sinh.p)/(cosh.p)^2 by Th31
    .=( (cosh.p)^2 - (sinh.p)*(sinh.p) )/(cosh.p)^2 by Th32
    .=1/(cosh.p)^2 by Th14;
  hence thesis by A1,A2,FDIFF_2:14;
end;

Lm23: tanh=sinh/cosh
proof
  not 0 in rng(cosh)
  proof
    assume 0 in rng(cosh);
    then ex p be Element of REAL st p in dom(cosh) & 0=cosh.p by PARTFUN1:3;
    hence contradiction by Th15;
  end;
  then
A1: dom (sinh/cosh) = dom sinh /\ (dom cosh \ cosh"{0}) & cosh"{0}={} by
FUNCT_1:72,RFUNCT_1:def 1;
  for p being Element of REAL st p in REAL holds tanh.p=(sinh/cosh).p
  proof
    let p;
    p in REAL by XREAL_0:def 1;
    then (sinh/cosh).p =sinh.p*(cosh.p)" by A1,Th30,RFUNCT_1:def 1
      .=(sinh.p)/(cosh.p) by XCMPLX_0:def 9
      .=tanh.p by Th17;
    hence thesis;
  end;
  hence thesis by A1,Th30,PARTFUN1:5;
end;

theorem
  tanh is_differentiable_in p & diff(tanh,p)=1/(cosh.p)^2 by Lm22,Lm23;

theorem Th34:
  sinh is_differentiable_on REAL & diff(sinh,p)=cosh.p
proof
A1: [#]REAL is open Subset of REAL & REAL c= dom(sinh) by FUNCT_2:def 1;
  for p st p in REAL holds sinh is_differentiable_in p by Th31;
  hence thesis by A1,Th31,FDIFF_1:9;
end;

theorem Th35:
  cosh is_differentiable_on REAL & diff(cosh,p)=sinh.p
proof
A1: [#]REAL is open Subset of REAL & REAL c=dom cosh by FUNCT_2:def 1;
  for r st r in REAL holds cosh is_differentiable_in r by Th32;
  hence thesis by A1,Th32,FDIFF_1:9;
end;

theorem Th36:
  tanh is_differentiable_on REAL & diff(tanh,p)=1/(cosh.p)^2
proof
  [#]REAL is open Subset of REAL & for p st p in REAL holds tanh
  is_differentiable_in p by Lm22,Lm23;
  hence thesis by Lm22,Lm23,Th30,FDIFF_1:9;
end;

Lm24: exp_R.p + exp_R.-p >= 2
proof
A1: exp_R.p >= 0 & exp_R.-p>= 0 by SIN_COS:54;
  2*sqrt((exp_R.p) * (exp_R.-p)) =2*sqrt(exp_R.(p+-p)) by Th12
    .=2*1 by SIN_COS:51,def 23,SQUARE_1:18;
  hence thesis by A1,Th1;
end;

theorem
  cosh.p >= 1
proof
  (exp_R.p + exp_R.-p)/2 >= 2/2 by Lm24,XREAL_1:72;
  hence thesis by Def3;
end;

theorem
  sinh is_continuous_in p by Th31,FDIFF_1:24;

theorem
  cosh is_continuous_in p by Th32,FDIFF_1:24;

theorem
  tanh is_continuous_in p by Lm22,Lm23,FDIFF_1:24;

theorem
  sinh|REAL is continuous by Th34,FDIFF_1:25;

theorem
  cosh|REAL is continuous by Th35,FDIFF_1:25;

theorem
  tanh|REAL is continuous by Th36,FDIFF_1:25;

theorem
  tanh.p<1 & tanh.p>-1
proof
A1: exp_R.p > 0 & exp_R.(-p) > 0 by SIN_COS:54;
  thus tanh.p<1
  proof
    assume tanh.p>=1;
    then
A2: (exp_R.p - exp_R.(-p))/(exp_R.p + exp_R.(-p))>=1 by Def5;
    exp_R.p > 0 & exp_R.(-p) > 0 by SIN_COS:54;
    then
A3: (exp_R.p - exp_R.(-p))/(exp_R.p + exp_R.(-p))*(exp_R.p + exp_R.(-p)) =
    exp_R.p - exp_R.(-p) by XCMPLX_1:87;
    (exp_R.p + exp_R.-p) >= 2 by Lm24;
    then
    (exp_R.p - exp_R.(-p))/(exp_R.p + exp_R.(-p))*(exp_R.p + exp_R.(-p ))
    >=1*(exp_R.p + exp_R.(-p)) by A2,XREAL_1:64;
    then (exp_R.p - exp_R.(-p))-exp_R.p >= (exp_R.p + exp_R.(-p))-exp_R.p by A3
,XREAL_1:9;
    then
A4: -exp_R.(-p)+exp_R.(-p) >= exp_R.(-p)+exp_R.(-p) by XREAL_1:6;
    exp_R.(-p) > 0 by SIN_COS:54;
    hence contradiction by A4;
  end;
  assume tanh.p<=-1;
  then
A5: (exp_R.p - exp_R.(-p))/(exp_R.p + exp_R.(-p))<=-1 by Def5;
  (exp_R.p + exp_R.-p) >= 2 by Lm24;
  then
  (exp_R.p - exp_R.(-p))/(exp_R.p + exp_R.(-p))*(exp_R.p + exp_R.(-p)) <=
  (-1)*(exp_R.p + exp_R.(-p)) by A5,XREAL_1:64;
  then exp_R.p - exp_R.(-p)<=-(exp_R.p + exp_R.(-p)) by A1,XCMPLX_1:87;
  then
A6: exp_R.p - exp_R.(-p)+exp_R.(-p) <=-exp_R.p + -exp_R.(-p)+exp_R.(-p) by
XREAL_1:6;
  exp_R.p > 0 by SIN_COS:54;
  hence contradiction by A6;
end;

subset_name: curated