\section{Conclusion}\label{sec-conc}
In this paper, we have described an algorithm, which given a timed automaton $A$, produces another timed automaton $A_4$ with the smallest number of clocks that is timed bisimilar to $A$.
Since $A_4$ is timed bisimilar to $A$, it also accepts the same timed language as $A$.
The problem appears in a different form
in \cite{LLW1} where it was left open.
In \cite{LLW1}, this problem has been described as
given a timed automaton $A$ and a set $C$ of clocks,
whether there exists a timed automaton that is
timed bisimilar to the original automaton $A$ and has $|C|$ clocks.
It is easy to see that a solution to our problem solves
the problem in \cite{LLW1} and vice versa.

For reducing the number of clocks of the timed automaton, we rely on a semantic representation of the timed automaton rather than its syntactic form as in \cite{DY1}. This helps us to reason about the behaviour of the timed automaton more effectively.
We could just as well have used a region graph construction,
the asymptotic complexity would remain the same.
However, the zone graph we use in our approach is usually much smaller in size than the region graph.
Besides, the size of the zone graph is independent of the constants used in the timed automaton.
There is an exponential increase in the number of locations while producing the TA with the minimal number of clocks.
However, there is no addition of new valuations to the underlying state space of the timed automaton in stage 2,
since the splitting of a location $l$,
described in stage 2, involves distributing the zones of $l$ across the locations $l$ is split into.

We note that pre-stability of the zones with respect to transitions is essential for the operations in stage 2.
Notice that regions are inherently pre-stable \cite{Sri1} and therefore the region graph construction would work here as well
though the number of locations obtained through splitting would be larger than in the case of zone graphs.

Our algorithm has a 2-EXPTIME upper bound,
however, it needs to be used only once to reduce the number of clocks.
Reachability analysis and often model checking
for timed automata like TCTL model checking are PSPACE-complete \cite{AD1,AL1,ACD1}.
Since reachability analysis and model checking the timed automata
are usually performed several times to verify different properties,
we may consider first reducing the number of clocks before doing the model checking.
With fewer clocks, the DBM representation of a zone will also be smaller in size
which may as well contribute towards reducing the computation time for model checking.
\chapter{Clock Reduction while Preserving Timed Bisimulation} \label{chap-clockred}

\section{Introducton} \label{sec-clockredIntro}
Timed automata \cite{AD1} is a formalism for modelling and analyzing real time systems. The complexity of model checking is dependent on the number of clocks of the timed automaton (TA) \cite{AD1,ACD1}.
Model checkers use the region graph or a zone graph construction
for analysing reachability and other properties in timed automata.
These graphs have sizes exponential in the number of clocks of the timed automaton.
The algorithms for model checking in turn depend on the sizes of these graphs.
Hence it is desirable to construct a timed automaton with the minimum number of clocks that preserves some property of interest
(such as timed language equivalence or timed bisimilarity).
Here we show that checking the existence of a timed automaton with fewer clocks
that is timed bisimilar to the original timed automaton is decidable.
Our method is constructive and we provide a 2-EXPTIME algorithm to construct the timed bisimilar automaton with the least possible number of clocks.
We also note that if the constructed TA has a smaller number of clocks,
then it implies that there exists an automaton with a smaller number of clocks accepting the same timed language.

\textbf{Related work:}
In \cite{DY1}, an algorithm has been provided to reduce the number of clocks of a given timed automaton.
It produces a new timed automaton that is timed bisimilar to the original one. The algorithm detects a set of \emph{active clocks} at every location and partitions these active clocks into classes such that all the clocks belonging to a class in the partition always have the same value.
However, this may not result in the minimum possible number of clocks.
The algorithm described in the paper \cite{DY1} works on the timed automaton directly rather than on its semantics.
For instance, if a constraint associated with clock $x$ implies a constraint
associated with clock $y$, and the conjunction of the two constraints appears on an edge,
then the constraint associated with clock $y$
may be eliminated.
However, the algorithm of \cite{DY1} does not capture such implication.
Also by considering
constraints on more than one outgoing edge from a location,
e.g.  $l_0 \to{a, x \le 3, \emptyset} l_1$ and $l_0 \to{a, x > 3, \emptyset} l_2$
collectively, we may sometimes eliminate the constraints
that may remove some clock. This too has not been accounted for by the algorithm of \cite{DY1}.

In a related line of work
\cite{Trip1}, it has been shown that no algorithm can decide the
following two things together.
\begin{enumerate}[(i)]
\item minimality of the number of clocks while preserving the timed language and
\item for the non-minimal case finding a timed language equivalent automaton with fewer clocks.
\end{enumerate}
Also as mentioned earlier, for a given timed automaton,
the problem of determining whether there exists another TA with fewer clocks accepting the same timed language is undecidable\index{clock reduction!preserving timed language} \cite{Fin1}.

Another result appearing in \cite{LLW1} which uses the region-graph construction is the following.
A $(C, M)$-automaton\index{(C, M)-automaton} is one with $\lvert C \lvert$ clocks and $M$ is the largest integer appearing in the timed automaton.
Given a timed automaton $A$, a set of clocks $C$ and an integer $M$, checking the existence of a ($C, M$)-automaton that is timed bisimilar to $A$ is shown to be decidable.
The method described there constructs a logical formula called the \emph{characteristic formula} and checks whether there exists a $(C,M)$-automaton that satisfies it.
Further, it is shown that a pair of automata satisfying the same characteristic formula are timed bisimilar.
The paper also considers the following problem $P$ and leaves it open.
\begin{quote}
Does there exist a timed automaton with $|C|$ clocks
that is timed bisimilar to $A$?
\end{quote}
We solve the following problem $Q$.
\begin{quote}
Given a timed automaton $A$, construct a timed bisimilar automaton $B$
with the least number of clocks.
\end{quote}
Solving the problem $Q$ also solves the problem $P$ and vice versa.
If $B$ has say $|C'|$ clocks and $|C| \ge |C'|$,
then the answer to the problem $P$ is in the affirmative.
For the problem $P$, if the answer is `yes' for $|C|$ clocks
but the answer is `no' for $|C|-1$ clocks,
then it is clear that the timed automaton
with the minimum possible number of clocks that is timed bisimilar to $A$
can have no less than $\lvert C \lvert$ clocks.

The rest of the paper is organized as follows: in Section \ref{sec-ta}, we describe timed automata and introduce several concepts that will be used in the paper. We also describe the construction of the zone graph
in Section \ref{sec-zoneGraph} used in reducing the number of clocks.
This construction of the zone graph has also been used in \cite{GKNA1}
for deciding several timed and time abstracted relations.
In Section \ref{sec-method}, we discuss our approach in detail along with a few examples. Section \ref{sec-conc} is the conclusion.

\section{Clock Reduction} \label{sec-method}\index{clock reduction}
Unlike the method described in \cite{DY1}, which works on the syntactic structure of the timed automaton, we use a semantic representation as given by the zone graph described in Section \ref{sec-ta} to capture the behaviour of the timed automaton. This helps us to reduce the number of clocks in a more effective way. For a given TA $A$, we first describe a sequence of stages
to construct a TA $A_4$ that is timed bisimilar to $A$.
Later we prove the minimality in terms of the number of clocks for the TA $A_4$.
The operations involved in our procedure use a \emph{difference bound matrix} \index{difference bound matrix (DBM)} representation of the zones.
The following  are important considerations in reducing the number of clocks.
\begin{itemize}
\item There may be some clock constraints on an edge of the TA that are never enabled.
Such edges and constraints may be removed. (Stage 1)
\item Splitting some locations may lead to a reduction in the number of the clocks. (Stage 2)
\item At some location, some clocks whose values may be expressed in terms of other clock values, may be removed. (Stage 2)
\item Two or more constraints on edges outgoing from a location
when considered collectively may lead to the removal of some constraints. (Stage 3)
\item An efficient way of renaming the clocks across all locations can reduce the total number of clocks further. (Stage 4)
\end{itemize}
Given a TA $A$, we apply the above operations in sequence to obtain the TA $A_4$. These operations are detailed below.

\input{stage1}
\input{stage2_new}
\input{stage3_new}
\input{stage4}
In the presence of an invariant condition, considering an edge $l \to{g, a, R} l'$, a zone $Z'$ of $l'$ is initially created
from zone $Z$ of $l$ such that $Z'=(Z \cap g)_{[R \leftarrow \overline{0}]}\!\!\uparrow \cap I(l')$, if $(Z \cap g)_{[R \leftarrow \overline{0}]}\!\!\uparrow \cap I(l') \neq \emptyset$, where $I(l')$ is the invariant on location $l'$.
The edge can be removed from the timed automaton if for all zones $Z$ of $l$, $(Z \cap g)_{[R \leftarrow \overline{0}]}\uparrow \cap I(l') = \emptyset$.
The invariant condition on location $l'$ can be entirely removed if for every incoming edge $l \to{g, a, R} l'$ and for each zone $Z$ of $l$, $(Z \cap g)_{[R \leftarrow \overline{0}]}\uparrow \cap I(l') = (Z \cap g)_{[R \leftarrow \overline{0}]}\uparrow$ holds.
Considering the clock relations in the zones of a location, an elementary constraint in the invariant too is removed if it is implied by the rest of the elementary constraints in the invariant.
In stage 4, for computing the set of active clocks $act(l)$ in location $l$, $clk(l)$ is defined as the
union of the set of clocks appearing in the constraints on the outgoing edges from $l$ and the set of clocks appearing in $I(l)$.

\input{proof3}

\section{Timed Automata}\label{sec-ta}
Formally, a \emph{timed automaton} (TA) \cite{AD1} is defined as a tuple $A=(L, Act, l_0, C, E)$ where
$L$ is a finite set of locations, $Act$ is a finite set
of visible actions, $l_0 \in L$ is the initial location,
$C$ is a finite set of clocks and $E \subseteq L \:\times\: \mathcal{B}(C) \:\times\: Act \:\times\: 2^C \:\times\: L$
is a finite set of \emph{edges}. The set of constraints or guards on the edges, denoted $\mathcal{B}(C)$, is
given by the grammar $g ::= \; x \bowtie k \:|\: g \wedge g$,
where $k \in \mathbb{N}$ and $x \in C$ and $\bowtie \: \in \: \{\le, <, =, >, \ge\}$.
Given two locations $l, l'$, a transition from $l$ to $l'$ is of the form
$(l, g, a, R, l')$ i.e. a transition from $l$ to $l'$  on action $a$ is possible if the constraints
specified by $g$ are satisfied; $R \subseteq C$ is a set of clocks which are reset to zero during
the transition.

The semantics of a timed automaton(TA) is described by a \emph{timed labelled
transition system} (TLTS) \cite{LAKJ1}. The timed labelled transition system $T(A)$ generated by $A$ is defined as
$T(A) = (Q, Lab, Q_0, \{\stackrel{\alpha}{\longrightarrow} | \alpha \in Lab\})$, where
$Q \: = \: \{(l,v)\:|\: l \in L, v \in {\mathbb{R}_{\ge 0}}^{|C|}\}$ is the set of  \emph{states}, each of which is of the form $(l,v)$, where $l$ is a location of the timed automaton
and $v$ is a valuation in $|C|$ dimensional real space
where each clock is mapped to a unique dimension in this space;
$Lab = Act \cup \mathbb{R}_{\ge 0}$ is the set of labels.
Let $v_0$ denote the valuation such that $v_0(x)=0$ for all
$x \in C$. $Q_0=(l_0, v_0)$ is the initial state of $T(A)$.
A transition may occur in one of the following ways:\\
(i) \emph{Delay transitions} : $(l,v)\stackrel{d}{\longrightarrow}(l, v+d)$. Here, $d \in \mathbb{R}_{\ge 0}$ and $v + d$ is the valuation
in which the value of every clock is incremented by $d$. \\
(ii) \emph{Discrete transitions} : $(l,v) \stackrel{a}{\longrightarrow} (l', v')$ if for an edge
$e=(l,g,a,R,l') \in \: E$,
$v$ satisfies the constraint $g$, (denoted $v \models g$),
$v' = v_{[R \leftarrow \overline{0}]}$,
where $v_{[R \leftarrow \overline{0}]}$ denotes that every clock in $R$ has been reset to 0, while the remaining clocks are unchanged.
From a state $(l,v)$, if $v \models g$, then there exists an $a$-transition to a state $(l', v')$;
after this,  the clocks in $R$ are reset while those in $C\backslash R$ remain unchanged.

For simplicity, we do not consider annotating
locations with clock constraints (known as \emph{invariant conditions} \cite{HNSY1}). Our results extend in a
straightforward manner to timed automata with invariant conditions.
In Section \ref{sec-method}, we provide the modifications to our method for dealing with location invariants.
We now define various concepts that will be used in the rest of the paper.
\begin{definition}
Let $A=(L, Act, l_0, E, C)$ be a timed automaton, and $T(A)$ be the TLTS corresponding to $A$.
\begin{enumerate}
\item \textbf{Timed trace}: A sequence of delays and visible actions $d_1 a_1 d_2 a_2 \dots d_n a_n$
is called a timed trace of the timed automaton $A$ iff there is a sequence of transitions
$p_0 \to{d_1} p_1 \to{a_1} p_1'\to{d_2}p_2\to{a_2}p_2'\cdots \to{d_n}p_n\to{a_n}p_n'$ in $T(A)$,
where each $p_i$ and $p_i'$, $0 \le i \le n$
are timed automata states and are of the form $(l, v)$.
$p_0$ is the initial state of the timed automaton.
Given a timed trace $tr = d_1 a_1 d_2 a_2 \dots d_n a_n$,
$untime(tr)$ denotes its projection on $Act$, i.e.
$untime(tr) = a_1 a_2 \dots a_n$.
\item \textbf{Zone and difference bound matrix (DBM)}: A zone $Z$ is a set of valuations $\{v \in \mathbb{R}_{\ge 0}^{|C|} \mid v \models \beta\}$, where
$\beta$ is of the form $\beta ::= \; x \bowtie k \:|\: x-y \bowtie k \:|\: \beta \wedge \beta$,
$k$ is an integer, $x,y \in C$ and $\bowtie \: \in \: \{\le, <, =, >, \ge\}$.
$Z \!\!\uparrow$ denotes the future of the zone $Z$.
$Z\!\!\uparrow=\{v+d \mid v \in Z, d \geq 0\}$ is the set of all valuations
reachable from $Z$ by delay transitions.
A zone is a convex set of clock valuations.
A difference bound matrix is a data structure used to
represent zones.
A DBM  \cite{BM1,Dill1} for a set $C = \{x_1, x_2, \dots, x_n\}$ of $n$ clocks is an $(n+1)$ square matrix $\Xi$
where an extra variable $x_0$ is introduced such that the value of $x_0$ is always $0$.
An element $\Xi_{ij}$ is of the form $(m_{ij}, \prec)$ where $\prec \in \{<, \le\}$
such that $x_i - x_j \prec m_{ij}$.
Thus an entry $m_{i0}$ denotes the constraint $x_i - x_0 \prec m_{i0}$ which is equivalent to $x_i \prec m_{i0}$.
\item \textbf{Canonical decomposition}: Let $\displaystyle{g =\bigwedge_{i=1}^n \gamma_i \in \mathcal{B}(C)}$, where
each $\gamma_i$ is an \emph{elementary constraint} of the form $x_i \bowtie k_i$,
such that $x_i \in C$ and $k_i$ is a non-negative integer.
A canonical decomposition of a zone $Z$ with respect to  $g$
is obtained by splitting $Z$ into a set of zones $Z_1, \dots, Z_m$ such that
for each $1 \leq j \leq m$, and $1 \leq i \leq n$, either $\forall v \in Z_j, v \models \gamma_i$ or $\forall v \in Z_j$, $v \not \models \gamma_i$.
For example, consider the zone $Z=x \geq 0 \wedge y \geq 0$ and the guard  $x \le 2 \wedge y > 1$.
$Z$ is split with respect to $ x \le 2$, and then with respect to $y > 1$, hence
into four zones :
$x \le 2 \wedge y \le 1$, $x > 2 \wedge y \le 1$, $x \le 2 \wedge y > 1$ and $x > 2 \wedge y > 1$.
An elementary constraint $x_i \bowtie k_i$ induces the hyperplane $x_i = k_i$ in a zone graph of the timed automaton.
\item \textbf{Zone graph}:
Given a timed automaton $A= (L, Act, l_0, C, E)$, a \emph{zone graph} $\mathcal{G}_A$ of $A$
is a transition system $(S, s_0, Lep, \rightarrow)$,
that is a finite representation of $T(A)$.
Here $Lep = Act \cup \{\varepsilon\}$.
$S \subseteq L \times \mathcal{Z}$ is the set of nodes of $\mathcal{G}_A$, $\mathcal{Z}$ being the set of zones.
The node $s_0 = (l_0, Z_0)$ is the initial node such that $v_0 \in Z_0$.
$(l_i, Z) \to{a} (l_j, Z')$ iff $l_i \to{g, a, R} l_j$ in $A$ and $Z' \subseteq ([Z \cap g]_{R \leftarrow \overline{0}})\!\!\uparrow$ obtained after canonical decomposition of $([Z \cap g]_{R \leftarrow \overline{0}})$.
For any $Z$ and $Z'$, $(l_i, Z) \to{\varepsilon}(l_i, Z')$ iff there exists a delay $d$ and a valuation $v$ such that $v\in Z$, $v+d \in Z'$ and $(l_i,v)\stackrel{d}{\longrightarrow}(l_i, v+d)$ is in $T(A)$. Here the zone $Z'$ is called a \emph{delay successor} of zone $Z$, while $Z$ is called the \emph{delay predecessor} of $Z'$.
Every node has an $\varepsilon$ transition to itself and the $\varepsilon$ transitions are also transitive.
The relation $\to{\varepsilon}$ is reflexive and transitive and so is the \emph{delay successor} relation.
We denote the set of delay successor zones of $Z$ by $ds(Z)$.
A zone $Z' \neq Z$  is called the \emph{immediate delay successor} of a zone $Z$ iff $Z' \in ds(Z)$ and
$\forall Z'' \in ds(Z): Z'' \neq Z$ and $Z'' \neq Z'$, $Z'' \in ds(Z')$.
We call a zone $Z$ corresponding to a location to be a \emph{base zone} if $Z$ does not have a delay predecessor
other than itself.
For both $a$ and $\varepsilon$ transitions, if $Z$ is a zone then $Z'$ is also a zone, i.e. $Z'$ is a convex set.
\item A hyperplane $x = k$ is said to \emph{bound} a zone $Z$ from above if
$\exists v \in Z \:[\:\forall d \in \mathbb{R}_{\ge 0} \:[\: (v+d)(x) \succ k \iff (v+d)(x) \notin Z \:]\:]$, where $\succ \in \{>, \ge \}$.
A zone, in general, can be bounded above by several hyperplanes.
A hyperplane $x = k$ is said to \emph{bound} a zone $Z$ fully from above if
$\forall v \in Z \:[\:\forall d \in \mathbb{R}_{\ge 0} \:[\: (v+d)(x) \succ k \iff (v+d)(x) \notin Z \:]\:]$.
Analogously, we can also say that a hyperplane $x = k$ bounds a zone from below if
$\exists v \in Z \:[\:\forall d \in \mathbb{R}_{\ge 0} \:[\: (v-d)(x) \prec k \iff (v-d)(x) \notin Z \:]\:]$, where $\prec \in \{<, \le \}$.
We can also define a hyperplane bounding a zone fully from below in a similar manner.

When not specified otherwise, in this paper, a hyperplane bounding a zone implies
that it bounds the zone from above.
A zone $Z$ is bounded above if it has an immediate delay successor zone.
\item \textbf{Pre-stability}: A zone $Z_1$ of location $l_1$ is pre-stable with respect to a transition $\eta \in Lep$
leading to a zone $Z_2$ of location $l_2$ if\\
i) for $\eta = a \in Act \: :$ $\forall v_1 \in Z_1, \exists v_2 \in Z_2$ such that $(l_1, v_1) \to{a} (l_2, v_2)$ or
$\forall v_1 \in Z_1, \not\!\exists v_2 \in Z_2$ such that $(l_1, v_1) \to{a} (l_2, v_2)$. \\
ii) for $\eta = \varepsilon \: :$ $\forall v_1 \in Z_1, \exists v_2 \in Z_2, \exists d \in \mathbb{R}_{\ge 0}$ such that $(l_1, v_1) \to{d} (l_2, v_2)$ or
$\forall v_1 \in Z_1$, $\not\!\exists v_2 \in Z_2$ such that $(l_1, v_1) \to{d} (l_2, v_2)$, for any $d \in \mathbb{R}_{\ge 0}$. \\
A zone graph is said to be pre-stable if the zones corresponding to each of the nodes
of the zone graph are pre-stable with respect to their outgoing transitions.
For example zone $Z_1$ in Figure \ref{fig-delayZone} is pre-stable with respect to zone $Z_4$,
however, $Z_1 \cup Z_2 \cup Z_3$ is not pre-stable with respect to zone $Z_4$
since $\exists v \in Z_3$ and hence $v \in Z_1 \cup Z_2 \cup Z_3$ such that for all $d \in \mathbb{R}_{\ge 0}$, $v + d \notin Z_4$.
\item \textbf{corner point and corner point traces}: A corner point\index{corner point} of a zone in the zone graph is a state
where corresponding to each clock $x$, a hyperplane of the form $x = c$ bounding a zone intersects
another hyperplane defining the zone.
Considering $M$ to be the largest constant appearing in the constraints or the location invariants of the timed automaton,
each coordinate of a corner point is of the form $n, n + \delta$ or $n - \delta$ where
$n \in \{0, 1, \dots, M\}$ and $\delta$ is a symbolic value for an infinitesimally small quantity.
In a zone $Z$, there can be two kinds of corner points, an \emph{entry}\index{corner point!entry} corner point and an \emph{exit}\index{corner point!exit} corner point. For an entry corner point, there is a unique exit corner point, that can be reached from the entry corner point by performing a delay $d$ such that any delay more than $d$ will lead to a state that is in a delay successor zone of $Z$.

For a zone, it is possible to have a pair of such entry and exit corner points that are not distinct, for example, in Figure \ref{fig-delayZone}, in zone $Z_4$, the corner point ($x = 5 + \delta, y = 7$) is both an entry as well as an exit corner point, while the corner point ($x = 5 + \delta, y = 5+\delta$) is another entry corner point for the same zone $Z_4$ and the corresponding exit corner point is $(x = 7, y = 7)$.
If the entry and the exit corner points are distinct,
then a non-infinitesimally small delay can be performed from the entry corner point
to reach the corresponding exit corner point whereas an infinitesimally small delay from the exit corner point causes it to evolve into an entry corner point of the immediate delay successor zone.
A zone that is not bounded above does not have any exit corner point.

For an entry corner point, the corresponding exit corner point is the \emph{next delay corner point} whereas for an exit corner point, the entry corner point in the immediate delay successor zone is the \emph{next delay corner point}.

Let $c_1$ and $c_2$ be constraints of the form $x \succ k_1$ and $y \prec k_2$ respectively,
where $k_1, k_2 \in \{0,1, \dots, M\}$
and $\succ \:\in \: \{>, \ge \}$ and $\prec \: \in \: \{<, \le \}$.

We consider delays of the form $d_{c_1, c_2}$ and $d_{c_2, c_1}$.
$d_{c_1, c_2}$ is the delay made from an entry corner point to reach an exit corner point
while $d_{c_2, c_1}$ is the delay
made from the exit corner point to reach an entry corner point of a delay successor zone.
For example, the delay $d_{x \ge 5, y < 7}$ denotes a delay
made from a state $p'$ with $p'(x) = 5$ and after making this delay
it evolves to state $p''$ such that $p''(y) = 7 - \delta$.
Delays of the form mentioned above are called \emph{corner point delays}\index{corner point delay} or \emph{cp-delays}.

A corner point trace\index{corner point trace} (cp-trace) is a timed trace ending at a corner point through delays and visible actions such that each delay takes the trace from one corner point to the next delay corner point.
All the delay actions present in a corner point trace are cp-delays.
\end{enumerate}
\end{definition}
We create a zone graph such that for any location $l$, the zones $Z$ and $Z'$ of
any two nodes $(l, Z)$ and $(l, Z')$ in the zone graph are disjoint and all zones of the zone graph are pre-stable.
This zone graph is constructed in two phases in time exponential in the number of the clocks.
The first phase performs a forward analysis of the timed automaton while the second phase
ensures pre-stability in the zone graph.
The forward analysis may cause a zone graph to become infinite \cite{BBFL1}. Several kinds of abstractions
have been proposed in the literature \cite{DT1,BBFL1,BBLP1} to make the zone graph finite.
We use \emph{location dependent maximal constants}
abstraction \cite{BBFL1} in our construction.
In phase 2 of the zone graph creation, the zones are further split to ensure that the resultant zone graph is pre-stable.

Some approaches for preserving convexity and implementing pre-stability have been discussed in \cite{TY1}.
As an example, consider the timed automaton in Figure \ref{fig-delayZone}. The pre-stable zones of location $l_1$ are shown in the right side of the figure.
In this paper, from now on, unless stated otherwise, \emph{zone graph} will refer to this form of the pre-stable zone graph that is described above.
An algorithmic procedure for the construction of the zone graph is given in \cite{GKNA1}.

A relation $\mathcal{R} \subseteq Q \times Q$ is a \emph{timed simulation} relation if the following conditions hold
for any two timed states $(p, q) \in \mathcal{R}$. \\
$\forall$ $a \in Act$, $p \to{a}p' \implies \exists q'$ : $q \to{a}q'$ and $(p', q') \in \mathcal{R}$ and \\
$\forall d \in \mathbb{R}_{\ge 0}$, $p \to{d}p' \implies \exists q'$ : $q \to{d}q'$ and $(p', q') \in \mathcal{R}$.\\
A \emph{timed bisimulation} relation is a symmetric timed simulation. Two timed automata are timed bisimilar if and only if their initial states are timed bisimilar. Using a product construction on region graphs,
timed bisimilarity for timed automata was shown to be decidable in EXPTIME \cite{Cer1}.
In \cite{WL1}, a product construction on zone graph was used
for deciding timed bisimulation relation for timed automata.
 \begin{figure}[t]
 \centering
 \begin{minipage}[b]{0.46\linewidth}
 \includegraphics[width=0.9\textwidth]{delayZone.pdf}
 \caption{\label{fig-delayZone} a timed automaton and the zones for location $l_1$}
 \end{minipage}
 \begin{minipage}[b]{0.49\linewidth}
 \includegraphics[width=0.95\textwidth]{zg.pdf}
 \caption{\label{fig-zoneGraph} zone graph for the TA in Figure \ref{fig-delayZone}}
 \end{minipage}
 \end{figure}

\section{Zone Graph}\label{sec-zoneGraph}\index{zone graph}
\begin{algorithm}
\caption{Construction of Zone Graph\newline
\emph{Input: Timed automaton $A$ \newline
Output: Zone  graph corresponding to $A$}
}\label{algo-zonegraph}
\begin{algorithmic}[1]
{\small
\State Calculate $max_x^l$ for each location $l \in L$ and each clock $x \in C$. \Comment{\textsf{\tiny This is required
for abstraction to ensure finite number of zones in the zone graph.}}
\State Initialize $Q$ to an empty queue.
\State $Enqueue(Q, \langle l_0, \emptyset \rangle)$. \Comment{\textsf{\tiny Every element is a pair consisting of a location and its parent}}
\State $successors\_added = false$. \Comment{\textsf{\tiny flag set to true whenever  successors of a location are added to $Q$}}
\While {Q not empty}
\State $\langle l, l_p \rangle = dequeue(Q)$
\If {$l_p\neq \emptyset$},
\State For the edge $l_p \to{g, a, R'} l$ in $A$, for each existing zone $Z_{l_p}$ of $l_p$, create
the zone $Z = (Z_{l_p} \uparrow \cap \:g_{[R' \leftarrow \overline{0}]})$ of $l$, when $Z \neq \emptyset$.
\State Abstract each of the newly created zones if necessary and for any newly created zone $Z$, for location $l$,
if $\exists Z_1$ of same location such that $Z \cap Z_1 \neq \emptyset$, then merge $Z$ and $Z_1$.
\State Update edges from zones of $l_p$ to zones of $l$ appropriately.
\State If {a new zone of $l$ is added or an existing zone of $l$ is modified} for all successors $l_j$ of $l$,
enqueue $\langle l_j, l \rangle$ to Q.
\State $successors\_added := true$.
\EndIf
\State $new\_zone\_l := true$. \Comment{\textsf{\tiny flag set to false when the canonical decomposition does not produce further zones}}
\While {$new\_zone\_l$}
\State Split the existing zones $Z$ of $l$ based on the canonical decomposition of the guards on the
outgoing edges of $l$ \Comment{\textsf{\tiny It is not always necessary for a split to happen.}}
\State For every zone $Z$ of $l$, consider $Z \uparrow$ and split it further based on the canonical
decomposition of the guards on the outgoing edges of $l$
\Comment{\textsf{\tiny Note that the zones created from this split are convex.}}
\State Abstract each of the newly created zones if necessary and update edges appropriately.
\If {new zones are not created} set $new\_zone\_l$ to $false$.
\EndIf
\EndWhile
\If {any new zones of $l$ are created or any existing zones of $l$ are modified due to the canonical
decomposition of the outgoing edges of $l$ and $successors\_added = false$}
\State for all the successor locations $l_j$ of $l$ to Q, enqueue $\langle l_j, l \rangle$ to Q.
\EndIf
\EndWhile
\State  /* \textsf{Phase 2 :  In this phase, pre-stability is enforced} */
\State $new\_zone = true$
\While {$new\_zone$}
\State $new\_zone  = false$
\ForAll {edges $l_i \to{g, a, R'} l_j$}
\ForAll {pairs of zones $Z_{lik}$, $Z_{ljm}$ such that $Z_{lik} \to{\alpha}Z_{ljm}$ is an edge in the
zone graph where $\alpha \in Act \cup \mathbb{R}_{\ge 0}$}
\If {$Z_{lik}$ is not pre-stable with respect to $Z_{ljm}$,}
Split $Z_{lik}$ to make it pre-stable with respect to $Z_{ljm}$.   \Comment{\textsf{\tiny Note that this split still maintains
convexity of $Z_{lik}$ since the zone is split entirely along an axis that is parallel to the diagonal in
the $\lvert C\lvert$-dimensional space.}}
\State $new\_zone := true$
\State Update the incoming and the outgoing edges for the newly created zones.
\EndIf
\EndFor
\EndFor
\EndWhile
}
\end{algorithmic}
\end{algorithm}

The detailed algorithm for creating the zone graph has been described in Algorithm \ref{algo-zonegraph}
and consists of two phases, the first one being a forward analysis of the timed automaton while the second phase
ensures pre-stability in the zone graph.
The presence of an edge to a location $l'$ in the timed automaton does not
imply that the the action of the edge can be performed or
location $l'$ is reachable.
For an edge $l \to{g, a, R} l'$, if for all zones $Z$ of $l$,
$Z \cap g = \emptyset$, then action $a$ cannot be performed at $l$
and $l'$ cannot be reached at $l$.
Such edges can be removed from the timed automaton.
For a location $l'$, if the actions corresponding to none of
the incoming edges can be performed, then the location $l'$
is unreachable and it can be removed along with its incoming and
outgoing transitions.

The set of valuations for every
location is initially split into zones based on the canonical decomposition of its outgoing transition\delete{as discussed in
\cite{TY1}}.

After phase 2, pre-stability ensures the following: For a node $(l, z)$ in the zone graph, with $v \in z$,
\delete{such that for a $v \models z$,}for a timed trace $tr$,
if $(l, v) \to{tr} (l'', v'')$, with $v'' \in z''$, then
$\forall v' \; \delete{\models}\mbox{in}\; z$, $\exists tr'.(l, v') \to{tr'} (l'', \tilde{v})$, with
$untime(tr') = untime(tr)$ and $\tilde{v} \in z''$.
\delete{Here $untime(tr)$ represents the sequence of visible actions in $tr$.}According to the construction given
in algorithm \ref{algo-zonegraph}, for a particular location of the timed automaton, the zones corresponding
to any two nodes are disjoint. Convexity of the zones and pre-stability property \delete{of the zones}together ensure
that a zone with elapse of time is intercepted by a single hyperplane of the form $x = h$,
where $x \in C$ and $h \in \mathbb{N}$\delete{in the future}.
This is stated formally in Lemma \ref{lem-pre-stability}.
Some approaches for preserving convexity and implementing pre-stability have been discussed \delete{in detail}in \cite{TY1}.
As an example consider the timed automaton in Figure \ref{fig-delayZone}. The zones corresponding to location $l_1$ as produced through algorithm \ref{algo-zonegraph} are shown in the right side of the figure.

For a state $q \in T(A)$, $\mathcal{N}(q)$ represents the node
 of the zone \delete{valuation}graph with the same location as that of $q$ such that the zone corresponding to $\mathcal{N}(q)$ includes
 the valuation of $q$. We often say that a state $q = (l, v)$ is in node $s = (l, Z)$ to indicate that $v \in Z$. For two zone \delete{valuation}graphs,
 $\mathcal{G}_{A_1} = (S_1,s_0,Lep,\rightarrow_1)$, $\mathcal{G}_{A_2} = (S_2,s_0',Lep,\rightarrow_2)$ and a relation
 $\mathcal{R}\subseteq S_1\times S_2$, $\mathcal{G}_{A_1}\:\mathcal{R}\:\mathcal{G}_{A_2}$ iff $(s_p,s_q)\in \mathcal{R}$.
An $\varepsilon$ transition represents a delay $d \in \mathbb{R}_{\ge 0}$.
\subsection{Creating pre-stable zones}

We describe below a procedure for splitting a zone $Z_1$ into a set $\mathcal{Z}_1$ of zones w.r.t a zone $Z_2$
such that the following properties are satisfied.
\begin{enumerate}
\item Each zone $Z_i$ in $\mathcal{Z}_1$ should either be a subset of $Z_2$ or disjoint from $Z_2$.
\item The DBMs in $\mathcal{Z}_1$ are mutually exclusive and disjoint.
\item The union of the DBMs in $\mathcal{Z}_1$ is $Z_1$.
\end{enumerate}

A naive approach would be to populate $\mathcal{Z}_1$ with two DBMs,
namely the (set) intersection and (set) difference of $Z_1$ with $Z_2$.
This ensures that the properties 1, 2 and 3 hold.
However, while the intersection of two DBMs is always a DBM, this cannot be said for their difference.
Thus, we seek to refine our process by replacing the difference of $Z_1$ and $Z_2$
with some DBMs which are mutually exclusive and which cover the difference of $Z_1$ and $Z_2$,
so that properties 2 and 3 continue to hold.

One naive way to achieve this is to initially populate $\mathcal{Z}_1$ with only $Z_1$,
and then iterate through all the elementary clock constraints $c_{ij}$ of $Z_2$.
The element $m_{ij}$ in the DBM is a constraint of the form $x_i - x_j < k$ or $x_i - x_j \le k$.
In each iteration, we replace each DBM in $\mathcal{Z}_1$ with two DBMs, namely its intersection and difference with $c_{ij}$.
Note that the complement of an elementary clock constraint is always an elementary clock constraint
and thus the difference can be calculated by simply taking the intersection
with the complement of the elementary constraint, which implies that
both the intersection and the difference considering a single elementary constraint as above are DBMs.
Of course, some of the DBMs thus produced may be empty, and we do not add these to $\mathcal{Z}_1$.
The theoretical bound on the number of DBMs in $\mathcal{Z}_1$ at the end is
$2^{\mbox{(number of elementary clock constraints in dbm2)}}$,
that is, $2^{O((n+1)^2)}$.
The term $(n+1)^2$ comes from the fact that there are $(n+1)^2$ entries in $Z_2$
each of which is an elementary clock constraint.
The above exponential bound is because of the fact that, in the worst case,
we double the number of DBMs in $\mathcal{Z}_1$ in each iteration.
This leaves us with a set of mutually exclusive zones which satisfies properties 1, 2 and 3.

However, this is inefficient in that it generates too many zones.
So, we refine the process of iterating through the elementary clock constraints $c_{ij}$.
We start with an empty $\mathcal{Z}_1$ and a pointer which initially points to $Z_1$.
At each iteration, we evaluate the intersection of $c_{ij}$ and the zone pointed to by the pointer
and also their difference which too is a zone since the difference is calculated
by taking the intersection with the complement of $c_{ij}$.
The pointer is then updated to point to the intersection DBM,
while the difference DBM is added to $\mathcal{Z}_1$.
This procedure continues until all the $c_{ij}$ have been seen.
In the end, if the pointer's target is non-empty, we add it to $\mathcal{Z}_1$.
At this point, the zones in $\mathcal{Z}_1$ still cover $Z_1$ and also satisfy properties 1 and 2.
This method leads to a maximum of one zone being added to $\mathcal{Z}_1$ for each $c_{ij}$,
and one zone being added at the end, which gives us a maximum of $O((n + 1)^2 + 1)$ zones in $\mathcal{Z}_1$.

Algorithm \ref{algo-splitZone} describes the procedure mentioned above formally.

\begin{algorithm}
\caption{This algorithm returns the intersection and the difference of two zones.\newline
\emph{Input: Two zones $Z_1$ and $Z_2$ \newline
Output: i) $Z$ which is $Z_1 \cap Z_2$ and ii) $\mathcal{Z}_1$, a set of zones whose union is $Z_1 - Z_2$}
}\label{algo-splitZone}
\begin{algorithmic}[1]
{\small
\State $Z := Z_1 \cap Z_2$
\State $\mathcal{Z}_1 := \emptyset$, $p := Z_1$ \Comment{\textsf{\tiny $p$ is a pointer pointing to a zone.}}
\ForAll {$i \in [0, n]$}
\ForAll {$j \in [0, n]$}
\State Let $c_{ij}$ be the constraint $x_i - x_j \prec m_{ij}$, where $(m_{ij}, \prec)$ is the corresponding entry in the DBM
of $Z_2$.
\State Let $\overline{c_{ij}}$ be the complement of the constraint $c_{ij}$.
\State $Z_{11} := p \cap c_{ij}$, $Z_{12} := p \cap \overline{c_{ij}}$
\If {$Z_{12} \neq \emptyset$}
\State $\mathcal{Z}_1 := \mathcal{Z}_1 \cup Z_{12}$
\EndIf
\If {$Z_{11} \neq \emptyset$}
\State $p := Z_{11}$
\EndIf
\EndFor
\EndFor
\If {$Z_{11} \neq \emptyset$}
\State $\mathcal{Z}_1 = \mathcal{Z}_1 \cup Z_{11}$
\EndIf
}
\end{algorithmic}
\end{algorithm}

Pre-stability is ensured corresponding to both \emph{delay} and \emph{discrete} transitions in the zone graph.
\begin{itemize}
\item {\bf Pre-stabilizing zones corresponding to delay transitions:}
Consider two nodes in the zone graph $(l, Z) \to{\varepsilon} (l, Z')$ such that $Z$ is not pre-stable w.r.t $Z'$.
The procedure involves splitting $Z$ w.r.t $Z'\!\!\downarrow$ such that each newly created zones
obtained by splitting $Z$ is pre-stable w.r.t. $Z'$.

We follow the procedure given in Algorithm \ref{algo-splitZone} to split a zone $Z$ w.r.t the zone $Z'$.
\item {\bf Pre-stabilizing zones corresponding to discrete transitions  in the zone graph:}
Consider two nodes in the zone graph $(l, Z) \to{a} (l', Z')$ such that $a \in Act$  and $Z$ is not pre-stable w.r.t $Z'$.
The method involving the pre-stabilization operation is similar to the $Pre_e(Z')$ operator used in backward analysis
of the timed automaton, where $e$ is an edge between say location $l$ and $l'$ of the timed automaton and $Z'$ is a zone of $l'$.
Formally for an edge $l \to{g, a, R} l'$ $Pre_e(Z')$ is defined as the following:
\begin{center}
$Pre_e(Z') = \{v \in \mathbb{R}_{\ge 0}^{\lvert C \lvert}\ \lvert \exists v' \in Z', \: \exists t \in \mathbb{R}_{\ge 0} \mbox{ such that }
v + t \models g \mbox{ and }$ \\
$[R\leftarrow 0] (v + t) = v'\}$
\end{center}
Let $\sembrack{R = 0}$ denotes the set of all clock valuations $v$
such that $v(x) = 0$ when $x \in R$ and $v(x) \in \mathbb{R}_{\ge 0}$ when $x \notin R$. $v(x)$ denotes the valuation of clock $x$.
In our case, we define the zone $\widehat{Z'}$ as follows:\\
$\widehat{Z'} = [R \leftarrow 0]^{-1}(Z' \cap \sembrack{R = 0})$.
Now $Z$ is split w.r.t $\widehat{Z'}$ following Algorithm \ref{algo-splitZone} that
makes all the newly created zones pre-stable w.r.t $Z'$.
\end{itemize}

The following lemma states an important property of the zone graph which will further be used for clock reduction.
\begin{lemma} \label{lem-pre-stability}
Pre-stability ensures that
if the zone $Z$ in any node $(l, Z)$ in the zone graph
is bounded above, then it is bounded fully from above
by a hyperplane $x = h$, where $x \in C$ and $h \in \mathbb{N}$.
\end{lemma}
\begin{proof}
We show the proof for two clocks. The same argument holds for arbitrary number of clocks.
Consider a $\lvert C \lvert$ dimensional zone $Z$.
Since $\lvert C \lvert = 2$, a zone $Z$ that is bounded above is of the form $k_{xy_2} \prec x - y \prec k_{xy_1}$,
$k_{x_1} \prec x \prec k_{x_2}$ and $k_{y_1} \prec y \prec k_{y_2}$,
where $x$ and $y$ are the two clocks, $\prec \in \{<, \le\}$.
Let $Z$ be bounded above by two hyperplanes $x = k_{x_2}$ and  $y = k_{y_2}$.
As can be seen from Figure \ref{fig-pre-proof}, there are two zones $Z_1$ and $Z_2$ such that
$Z_1$ is defined by the inequations $x-y \prec k_{xy_1}$, $k_{x_2} \prec x$ and $y \prec k_{y_2}$,
while $Z_2$ is defined as $k_{xy_2} \prec x-y$, $x \prec k_{x_2}$ and $k_{y_2} \prec y$.
Making $Z$ pre-stable will divide it into two parts:
\begin{itemize}
\item $Z'$ defined by the inequations
$k_{x_2} -  k_{y_2} \prec x - y \prec k_{xy_1}$, $k_{x_1} \prec x \prec k_{x_2}$ and $k_{y_1} \prec y$ and
\item $\widehat{Z}$ define by the inequations
$k_{xy_2} \prec x - y \prec k_{x_2} -  k_{y_2}$, $k_{x_1} \prec x$ and $y \prec k_{y_2}$
\end{itemize}
$Z'$ is bounded fully from above by the hyperplane $k_{x_2}$ while $\widehat{Z}$ is bounded fully from above by the hyperplane $k_{y_2}$.
Note that the inequations defining the zones $Z'$ and $\widehat{Z}$ may vary depending on
the relation among the various constants $k_{xy_1}, k_{xy_2}, k_{x_1}, k_{x_2}, k_{y_1}$ and $k_{y_2}$.
In all cases, pre-stability will ensure that there exists a single hyperplane that fully bounds a zone from above.
\end{proof}
\begin{figure}
\centering
\vspace{-10pt}
\includegraphics[width=0.3\textwidth]{pre-stability.pdf}
\caption{\label{fig-pre-proof}For every zone, there exists a hyperplane that bounds it fully from above in a pre-stable zone graph}
\vspace{-10pt}
\end{figure}


\section{Introduction}
To study the manifolds with negative curvature Bishop and O'Neill introduced the idea of warped products [2]. Afterwards, this idea was used to model the standard space time, especially in the neighborhood of massive stars and black holes [11]. However, warped product spaces were developed in Riemannian manifolds enormously [8, 9]. The geometry of warped product submanifolds is intensified after B.Y. Chen's work on warped product CR-submanifolds of Kaehler manifolds [8]. Motivated by Chen's papers, many geometers studied warped product submanifolds in almost Hermitian as well as almost contact metric manifolds [1, 9, 13, 15].

\parindent=8mm
On the other hand, almost contact manifolds with Killing structures tensors were defined in [3] as nearly cosymplectic manifolds.  Later on, Blair and Showers [5] studied nearly cosymplectic structure $(\phi, \xi, \eta, g)$ on a Riemannian manifold  $\bar M$ with $\eta$ closed from the topological viewpoint. An almost contact metric structure $(\phi, \xi, \eta, g)$ satisfying  $(\bar\nabla_X\phi)X=0$ is called a nearly cosymplectic structure. If we consider $S^5$ as a totally geodesic hypersurface of $S^6$; then it is known that $S^5$ has a non cosymplectic nearly cosymplectic structure. It was shown that the normal nearly cosymplectic manifolds are cosymplectic (see [4]).

\parindent=8mm Next, the slant submanifolds of an almost contact metric manifold were defined and studied by J.L. Cabrerizo et.al [7]. The notion of semi-slant submanifolds of almost Hermitian manifolds was introduced by N. Papaghuic [12]. In fact, semi-slant submanifolds in almost Hermitian manifolds are defined on the line of CR-submanifolds. In the setting of almost contact metric manifolds, semi-slant submanifolds are defined and investigated by Cabrerizo et. al [6].

\parindent=8mm Recently, we studied warped product semi-slant submanifolds of cosymplectic manifolds [10]. We have seen that there do not exist warped product semi-slant submanifolds in cosymplectic manifolds. As the nearly cosymplectic structure is generalizes the cosymplectic ones in the same sense as nearly Kaehler generalizes Kaehler structure. Therefore, the geometric study of warped product semi-slant submanifolds in nearly cosymplectic is interesting. Such type of warped products exist in nearly cosymplectic case while in case of cosymplectic do not exist. In this paper, we prove that the warped product semi-slant submanifolds of the type $N_\theta\times{_{f}N_T}$ do not exist in a nearly cosymplectic manifold $\bar M$. However, we obtain some interesting results on the existence of the warped product submanifolds of the type $N_T\times{_{f}N_\theta}$ of a nearly cosymplectic manifold $\bar M$, where $N_T$ and $N_\theta$ are the invariant and proper slant submanifolds of $\bar M$, respectively. We also establish a general sharp inequality for squared norm of the second fundamental form in terms of the warping function and the slant angle for the warped product semi-slant submanifolds in the form $N_T\times{_{f}N_\theta}$ in a nearly cosymplectic manifold $\bar M$.

\section{Preliminaries}
A $(2n+1)-$dimensional $C^\infty$ manifold $\bar M$ is said to have an {\it{almost contact structure}} if there exist on $\bar M$ a tensor field $\phi$ of type $(1,1)$, a vector field $\xi$ and a $1-$form $\eta$ satisfying [5]
$$\phi^2=-I+\eta\otimes\xi,~~\phi\xi=0,~~\eta\circ\phi=0,~~\eta(\xi) = 1.\eqno(2.1)$$
There always exists a Riemannian metric $g$ on an almost contact manifold $\bar M$ satisfying the following compatibility condition
$$\eta(X)=g(X, \xi),~~~g(\phi X, \phi Y)=g(X, Y)-\eta(X)\eta(Y)\eqno(2.2)$$
where $X$ and $Y$ are vector fields on $\bar M$ [5].

\parindent=8mm
An almost contact structure $(\phi, \xi, \eta)$ is said to be {\it normal} if the almost complex structure $J$ on the product manifold $\bar M\times\mathbb R$ given by
$$J(X, f\frac{d}{dt})=(\phi X-f\xi,~\eta(X)\frac{d}{dt}),$$
where $f$ is a $C^\infty-$function on $\bar M\times\mathbb R$ has no torsion i.e., $J$ is integrable, the condition
for normality in terms of $\phi,~\xi$ and $\eta$ is $[\phi, \phi]+2d\eta\otimes\xi=0$ on $\bar M$, where $[\phi, \phi]$ is the Nijenhuis tensor of $\phi$. Finally the {\it{fundamental 2-form}} $\Phi$ is defined by $\Phi(X, Y)=g(X, \phi Y)$.

\parindent=8mm
An almost contact metric structure $(\phi, \xi, \eta, g)$ is said to be {\it{cosymplectic}}, if
it is normal and both $\Phi$ and $\eta$ are closed [5]. The structure is said to be {\it{nearly cosymplectic}} if $\phi$ is Killing, i.e., if
$$(\bar\nabla_X\phi)Y+(\bar\nabla_Y\phi)X=0,\eqno(2.3)$$
for any $X,~Y\in T\bar M$, where $T\bar M$ is the tangent bundle of $\bar M$ and $\bar\nabla$ denotes the Riemannian connection of the metric $g$. Equation (2.3) is equivalent to $(\bar\nabla_X\phi)X=0$, for each $X\in T\bar M$. The structure is said to be {\it{closely cosymplectic}} if $\phi$ is Killing and $\eta$ is closed. It is well known that an almost contact metric manifold is {\it{cosymplectic}} if and only if $\bar\nabla\phi$ vanishes identically, i.e., $(\bar\nabla_X\phi)Y=0$ and $\bar\nabla_X\xi=0$.\\

\noindent{\bf{Proposition 2.1 [5].}} {\it{On a nearly cosymplectic manifold the vector field $\xi$ is
Killing.}}

\parindent=8mm
From the above proposition, we have $g(\bar\nabla_X\xi, X)=0$, for any vector field $X$ tangent to $\bar M$, where $\bar M$ is a nearly cosymplectic manifold.

\parindent=8mm
Let $M$ be submanifold of an almost contact metric manifold $\bar M$ with induced metric $g$ and let $\nabla$ and $\nabla^{\perp}$ be the induced connections on the tangent bundle $TM$ and the normal bundle $T^{\perp}M$ of $M$, respectively. Denote by ${\cal{F}}(M)$ the algebra of smooth functions on $M$ and by $\Gamma(TM)$ the ${\cal{F}}(M)$-module of smooth sections of $TM$ over $M$. Then the Gauss and Weingarten formulas are given by
$$\bar\nabla_XY=\nabla_XY+h(X, Y)\eqno(2.4)$$
$$\bar\nabla_XN=-A_NX+\nabla^{\perp}_XN,\eqno(2.5)$$
for each $X,~Y\in\Gamma(TM)$ and $N\in \Gamma(T^\perp M)$, where $h$ and $A_N$ are the second
fundamental form and the shape operator (corresponding to the normal vector
field $N$) respectively for the immersion of $M$ into $\bar M$. They are related as
$$g(h(X, Y), N)=g(A_NX, Y),\eqno(2.6)$$
where $g$ denotes the Riemannian metric on $\bar M$ as well as induced on $M$.

\parindent=8mm
For any $p\in M$, let $\{ e_1,\cdots, e_m,\cdots, e_{2n+1}\}$ be an orthonormal frame for the tangent space $T_p\bar M$, such that $e_1,\cdots, e_m$ are tangent to $M$ at $p$. We denote by $H$ the mean curvature vector, that is
$$H(p)=\frac {1}{m}\sum_{i=1}^mh(e_i, e_i).$$
Also, we set
$$h_{ij}^r = g( h(e_i, e_j), e_r), ~~~ i, j \in \{1,\cdots, m\}, ~r\in \{m+1,\cdots, 2n+1\}$$
and
$$||h||^2 = \sum _{i, j = 1}^m g( h(e_i, e_j), h(e_i, e_j)).$$

\parindent=8mm
For any $X\in \Gamma(TM)$, we write
$$\phi X=PX+FX,\eqno(2.7)$$
where $PX$ is the tangential component and $FX$ is the normal component of $\phi X$.

\parindent=8mm
Similarly for any $N\in \Gamma(T^\perp M)$, we write
$$\phi N=BN+CN,\eqno(2.8)$$
where $BN$ is the tangential component and $CN$ is the normal component of $\phi N$. The covariant derivative of the tensor field $\phi$ is defined as
$$(\bar\nabla_X\phi)Y=\bar\nabla_X\phi Y-\phi\bar\nabla_XY,\eqno(2.9)$$
for all $X,~Y\in \Gamma(T\bar M)$.

\parindent=8mm
Now, denote by ${\cal P}_XY$ and ${\cal Q}_XY$ the tangential and normal
parts of $(\bar\nabla_X\phi)Y$, i.e.,
$$(\bar\nabla_X\phi)Y={\cal P}_XY+{\cal Q}_XY\eqno(2.10)$$
for all $X, Y\in \Gamma(TM)$. Making use of (2.7)-(2.10) and the Gauss and Weingarten formulae, the following
equations may easily be obtained
$${\cal P}_XY=(\bar\nabla_XP)Y-A_{FY}X-Bh(X,Y)\eqno(2.11)$$
$${\cal Q}_XY=(\bar\nabla_XF)Y+h(X, PY)-Ch(X, Y)\eqno(2.12)$$
where the covariant derivative of $P$ and $F$ are defined by
$$(\bar\nabla_X P)Y=\nabla_XPY-P\nabla_XY\eqno(2.13)$$
$$(\bar\nabla_XF)Y=\nabla^\perp_XFY-F\nabla_XY\eqno(2.14)$$
for all $X,Y\in \Gamma(TM)$.

\parindent=8mm
Similarly, for any $X\in\Gamma(TM)$ and $N\in \Gamma(T^\perp M)$, denoting the tangential and
normal parts of $(\bar\nabla_X\phi)N$ by ${\cal P}_XN$ and ${\cal Q}_XN$
respectively, we obtain
$${\cal P}_XN=(\bar\nabla_XB)N+PA_NX-A_{CN}X \eqno(2.15)$$
$${\cal Q}_XN=(\bar\nabla_XC)N+h(BN, X)+FA_NX \eqno(2.16)$$
where the covariant derivative of $B$ and $C$ are defined
by
$$(\bar\nabla_XB)N=\nabla_XBN-B\nabla^\perp_XN\eqno(2.17)$$
$$(\bar\nabla_XC)N=\nabla^\perp_XCN-C\nabla^\perp_XN.\eqno(2.18)$$

\parindent=8mm
It is straightforward to verify  the following properties of ${\cal P}$ and ${\cal Q}$, which we enlist here for later use:
$$(p_1)~~(i)~~~~{\cal P}_{X + Y}W={\cal P}_XW + {\cal P}_YW,
~~~~~~(ii)~~~{\cal Q}_{X + Y}W={\cal Q}_XW + {\cal Q}_YW,~~~$$
$$~(p_2)~~~(i)~~~~{\cal P}_X(Y+W)={\cal P}_XY + {\cal P}_XW,~~~~(ii)~~~~{\cal Q}_X(Y + W)={\cal Q}_XY + {\cal Q}_XW,$$
$$(p_3)~~(i)~~~~g({\cal P}_XY,~ W)=- g(Y, {\cal P}_XW),~~~(ii)~~~~g({\cal Q}_XY,~ N)=-g(Y, {\cal P}_XN),$$
for all $X, Y, W\in\Gamma(TM)$ and $N\in\Gamma(T^\perp M)$.

\parindent=8mm
On a submanifold $M$ of a nearly cosymplectic  manifold, by equations (2.3) and (2.10) we have
$$(a)~~{\cal P}_XY + {\cal P}_YX=0,~~(b)~~{\cal Q}_XY + {\cal Q}_YX= 0\eqno(2.19)$$
for any $X,~Y\in\Gamma(TM)$.

\parindent=8mm
A submanifold $M$ of an almost contact metric manifold $\bar M$ is said to be {\it{invariant}} if $F$ is identically zero, that is, $\phi X\in \Gamma(TM)$ and {\it{anti-invariant}} if $P$ is identically zero, that is, $\phi X\in \Gamma(T^\perp M)$, for any $X\in \Gamma(TM)$.

\parindent=8mm
We shall always consider $\xi$ to be tangent to the submanifold $M$. There is another class of submanifolds that is called the slant submanifold. For each non zero vector $X$ tangent to $M$ at $x$, such that $X$ is not proportional to $\xi_x$, we denote by $0\leq\theta (X) \leq\pi /2$, the angle between $\phi X$ and $T_x M$ is called the slant angle. If the slant angle $\theta (X)$ is constant for all $X\in T_x M - \langle\xi_x\rangle$ and $\ x\in M$, then M is said to be {\it{slant submanifold}} [7]. Obviously if $\theta = 0,~M$ is invariant and if $\theta$ = $\pi /2,~M$ is an anti-invariant submanifold. A slant submanifold is said to be {\it{proper slant}} if it is neither invariant nor anti-invariant submanifold.

\parindent=8mm
We recall the following result for a slant submanifold.\\

\noindent
{\bf {Theorem 2.1 [7]}} {\it Let M be a submanifold of an almost contact metric manifold $\bar M$, such that $\xi$ is tangent to $M$. Then $M$ is slant if and only if there exists a constant  $\lambda\in [0, 1]$ such that
$$P^2=\lambda (- I + \eta \otimes \xi)\eqno (2.20)$$
Furthermore, if $\theta$ is slant angle of $M$, then $\lambda = \cos ^2\theta.$}

\parindent=8mm
Following relations are straightforward consequences of relation (2.20)
$$g(PX,PY)=\cos^2\theta ( g(X,Y)- \eta (Y)\eta (X))\eqno(2.21)$$
$$g(FX,FY)=\sin^2\theta ( g(X,Y)-\eta (Y)\eta (X))\eqno(2.22)$$
for all $X,Y\in \Gamma(TM)$.

\parindent=8mm
A submanifold $M$ of an almost contact metric manifold $\bar M$ is said to be a {\it semi-slant} if there exist two orthogonal distributions
$D_1$ and $D_2$ satisfying:
\begin{enumerate}
\item[{(i)}] $TM={\cal D}_1\oplus {\cal D}_2\oplus \langle\xi\rangle$
\item[{(ii)}] ${\cal D}_1$ is an invariant  i.e., $\phi {\cal D}_1\subseteq {\cal D}_1$
\item[{(iii)}] ${\cal D}_2$ is a slant distribution with slant angle $\theta\neq \frac{\pi}{2}.$
\end{enumerate}

\parindent=8mm
A semi-slant submanifold $M$ of an almost contact  manifold $\bar M$ is {\it mixed geodesic} if
$$h(X, Z)=0\eqno(2.23)$$
for any $X\in {\cal D}_1$ and $Z\in {\cal D}_2$. Moreover, if $\mu$ is the $\phi-$invariant subspace of the normal bundle $T^\perp M$, then in case of semi-slant submanifold, the normal bundle $T^\perp M$ can be decomposed as
$$T^\perp M= F{\cal D}_2\oplus \mu.\eqno(2.24)$$

\section{Warped product semi-slant submanifolds}
Bishop and O'Neill [2] introduced the notion of warped product manifolds. These manifolds are natural generalizations of Riemannian product manifolds. They defined these manifolds as: Let $(N_1, g_1)$ and $(N_2, g_2)$ be two Riemannian manifolds and $f>0$ a differentiable function on $N_1$. Consider the product manifold $N_1\times N_2$ with its projections $\pi_1:N_1\times N_2\to N_1$ and $\pi_2:N_1\times N_2\to N_2$. Then the warped product of $N_1$ and $N_2$ denoted by $M=N_1\times{_f}N_2$ is a Riemannian manifold $N_1\times N_2$ equipped with the Riemannian structure such that
$$g(X, Y)=g_1({\pi_1}_\star X, {\pi_1}_\star Y)+(f\circ\pi)^2g_2({\pi_2}_\star X, {\pi_2}_\star Y)$$
for each $X, Y\in\Gamma(TM)$ and $\star$ is a symbol for the tangent map. Thus we have
$$g=g_1+f^2g_2.\eqno(3.1)$$
The function $f$ is called the {\it{warping function}} of the warped product [2, 11]. A warped product manifold $N_1\times{_f}N_2$ is said to be {\it{trivial}} if the warping function $f$ is constant. We recall the following general result for a warped product manifold for later use.\\

\noindent
{\bf {Lemma 3.1 [2]}} {\it Let $M = N_1\times_{f}N_2$ be a warped product manifold. Then}
\begin{enumerate}
\item[{(i)}] {\it$\nabla_XY\in TN_1$ is the lift of $\nabla_XY$ on $N_1$}
\item[{(ii)}] {\it $\nabla_XZ=\nabla_ZX=(X\ln f)Z$}
\item[{(iii)}] {\it $\nabla_ZW = \nabla^{N_2}_ZW - g(Z, W) \nabla\ln f$}
\end{enumerate}
{\it{for each X, Y$\in\Gamma(TN_1)$ and $Z, W\in \Gamma(TN_2)$ where $\nabla$ and $\nabla ^{N_2}$ denote the Levi-Civita connections on $M$ and $N_2$, respectively, and $\nabla\ln f$ is the gradient of $\ln f$.}}

\parindent=8mm Let $M$ be a Riemannian manifold of dimension $k$ with the inner product $g$ and $\{e_1,\cdots,e_k\}$ be an orthonormal frame on $M$. Then for a differentiable function $f$ on $M$, the gradient $\nabla f$ of a function $f$ on $M$ is defined by
$$g(\nabla f, U)=Uf,\eqno(3.2)$$
for any $U\in \Gamma(TM)$. As a consequence, we have
$$\|\nabla f\|^2=\sum\limits_{i=1}\limits^k (e_i(f))^2\eqno(3.3)$$
where $\nabla f$ is the gradient of the function $f$ on $M$.

\parindent=8mm
Now, we consider the warped product semi-slant submanifolds tangent to the structure vector field $\xi$ which are either in the form $M=N_T\times{_{f}N_\theta}$ or $M=N_\theta\times{_{f}N_T}$ in a nearly cosymplectic  manifold $\bar M$, where $N_T$ and $N_\theta$ are invariant and proper slant submanifolds of a nearly cosymplectic manifold $\bar M$, respectively. On a warped product submanifold $M=N_1\times{_{f}N_2}$ of a nearly cosymplectic manifold $\bar M$, we have the following result.\\

\noindent
{\bf{Theorem 3.1. [15]}} {\it A warped product submanifold $M=N_1\times{_{f} N_2}$ of a nearly cosymplectic manifold $\bar M$ is an usual Riemannian product if the structure vector field $\xi$ is tangent to $N_2$, where $N_1$ and $N_2$ are the Riemannian submanifolds of $\bar M$.}

\parindent=8mm
From the above theorem for the existence of warped products we always consider the structure vector field $\xi$ is tangent to the base. Now, we start with the warped product semi-slant submanifolds of the type $M=N_\theta\times{_{f}N_T}$ of a nearly cosymplectic manifold $\bar M$.\\

\noindent
{\bf{Theorem 3.2}} {\it There do not exist the warped product semi-slant submanifolds M=$N_\theta \times{_{f}N_T}$ of a nearly cosymplectic manifold $\bar M$, where $N_\theta$ and $N_T$ are proper slant and invariant submanifolds of $\bar M$, respectively.}\\

\noindent
{\it Proof.} The proof is similar to the semi-invariant case, which we have proved in [15].$~\blacksquare$

\parindent=8mm
Now, we discuss the other case and all results before the geometric inequality are preparatory and we can not skip all these. We shall use these preparatory results in our main theorem.\\

\noindent
{\bf{Lemma 3.2}} {\it Let $M=N_T\times{_{f}N_\theta}$ be a warped product semi-slant submanifold of a nearly cosymplectic manifold $\bar M$, where $N_T$ and $N_\theta$ are invariant and proper slant submanifolds of $\bar M$, respectively. Then}
\begin{enumerate}
\item[{(i)}] $\xi\ln f=0$
\item[{(ii)}] $g(h(X, Y), FZ)=0$
\item[{(iii)}] $g(h(PX, Z), FZ)=(X\ln f)\|Z\|^2$
\item[{(iv)}] $g(h(X, Z), FPZ)=-g(h(X, PZ), FZ)$
\item[{(v)}]  $g({\cal P}_XZ, PZ)= 2  g(h(X, Z), FPZ)$
\end{enumerate}

\noindent {\it{for any $X, Y\in\Gamma(TN_T)$ and $Z\in\Gamma(TN_\theta)$.}}\\

\noindent{\it Proof.} Let $M=N_T\times{_{f}N_\theta}$ be a warped product semi-slant submanifold of a nearly cosymplectic manifold $\bar M.$ We assume that the structure vector field $\xi$ is tangent to $N_T$, then for any $Z\in\Gamma(TN_\perp)$, we have
$$\bar\nabla_Z\xi=\nabla_Z\xi+h(Z, \xi).$$
Taking the inner product with $Z$ and using Proposition 2.1 and Lemma 3.1 $(ii)$, we obtain
$$(\xi\ln f)\|Z\|^2=0.$$
This means that either $M$ is invariant or $\xi\ln f=0$, which proves $(i)$. Now, we consider $X, Y\in\Gamma(TN_T)$ and $Z\in\Gamma(TN_\theta)$, the we have
$$g(h(X, Y), FZ)=g(\bar\nabla_XY, FZ)=-g(Y, \bar\nabla_XFZ).$$
Using (2.7) and then (2.9), we obtain
$$g(h(X, Y), FZ)=-g(Y, (\bar\nabla_X\phi)Z)-g(Y, \phi\bar\nabla_XZ)+g(Y, \bar\nabla_XPZ).$$
Then from (2.2), (2.4) and Lemma 3.1 (ii), the second and last terms of right hand side vanish identically and hence by (2.10), we derive
$$g(h(X, Y), FZ)=-g(Y, {\cal P}_XZ).$$
Thus, on using the property $p_3~(i)$, we get
$$g(h(X, Y), FZ)=g({\cal P}_XY, Z).$$
Hence, by skew-symmetry of ${\cal P}_XY$ and symmetry of $h(X, Y)$, we get the second part of the lemma. For the third part, consider for any $X\in\Gamma(TN_T)$ and $Z\in\Gamma(TN_\theta)$, we have
$$g(h(PX, Z), FZ)=g(\bar\nabla_{Z}PX, \phi Z-PZ)~~~~~~~~~~~~~~$$
$$~~~~~~~~~~~~~~~~~~~~~~~~~=-g(PX, \bar\nabla_{Z}\phi Z)-g(\bar\nabla_{Z}PX, PZ).$$
From (2.4), (2.9) and Lemma 3.1 $(ii)$, the above equation reduced to
$$g(h(PX, Z), FZ)=-g(PX, (\bar\nabla_{Z}\phi) Z)-g(PX, \phi\bar\nabla_{Z}Z)-(PX\ln f)g(Z, PZ).$$
On using the structure equation of nearly cosymplectic and the fact that $Z$ and $PZ$ are orthogonal vector fields, the first and last terms of right hand side are identically zero. Then from (2.2), we derive
$$g(h(PX, Z), FZ)=g(\phi^2X, \bar\nabla_{Z}Z).$$
Using (2.1), we get
$$g(h(PX, Z), FZ)=-g(X, \bar\nabla_{Z}Z)+\eta(X)g(\xi, \bar\nabla_{Z}Z).$$
By the property of Riemannian connection $\bar\nabla$, the above equation takes the form
$$g(h(PX, Z), FZ)=g(\bar\nabla_{Z}X, Z)-\eta(X)g(\bar\nabla_{Z}\xi, Z).$$
Then from (2.4), Proposition 2.1 and Lemma 3.1 $(ii)$, we obtain
$$g(h(PX, Z), FZ)=(X\ln f)\|Z\|^2,$$
which is third part of the lemma. For the other parts, consider
$$g(\nabla_{PZ}\phi X, Z)=g(\bar\nabla_{PZ}\phi X, Z)$$
for any $X\in\Gamma(TN_T)$ and $Z\in\Gamma(TN_\theta)$. Using the property of Riemannian connection $\bar\nabla$ and Lemma 3.1 $(ii)$, we get
$$(\phi X\ln f)g(Z, PZ)=-g(\phi X, \bar\nabla_{PZ}Z).$$
Using the fact that $Z$ and $PZ$ are orthogonal vector fields, the above equation reduces to
$$0=g(X, \phi\bar\nabla_{PZ}Z).$$
Then form (2.9), we derive
$$0=g(X, \bar\nabla_{PZ}\phi Z)-g(X, (\bar\nabla_{PZ}\phi)Z).$$
By (2.7) and (2.10), we obtain
$$0=g(X, \bar\nabla_{PZ}PZ)+g(X, \bar\nabla_{PZ}FZ)-g(X, {\cal{P}}_{PZ}Z).$$
Using (2.4), (2.5) and (2.19) $(a)$, we get
$$0=-g(\nabla_{PZ}X, PZ)-g(X, A_{FZ}PZ)+g(X, {\cal{P}}_{Z}PZ).$$
Then from the property $p_3~(i)$ and Lemma 3.1 $(ii)$, we obtain
$$0=-(X\ln f)g(PZ, PZ)-g(h(X, PZ), FZ)-g({\cal{P}}_{Z}X, PZ).$$
Again using (2.19) $(a)$, (2.21) and the fact that $\xi$ is tangent to $N_T$, we derive
$$g({\cal{P}}_{X}Z, PZ)=(X\ln f)\cos^2\theta\|Z\|^2+g(h(X, PZ), FZ).\eqno(3.4)$$
Interchanging $Z$ by $PZ$ and then using (2.20), (2.21) and the fact that $\xi$ is tangent to $N_T$, we obtain
$$-\cos^2\theta g({\cal{P}}_{X}PZ, Z)=(X\ln f)\cos^4\theta\|Z\|^2-\cos^2\theta g(h(X, Z), FPZ).$$
By the property $p_3~(i)$, the above equation will be
$$g({\cal{P}}_{X}Z, PZ)=(X\ln f)\cos^2\theta\|Z\|^2-g(h(X, Z), FPZ).\eqno(3.5)$$
Thus, the fourth part of the lemma follows from (3.4) and (3.5). Now, for the part $(v)$, we consider
$$g(h(X, PZ), FZ)=g(\bar\nabla_XPZ, FZ)$$
for any $X\in\Gamma(TN_T)$ and $Z\in\Gamma(TN_\theta)$. Using the property of Riemannian connection $\bar\nabla$ and then using (2.7), we have
$$g(h(X, PZ), FZ)=-g(PZ, \bar\nabla_X\phi Z)+g(PZ, \bar\nabla_XPZ).$$
Using (2.9), Lemma 3.1 $(ii)$, (2.21) and the fact that $\xi$ is tangent to $N_T$, we obtain
$$g(h(X, PZ), FZ)=-g(PZ, \phi\bar\nabla_XZ)-g(PZ, (\bar\nabla_X\phi)Z)+(X\ln f)\cos^2\theta\|Z\|^2.$$
Then from (2.2) and (2.10), we get
$$g(h(X, PZ), FZ)=g(\phi PZ, \bar\nabla_XZ)-g(PZ, {\cal{P}}_XZ)+(X\ln f)\cos^2\theta\|Z\|^2.$$
Using (2.4) and (2.7), we derive
$$g(h(X, PZ), FZ)=g(P^2Z, \nabla_XZ)+g(h(X, Z), FPZ)~~~~~~~~~~~~~~~$$
$$~~~~~~~~~~~~~~~-g(PZ, {\cal{P}}_XZ)+(X\ln f)\cos^2\theta\|Z\|^2.$$
Again, from the fact that $\xi$ is tangent to $N_T$ and using (2.20), the above equation reduces to
$$g({\cal{P}}_XZ, PZ)=g(h(X, Z), FPZ)-g(h(X, PZ), FZ).\eqno(3.6)$$
Thus, the fifth part of the lemma follows from (3.6) and the fourth part of this lemma. This proves the lemma completely.$~\blacksquare$\\

\noindent
{\bf{Theorem 3.3}} {\it Let $M=N_T\times{_{f}N_\theta}$ be a warped product semi-slant submanifold in a nearly cosymplectic manifold $\bar M$, where $N_T$ and $N_\theta$ are invariant and proper slant submanifolds of $\bar M$, respectively. Then,
$$g({\cal P}_XZ, PZ)=\frac{2}{3}(X\ln f)\cos^2\theta\|Z\|^2$$
for any $X\in\Gamma(TN_T)$ and $Z\in\Gamma(TN_\theta)$.}\\

\noindent
{\it Proof.} From (3.4) and (3.5), we have
$$2g({\cal{P}}_{X}Z, PZ)=2(X\ln f)\cos^2\theta\|Z\|^2+g(h(X, PZ), FZ)$$
$$-g(h(X, Z), FPZ).\eqno(3.7)$$
Then, by Lemma 3.2 $(iv)$, we obtain
$$2g({\cal{P}}_{X}Z, PZ)=2(X\ln f)\cos^2\theta\|Z\|^2-2g(h(X, Z), FPZ).\eqno(3.8)$$
Thus, from Lemma 3.2 $(v)$ and (3.8), we obtain the desire result. $~\blacksquare$

\parindent=8mm From the above theorem we have the following consequence.\\

\noindent
{\bf{Corollary 3.1}} {\it{ A warped product semi-slant submanifold $M=N_T\times{_{f}N_\theta}$ of a nearly cosymplectic manifold $\bar M$ is simply a Riemannian product of $N_T$ and $N_\theta$ if and only if ${\cal P}_XZ\in\Gamma(TN_T)$, for any X$\in\Gamma(TN_T)$ and Z$\in\Gamma(TN_\theta)$.}}\\

\noindent
{\bf{Theorem 3.4}} {\it Let $M=N_T\times{_{f}N_\theta}$ be a warped product semi-slant submanifold in a nearly cosymplectic manifold $\bar M$, where $N_T$ and $N_\theta$ are invariant and proper slant submanifolds of $\bar M$, respectively. Then,
$$g(h(X, Z), FPZ)=-g(h(X, PZ), FZ)=\frac{1}{3}(X\ln f)\cos^2\theta\|Z\|^2,$$
for any $X\in\Gamma(TN_T)$ and $Z\in\Gamma(TN_\theta)$.}\\

\noindent
{\it Proof.} The first equality is nothing but Lemma 3.2 (iv) and the second equality is directly followed by the equation (3.8) and Lemma 3.2 $(v).~\blacksquare$

\parindent=8mm From the above theorem we have the following corollaries.\\

\noindent
{\bf{Corollary 3.2}} {\it{A semi-slant warped product submanifold $M=N_T\times{_{f}N_\theta}$ of a nearly cosymplectic manifold is simply a Riemannian product of $N_T$ and $N_\theta$ if and only if $h(X, Z)\in\Gamma(\mu)$, for all  X$\in\Gamma(TN_T)$ and Z$\in\Gamma(TN_\theta)$, where $\mu$ is the invariant normal subbundle of $T^\perp M$.}}\\

\noindent
{\bf{Corollary 3.3}} {\it{There does not exist a mixed geodesic warped product semi-slant submanifold of a nearly cosymplectic manifold.}}

\parindent=8mm From Lemma 3.2 $(i),~(iv)$ and Theorem 3.4, we obtain
$$g(h(\xi, Z), FPZ)=g(h(\xi, PZ), FZ)=0,\eqno(3.9)$$
for any $Z\in\Gamma(TN_\theta)$.

\parindent=8mm Using the previous results, we derive the following geometric inequality for the length of second fundamental form.\\

{\bf{Theorem 3.5}} {\it Let M=$N_T\times {_{f} N_\theta}$ be a warped product semi-slant submanifold of a nearly cosymplectic manifold $\bar M$, where $N_T$ and $N_\theta$ are invariant and proper slant submanifolds of $\bar M$, respectively. Then}
\begin{enumerate}
\item[{(i)}] {\it{The squared norm of the second fundamental form $h$ satisfies}}
$$\|h\|^2 \geq 4q\csc^2\theta\{1+\frac{1}{9}\cos^2\theta\}\|\nabla \ln f\|^2\eqno(3.10)$$
{\it{where $\nabla \ln f$ is the gradient of $\ln f$ and $2q$ is the dimension of $N_\theta.$}}
\item[{(ii)}] {\it{If the equality sign of (3.10) holds, then $N_T$ is totally geodesic in $\bar M$ and $N_\theta$ is a totally umbilical submanifold of $\bar M$. Moreover, $M$ is a minimal submanifold of $\bar M$.}}
\end{enumerate}

\noindent
{\it {Proof.}} Let $\bar M$ be a $(2n+1)$-dimensional nearly cosymplectic manifold and let $M=N_T\times {_{f} N_\theta}$ be a warped product semi-slant submanifold of $\bar M$. Let us denote by ${\cal D}$ and ${\cal D}_\theta$ the tangent bundles on $N_T$ and $N_\theta$, respectively and let $dim N_T=2p+1$ and $dim N_\theta=2q$, then $m=2p+2q+1$. Let $\{X_1,\cdots, X_p, X_{p+1}=\phi X_1,\cdots, X_{2p}=\phi X_p, X_{2p+1}=\xi\}$ and $\{ Z_1,\cdots, Z_{q} , Z_{q+1}=\sec\theta PZ_1,\cdots, Z_{2q}=\sec\theta PZ_{q}\}$ be the local orthonormal frames of ${\cal D}$ and ${\cal D}_\theta$, respectively. Then, the orthonormal frames of $F{\cal D}_\theta$ and $\mu$ are $\{Z_1^\star=\csc\theta FZ_1,\cdots, Z_q^\star=\csc\theta FZ_{q}, Z_{q+1}^\star=\csc\theta\sec\theta FPZ_1,\cdots, Z_{2q}^\star=\csc\theta\sec\theta FPZ_{q}\}$ and $\{N_{m+2q+1}^\star\cdots, N_{2n+1}^\star\}$, respectively, where $N$ is the normal vector in the invariant normal subbundle of $T^\perp M$. The dimensions of $F{\cal D}_\theta$ and $\mu$ will be $2q$ and $2n+1-m-2q$, respectively. The length of second fundamental form $h$ is defined as
$$\|h\|^2=\sum_{r=m+1}^{2n+1}\sum_{i, j=1}^{m} g(h(U_i, V_j), N_r)^2\eqno(3.11)$$
for any vector fields $U_i, V_j$ tangent to $M$ and $N_r$ normal to $M$. Now, for the assumed frames of $F{\cal D}_\theta$ and $\mu$, the above equation can be written as
$$\|h\|^2=\sum_{r=m+1}^{m+2q}\sum_{i, j=1}^{m} g(h(U_i, V_j), Z_r^\star)^2+\sum_{r=m+2q+1}^{2n+1}\sum_{i, j=1}^{m}g(h(U_i, V_j), N_r^\star)^2.\eqno(3.12)$$
The first term in the right hand side of the above equation is the $F{\cal D}_\theta$-component and the second term is $\mu$-component. Here, we equate only $F{\cal D}_\theta$-component term, thus we have
$$\|h\|^2\geq\sum_{r=m+1}^{m+2q}\sum_{i, j=1}^{m} g(h(U_i, V_j), Z_r^\star)^2.\eqno(3.13)$$
Thus, for the given frames of ${\cal D}$ and ${\cal D}_\theta$, the above equation will be
$$\|h\|^2\geq\sum_{r=1}^{2q} \sum_{i, j=1}^{2p+1} g(h(X_i, X_j), Z_r^\star)^2+2\sum_{r=1}^{2q}\sum_{i=1}^{2p+1}\sum_{ j=1}^{2q} g(h(X_i, Z_j), Z_r^\star)^2$$
$$~~~~+\sum_{r=1}^{2q} \sum_{i, j=1}^{2q} g(h(Z_i, Z_j), Z_r^\star)^2.\eqno(3.14)$$
By Lemma 3.2 (ii), the first term of the right hand side is identically zero and we shall compute the next term and leave the third term
$$\|h\|^2\geq2\sum_{r=1}^{2q}\sum_{i=1}^{2p+1}\sum_{j=1}^{2q}g(h(X_i, Z_j), Z_r^\star)^2.$$
As $j, r=1,\cdots, 2q$, then the above equation can be written for one summation as
$$\|h\|^2\geq2\sum_{i=1}^{2p+1}\sum_{j=1}^{2q}g(h(X_i, Z_j), Z_j^\star)^2.$$
Separating the $h(\xi, Z)$-components, the above inequality will be
$$\|h\|^2\geq2\sum_{i=1}^{2p}\sum_{j=1}^{2q}g(h(X_i, Z_j), Z_j^\star)^2+2\sum_{j=1}^{2q}g(h(\xi, Z_j), Z_j^\star)^2.\eqno(3.15)$$
Now, we solve the second term of right hand side of (3.15) as follows
$$\sum_{j=1}^{2q}g(h(\xi, Z_j), Z_j^\star)^2=\csc^2\theta\sum_{j=1}^{q}g(h(\xi, Z_j), FZ_j)^2$$
$$~~~~~~~~~~~~+\csc^2\theta\sec^4\theta\sum_{j=1}^{q}g(h(\xi, PZ_j), FPZ_j)^2$$
$$~~~~~~~~~~+\csc^2\theta\sec^2\theta\sum_{j=1}^{q}g(h(\xi, Z_j), FPZ_j)^2$$
$$~~~~~~~~~~~+\csc^2\theta\sec^2\theta\sum_{j=1}^{q}g(h(\xi, PZ_j), FZ_j)^2.\eqno(3.16)$$
From (3.9), the last two terms of right hand side of (3.16) are identically zero and we will compute the first two terms as follows. We know that
$$g(h(\xi, Z), FZ)=g(\bar\nabla_Z\xi, FZ)=-g(\xi, \bar\nabla_ZFZ).$$
Using (2.7) and then (2.9), we obtain
$$g(h(\xi, Z), FZ)=-g(\xi, (\bar\nabla_Z\phi)Z)-g(\xi, \phi\bar\nabla_ZZ)+g(\xi, \bar\nabla_ZPZ).$$
Using the nearly cosymplectic character the first term of right hand side is zero, second and last terms are also zero by using (2.1), property of Riemannian connection and either orthogonality of vectors $Z$ and $PZ$ or the fact that $\xi\ln f=0$ or both, hence
$$g(h(\xi, Z), FZ)=0.\eqno(3.17)$$
If we interchange $Z$ by $PZ$, then $g(h(\xi, PZ), FPZ)=0.$ Put all these values in (3.16), we obtain
$$\sum_{j=1}^{2q}g(h(\xi, Z_j), FZ_j)^2=0.\eqno(3.18)$$
Thus, from (3.15) and (3.18), we derive
$$\|h\|^2\geq2\sum_{i=1}^{2p}\sum_{j=1}^{2q}g(h(X_i, Z_j), Z_j^\star)^2.$$
Using the frame of $F{\cal D}_\theta$, the above inequality can be written as
$$\|h\|^2\geq2\csc^2\theta\sum_{i=1}^{p}\sum_{j=1}^{q}g(h(X_i, Z_j), FZ_j)^2~~~~~~~~~~$$
$$~~~~~~~~~~~+2\csc^2\theta\sec^4\theta\sum_{i=1}^{p}\sum_{j=1}^{q}g(h(X_i, PZ_j), FPZ_j)^2$$
$$+2\csc^2\theta\sum_{i=1}^{p}\sum_{j=1}^{q}g(h(\phi X_i, Z_j), FZ_j)^2~$$
$$~~~~~~~~~~~~~~+2\csc^2\theta\sec^4\theta\sum_{i=1}^{p}\sum_{j=1}^{q}g(h(\phi X_i, PZ_j), FPZ_j)^2~~$$
$$~~~~~~~~~~~+2\csc^2\theta\sec^2\theta\sum_{i=1}^{p}\sum_{j=1}^{q}g(h(X_i, PZ_j), FZ_j)^2$$
$$~~~~~~~~~~~+2\csc^2\theta\sec^2\theta\sum_{i=1}^{p}\sum_{j=1}^{q}g(h(X_i, Z_j), FPZ_j)^2~~~~~~~$$
$$~~~~~~~~~~~+2\csc^2\theta\sec^2\theta\sum_{i=1}^{p}\sum_{j=1}^{q}g(h(\phi X_i, PZ_j), FZ_j)^2$$
$$~~~~~~~~~~~+2\csc^2\theta\sec^2\theta\sum_{i=1}^{p}\sum_{j=1}^{q}g(h(\phi X_i, Z_j), FPZ_j)^2.\eqno(3.19)$$
The first four terms of above inequality will be solved as follows. From Lemma 3.2 (iii), we have
$$g(h(\phi X, Z), FZ)=(X\ln f)g(Z, Z).\eqno(3.20)$$
Interchanging $X$ by $\phi X$ and then using (2.1), we derive
$$-g(h(X, Z), FZ)+\eta(X)g(h(\xi, Z), FZ)=(\phi X\ln f)g(Z, Z).$$
But from (3.17), the second term of left hand side of above equation vanishes identically, thus we obtain
$$g(h(X, Z), FZ)=-(\phi X\ln f)g(Z, Z).\eqno(3.21)$$
Interchanging $Z$ by $PZ$ in (3.20) and (3.21) and using (2.21) and the fact that $\xi$ is tangent to $N_T$, thus we obtain the following equations, respectively
$$g(h(\phi X, PZ), FPZ)=(X\ln f)\cos^2\theta g(Z, Z)\eqno(3.22)$$
and
$$g(h(X, PZ), FPZ)=-(\phi X\ln f)\cos^2\theta g(Z, Z).\eqno(3.23)$$
The last four terms of (3.19) will be solved by Theorem 3.4 as follows
$$g(h(X, Z), FPZ)=-g(h(X, PZ), FZ)=\frac{1}{3}(X\ln f)\cos^2\theta g(Z, Z).\eqno(3.24)$$
Interchanging $X$ by $\phi X$ in (3.24), we obtain
$$g(h(\phi X, Z), FPZ)=-g(h(\phi X, PZ), FZ)=\frac{1}{3}(\phi X\ln f)\cos^2\theta g(Z, Z).\eqno(3.25)$$
Put all these values in (3.19), we derive
$$\|h\|^2\geq4\csc^2\theta\sum_{i=1}^{p}\sum_{j=1}^{q}\{1+\frac{1}{9}\cos^2\theta\}(X_i\ln f)^2g(Z_j, Z_j)^2$$
$$~~~~~~~~+4\csc^2\theta\sum_{i=1}^{p}\sum_{j=1}^{q}\{1+\frac{1}{9}\cos^2\theta\}(\phi X_i\ln f)^2g(Z_j, Z_j)^2.$$
Hence, from (3.3) the above expression will be
$$\|h\|^2\geq4\csc^2\theta\sum_{j=1}^{q}\{1+\frac{1}{9}\cos^2\theta\}\|\nabla\ln f\|^2g(Z_j, Z_j)^2$$
$$~~~~~~~~~~~=4q\csc^2\theta\{1+\frac{1}{9}\cos^2\theta\}\|\nabla\ln f\|^2,$$
which is the inequality (3.10). If the equality holds in (3.10), then by (3.12) and (3.14), we obtain
$$h({\cal D},{\cal D})=0,~~~h({\cal D}_\theta, {\cal D}_\theta)=0,~~h({\cal D}, {\cal D}_\theta)\subset F{\cal D}_\theta.\eqno(3.26)$$
Now, for any $Z, W\in\Gamma({\cal D}_\theta)$ and $X\in\Gamma({\cal D})$, we have
$$g(h^\star(Z, W), X)=g(\nabla_ZW, X)=-g(W, \nabla_ZX)=-(X\ln f)g(Z, W)$$
where $h^\star$ is the second fundamental form of $N_\theta$ in $M$. Using (3.2), we obtain
$$g(h^\star(Z, W), X)=-g(Z, W)g(\nabla\ln f, X),$$
where $\nabla\ln f$ is the gradient of $\ln f$. Thus from the last relation, we get
$$h^\star(Z, W)=-\nabla\ln fg(Z, W).\eqno(3.27)$$
Since $N_T$ is totally geodesic submanifold in $M$ (by Lemma 3.1 (i)), using this fact with the first condition of (3.26), we get $N_T$ is totally geodesic in $\bar M$. Also, the second condition of (3.26) with (3.27) implies that $N_\theta$ is totally umbilical in $\bar M$. Moreover all conditions of (3.26) imply that $M$ is a minimal submanifold of $\bar M$. This completes the proof of the theorem.~$\blacksquare$\\

\parindent=8mm If we consider $\theta=\frac{\pi}{2}$, then the inequality (3.10) generalizes the inequality which we have obtained for contact CR-warped products in [16].


\section{Introduction}

Energy-momentum localization plays an important role among the open issues  appeared over the years in general relativity. The difficulty which arises consists in constructing a properly defined  energy density of
gravitating systems. As a result, up today there is no generally accepted satisfactory description
for gravitational energy.

A number of researchers have used different  methods
for the energy-momentum localization. Standard research methods include the use
of different tools such as superenergy tensors \cite{Bel}, quasi-local expressions
\cite{Brown} and the well-known energy-momentum complexes of Einstein \cite{Einstein},
Landau-Lifshitz \cite{Landau}, Papapetrou \cite{Papapetrou}, Bergmann-Thomson \cite{Bergmann}, M\o ller \cite{Moller}, Weinberg \cite{Weinberg}, and Qadir-Sharif \cite{Qadir}. The coordinate system dependence of these computational tools remains the main problem encountered. An alternative to avoid the dependence on coordinates is the tele-parallel theory of gravitation \cite{Moller_1} which has been used in many studies on the calculation of the energy and momentum distribution.

Concerning the pseudotensorial prescriptions, only the M\o ller energy-momentum complex is coordinate independent. Schwarzschild Cartesian coordinates and Kerr-Schild Cartesian coordinates are used to
compute the energy-momentum in the case of the other pseudotensorial
prescriptions. In the last decades, despite the critique on energy-momentum
complexes, their application has yielded physically reasonable results in
the case of many spacetime geometries, in particular for geometries in $(3+1)$, $(2+1)$ and $(1+1)$ dimensions \cite{Virbhadra}. It is worth noticing that there is an agreement between on one hand the Einstein, Landau-Lifshitz, Papapetrou, Bergmann-Thomson, Weinberg and M\o ller prescriptions and, on the other hand,
the quasi-local mass definition introduced by Penrose \cite{Penrose} and developed by
Tod \cite{Tod} for some gravitational backgrounds. It should be stressed that
several pseudotensorial definitions yield the same
results for any metric of the Kerr-Schild class as well as  for
solutions more general than the Kerr-Schild class (see, e.g., \cite{Aguirregabiria}).
Further, one can underline the similarity of some of the
aforementioned results with those obtained by using the tele-parallel theory
of gravitation \cite{Gamal}. Finally, one can also point out the ongoing elaboration of the definition and application of energy-momentum complexes, as well as the attempt for their rehabilitation \cite{Chang}.

The present work has the following structure: in Section 2 we describe the new four-dimensional charged regular black hole solution \cite{Balart} under study. Section 3  contains the presentation of the Einstein, Landau-Lifshitz, Weinberg and M\o ller energy-momentum complexes utilized for the calculations. Section 4 contains  the computations of the energy and momentum distributions. Finally, in the Discussion given in Section 5, we present our results and examine some limiting and particular cases.
Throughout the paper we use geometrized units $(c=G=1)$, while the signature
chosen is $(+,-,-,-)$. The calculations for the Einstein, Landau-Lifshitz and Weinberg energy-momentum complexes are performed in Schwarzschild Cartesian coordinates. Greek indices range from $0$ to $3$,
while Latin indices run from $1$ to $3$.

\section{A New Regular Black Hole Solution in Einstein-nonlinear
Electrodynamics}

Attempts to avoid the central singularity in the presence of a horizon for a classical black hole solution have started already in the 1960's with Bardeen's regular black hole solution \cite{Bardeen} followed up by a number of works in the subject (see, \cite{Ansoldi} for a review of existing models and many relevant references). The interest in this research area was revived after E. Ay\'on-Beato and A. Garc\'{\i}a \cite{Ayon} studied a regular black hole solution coupled to nonlinear electrodynamics, whereby the Einstein-Hilbert action of general relativity is augmented by the addition to the Lagrangian of a term depending nonlinearly on the electromagnetic field tensor. Since then, a number of very interesting works has been published on this subject \cite{Beato}. In fact, it is worth noticing that the aforesaid authors were the first to establish a connection of Bardeen's regular black hole solution with a nonlinear electrodynamics leading to the current interpretation o!
 f an exact solution supported by a nonlinear magnetic monopole \cite{Garcia}. Recently, L. Balart and E.C. Vagenas \cite{Balart} constructed a family of new spherical, static and charged black hole solutions in the context of Einstein-nonlinear electrodynamics by imposing the following conditions: (i) the weak energy condition has to be
satisfied, (ii)  the energy-momentum tensor has to respect the symmetry $ T_{0}^{0}=T_{1}^{1}$ and, (iii)
the metric must have an asymptotical
behaviour as the Reissner-Nordstr\"{o}m black hole solution.

Starting with a general static and spherically symmetric line element:
\begin{equation}\label{line_element}
ds^2=B(r) dt^2- B^{-1}(r) dr^2-r^2(d\theta^2+\sin^2\theta d\phi^2),
\end{equation}
with $B(r)=1-\frac{2m(r)}{r}$ and by relaxing the third condition, the authors have also constructed a more
general family of charged regular black hole solutions  which do not behave asymptotically as the
Reissner-Nordstr\"{o}m black hole metric. In this case, the mass function
depends on two  integer constant parameters $\mu \geqslant 4$ and $\alpha\geqslant 1$:
\begin{equation}\label{mass_function}
\begin{split}
m(r) &=\frac{r^{3}q^{2}}{6}\frac{\Gamma(\frac{1}{\alpha})\Gamma
(\frac{\mu }{\alpha })}{\Gamma (\frac{4}{\alpha })\Gamma
(\frac{\mu -3}{\alpha })}
\left(\frac{6\Gamma (\frac{4}{\alpha })}{\Gamma  (\frac{1}{\alpha})
\Gamma (\frac{\alpha +3}{\alpha })}\frac{M\,}{q^{2}}\right)^{4}\times\\
&\\
&\,_{2}F_{1}\left(\frac{3}{\alpha },\frac{\mu }{\alpha },\frac{3+\alpha }{\alpha }; -(\frac{6\,\Gamma (\frac{4}{\alpha })}{\Gamma (\frac{1}{\alpha })\Gamma (\frac{\alpha +3}{\alpha })}\frac{M\,}{q^{2}}r)^{\alpha}\right),
\end{split}
\end{equation}
where $q$ is the electric charge, $\Gamma$ is the Gamma function, $M$ is the black hole mass and $\,_{2}F_{1}$ is the Gauss hypergeometric function.
Particular cases such as $\mu =5$, $\alpha =2$ and $\mu =6$, $\alpha =3$ lead to the metrics obtained by Bardeen \cite{Bardeen} and Hayward \cite{Hayward}.

The line element (\ref{line_element}) with the mass function (\ref{mass_function}) does not behave asymptotically  as the Reissner-Nordstr\"{o}m black hole solution unless $\mu =4$.  In fact, for $\mu =4$ and $\alpha =3$ in Eq.(\ref{mass_function}) Balart and Vagenas obtained the
asymptotically Reissner-Nordstr\"{o}m solution
\begin{equation}\label{RN_solution}
B(r)=1-\frac{2M}{r}\left[1-\frac{1}{[1+(\frac{2Mr}{q^2 })^3]^{1/3}}\right].
\end{equation}

In the present paper, we will study the, non-asymprotically Reissner-Nordstr\"{o}m, gravitational background obtained for the particular
case $\mu =3$ and $\alpha$ arbitrary. Consequently the metric for the new charged regular black hole solution to be considered in what follows has the line element (\ref{line_element}) with
\begin{equation}\label{metric_coef_B}
B(r)= 1-\frac{2M}{r}\left(1-\frac{1}{[1+(\frac{2Mr}{q^{2}})^{3}]^{(\alpha -3)/3}}\right).
\end{equation}

\begin{figure}[b!]
\includegraphics[width=84mm]{gtt.png}
\caption{Two horizons exist for the charged regular black hole at the points where $g_{tt}$ meets the $r$-axis.}
\label{fig7}
\end{figure}

In fact, this regular black hole solution exhibits two event horizons. By choosing, for example, the numerical values $M=1$ and $q=0.9$, we get  Figure 1 showing the two horizons. One can see that  the position of the inner horizon is shifted towards bigger values of $r$ as the parameter $\alpha$ decreases, while the position of the outer horizon remains unaffected.

\section{Einstein, Landau-Lifshitz, Weinberg and M\o ller Energy-Momentum
Complexes}

The Einstein energy-momentum complex \cite{Einstein} for a $(3+1)$ dimensional
gravitational background is given by the expression
\begin{equation}
\theta _{\nu }^{\mu }=\frac{1}{16\pi }h_{\nu ,\,\lambda }^{\mu \lambda },
\end{equation}
where the superpotentials $h_{\nu }^{\mu \lambda }$ are given as
\begin{equation}
h_{\nu }^{\mu \lambda }=\frac{1}{\sqrt{-g}}g_{\nu \sigma }[-g(g^{\mu \sigma
}g^{\lambda \kappa }-g^{\lambda \sigma }g^{\mu \kappa })]_{,\kappa }
\end{equation}
and satisfy the corresponding antisymmetric property
\begin{equation}
h_{\nu }^{\mu \lambda }=-h_{\nu }^{\lambda \mu }.
\end{equation}
The components $\theta _{0}^{0}$ and $\theta _{i}^{0}$ represent the energy and
the momentum density, respectively. In this prescription the
local conservation law is respected:
\begin{equation}
\theta _{\nu,\mu }^{\mu }=0.
\end{equation}
Thus, the energy and momentum can be calculated by
\begin{equation}
P_{\mu }=\int \int \int \theta _{\mu }^{0}\,dx^{1}dx^{2}dx^{3}.
\end{equation}
Applying Gauss' theorem, the energy-momentum becomes
\begin{equation}\label{Einstein_P}
P_{\mu }=\frac{1}{16\pi }\int\int h_{\mu }^{0i}n_{i}dS,
\end{equation}
with $n_{i}$ the outward unit normal vector over the surface $dS$.

In the Landau-Lifshitz prescription, the corresponding energy-momentum
complex \cite{Landau} is defined as
\begin{equation}
L^{\mu \nu }=\frac{1}{16\pi }S_{,\,\rho \sigma }^{\mu \nu \rho \sigma },
\end{equation}
with the Landau-Lifshitz superpotentials:
\begin{equation}
S^{\mu \nu \rho \sigma }=-g(g^{\mu \nu }g^{\rho \sigma }-g^{\mu \rho }g^{\nu
\sigma }).
\end{equation}
The $L^{00}$ and $L^{0i}$ components are the energy and the momentum density,
respectively. The local conservation law holds:
\begin{equation}
L_{,\,\nu }^{\mu \nu }=0.
\end{equation}
By integrating $L^{\mu \nu }$ over the 3-space one gets for the
energy-momentum:
\begin{equation}
P^{\mu }=\int \int \int L^{\mu 0}\,dx^{1}dx^{2}dx^{3}.
\end{equation}
By using Gauss' theorem we get
\begin{equation}\label{LL_P}
P^{\mu }=\frac{1}{16\pi }\int \int S_{,\nu }^{\mu 0i\nu }n_{i}dS=\frac{1}{16\pi }\int \int U^{\mu 0i}n_{i}dS.
\end{equation}

The Weinberg energy-momentum complex \cite{Weinberg} is given as
\begin{equation}
W^{\mu \nu }=\frac{1}{16\pi }D_{,\,\lambda }^{\lambda \mu \nu },
\end{equation}
where $D^{\lambda \mu \nu }$ represents the corresponding superpotentials:
\begin{equation}
D^{\lambda \mu \nu }=\frac{\partial h_{\kappa }^{\kappa }}{\partial
x_{\lambda }}\eta ^{\mu \nu }-\frac{\partial h_{\kappa }^{\kappa }}{\partial
x_{\mu }}\eta ^{\lambda \nu }-\frac{\partial h^{\kappa \lambda }}{\partial
x^{\kappa }}\eta ^{\mu \nu }+\frac{\partial h^{\kappa \mu }}{\partial
x^{\kappa }}\eta ^{\lambda \nu }+\frac{\partial h^{\lambda \nu }}{\partial
x_{\mu }}-\frac{\partial h^{\mu \nu }}{\partial x_{\lambda }},
\end{equation}
with
\begin{equation}
h_{\mu \nu }=g_{\mu \nu }-\eta _{\mu \nu }.
\end{equation}
Here the $W^{00}$ and $W^{0i}$ components represent the energy and
the momentum density, respectively. In the Weinberg prescription the local
conservation law is also respected:
\begin{equation}
W_{,\,\nu }^{\mu \nu }=0.
\end{equation}
The integration of $W^{\mu \nu }$ over the 3-space gives the energy-momentum:
\begin{equation}
P^{\mu }=\int \int \int W^{\mu 0}\,dx^{1}dx^{2}dx^{3}.
\end{equation}
Applying Gauss' theorem and integrating over the surface of a
sphere of radius $r$, we get the following expression for the energy-momentum distribution:
\begin{equation}\label{W_P}
P^{\mu }=\frac{1}{16\pi }\int \int D^{i0\mu }n_{i}dS.
\end{equation}

The M{\o }ller energy-momentum complex \cite{Moller} is
\begin{equation}
\mathcal{J}_{\nu }^{\mu }=\frac{1}{8\pi }M_{\nu \,\,,\,\lambda }^{\mu
\lambda },
\end{equation}
with the M{\o }ller superpotentials $M_{\nu }^{\mu \lambda }$
\begin{equation}
M_{\nu }^{\mu \lambda }=\sqrt{-g}\left( \frac{\partial g_{\nu \sigma }}{\partial x^{\kappa }}-\frac{\partial g_{\nu \kappa }}{\partial x^{\sigma }} \right) g^{\mu \kappa }g^{\lambda \sigma },
\end{equation}
being  antisymmetric:
\begin{equation}
M_{\nu }^{\mu \lambda }=-M_{\nu }^{\lambda \mu }.
\end{equation}
M{\o }ller's energy-momentum complex also satisfies the local conservation
law
\begin{equation}
\frac{\partial \mathcal{J}_{\nu }^{\mu }}{\partial x^{\mu }}=0,
\end{equation}
while $\mathcal{J}_{0}^{0}$  and $\mathcal{J}_{i}^{0} $ are the energy and the momentum density respectively. Thus, the energy-momentum is given by
\begin{equation}
P_{\mu }=\int \int \int \mathcal{J}_{\mu }^{0}dx^{1}dx^{2}dx^{3}.
\end{equation}
With the use of Gauss' theorem one gets
\begin{equation}\label{M_P}
P_{\mu }=\frac{1}{8\pi }\int\int M_{\mu }^{0i}n_{i}dS.
\end{equation}

\section{Energy and Momentum Distributions for the Charged Regular Black Hole}
In order to perform the calculations using the Einstein energy-momentum
complex, it is necessary to transform the metric given by the line element
(\ref{line_element}) in Schwarzschild Cartesian coordinates. Thus, we obtain the line element
in the form:
\begin{equation}\label{line_element_SC}
ds^{2} =B(r)dt^{2}-(dx^{2}+dy^{2}+dz^{2})-\frac{B^{-1}(r)-1}{r^{2}}(xdx+ydy+zdz)^{2}\text{.}
\end{equation}
Now, after computing the superpotentials in the Einstein prescription we get
the following non-vanishing components:
\begin{equation}\label{h001}
h_{0}^{01}=\frac{2x}{r^{2}}\frac{2M}{r}\left\{1-\frac{1}{[1+(\frac{2Mr}{q^{2}})^{3}]^{(\alpha -3)/3}}\right\},
\end{equation}
\begin{equation}\label{h002}
h_{0}^{02}=\frac{2y}{r^{2}}\frac{2M}{r}\left\{1-\frac{1}{[1+(\frac{2Mr}{q^{2}})^{3}]^{(\alpha -3)/3}}\right\},
\end{equation}
\begin{equation}\label{h003}
h_{0}^{03}=\frac{2z}{r^{2}}\frac{2M}{r}\left\{1-\frac{1}{[1+(\frac{2Mr}{q^{2}})^{3}]^{(\alpha -3)/3}}\right\}.
\end{equation}
Using the line element (\ref{line_element}), with the metric coefficient (\ref{metric_coef_B}), the equation (\ref{Einstein_P}) for the energy and the expressions (\ref{h001})-(\ref{h003}) for the superpotentials, we obtain the energy distribution for the  regular charged black hole in the Einstein prescription:
\begin{equation}\label{Einstein_Energy}
E_{E}=M\left\{1-\frac{1}{[1+(\frac{2Mr}{q^{2}})^{3}]^{(\alpha -3)/3}}\right\}.
\end{equation}
Further, all the momenta are zero.

In Figures 2, 3 we can see the graph of the energy distribution as a function of the distance for different values of $\alpha$, while Figure 4 presents the energy distribution as a function of charge and mass.

\begin{figure}[h!]
\includegraphics[width=84mm]{e1.png}
\caption{Energy distribution calculated by the Einstein prescription between the inner and outer horizons.}
\label{fig5}
\end{figure}

\begin{figure}[h!]
\includegraphics[width=84mm]{e2.png}
\caption{Energy distribution calculated by the Einstein prescription outside the outer horizon.}
\label{fig6}
\end{figure}

\begin{figure}[h!]
\includegraphics[width=84mm]{3de.png}
\caption{Variation of the energy distribution calculated by the Einstein prescription with respect to charge and mass outside the outer horizon.}
\label{fig1}
\end{figure}

In order to apply the Landau-Lifshitz prescription, we use the $U^{\mu 0i}$ quantities defined in Eq.(\ref{LL_P}) to compute the energy--momentum distribution.
The nonvanishing components of the $U^{\mu 0i}$ quantities are:
\begin{equation}\label{Uttx}
U^{ttx}=\frac{2x}{r^{2}}
\frac{\frac{2M}{r}\left[1-\frac{1}{[1+(\frac{2Mr}{q^{2}})^{3}]^{(\alpha -3)/3}}\right]}
{1-\frac{2M}{r}\left[1-\frac{1}{[1+(\frac{2Mr}{q^{2}})^{3}]^{(\alpha -3)/3}}\right]},
\end{equation}
\begin{equation}\label{Utty}
U^{tty}=\frac{2y}{r^{2}}
\frac{\frac{2M}{r}\left[1-\frac{1}{[1+(\frac{2Mr}{q^{2}})^{3}]^{(\alpha -3)/3}}\right]}
{1-\frac{2M}{r}\left[1-\frac{1}{[1+(\frac{2Mr}{q^{2}})^{3}]^{(\alpha -3)/3}}\right]},
\end{equation}
\begin{equation}\label{Uttz}
U^{ttz}=\frac{2z}{r^{2}}
\frac{\frac{2M}{r}\left[1-\frac{1}{(1+(\frac{2Mr}{q^{2}})^{3})^{(\alpha -3)/3}}\right]}
{1-\frac{2M}{r}\left[1-\frac{1}{[1+(\frac{2Mr}{q^{2}})^{3}]^{(\alpha -3)/3}}\right]}.
\end{equation}
Substituting (\ref{Uttx})-(\ref{Uttz}) in (\ref{LL_P}), we get the energy distribution:
\begin{equation}\label{LL_Energy}
E_{LL}=\frac{M\left[1-\frac{1}{[1+(\frac{2Mr}{q^{2}})^{3}]^{(\alpha -3)/3}}\right]}
{1-\frac{2M}{r}\left[1-\frac{1}{[1+(\frac{2Mr}{q^{2}})^{3}]^{(\alpha-3)/3}}\right]}.
\end{equation}
Also, for this prescription  all the momenta are found to vanish.

In the case of the Weinberg prescription we calculate the non-vanishing
components:
\begin{equation}\label{Dxtt}
D^{xtt}=\frac{2x}{r^{2}}
\frac{\frac{2M}{r}\left[1-\frac{1}{[1+(\frac{2Mr}{ q^{2}})^{3}]^{(\alpha -3)/3}}\right]}
{1-\frac{2M}{r}\left[1-\frac{1}{[1+(\frac{2Mr}{q^{2}})^{3}]^{(\alpha -3)/3}}\right]},
\end{equation}
\begin{equation}\label{Dytt}
D^{ytt}=\frac{2y}{r^{2}}
\frac{\frac{2M}{r}\left[1-\frac{1}{[1+(\frac{2Mr}{q^{2}})^{3}]^{(\alpha -3)/3}}\right]}
{1-\frac{2M}{r}\left[1-\frac{1}{[1+(\frac{2Mr}{q^{2}})^{3}]^{(\alpha -3)/3}}\right]},
\end{equation}
\begin{equation}\label{Dztt}
D^{ztt}=\frac{2z}{r^{2}}
\frac{\frac{2M}{r}\left[1-\frac{1}{[1+(\frac{2Mr}{q^{2}})^{3}]^{(\alpha -3)/3}}\right]}
{1-\frac{2M}{r}\left[1-\frac{1}{[1+(\frac{2Mr}{q^{2}})^{3})^{(\alpha -3)/3}}\right]}.
\end{equation}
Inserting the expressions obrained in (\ref{Dxtt})-(\ref{Dztt}) into (\ref{W_P}) we obtain for the
energy distribution inside a $2$-sphere of radius $r$:
\begin{equation}\label{W_Energy}
E_{W}=M\frac{\left[1-\frac{1}{[1+(\frac{2Mr}{q^{2}})^{3}]^{(\alpha -3)/3}}\right]}
{1-\frac{2M}{r}\left[1-\frac{1}{[1+(\frac{2Mr}{q^{2}})^{3}]^{(\alpha -3)/3}}\right]},
\end{equation}
while all the momenta vanish. As it can be seen, the energy in the Weinberg prescription (\ref{W_Energy}) is identical with the energy in the Landau-Lifshitz prescription (\ref{LL_Energy}).

In Figures 5, 6 we can see the graph of the energy distribution as a function of the distance for different values of $\alpha$, while Figure 7 presents the energy distribution as a function of charge and mass.

\begin{figure}[h!]
\includegraphics[width=84mm]{l1.png}
\caption{Energy distribution calculated by the Landau-Lifshitz prescription between the inner and outer horizons.}
\label{fig8}
\end{figure}

\begin{figure}[h!]
\includegraphics[width=84mm]{l2.png}
\caption{Energy distribution calculated by the Landau-Lifshitz prescription outside the outer horizon.}
\label{fig9}
\end{figure}

\begin{figure}[h!]
\includegraphics[width=84mm]{3dl.png}
\caption{Variation of the energy distribution calculated by the Landau-Lifshitz prescription with respect to charge and mass outside the outer horizon.}
\label{fig2}
\end{figure}

In the M\o ller prescription, the only non-vanishing superpotential is:
\begin{equation}\label{M001}
M_{0}^{01}=\left\{\frac{2M\left[1-\frac{1}{(1+\frac{8M^{3}r^{3}}{q^{6}})^{(\frac{\alpha }{3}-1)}}\right]}
{r^{2}}-\frac{48M^{4}(\frac{\alpha }{3}-1)r}{\left[1+\frac{8M^{3}r^{3}}{q^{6}}\right]^{(\frac{\alpha }{3}-1)}q^{6}\left[1+\frac{8M^{3}r^{3}}{q^{6}}\right]}\right\}r^{2}\sin \theta.
\end{equation}

Applying the last result for the line element (\ref{line_element}) with the metric coefficient (\ref{metric_coef_B}) and using the expression (\ref{M_P}), we get the energy distribution in the M\o ller prescription:
\begin{equation}\label{Moller_Energy}
E_{M}=\frac{r^{2}}{2}
\left\{\frac{2M\left[1-\frac{1}{(1+\frac{8M^{3}r^{3}}{q^{6}})^{(\frac{\alpha }{3}-1)}}\right]}{r^{2}}-
\frac{48M^{4}(\frac{\alpha }{3}-1)r}{\left[1+\frac{8M^{3}r^{3}}{q^{6}}\right]^{(\frac{\alpha }{3}-1)}q^{6}\left[1+\frac{8M^{3}r^{3}}{q^{6}}\right]}\right\}.
\end{equation}
Also, in this prescription it is found that  all the momenta vanish.

In Figures 8, 9 we can see the graph of the energy distribution as a function of the distance for different values of $\alpha$, while Figure 10 presents the energy distribution as a function of charge and mass.

\begin{figure}[t!]
\includegraphics[width=84mm]{m1.png}
\caption{Energy distribution calculated by the M\o ller prescription between the inner and outer horizons.}
\label{fig10}
\end{figure}

\begin{figure}[t!]
\includegraphics[width=84mm]{m2.png}
\caption{Energy distribution calculated by the M\o ller prescription outside the outer horizon.}
\label{fig11}
\end{figure}

\begin{figure}[t!]
\includegraphics[width=84mm]{3dm.png}
\caption{Variation of the energy distribution calculated by the M\o ller prescription with respect to charge and mass outside the outer horizon.}
\label{fig3}
\end{figure}

It is pointed out, according to Figure 8, that for $\alpha \leqslant  6$ the M\o ller energy of the regular black hole is negative and decreases monotonically for the whole range of the values of $r$ considered.

Finally, in Figure 11 we give a comparison of the energy distributions as obtained by the different prescriptions for $\alpha=8$.
\begin{figure}[h!]
\includegraphics[width=84mm]{3E.png}
\caption{Comparison of the energy distributions calculated by the Einstein, Landau-Lifshitz and M\o ller prescriptions outside the outer horizon.}
\label{fig4}
\end{figure}

\section{Discussion}

This paper is focused on the study of the energy-momentum for the gravitational field of a new four-dimensional, spherically symmetric, static and charged regular black hole solution developed
in the context of Einstein-nonlinear electrodynamics coupling, by applying the Einstein,
Landau-Lifshitz, Weinberg and M\o ller energy-momentum complexes. The metric considered, which has two horizons, does not asymptotically behave as the Reissner-Nordstr\"{o}m black hole
solution unless $\mu =4$, where $\mu $ is a positive integer parameter.
The expressions for the energy are well-defined in all the aforesaid prescriptions and
depend on the mass $M$ of the black hole, its charge $q$, a positive integer $\alpha $ and the radial coordinate $r$. Both the Landau-Lifshitz and Weinberg
prescriptions give the same result for the energy distribution, while all the momenta are zero in all four prescriptions used.

In order to examine the physical meaning of the results obtained, we study the limiting
behavior of the energy for $r\rightarrow \infty$, $r=0$ and $q=0$ in three
cases: (i) $1\leqslant \alpha <3$, (ii) $\alpha =3$ and (iii) $\alpha >3$. However, since the cases (i) and (ii) do not yield any physically meaningful result, we present only the results for the case (iii) in Table 1.

\begin{table}[t!]
\centering
\begin{tabular}{|l|c|c|c|}
\hline
Prescription & $r\rightarrow \infty$ & $r=0$ & $q=0$ \\
\hline
&&&\\
Einstein & $M$ & $0$ & $M$ \\
\hline
&&&\\
Landau-Lifshitz & $M$ & $0$ & $M(1-\frac{2M}{r})^{-1}$ \\
\hline
&&&\\
Weinberg & $M$ & $0$ & $M(1-\frac{2M}{r})^{-1}$ \\
\hline
&&&\\
M\o ller & $M$ & $0$ & $M$\\
\hline
\end{tabular}
\caption{Limiting cases for $\alpha>3$}
\end{table}

For $\mu =4$ and $\alpha=3$ we calculate the energy-momentum distribution for the asymprotically Reissner-Nordstr\"{o}m metric with the metric coefficient (\ref{RN_solution}).

For this black hole solution, the energy distribution obtained by using the Einstein, Landau-Lifshitz, Weinberg and M\o ller energy-momentum complexes is given in the following Table 2, along with the limiting cases for $r\rightarrow \infty $, $r=0$ and $q=0$.

\begin{table}[t!]
\centering
\begin{tabular}{|l|c|c|c|c|}
\hline
\textbf{Prescription} & \textbf{Energy} &$r\rightarrow \infty $ & $r=0$ & $q=0$ \\
\hline
&&&&\\
Einstein & $M\left\{1-\frac{1}{[1+(\frac{2Mr}{q^{2}})^{3}]^{1/3}}\right\}$ & $M$ & $0$ & $M$ \\
&&&&\\
\hline
&&&&\\
Landau-Lifshitz & $\frac{M\left[(1-\frac{1}{[1+(\frac{2Mr}{q^{2}})^{3}]^{1/3}}\right]}
{1-\frac{2M}{r}\left[1-\frac{1}{[1+(\frac{2Mr}{q^{2}})^{3}]^{1/3}}\right]}$ & $M$ & $0$ & $M((1-\frac{2\,M}{r})^{-1}$\\
&&&&\\
\hline
&&&&\\
Weinberg & $\frac{M\left[1-\frac{1}{[1+(\frac{2Mr}{q^{2}})^{3}]^{1/3}}\right]}
{1-\frac{2M}{r}\left[1-\frac{1}{[1+(\frac{2Mr}{q^{2}})^{3}]^{1/3}}\right]}$ & $M$ & $0$ & $M((1-\frac{2\,M}{r})^{-1}$\\
&&&&\\
\hline
&&&&\\
M\o ller & $\frac{r^{2}}{2}
\left\{-\frac{2M\left[(\frac{8\,M^{3}r^{3}}{q^{6}}+1)^{-\frac{1}{3}}-1\right]}{r^{2}}-
\frac{16M^{4}r}{q^{6}\left[\frac{8M^{3}r^{3}}{q^{6}}+1\right]^{\frac{4}{3}}}\right\}$& $M$ & $0$ & $M$\\
&&&&\\
\hline
\end{tabular}
\caption{$\mu=4$ and $\alpha=3$}
\end{table}

As we can see, the limiting behavior of the energy for the asymptotically Reissner-Nordstr\"{o}m solution
($\mu=4$, $\alpha=3$) is the same as in the case of the non-asymptotically Reissner-Nordstr\"{o}m solution given in Table 1.

Thus, by comparing the results obtained for the energy with all the energy-momentum complexes applied, we see that, for $q=0$ and at infinity, the Einstein prescription and the M\o ller prescription give the same results (namely, the ADM mass $M$) with that obtained for the Schwarzschild black hole solution. On the other hand, the Landau-Lifshitz and the Weinberg prescriptions, while at infinity lead to the same result (the ADM mass $M$), in the chargeless case $q=0$ give the energy distribution  $E_{LL}=E_{W}=M(1-\frac{2M}{r})^{-1}$ which is in total agreement with the expression obtained for the Schwarzschild black hole using Schwarzschild Cartesian coordinates (see, K.S. Virbhadra (1999), in \cite{Aguirregabiria}).

Finally,  all the four prescriptions used yield a vanishing energy at $r=0$.

Based on the present study, as well as on the aforementioned agreement with previous results
obtained in the case of the Schwarzschild metric, we conclude that the Einstein and M\o ller
energy-momentum complexes give physically meaningful results and, thus,
can be considered as the most reliable tools among the energy-momentum complexes for the study of the energy-momentum localization of gravitational systems.

An interesting future perspective lies in the calculation of the
energy-momentum of the new regular black hole solution considered in the present work by applying  other energy-momentum complexes and/or the tele-parallel equivalent to general relativity (TEGR).