\documentclass[reqno]{amsart} 
\usepackage{psfig}
\begin{document}
{\noindent\small {\em Electronic Journal of Differential Equations},
Vol.~2000(2000), No.~58, pp.~1--32.\newline
ISSN: 1072-6691. URL: http://ejde.math.swt.edu or http://ejde.math.unt.edu
\newline ftp ejde.math.swt.edu \quad ejde.math.unt.edu (login: ftp)}
\thanks{\copyright 2000 Southwest Texas State University  and 
University of North Texas.} 
\vspace{1cm}

\title[\hfilneg EJDE--2000/58\hfil Three-dimensional Kolmogorov problem]
{Steady-state bifurcations of the three-dimensional Kolmogorov problem} 

\author[Z. M. Chen \& S. Wang\hfil EJDE--2000/58\hfilneg]
{Zhi-Min Chen \&  Shouhong Wang}


\address{Zhi-Min Chen \hfill\break 
Department of Ship Science, Southampton University,
Southampton SO17 1BJ, UK \hfill\break and:
Department of Mathematics, Tianjin University,
 China}
 \email{zhimin@ship.soton.ac.uk}
 
 \address{Shouhong Wang \hfill\break 
Department of Mathematics, Indiana University,
Bloomington, IN 47405, USA}
\email{showang@indiana.edu}

 
\date{}
\thanks{Submitted February 12, 1999. Revised May 26, 2000. Published August 30, 2000.}
\subjclass{35Q30, 76D05, 58J55, 35B32}
\keywords{3D Navier-Stokes equations, Kolmogorov flow, multiple steady states,
\hfill\break\indent supercritical pitchfork bifurcation, continuous fractions}


\begin{abstract}
  This paper studies the spatially periodic incompressible fluid motion in
  $\mathbb R^3$ excited by the external force $k^2(\sin kz, 0,0)$ with 
  $k\geq 2$ an integer.  This driving force gives rise to the existence of the
  unidirectional basic steady flow $u_0=(\sin kz,0, 0)$ for any Reynolds number.
  It is shown in Theorem 1.1 that there exist a number of critical Reynolds
  numbers such that $u_0$ bifurcates into either 4 or 8 or 16 different steady
  states, when the Reynolds number increases across each of such numbers.

  Thanks to the Rabinowitz global bifurcation theorem, all of the bifurcation
  solutions are extended to global branches for $\lambda \in (0, \infty)$.
  Moreover we prove that when $\lambda$ passes each critical value, a) all the
  corresponding global branches do not intersect with the trivial branch 
  $(u_0,\lambda)$, and b) some of them never intersect each other; 
  see Theorem 1.2.
\end{abstract} 

\maketitle


\renewcommand{\theequation}{\arabic{section}.\arabic{equation}}
\renewcommand\thefigure{\arabic{section}.\arabic{figure}}
\newtheorem{Theorem}{Theorem}[section]
\newtheorem{Lemma}{Lemma}[section]
\newtheorem{Remark}{Remark}[section]

\def\dbH{\dot{\mathbb H}}
\def\tilde{\widetilde}


\section{Introduction}
\setcounter{equation}{0}

The three-dimensional  Kolmogorov  problem, which was first formulated
by Kolmogorov (see \cite{AM}),  refers  to  the Navier-Stokes equations
defining the spatially periodic  fluid motion in the following form:
\begin{gather}
\frac{\partial u}{\partial t} -\Delta u+\lambda (u\!\cdot\!\nabla) u + \lambda\nabla p
=k^2(\sin kz,0,0), \label{ns1}\\
 \text{div }u = 0, \label{ns2}\\
 u(t,x,y,z) =u(t,x+2\pi,y,z)=u(t,x,y+2\pi,z)=u(t,x,y,z+2\pi),  \label{ns3}\\
 \int_{{\mathbb T}^3}  u\,dx\,dy\,dz= 0\,. \label{ns4}
\end{gather}
Here $u=(u_1,u_2,u_3)$ is the velocity field, $p$ the pressure,
$\lambda >0$ the Reynolds number in this dimensionless formulation,
$k$ a positive integer, and
${\mathbb T}^n={{\mathbb R}}^n/(2 \pi {\mathbb Z})^n$  ($n=2, 3$) the n-dimensional flat
torus.  In particular $u_0=(\sin kz,0,0)$ solves this problem for all
$\lambda$.


There have been extensive mathematical and physical studies for the
Kolmogorov problem as well as for  general fluid equations; see for
instance \cite{CP1,FSV, Green,Iu,Kir, L1, Mar, MS,Obu,Oka,Sirov,Rab, Rab68,
Rab73,serrin, velte}, and the references therein.
From the bifurcation point of view, most of the literatures are devoted to
the existence of secondary steady-states of the Navier-Stokes
equations, and in particular the detailed rigorous mathematical analysis  on
the bifurcation of the solutions of the Kolmogorov problem
as well as the general
Navier-Stokes equations is still in its early stage,
due partially to the difficulty in estimating the eigenvalues
of the associated  linearized operators.

The main objective of this paper is to establish the existence
of multi-branches  of steady-state solutions bifurcated respectively
from $u_0$ at some critical Reynolds numbers.
The multi-branches are obtained in some flow invariant subspaces,
which are defined using special Fourier modes 
in terms of some $n$--tuple integer vectors.

Let $\mathbb N$ be the set of all positive integers, and
$\mathbb Z$ be the integer set.
The following  condition
for an integer vector $(l,j,i,k)$
will be used often to define some
invariant function spaces:
\begin{equation} \label{4-tuple}
\gathered
 (l^2+j^2+i^2-k^2)(l^2+j^2+(k-i)^2-k^2)<0, \\
 0\leq l<k,\,\, j\in {\mathbb Z},\,\,0\leq  i <k\,.
\endgathered
\end{equation}
For such an integer vector $(l,j,i,k)$, we set
\begin{equation} \label{n-zero}
n_0=\max\{n\in \mathbb N  \, |\,\,n^2l^2+n^2j^2+\min\{\{ni\}^2,\,
(k-\{ni\})^2\}-k^2<0\},
\end{equation}
and
\begin{equation} \label{i-zero}
i_0=\left\{
\begin{array}{ll}
4,& \mbox{  when } j=\{n_0i\}=0,\\
8,& \mbox{ when } |j|+\{n_0i\}> 0 \mbox{ and }j\{n_0i\}=0,\\
16,& \mbox{ when } j\{n_0i\}\neq 0.
\end{array}\right.
\end{equation}
Here $\{n\}$  is defined by
\begin{equation} \label{mod-k}
\{n\} \equiv n \mbox{\rm \quad (mod  $k$)}.
\end{equation}


Let $H^m({\mathbb T}^3)$ be the usual Sobolev space of scalar and periodic
functions
endowed with the usual $H^m$
norm $|| \cdot ||_{H^m}$, and let ${\mathbb H}^m({\mathbb T}^3)=H^m({\mathbb T}^3)^3$
be the vectorial
Sobolev space. For $m\ge 1$,  we use the following
function spaces of divergence-free vector fields:
$${\mathbb H}^m_{\sigma}=
\left\{u\in {\mathbb H}^{m}({\mathbb T}^3) \, |\,\text{ div } u=0
 \right\}.$$
Let
$\dot{H}^m({\mathbb T}^3)=H^m({\mathbb T}^3)/{\mathbb R}$, $\dbH^m({\mathbb T}^3) = {\mathbb H}^m({\mathbb T}^3)/ {\mathbb R}$,
and $\dbH^2_{\sigma}= {\mathbb H}^2_{\sigma}/{\mathbb R}$.
When $m=2$, we can use $\|u\|_{H^2_\sigma}=\|\Delta u\|_{L^2}$  as the norm
of $\dbH^2_{\sigma}$.

\begin{figure}
$$\psfig{figure=fig1.eps,width=10cm}$$
\caption{The global branch $u_{l,j,i,k,n,\lambda}$ in $\mathcal C_1$.
Here the dotted horizontal lines represent $ \| u\|_{H^1} = c_2 k^2 $, and
dotted vertical line is $\displaystyle \lambda = \frac{c_1}{k^2}$.}
\end{figure}

\begin{figure}
$$\psfig{figure=fig2.eps,width=10cm}$$
\caption{The global branch $u_{l,j,i,k,n,\lambda}$ in $\mathcal C_2$.
Here the dotted vertical line is $\displaystyle \lambda = \frac{c_1}{k^2}$,
while the dotted curves represent
$\| u\|_{H^2} = c_3 k^2 (\lambda^2 k^4 + 1)$.}
\end{figure}
     


As we mentioned before, the main objective of this article is to
find and possibly classify steady-state bifurcations of
the Kolmogorov problem in some subspaces  of $\dbH^2_{\sigma}$ invariant to
the Navier-Stokes equations. Now we state two main theorems of this article
without specify invariant subspaces, which will be explicitly given in
Sections 4--6.

Thanks to the Rabinowitz global bifurcation theorem,  we obtain 
in these two theorems  $i_0$ 
different global branches of steady-state solutions bifurcating from $u_0$ at 
a single critical value. Each of those branches  undergoes 
local supercritical bifurcation  around the bifurcation point,  
whereas  half ($i_0/2$) of those branches  never touch each other 
 away from the bifurcation point.


\begin{Theorem}
Let  $(l,j,i,k)$ be an integer vector satisfying (\ref{4-tuple}).
Then  there exists  a critical Reynolds number
$\lambda_{l,j,i,k}>0$ such that     (\ref{ns1}-\ref{ns4}) with 
$0<\lambda <\infty $ admit
$i_0$  steady-state solutions
$u_{l,j,i,k,n,\lambda}$ $(n=1,...,i_0)$   branching  off 
$(\lambda_{l,j,i,k},u_0)$ continuously 
such that:
      \footnote{See (\ref{5.1}--\ref{5.5}) for specific spaces
       where these bifurcation solutions are located.}
\begin{enumerate}
\item[a).] if $\lambda_{l,j,i,k}<\lambda$, then 
$$   u_{l,j,i,k,n,\lambda}  =u_0,$$
\item[b).] if $\lambda_{l,j,i,k}<\lambda $  and 
$\|\Delta(u_{l,j,i,k,n,\lambda}-u_0)\|_{L^2} + |\lambda-\lambda_{l,j,i,k}|
<\epsilon $ for some small $\epsilon> 0 $, then 
\begin{equation}\label{1.3}
u_0\neq u_{l,j,i,k,n,\lambda}   \neq u_{l,j,i,k,n',\lambda}\neq u_0.
\end{equation}

\end{enumerate}

\end{Theorem}

We remark that for one value of
$\lambda \in (0, \infty)$, there might be more than one steady-state
solutions on the global bifurcation branch as shown in Figures 1.1  and 1.2.
It is easy to see that all steady-state solutions of the Kolmogorov problem
lie in both the following subsets of $\dbH^1_\sigma\times (0, \infty)$
and $\dbH^2_\sigma\times (0, \infty)$:
\begin{align}
 \mathcal C_1= &
 \left\{ (u, \lambda) \in \dbH^1_\sigma \times (0, \infty) \quad | \quad
      \| u\|_{H^1} \le c_2 k^2, \quad \lambda \ge \frac{c_1}{k^2} \right\} \\
    &  \cup \left\{ (0, \lambda) \quad |
           \quad 0 < \lambda \le \frac{c_1}{k^2} \right\},
                \nonumber \\
 \mathcal C_2= &
 \left\{ (u, \lambda) \in \dbH^2_\sigma \times (0, \infty) \quad | \quad
          \| u\|_{H^2} \le c_3 k^2 (\lambda^2 k^4 + 1),
\quad \lambda \ge \frac{c_1}{k^2} \right\} \\
   & \cup \left\{ (0, \lambda)\quad | \quad
                  0 < \lambda \le \frac{c_1}{k^2} \right\}.
              \nonumber
\end{align}
Here the absolute constants $c_1$, $c_2$  and $c_3$
are given in Lemma 3.2 in Section 3.
The set $\mathcal C_1$   is as shown in Figure 1.1.
Since stronger  $H^2$ norm is used in $\mathcal C_2$, the schematic picture
of  $\mathcal C_2$ shown in Figure 1.2
is similar to  $\mathcal C_1$ depicted in Figure 1.1
but with  the dotted horizontal lines replaced by quadratic line of
$\lambda$ given by
$ \| u\|_{H^2} = c_3 k^2 (\lambda^2 k^4 + 1).$

\begin{Theorem} Let  $(l,j,i,k)$ be an integer vector
satisfying (\ref{4-tuple}). Then all bifurcation solutions
$u_{l,j,i,k,n,\lambda}$ given by Theorem 1.1  
satisfying the following properties:

\begin{enumerate}
\item each branch $u_{l,j,i,k,n,\lambda}$ extends to $\lambda=\infty$ in
$\mathcal C_1 \cap  \mathcal C_2$ as shown in  Figures 1.1  and 1.2;

\item each branch $u_{l,j,i,k,n,\lambda}$ intersects with the $\lambda$--axis
only at the line segment $pq$. In particular, $u_{l,j,i,k,n,\lambda}$
never touches the $\lambda$--axis for $\lambda > \lambda_{l,j,i,k}$;

\item for $ j \neq 0$   and for $i_1=0$ or $i_0/2$,
\begin{equation}
u_{l,j,i,k,n,\lambda}\neq u_{l,j,i,k,n',\lambda},\label{1.4}
\end{equation}
if 
$$\lambda_{l,j,i,k}<\lambda<\infty,\,
i_1+1\leq n\leq i_1+\frac{i_0}{4}<n'\leq i_1+\frac{i_0}{2}.
$$

\end{enumerate}

\end{Theorem}


\bigskip

The above two main theorems are  the first ones showing
the existence of either 4 or 8 or 16 global branches of steady-state
solutions undergoing supercritical pitchfork bifurcations
from a single bifurcation point for
a Navier-Stokes problem.
The method we employ in the proof of these theorems is the
Rabinowitz global bifurcation theorem combined with
continuous fraction  method
first introduced by Meshalkin and Sinai \cite{MS}
and with  Fourier analysis.

The central gravity of the proof relies
on estimating the eigenvalues of the  linearized Navier-Stokes  operator
associated with   (\ref{ns1}-\ref{ns4}). One of the main
difficulties is that  the  eigenvalues of the  linearized
Navier-Stokes operator  always have
the multiplicity $2m$ for some integer $m$.
Thus  careful examination  is necessary for the bifurcation analysis.
To overcome this difficulty, we reduce the problem in some flow
invariant subspaces.   More specifically, we find that for  such an eigenvalue
there exist   exactly $2m$ flow invariant subspaces, in each of which
the eigenvalue is   simple in the algebraic sense.
Then the global bifurcation theorem introduced by
P. Rabinowitz \cite{Rab} can be used to complete the
bifurcation analysis.

Furthermore, restricting to  each flow invariant subspace, 
the Navier-Stokes equations undergo local supercritical bifurcation
around the bifurcation point. For the first time,  we obtain in Theorem 5.2 
the precise local asymptotic expansions of 
the bifurcation branches in terms of the Reynolds number $\lambda$ near the 
bifurcation point; see (\ref{expansion}). Technically speaking, 
the local asymptotic expansion of each bifurcation branch is obtained 
by studying the projection of Navier-Stokes equations to the 
unstable direction, and by carefully examining the nonlinear 
interaction of the corresponding unstable mode given by $u^*$ as in 
(\ref{u*})).

The two-dimensional (2D)  Kolmogorov problem 
is the 2D Navier-Stokes equations with 
the special zonal forcing $k^2(\sin ky,0)$.
In this case, the first study of this problem is due to
Meshalkin and Sinai \cite{MS} using the continuous fraction method.
They proved  the
absence of instability phenomena
of  the 2D Kolmogorov problem with $k=1$ regardless  the magnitude of
the Reynolds number.
This result was also confirmed in Marchioro \cite{Mar} with
an alternative approach.
The bifurcation analysis to the 2D Kolmogorov problem was
first examined by Iudovich \cite{Iu}
in a two-dimensional elongated domain
$[0, 2\pi/a] \times [0, 2\pi]$   with $0<a<1$, and 
the supercritical pitchfork bifurcation phenomenon 
with respect to some critical values $\lambda_{l,k}$ with $l=0,1,...,k-1$
was examined in \cite{CP1} based on the analysis and numerical computations.

Back to the three-dimensional case, the
governing linearized Navier-Stokes system may be reduced to a
two-dimensional
one via the Squire transformation.
Thus the three-dimensional linear instability phenomenon  may
be observed
from a two-dimensional problem; see  \cite{Lin} for details.
We note, however, that
the Squire transformation is not invertible, and thus is not
suitable for   bifurcation analysis.

It is not difficult to see that the steady-state  solutions
not branching off   the first bifurcation
point are unstable. One may think that unstable steady-state solutions
are not of physical interest, but they are
necessary in the  understanding of the transition to turbulence.
Moreover it is normally difficult to give a
numerical scheme ensuring the time-dependent discretized solution
converging  to  an unstable steady-state solution as $t\to \infty$;
see, for example,  \cite{Con}. In this article a bifurcated steady-state
solution is obtained  in a flow invariant space, in which the solution is
locally  stable at  least for  $\lambda $ close to the associated
bifurcation value.  Therefore one may derive  a convergent numerical
scheme with respect to every bifurcated
steady-state solution above.
In addition, noticing  (see \cite{CF,Tem})
that steady-state solutions to a
three-dimensional Navier-Stokes system are regular.
Thus the local stability of a steady-state solution  in
a flow invariant subspace
shows the  global existence
of time-dependent regular solutions  starting from
large initial data near the  steady-state solution,
although  the global existence of a  regular time-dependent  solution
to a  three-dimensional Navier-Stokes system remains to be
a fundamental open question, and only partial regularity
for a solution are available; see e.g. \cite{Nir, JLL69, PLL96, TX}.

\setcounter{equation}{0}
\section{Spectrum of the Linearized Problem}

From now on, we consider only the steady-state problem of the
Navier-Stokes equations (\ref{ns1}-\ref{ns4}), which can be rewritten
as
\begin{equation}\label{2.1}
-\Delta u+\lambda B(u,u)=k^2(\sin kz,0,0),    \,\,
u \in  \dbH^2_\sigma,
\end{equation}
where $B(u,u)=P[u \cdot \nabla u]$  is the bilinear term, and
$P: \dbH^2({\mathbb T}^3) \to \dbH^2_{\sigma}$ is the Leray
projection operator.


Linearizing  (\ref{2.1}) around the steady-state
$u_0=(\sin kz, 0, 0) $, we derive
$- \Delta u+ \lambda A u =0 $. Here the linearized operator $A$
is defined by
\begin{equation}
 A u= P[(u_0\cdot\nabla) u +
(u\cdot \nabla u_0)]=P\left[\sin kz \frac{\partial u}{\partial x} 
+ u_3 k \cos kz (1, 0, 0)\right],
\end{equation}
for any $u=(u_1, u_2, u_3)\in \dbH^2_\sigma$.

To obtain the steady-state
bifurcation result, we are interested in the
real eigenvalues of the linear  operator  in $\dbH^2_\sigma$. Namely we study
the nontrivial solutions of
\begin{equation}
 \Delta u -\lambda A u = \rho u. \label{2.2}
\end{equation}

For reasons which will be clear later we first introduce some
invariant subspaces of $\Delta -   \lambda A $.
For an integer vector $(l,j,i,k)$ with  $l,k \geq 1$, $0\leq i<k$ and
$j\in Z$, we define the following
subspaces   of $\dbH^2_\sigma$:
\begin{equation} \nonumber
\begin{split}
&  {\mathbb E}_{l,j,i,k} = \left\{u\in \dbH^2_\sigma \, |\,u= \sum_{n\in \mathbb Z}
        (\xi_n,\eta_n,\zeta_n)\sin (lx+j y+iz+nkz)\right\}, \\
& \tilde{{\mathbb E}}_{l,j,i,k}= \left\{u\in \dbH^2_\sigma  \,  |\,u=
 \sum_{n\in \mathbb Z} (\xi_n,\eta_n,\zeta_n)\cos (lx+j y+iz+nkz)\right\}.
\end{split}
\end{equation}
The spectrum of  this operator
in these subspaces is examined in the following theorem.

\begin{Theorem} Let $(l,j,i,k)$ be an integer vector subject to
the condition  (\ref{4-tuple}).
Then there exists  an eigenvalue
$\rho_{l,j,i,k}:\, (0,\infty)\rightarrow  {\mathbb R} $ of
(\ref{2.2}) such that for  any $\lambda >0$
$$\gathered
 \frac{d\rho_{l,j,i,k}(\lambda)}{d\lambda} >0,   \\
 \lim_{\lambda\rightarrow \infty}\rho_{l,j,i,k}(\lambda)=\infty,\\
 \lim_{\lambda\rightarrow 0^+}\rho_{l,j,i,k}(\lambda)=-l^2-j^2-\min\{i^2,
(k-i)^2\}.
\endgathered $$
Furthermore for any $\lambda >0$   and
 $\rho > -(l^2+j^2+\min\{i^2,(k-i)^2\})$,
\begin{equation}\nonumber
\begin{split}
& \mbox{\rm dim}\bigcup_{m\in {\mathbb N}}
\{u\in {\mathbb E}_{l,j,i,k}|\,\,(\Delta  -\lambda A -\rho)^m u=0\}\leq 1 ,\\
& \mbox{\rm dim} \bigcup_{m\in {\mathbb N}}
\{u\in \tilde{{\mathbb E}}_{l,j,i,k}|\,\,
(\Delta  -\lambda A -\rho)^m u=0\}
\leq 1, \\
& \mbox{\rm dim}\bigcup_{m\in {\mathbb N}}
\{u\in {\mathbb E}_{l,j,k-i,k}|\,\,
(\Delta -\lambda A-\rho )^mu=0\}\leq 1 \quad  \mbox{\rm if }i\neq 0, \\
& \mbox{\rm dim}\bigcup_{m\in {\mathbb N}}
\{u\in \tilde{{\mathbb E}}_{l,j,k-i,k}|\,\,
(\Delta -\lambda A -\rho)^m u=0\}
\leq 1 \quad \mbox{ if }i\neq 0,
\end{split}
\end{equation}
where the equalities hold if and only if $\rho=\rho_{l,j,i,k}(\lambda)$.
\end{Theorem}

The proof of Theorem 2.1 will be accomplished by converting the
eigenvalue problem $ \Delta u -\lambda A u = \rho u$ to
coupled continuous fraction equations. This approach was first used by
Meshalkin and Sinai \cite{MS} in a stability problem of fluid flows.

To proceed, we first recall a theorem on continuous fraction \cite{MS}.
Consider three term recurrence equations:
\begin{equation}
d_n e_n +e_{n-1} -e_{n+1}=0, \quad n \in {\mathbb Z},
\label{2.5}
\end{equation}
where $d_n$  and $e_n$ are complex numbers. Then we have

\begin{Theorem} Assume
\begin{equation}\label{2.6}
\gathered
 \text{Re} \, d_n > 0, for \,n\not= 0, \\
 \lim_{|n| \to \infty} \text{Re} \, d_n = \infty.
\endgathered
\end{equation}
Then the following two assertions are equivalent:
\begin{enumerate}
\item[(i)] There exists a nontrivial solution $\{ e_n\}_{n\in {\mathbb Z}}$ of
(\ref{2.5}) such that
\begin{equation}
\sum_{n\in {\mathbb Z}}|e_n|^2 < \infty. \label{2.7}
\end{equation}

\item[(ii)] The following equation holds true:
\begin{equation} \label{2.8}
d_0 + \displaystyle\frac{1}
{ \displaystyle d_{-1} +\displaystyle\frac{1}{\displaystyle d_{-2}+
\displaystyle\frac{1}{\ddots}}}
= \displaystyle\frac{- 1}
{ \displaystyle d_{1} +\displaystyle\frac{1}{\displaystyle d_{2}+
\displaystyle\frac{1}{\ddots}}}.
\end{equation}
\end{enumerate}
Furthermore, the solution obtained in (i) or (ii) is
unique up to a constant factor
and satisfies that
\begin{equation}
\text{ either $e_n = 0$  for all $n\in {\mathbb Z}$, or
 $e_n \not= 0$  for all  $n \in {\mathbb Z}$.}
\label{2.9}
\end{equation}
\end{Theorem}
\begin{proof}
1. Let $\{e_n\}_{n\in {\mathbb Z}}$ be a nontrivial solution of (\ref{2.5}).
If $e_{n_0}=0$ for some $n_0\geq 0$, multiplying the $n$-th equation of
(\ref{2.5}) by $\bar e_n$
and summing up the resultant equations for $n\geq n_0+1$
yields
$$\sum_{n\geq n_0+1}\mbox{Re} d_n|e_n|^2=0.$$
This together with (\ref{2.6}) gives $e_{n}=0$ for $n\geq n_0+1$, and
thus (\ref{2.5}) shows $e_n=0$ for all $n\in {\mathbb Z}$.
In the same way, the assumption
$e_{n_0}=0$ for some $n_0\leq -1$   implies
$e_n=0$ for $n\leq n_0-1$ and then $e_n=0$ for all $n\in {\mathbb Z}$.
Namely (\ref{2.9})  holds true.


2.  Dividing the $n$th and $-n$th equation of (\ref{2.5}) by
$e_n$ and $e_{-n}$ respectively yields that for any $n\ge 0$,
\begin{align}
&  \frac{e_{n}}{e_{n-1}}=\frac{-1}{d_n-\displaystyle\frac{e_{n+1}}{e_n} },
            \nonumber \\
& \frac{e_{-n}}{e_{-n+1}}=\frac{1}
{d_{-n}+\displaystyle\frac{e_{-n-1}}{e_{-n}}}. \nonumber
\end{align}
Namely, for any $n\ge 0$,
$$
\frac{e_{\pm n}}{e_{\pm (n-1)}}=\frac{\mp 1}
{d_{\pm n}+\displaystyle\frac{1}{d_{\pm(n+1)}
+\displaystyle\frac{1}{\ddots +\displaystyle\frac{1}{d_{\pm (n+m)}
\mp\displaystyle\frac{e_{\pm( n+m+1)}}{e_{\pm (n+m)}} }  }}  }.$$
This together with (\ref{2.6}-\ref{2.7}) implies that for any $n\ge 0$,
\begin{equation}\label{2.10}
\frac{e_{\pm n}}{e_{\pm (n-1)}}
=\gamma_{\pm n}\overset{ \text{def} }{=}
\frac{\mp 1}
{d_{\pm n}+\displaystyle\frac{1}{d_{\pm(n+1)}
+ \displaystyle\frac{1}{d_{\pm (n+2)}
+\displaystyle\frac{1}{\ddots}}}}.
\end{equation}
In particular, (\ref{2.8}) follows from
$(e_0/e_{-1})^{-1}=e_{-1}/e_0$.

3. Furthermore we infer from (\ref{2.10}) that for any $n\geq 1$,
\begin{equation}\nonumber
\begin{split}
&e_n=e_0\gamma_1\cdots\gamma_{n}, \\
 &e_{-n}=e_0\gamma_{-1}\cdots \gamma_{-n},
\end{split}
\end{equation}
which are uniquely determined up to constants, i.e. an arbitrary choice
of $e_0$.

4. Moreover, if (\ref{2.8}) is valid, we can also define
$\{\gamma_n\}_{n\in {\mathbb Z}}$  and a nontrivial
solution $\{e_n\}_{n\in {\mathbb Z}}$ as above.
\end{proof}

\begin{proof}[Proof  of Theorem 2.1] We divide the proof into several steps.

{\sl Step 1.}
We observe that the condition (\ref{4-tuple}) is equivalent to either
\begin{equation}\label{4-tuple1}
\gathered
 l^2+j^2+i^2 < k^2 , \\
 l^2+j^2+(k-i)^2 > k^2 , \\
 0\leq i < k/2, \\
 j \in {\mathbb Z}, \\
 0\leq l <  k,
\endgathered
\end{equation}
or
\begin{equation}\label{4-tuple2}
\gathered
 l^2+j^2+i^2 > k^2 , \\
 l^2+j^2+(k-i)^2 < k^2 , \\
 k/2 <  i < k, \\
 j \in {\mathbb Z}, \\
 0\leq l < k.
\endgathered
\end{equation}

We first consider the case where (\ref{4-tuple1}) holds true.

{\sl Step 2. An equivalent form of the non homogeneous problem
$(-\Delta +\rho) u +\lambda A u=(-\Delta +\rho)f $.}

For either  $u,\,f \in {\mathbb E}_{l,j,i,k}$ or
$u,\,f\in \tilde{{\mathbb E}}_{l,j,i,k}$, 
we write 
$$
u=\sum_{n\in {\mathbb Z}} (\xi_n,\eta_n,\zeta_n)\phi_n,\,\,\,
f=\sum_{n\in {\mathbb Z}} (a_n,b_n,c_n)\phi_n, $$
where 
$$\phi_n =\sin(lx+jy+iz+nkz)\,\,\mbox{ if  } \,\, u,\,f 
\in {\mathbb E}_{l,j,i,k},$$
$$\phi_n =\cos(lx+jy+iz+nkz)
    \,\,\mbox{ if } \,\,
u,\,f\in \tilde{{\mathbb E}}_{l,j,i,k}.$$
Here 
$$ (\xi_n,\eta_n,\zeta_n)=(\xi_n(l,j,i,k,\lambda),
\eta_n(l,j,i,k,\lambda),\zeta_n
(l,j,i,k,\lambda)),$$
$$(a_n,b_n,c_n)=(a_n(l,j,i,k,\lambda),
b_n(l,j,i,k,\lambda),c_n
(l,j,i,k,\lambda)).$$

Let 
\begin{equation} \label{beta-n}
\beta_n=\beta_n(l,j,i,k)=l^2+j^2+(nk+i)^2.
\end{equation}
Then direct computation yields that
$(-\Delta +\rho) u +\lambda A u=(-\Delta +\rho)f $
can be specified as follows:
\begin{gather}
 \sum_{n\in {\mathbb Z}}\!
\left((\beta_n+\rho)(\xi_n-a_n)
+\!\frac{\lambda l}{2}(\xi_{n-1}-\xi_{n+1})+\!\frac{\lambda k}{2}
(\zeta_{n-1}+\zeta_{n+1})\right)\phi_n = -\partial_{x}p, \label{2.11}\\
 \sum_{n\in {\mathbb Z}}\!
\left((\beta_n+\rho)(\eta_n -b_n)
+\!\frac{\lambda l}{2}(\eta_{n-1}-\eta_{n+1})\right)\phi_n
 = -\partial_{y}p,  \label{2.12}\\
\sum_{n\in {\mathbb Z}}\!
\left(
(\beta_n+\rho)(\zeta_n-c_n)
+\!\frac{\lambda l}{2}(\zeta_{n-1}-\zeta_{n+1})\right)\phi_n
 = -\partial_{z}p, \label{2.13}\\
 \sum_{n\in {\mathbb Z}}(l\xi_n+j\eta_n+(i+nk)\zeta_n)\phi_n = 0, \label{2.14} \\
\sum_{n\in {\mathbb Z}}(la_n+jb_n+(i+nk)c_n)\phi_n = 0. \label{n2.14}
\end{gather}

  It follows from  (\ref{2.11}-\ref{2.12})  that
\begin{eqnarray*}
\lefteqn{ \sum_{n\in {\mathbb Z}}\!
\left((\beta_n+\rho)(l\xi_n\!+\!j\eta_n-la_n-jb_n)
\!+\!\frac{\lambda l}{2}(l\xi_{n-1}\!+\!j\eta_{n-1}
-l\xi_{n+1}-j\eta_{n+1}) \right)\phi_n }\\
\lefteqn{+ \sum_{n\in {\mathbb Z}}\frac{\lambda lk}{2}
(\zeta_{n-1}\!+\!\zeta_{n+1})\phi_n }\\
&=&
-l\partial_{x}p-j\partial_{y}p. \hspace{8cm}
\end{eqnarray*}
Applying the operator $\partial_{z}$ to this equation and the operator
$-l\partial_{x}-j\partial_{y} $ to  (\ref{2.13}) respectively, and
summing the resultant equations, we have
\begin{eqnarray*}
\lefteqn{\left((\beta_n\!+\!\rho)(l\xi_n\!+\!j\eta_n-la_n-jb_n)
\!+\!\frac{\lambda
l}{2}(l\xi_{n-1}\!+\!j\eta_{n-1}\!-\!l\xi_{n+1}\!-\!j\eta_{n+1})
\right)(nk+i)   }\\
\lefteqn{ +\frac{\lambda lk}{2}
(\zeta_{n-1}\!+\!\zeta_{n+1})(nk+i) }\\
&=&
\left((\beta_n+\rho)(\zeta_n -c_n)
+\!\frac{\lambda l}{2}(\zeta_{n-1}\!-\!\zeta_{n+1})\right)
(l^2+j^2),\; n\in {\mathbb Z}. \hspace{2.5cm}
\end{eqnarray*}
Moreover, eliminating the pressure function $p$ 
from (\ref{2.11}-\ref{2.12}) yields
\begin{align*}
& \left((\beta_n+\rho)(\xi_n -a_n)
+\!\frac{\lambda l}{2}(\xi_{n-1}\!-\!\xi_{n+1})+\!\frac{\lambda k}{2}
(\zeta_{n-1}+\zeta_{n+1})\right)j\\
& =\left((\beta_n+\rho)(\eta_n -b_n)
+\!\frac{\lambda l}{2}(\eta_{n-1}\!-\!\eta_{n+1})\right)l,\,\,n\in {\mathbb Z}.
\end{align*}
Thus by (\ref{2.14}-\ref{n2.14}), we
see that the non homogeneous  problem
$$(-\Delta+\rho) u +\lambda A u =(-\Delta +\rho)f,\,\,\,\,
u,\,f \in {\mathbb E}_{l,j,i,k} \mbox{  or } u,\,f \in \tilde{{\mathbb E}}_{l,j,i,k}
$$
is
equivalent  to the following set of coupled algebraic equations
\begin{align}
& \frac{2\beta_n(\beta_n+\rho)}{\lambda l}\zeta_n
+(\beta_{n-1}-k^2)\zeta_{n-1}-(\beta_{n+1}-k^2)\zeta_{n+1}
\label{2.15}   \\ 
& \qquad = \frac{2\beta_n(\beta_n+\rho)}{\lambda l}c_n, \nonumber  \\
&  \frac{2(\beta_n+\rho)}{\lambda }(j\xi_n-l\eta_n)
+l(j\xi_{n-1}-l\eta_{n-1})-l(j\xi_{n+1}-l\eta_{n+1})  \label{2.16}\\
& \qquad   =-kj
(\zeta_{n-1}+\zeta_{n+1}) +
\frac{2(\beta_n+\rho)}{\lambda }(ja_n-lb_n), \nonumber \\
& l\xi_n+j\eta_n=-(i+nk)\zeta_n, \label{2.17} \\
& l a_n+jb_n=-(i+nk)c_n,   \label{n2.17}
\end{align}
for any $n\in {\mathbb Z}$.

{\sl Step 3.  Claim:  $\{\xi_n\}_{n\in {\mathbb Z}}$ and
$\{\eta_n\}_{n\in {\mathbb Z}}$ are uniquely   determined by
$\{ja_n -lb_n\}_{n\in {\mathbb Z}}$  and
$\{\zeta_n\}_{n\in {\mathbb Z}}$
when $\rho>-(l^2+j^2+i^2)$.
}

Let $\mathcal S$  be
\begin{equation}
\label{ss}
\mathcal S=\left\{\{\tau_n\}_{n\in {\mathbb Z}}  \, |\,
       \sum_{n\in {\mathbb Z}}|\tau_n|^2<\infty\right\},
       \end{equation}
and  $L:\mathcal S \to  \mathcal S$  be a linear operator defined by
$$L\{\tau_n\}_{n\in {\mathbb Z}} =\left\{\frac{\tau_{n-1}-\tau_{n+1}}
{\beta_n+\rho}\right\}_{n\in {\mathbb Z}}.$$
Then  (\ref{2.16}) can be rewritten as
$$
\left( \frac{2}{\lambda l }+L\right)
\{j\xi_n-l\eta_n\}_{n\in {\mathbb Z}}
=-\left\{\frac{kj (\zeta_{n-1}+\zeta_{n+1})}{(\beta_n +\rho)l}
         \right\}_{n\in {\mathbb Z}}+
         \left\{ \frac{2}{\lambda l} (ja_n-lb_n)\right\}_{n\in {\mathbb Z}}.
$$
It is easy to see  that  $L:\mathcal S \to  \mathcal S$
is a compact operator for $\rho>-\beta_0$.

We now  show
that $-2/\lambda l$ is not an eigenvalue of
$L$ or  the coupled algebraic equations
\begin{equation}
\label{2.18}
\frac{2(\beta_n+\rho)}{\lambda l}\tau_n
+\tau_{n-1}-\tau_{n+1}=0,\,\,n\in {\mathbb Z}
\end{equation}
has no nontrivial solution $\{\tau_n\}_{n\in {\mathbb Z}}  \in \mathcal S$.
Otherwise, let $\{ \tau_n\}_{n\in {\mathbb Z}}$  be a  nontrivial solution, then
by Theorem 2.2,
\begin{equation}\label{2.19}           \nonumber
-\displaystyle\frac{2(\beta_0+\rho)}{\lambda l}
= \displaystyle\frac{1}
{\displaystyle\frac{2(\beta_1+\rho)}{\lambda l}
+\displaystyle\frac{1}{\displaystyle\frac{2(\beta_{2}+\rho)}{\lambda l}+
\displaystyle\frac{1}{\ddots}}}
+
\displaystyle\frac{1}
{\displaystyle\frac{2(\beta_{-1}+\rho)}{\lambda l}
+\displaystyle\frac{1}{\displaystyle\frac{2(\beta_{-2}+\rho)}{\lambda l}+
\displaystyle\frac{1}{\ddots}}}.
\end{equation}
Since $-(\beta_0+\rho)/(\lambda l)<0$
while  the right-hand side of this equation is positive,
this leads to a contradiction.
Hence  (\ref{2.18}) has no nontrivial solution.


By the Riesz--Schauder theory,
it is easy to see that  $(2/\lambda l+L)^{-1}:\mathcal S \to \mathcal S$
is a bounded  operator, and
$$\{j\xi_n-l\eta_n\}_{n\in {\mathbb Z}}=- \big(\frac{2}{\lambda l}+
L\big)^{-1}
\big(\big\{\frac{kj(\zeta_{n-1}
+\zeta_{n+1})}{(\beta_n +\rho)l}\big\}_{n\in {\mathbb Z}}+
\big\{\frac{2}{\lambda l}(ja_n-lb_n)\big\}_{n\in bZ}\big).
$$
This together with  (\ref{2.17}-\ref{n2.17}) implies the desired assertion.




{\sl Step 4. Existence and uniqueness of
$\rho=\rho_{l,j,i,k}(\lambda)$ for the eigenvalue problem
\begin{equation}       \label{eig}
\Delta u -\lambda A u =\rho u,\,\,\, u\in {\mathbb E}_{l,j,i,k}\cup
\tilde{{\mathbb E}}_{l,j,i,k}.
\end{equation}
}

It follows from
(\ref{4-tuple1})  and (\ref{beta-n}) that
$\beta_0-k^2<0$ and $\beta_n-k^2>0$ for $n\neq 0$.
By Steps 2 and 3, (\ref{eig}) is uniquely determined by
(\ref{2.15}) with $c_n=0$ or set of equations
\begin{equation}
\frac{2\beta_n(\beta_n+\rho)}{\lambda l}\zeta_n
+(\beta_{n-1}-k^2)\zeta_{n-1}-(\beta_{n+1}-k^2)\zeta_{n+1}
=0,\label{3-term}
\end{equation}
which is in the form of (\ref{2.5}) with $\tau_n =(\beta_n-k^2)\zeta_n$.
By Theorem 2.2, it becomes to show the existence and uniqueness of
$\rho=\rho_{l,j,i,k}$ in the following algebraic equation
\begin{eqnarray}\label{2.20}
-\displaystyle\frac{2\beta_0(\beta_0+\rho)}{\lambda l(\beta_0-k^2)}
&=&
\displaystyle\frac{1}
{\displaystyle\frac{2\beta_1(\beta_1+\rho)}{\lambda l(\beta_1-k^2)}
+\displaystyle\frac{1}{\displaystyle\frac{2\beta_2(\beta_{2}+\rho)}
{\lambda l(\beta_2-k^2)}+
\displaystyle\frac{1}{\ddots}}} \\
&&+
\displaystyle\frac{1}
{\displaystyle\frac{2\beta_{-1}(\beta_{-1}+\rho)}{\lambda
l(\beta_{-1}-k^2)}
+\displaystyle\frac{1}{\displaystyle\frac{2\beta_{-2}
(\beta_{-2}+\rho)}{\lambda
l
(\beta_{-2}-k^2)}+
\displaystyle\frac{1}{\ddots}}}. \nonumber
\end{eqnarray}
Multiplying  this equation  by $-l(\beta_0-k^2)(\beta_0(\beta_0+\rho))^{-1}$
gives
\begin{equation}\label{2.22}
\frac{2}{\lambda }=
\displaystyle\frac{1}{\displaystyle\frac{g_1(\rho)}{\lambda}
+ \displaystyle\frac{1}{\displaystyle\frac{g_2(\rho)}{\lambda}
+ \displaystyle\frac{1}{\ddots}}}
+\displaystyle\frac{1}{\displaystyle\frac{g_{-1}(\rho)}{\lambda}
+ \displaystyle\frac{1}{\displaystyle\frac{g_{-2}(\rho)}{\lambda}
+ \displaystyle\frac{1}{\ddots}}},
\end{equation}
where
\begin{equation} \nonumber
g_{\pm n}(\rho)=
\left\{
    \begin{split}
       & \frac{2\beta_{\pm n}\beta_0(\beta_{\pm n}+\rho)(\beta_0+\rho)}
  {l^2(k^2-\beta_0) (\beta_{\pm n}-k^2)},\,\,\mbox{  if  $n$ is odd},\\
       & \frac{2\beta_{\pm n}(\beta_{\pm n}+\rho)(k^2-\beta_0)}
{\beta_0(\beta_{\pm n}-k^2)(\beta_{0}+\rho)},\,\,\mbox{ if $n$ is even}.
       \end{split}  \right.
\end{equation}
Let $G(\lambda,\rho)$ be the right-hand side of  (\ref{2.22}), and
for $n\geq 1$ we set
$$G_{\pm n}(\lambda,\rho)=
\displaystyle\frac{1}{\displaystyle\frac{g_{\pm 1}(\rho)}{\lambda}
+ \displaystyle\frac{1}{\displaystyle\frac{g_{\pm 2}(\rho)}{\lambda}
+ \displaystyle\frac{1}{\ddots+
\displaystyle\frac{1}{\displaystyle\frac{g_{\pm n}(\rho)}{\lambda}}}}}.$$
We have
$$G(\lambda, \rho)=\lim_{n\to\infty}G_{n}(\lambda,\rho)+
\lim_{n\to\infty}G_{-n}(\lambda,\rho) .
$$
For any $\lambda >0$, $\epsilon >0$  and $M > 0$ arbitrarily large,
we see that
$\{G_{n}(\lambda, \cdot)\}_{n\geq 1}$ and
$\{G_{-n}(\lambda, \cdot)\}_{n\geq 1}$ are Cauchy sequences
of the Banach space $C^1([-\beta_0+\epsilon, M])$, and thus
$G(\lambda,\cdot)    \in C^1([-\beta+\epsilon, M])$,
\begin{eqnarray*}
\frac{\partial G(\lambda,\rho)}{\partial \rho}
&=&\lim_{n\to \infty}\frac{\partial G_{n}(\lambda,\rho)}{\partial \rho}
       +
\lim_{n\to \infty}\frac{\partial G_{-n}(\lambda,\rho)}
{\partial \rho}  \\
&=& \sum_{n=1}^\infty  \frac{(-1)^n}{\lambda}
      \frac{\mbox{d}g_{n}(\rho)}{\mbox{d}\rho}
\hat{G}_1^2(\lambda,\rho)  \cdots
 \hat{G}_n^2(\lambda,\rho)  \\
&  &  +
\sum_{n=1}^\infty    \frac{(-1)^n}{\lambda}
   \frac{\mbox{d}g_{-n}(\rho)}{\mbox{d}\rho}
 \hat{G}_{-1}^2(\lambda,\rho)  \cdots
 \hat{G}_{-n}^2(\lambda,\rho),
\end{eqnarray*}
with
$$\hat{G}_{\pm n}(\lambda,\rho)=
\displaystyle\frac{1}{\displaystyle\frac{g_{\pm n}(\rho)}{\lambda}
+ \displaystyle\frac{1}{\displaystyle\frac{g_{\pm (n+1)}(\rho)}{\lambda}
+ \displaystyle\frac{1}{\ddots}}}.
$$
By direct calculation, we infer from (\ref{4-tuple1})
that  for any $n \ge 1$,
${\mbox{d}g_{\pm n}(\rho)}/{\mbox{d}\rho} >0$  when  $n$ is odd,
and ${\mbox{d}g_{\pm n} (\rho)}/ {\mbox{d}\rho}<0$  when  $n$ is even.
Therefore when $\rho >- \beta_0,$
we obtain
\begin{equation}\label{2.23}
\frac{\partial G(\lambda,\rho)}{\partial\rho}  <0.
\end{equation}
Hence the observation
$$\lim_{\rho\searrow -\beta_0}G(\lambda,\rho)= \infty\mbox{ and }
\lim_{\rho\rightarrow \infty}G(\lambda,\rho)= 0$$
implies the existence and uniqueness  of
$\rho=\rho_{l,j,i,k}(\lambda)>-
\beta_0$
satisfying
\begin{equation}
\frac{2}{\lambda }=G(\lambda,
\rho_{l,j,i,k}(\lambda)).\label{2.24}
\end{equation}


{\sl Step 5.}  Notice that
$$
\lambda G(\lambda,\rho)=
\displaystyle\frac{1}{
\displaystyle\frac{g_1(\rho)}{\lambda^2}
+ \displaystyle\frac{1}{g_2(\rho)
+\displaystyle\frac{1}{\displaystyle\frac{g_3(\rho)}{\lambda^2}
+ \displaystyle\frac{1}{\ddots}}} }
+\displaystyle\frac{1}{
\displaystyle\frac{g_{-1}(\rho)}{\lambda^2}
+ \displaystyle\frac{1}{g_{-2}(\rho)
+\displaystyle\frac{1}{\displaystyle\frac{g_{-3}(\rho)}{\lambda^2}
+ \displaystyle\frac{1}{\ddots}}} }.
$$
Arguing as in the derivation of (\ref{2.23}), we have
$$\frac{\partial (\lambda G(\lambda,\rho))}{\partial \lambda }
>0\,\ \mbox{ for } \lambda >0 \,\mbox{ and }   \rho>-\beta_0.$$
 (\ref{2.24}) gives for
$\rho=\rho_{l,j,i,k}(\lambda)$,
\begin{eqnarray*}
0&=&\frac{\mbox{d}(\lambda G(\lambda,\rho))}
{\mbox{d}\lambda}\\
&=&\frac{\partial (\lambda G(\lambda,\rho))}{\partial\lambda}+
\lambda\frac{\partial G(\lambda,\rho)}{\partial\rho}
\frac{\mbox{d}\rho}{\mbox{d}\lambda}  \\
&>&\lambda\frac{\partial G(\lambda,\rho)}
{\partial\rho}\frac{\mbox{d}\rho}{\mbox{d}\lambda}.
\end{eqnarray*}
This shows, by  (\ref{2.23}), that
$\mbox{d}\rho_{l,j,i,k}/\mbox{d}\lambda>0.$

Furthermore, multiplying  (\ref{2.20}) by $\lambda$ gives
\begin{eqnarray}
\lefteqn{ -\frac{2\beta_0(\beta_0+\rho_{l,j,i,k}(\lambda))}
{ l(\beta_0-k^2)} }\\
&=&
\displaystyle\frac{1}
{\displaystyle\frac{2\beta_1(\beta_1+\rho_{l,j,i,k}(\lambda))}
{\lambda^2l(\beta_1-k^2)}
+\displaystyle\frac{1}{\displaystyle\frac{2\beta_2(\beta_{2}+
\rho_{l,j,i,k}(\lambda))}
{l(\beta_2-k^2)}+
\displaystyle\frac{1}{\ddots}}}\nonumber \\
&&+
\displaystyle\frac{1}
{\displaystyle\frac{2\beta_{-1}(\beta_{-1}+\rho_{l,j,i,k}(\lambda))}
{\lambda^2 l(\beta_{-1}-k^2)}
+\displaystyle\frac{1}{\displaystyle\frac{2\beta_{-2}(\beta_{-2}
+\rho_{l,j,i,k}(\lambda))}{l
(\beta_{-2}-k^2)}+
\displaystyle\frac{1}{\ddots}}}. \nonumber
\end{eqnarray}
Passing to the limits  $\lambda\rightarrow 0^+$ and
$\lambda\rightarrow \infty$
respectively, we obtain immediately
\begin{equation}
\nonumber
\begin{split}
 &  \lim_{\lambda\rightarrow 0^+}\rho_{l,j,i,k}(\lambda)=-l^2-j^2-i^2, \\
 & \lim_{\lambda\rightarrow \infty}\rho_{l,j,i,k}(\lambda)=\infty.
\end{split}
\end{equation}

In addition, it follows from
(\ref{2.15}) with $\rho=\rho_{l,j,i,k}(\lambda)$ and Theorem 2.2  that
$$\gamma_{\pm n} =\frac{(\beta_{\pm n}-k^2)\zeta_{\pm n}}
{(\beta_{\pm (n-1)}-k^2)\zeta_{\pm(n-1)}}
\,\, \mbox{ for } n\in N$$
with  $\zeta_m$ subject to the following condition
$$
\begin{array}{lcll}
\zeta_0&=&c,&\vspace{2mm}\\
\zeta_{\pm n}&=&c\displaystyle\frac{(\beta_0-k^2)
\gamma_{\pm 1}\cdots\gamma_{\pm n}}{\beta_{\pm n}-k^2},
&n\in N.
\end{array}
$$
Here $c$ is an arbitrary constant. Notice that
$\gamma_{\pm n}$ are uniquely determined by
$\rho_{l, i, j, k}(\lambda)$  and $\beta_0, \beta_{\pm 1}, \cdots$.
Thus all the eigenfunctions of the
spectral problem
$\Delta u -\lambda A u=\rho_{l,j,i,k}u$ form
a one-dimensional subspace in ${\mathbb E}_{l,j,i,k}$ and
$\tilde{{\mathbb E}}_{l,j,i,k}$
respectively.

{\sl Step 6. Claim: ker$(\Delta-\lambda A-\rho)^m
=\mbox{ker}(\Delta-\lambda A-\rho)$ $(m>1,\,\, \rho >-\beta_0).$}

It suffices to consider the operator with
 $\rho=\rho_{l,j,i,k}$ in the space ${\mathbb E}_{l,j,i,k}$.
We start with the case $m=2$. Let 
$$u=\sum_{n\in {\mathbb Z}}(\xi_n,\eta_n,\zeta_n)
\sin(lx+jy+iz+nkz)\in {\mathbb E}_{l,j,i,k}$$
such that
$$(\Delta-\lambda A -\rho)^2u=0,$$ 
or equivalently, 
$$(1 -\lambda (\Delta-\rho)^{-1}A )^2u=0.$$
This implies the existence of
$$u'=\sum_{n\in {\mathbb Z}}(\xi'_n,\eta'_n,\zeta'_n)
\sin(lx+jy+iz+nkz)\in {\mathbb E}_{l,j,i,k}$$
 such that
\begin{equation}\label{ge}
\gathered 
 (-\Delta +\rho)u +\lambda A u = (-\Delta +\rho) u', \\
 (-\Delta +\rho)u'+\lambda Au' = 0.
 \endgathered 
\end{equation}
Thanks to (\ref{2.15}-\ref{n2.17}), the above system (\ref{ge}) 
is equivalent to 
the following set of coupled algebraic equations:
\begin{align}
&  \frac{2\beta_n(\beta_n+\rho)}{\lambda l}\zeta_n
+(\beta_{n-1}-k^2)\zeta_{n-1}-(\beta_{n+1}-k^2)\zeta_{n+1} \label{n2.15} \\
& \qquad = \frac{2\beta_n(\beta_n+\rho)}{\lambda l}\zeta'_n,\nonumber \\
& \frac{2(\beta_n+\rho)}{\lambda }(j\xi_n-l\eta_n)
+l(j\xi_{n-1}-l\eta_{n-1})-l(j\xi_{n+1}-l\eta_{n+1}) \\
& \qquad =-kj(\zeta_{n-1}+\zeta_{n+1})
+\frac{2(\beta_n+\rho)}{\lambda }(j\xi'_n-l\eta'_n) , \nonumber \\
& l\xi_n+j\eta_n=-(i+nk)\zeta_n, \\
& \frac{2\beta_n(\beta_n+\rho)}{\lambda l}\zeta'_n
+(\beta_{n-1}-k^2)\zeta'_{n-1}-(\beta_{n+1}-k^2)\zeta'_{n+1}
= 0,\label{nn2.15} \\
&  \label{nn2.16}
\frac{2(\beta_n+\rho)}{\lambda }(j\xi'_n-l\eta'_n)
+l(j\xi'_{n-1}-l\eta'_{n-1})-l(j\xi'_{n+1}-l\eta'_{n+1})
\\
& \qquad =-kj(\zeta'_{n-1}+\zeta'_{n+1}), \nonumber  \\
&  \label{nn2.17}
l\xi_n'+j\eta_n'=-(i+nk)\zeta_n',
\end{align}
for any $n\in {\mathbb Z}$.


Thus it amounts to showing that $\zeta'_n=\eta'_n=\xi'_n =0$ for any 
$n\in {\mathbb Z}$.
Indeed, (\ref{n2.15}) gives the following equation in
$\mathcal S$ defined by (\ref{ss})
   \begin{equation}\label{ttt}
   (1+T)\{(\beta_n-k^2)\zeta_n\}_{n\in {\mathbb Z}}
   =\{(\beta_n-k^2)\zeta'_n\}_{n\in {\mathbb Z}},
  \end{equation}
with the operator $T$ in the form
$$ T\{\tau _n\}_{n\in {\mathbb Z}} =\left\{
\frac{\lambda l(\beta_n-k^2)
}{2\beta_n(\beta_n+\rho)}(\tau_{n-1}-\tau_{n+1})\right\}_
{n\in {\mathbb Z}}.$$
Since $T: \mathcal S\mapsto \mathcal S$ is compact, by the Riesz-Schauder
theory, (\ref{ttt}) is solvable if and only if
$$  \sum_{n\in {\mathbb Z}}(\beta_n-k^2)\zeta'_n \tau_n=0\mbox{ for  all }
 \{\tau_n\}_{n\in{\mathbb Z}}\in \mathcal S \mbox{ with }
 (1+T^*)\{\tau_n\}_{n\in {\mathbb Z}}=0,
 $$
 where $T^*$ denotes the dual operator of $T$.

 To specify the kernel of $1+T^*$, we note that
 $$T^*\{\tau _n\}_{n\in {\mathbb Z}} =\left\{
-\frac{\lambda l(\beta_{n-1}-k^2)}
{2\beta_{n-1}(\beta_{n-1}+\rho)}\tau_{n-1}
+\frac{\lambda l(\beta_{n+1}-k^2)}
{2\beta_{n+1}(\beta_{n+1}+\rho)}\tau_{n+1}
\right\}_
{n\in {\mathbb Z}}.$$
 Thus $(1+T^*)\{\tau_n\}_{n\in {\mathbb Z}}=0$ becomes, for $n\in {\mathbb Z}$,
 $$\frac{2\beta_n(\beta_n+\rho)}{\lambda l(\beta_n-k^2)}
 \sigma_n +\sigma_{n-1}-\sigma_{n+1}=0,\,\,\,\,
 \sigma_m=
(-1)^m\frac{\lambda l(\beta_{m}-k^2)}
{2\beta_{m}(\beta_{m}+\rho)}\tau_{m}.
$$
or $\{\sigma_n\}_{n\in{\mathbb Z}}\in \mbox{ker}(1+T)$.
If $\{\zeta'_n\}_{n\in{\mathbb Z}}\neq 0$,
applying Theorem 2.2 to (\ref{nn2.15})  gives
$\zeta'_n\neq 0$ for all $n\in {\mathbb Z}$ and
$$\{\sigma_n\}_{n\in {\mathbb Z}}=
\left\{ (-1)^n\frac{\lambda l(\beta_{n}-k^2)}
{2\beta_{n}(\beta_{n}+\rho)}\tau_{n}\right\}_{n\in {\mathbb Z}}
=c\{(\beta_n-k^2)\zeta'_n\}_{n\in{\mathbb Z}}
$$
for some constant $c\neq 0$, whenever $\{\tau_n\}_{n\in {\mathbb Z}}\neq 0$.
This implies
\begin{equation}\label{ab}
\sum_{n\in{\mathbb Z}}(\beta_n-k^2)\zeta'_n\tau_n
=c\sum_{n\in{\mathbb Z}}(-1)^n
\frac{2\beta_{n}(\beta_{n}+\rho)(\beta_n-k^2)}{\lambda l}{\zeta'}^2_n.
\end{equation}
On the other hand, multiplying by $(\beta_n-k^2)\zeta'_n$
the $n$th equation    of (\ref{nn2.15})
  and summing up the resultant equations yield
$$\sum_{n\in{\mathbb Z}}
\frac{2\beta_{n}(\beta_{n}+\rho)(\beta_n-k^2)}{\lambda l}{\zeta'}^2_n=0.  $$
        This together with (\ref{ab}) implies
$$\sum_{n\in{\mathbb Z}}(\beta_n-k^2)\zeta'_n\tau_n
=-c\sum_{n\in{\mathbb Z}}
\frac{2\beta_{2n+1}(\beta_{2n+1}+\rho)(\beta_{2n+1}-k^2)
}{\lambda l}{\zeta'}^2_{2n+1} \neq 0.
$$
This leads to a contradiction, and thus
 $\{\zeta'_n\}_{n\in{\mathbb Z}}=0$. From  (\ref{nn2.16}-\ref{nn2.17})
 and  Step 3, we have
 $\{\eta'_n\}_{n\in{\mathbb Z}}=0$ and $\{\xi'_n\}_{n\in{\mathbb Z}}=0$.
Hence we have $u'=0$ or
$$\mbox{ker}(\Delta-\lambda A -\rho)^2
=\mbox{ker}(\Delta-\lambda A-\rho).$$

To reach the case where $m>2$,
suppose $(1+\lambda(-\Delta+\rho)^{-1} A)^mu=0$ for some $u\in {\mathbb E}_{l,j,i,k}$.
Then there exits a $u'\in {\mathbb E}_{l,j,i,k}$ so that
\begin{eqnarray*}
(-\Delta +\rho)v +\lambda A v&=&(-\Delta +\rho)u',\\
(\Delta -\lambda A -\rho)u'&=&0,       \\
(1+\lambda (-\Delta+\rho)^{-1}A )^{m-2}u&=&v.
\end{eqnarray*}
It follows from the argument for the case  $m=2$ that $u'=0$.
Hence, we obtain the case for $m>2$ by induction.

{\sl Step 7.}  We now consider the case where (\ref{4-tuple2}) holds true.

Obviously, it follows from the Steps 2 and 3 that
spectral problem
$\Delta u -\lambda A u=\rho u$ with
$u\in {\mathbb E}_{l,j,i,k}\cup\tilde{{\mathbb E}}_{l,j,i,k}$
is uniquely determined by  the system
$$
\frac{2\beta_{-n-1}(\beta_{-n-1}+\rho)}{\lambda l}\zeta_{-n-1}
+(\beta_{-n-2}-k^2)\zeta_{-n-2}-(\beta_{-n}-k^2)\zeta_{-n}
=0, \,\,\,\mbox{ for } n\in {\mathbb Z}$$
with $\beta_n=\beta_n(l,j,i,k)$.
Setting $\beta'_n=\beta_{n}(l,j,k-i,k)$ and observing that \\
 $\beta_{-n-1}(l,j,i,k)=\beta_{n}(l,j,k-i,k)$, we see that
$$
\frac{2\beta'_{n}(\beta'_{n}+\rho)}{\lambda l}\zeta_{-n-1}
+(\beta'_{n+1}-k^2)\zeta_{-n-2}-(\beta'_{n-1}-k^2)\zeta_{-n}
=0, \,\,\,\mbox{ for } n\in {\mathbb Z}.$$
By setting $\zeta'_n =(-1)^n\zeta_{-n-1}$, we have
$$\frac{2\beta'_{n}(\beta'_{n}+\rho)}{\lambda l}\zeta'_{n}
+(\beta'_{n-1}-k^2)\zeta'_{n-1}-(\beta'_{n+1}-k^2)\zeta'_{n+1}
=0, \,\,\,\mbox{ for } n\in {\mathbb Z}.$$
Thus the above argument  implies the unique existence of the eigenvalue
$\rho_{l,j,i,k}(\lambda)$ so that
$\Delta u -\lambda A u=\rho_{l,j,i,k}u$ with $k/2<i<k$ also  form
a one-dimensional subspace in ${\mathbb E}_{l,j,i,k}$ and
$\tilde{{\mathbb E}}_{l,j,i,k}$
respectively.
In particular,
$$\lim_{\lambda \to 0^{+}}\rho_{l,j,i,k}(\lambda)=-
\beta'_0=-(l^2+j^2+(k-i)^2).$$

Likewise, we obtain the assertion on the spaces $ {\mathbb E}_{l,j,k-i,k}$ and
$\tilde{{\mathbb E}}_{l,j,k-i,k}$.
Obviously, we have $\rho_{l,j,i,k}(\lambda)=\rho_{l,j,k-i,k}(\lambda)$.

The proof of Theorem 2.1 is complete.
\end{proof}


For reader's convenience, we state an elementary lemma.

\begin{Lemma}
Let  $(l,j,i,k)$ be an integer vector such that  $l,k\geq 1$, $0\leq i <k$
and
$j\in {\mathbb Z}$. Then
$$l^2+j^2+i^2\geq k^2\quad \mbox{\rm   implies } \quad
n^2(l^2+j^2)+\{ni\}^2> k^2, \quad \forall n\geq 2,$$
and
$$l^2+j^2+(k-i)^2\geq k^2
\quad \mbox{\rm    implies }\quad
n^2(l^2+j^2)+(k-\{ni\})^2> k^2,
\quad \forall n\geq 2.$$
\end{Lemma}


\begin{proof}
    Notice that the desired conclusion holds true if
    $n^2(l^2+j^2)>k^2$. It suffices to consider the integers
     $n$ subject to
    the condition
    $  n^2(l^2+j^2)\leq k^2$.

In the case where  $l^2+j^2+i^2\geq k^2$,  we have
$$i^2\geq k^2-l^2-j^2\geq \frac{(n^2-1)k^2}{n^2}>
\frac{(n-1)^2k^2}{n^2},$$
which yields
$$ nk> ni =(n-1)k + \{ni\} > (n-1)k.$$
Namely
$\{ni\}=ni-(n-1)k>0$. Hence
\begin{eqnarray*}
n^2(l^2+j^2)+\{ni\}^2&=&n^2(l^2+j^2)+(ni-(n-1)k)^2\\
&=& n^2(l^2+j^2+i^2)+n^2k^2-2nk^2-2n(n-1)ki+k^2\\
&\geq & 2n(n-1)k(k-i)+k^2>k^2.
\end{eqnarray*}

For the second case where  $l^2+j^2+(k-i)^2\geq k^2$, we have
$$(k-i)^2\geq k^2-l^2-j^2\geq
\frac{(n^2-1)k^2}{n^2}>\frac{(n-1)^2k^2}{n^2},$$
which shows  $\{ni\}=ni <k$,
and the second part of the lemma can be proved in the same fashion.
\end{proof}

We are now in  position to state the main result of this section.

\begin{Theorem}
Let an  integer vector $(l,j,i,k)$   satisfy  (\ref{4-tuple}),
$n_0$  be given by (\ref{n-zero}), and $\{ni\}$  be defined as (\ref{mod-k}).
Then the following assertions hold true:

(i). For $\lambda>0$, $n\in {\mathbb N}$  with $n>n_0$
and $\rho >-(l^2+j^2+\min\{i^2,(k-i)^2\})$,
\begin{equation} \nonumber
\begin{split}
&  \mbox{\rm dim}\{u\in {\mathbb E}_{nl,nj,\{ni\},k}
\cup \tilde{{\mathbb E}}_{nl,nj,\{ni\},k}
|\,\,\Delta u -\lambda A u-\rho u=0\}=0,\\
& \mbox{\rm dim}\{u\in {\mathbb E}_{nl,nj,k- \{ni\},k}
\cup\tilde{{\mathbb E}}_{nl,nj,k- \{ni\},k}
|\,\,\Delta u -\lambda A u-\rho u=0\}=0, \text{ if } \{ni\}\neq 0.
\end{split}
\end{equation}

(ii). There exist $n_0 $ eigenvalues
$\rho_{nl,nj,\{ni\},k}:(0,\infty)\rightarrow {\mathbb R}$ for
$n=1,...,n_0$ such that
\begin{equation}\nonumber
\begin{split}
& \frac{d\rho_{nl,nj,\{ni\},k}(\lambda)}{d\lambda} >0, \\
& \lim_{\lambda\rightarrow \infty}\rho_{nl,nj,\{ni\},k}
(\lambda)=\infty, \\
& \lim_{\lambda\rightarrow 0^+}\rho_{nl,nj,\{ni\},k}(\lambda)=
-n^2(l^2+j^2)-\min\{\{ni\}^2,(k-\{ni\})^2\},
\end{split}
\end{equation}
and
\begin{equation}\nonumber
\begin{split}
& \mbox{\rm dim}\bigcup_{m\in {\mathbb N}}
\{u\in {\mathbb E}_{nl,nj,\{ni\},k}
|\,\,(\Delta  -\lambda A -\rho)^m u=0\}\leq 1, \\
& \mbox{\rm dim}\bigcup_{m\in{\mathbb N}}
\{u\in \tilde{{\mathbb E}}_{nl,nj,\{ni\},k}|\,\,(\Delta
-\lambda A -
\rho )^mu=0\}\leq 1, \\
&  \mbox{\rm dim}\bigcup_{m\in{\mathbb N}}
\{u\in {\mathbb E}_{nl,nj,k-\{ni\},k}
|\,\,(\Delta  -\lambda A -\rho)^m u=0\}\leq 1 \,\mbox{\rm  if } \,\{ni\}
\neq 0, \\
&
\mbox{\rm dim}\bigcup_{m\in{\mathbb N}}
\{u\in \tilde{{\mathbb E}}_{nl,nj,k-\{ni\},k}|\,\,
(\Delta -\lambda A - \rho )^mu=0\}\leq 1  \,\,\mbox{\rm  if } \,\{ni\}\neq 0,
\end{split}
\end{equation}
where the equalities hold if and only if
$\rho=\rho_{nl,nj,\{ni\},k}(\lambda)$ for $\lambda >0$.
\end{Theorem}

\begin{proof}

(i). Let  $n\geq  n_0+1$.
Without loss of generality, we may assume that
$\{ni\}<k-\{ni\}$.
By the definition (\ref{n-zero}) of $n_0$, we observe that
$$n^2(l^2+j^2)+(k-\{ni\})^2-k^2 > n^2(l^2+ j^2)+\{ni\}^2-k^2 \ge 0. $$

As in  the proof of Theorem 2.1, for
$u\in {\mathbb E}_{nl,nj,\{ni\},k}\cup\tilde{{\mathbb E}}_{nl,nj,\{ni\},k},$
the spectral problem $\Delta u -\lambda A u-\rho u=0 $,
is uniquely determined  by  the following analog of (\ref{2.15}):
\begin{equation}
\label{2.26}
\frac{2\beta_m(\beta_m+\rho)}{\lambda n l}\zeta_m
+(\beta_{m-1}-k^2)\zeta_{m-1}-(\beta_{m+1}-k^2)\zeta_{m+1}
=0,\,\,m\in {\mathbb Z},
\end{equation}
where  $\beta_m$ now equals
$\beta_m(nl,nj,\{ni\},k)=n^2(l^2+j^2)+(mk+\{ni\})^2$.

Multiplying the $m$--th equation of (\ref{2.26})
by $(\beta_m-k^2)\zeta_m$
and summing the resultant equations yield
$$
\sum_{m\in {\mathbb Z}}
\frac{2\beta_m(\beta_m+\rho)(\beta_m-k^2)}{\lambda n l}|\zeta_m|^2 =0 .
$$
Notice  that $\beta_m-k^2>0$ for $m\neq 0$.
Thus  $\zeta_m\equiv 0$ for $m\neq 0$, and $\zeta_0=0$
as well in view of  (\ref{2.26}).

In the case where $\{ni\}\neq 0$   and
$u\in{\mathbb E}_{nl,nj,k-\{ni\},k}\cup\tilde{{\mathbb E}}_{nl,nj,k-\{ni\},k},$
the spectral problem $\Delta u -\lambda A u-\rho u=0 $
is uniquely determined  by
\begin{equation}
\nonumber
\frac{2\beta_m(\beta_m+\rho)}{\lambda n l}\zeta_m
+(\beta_{m-1}-k^2)\zeta_{m-1}-(\beta_{m+1}-k^2)\zeta_{m+1}
=0,\,\,m\in {\mathbb Z},
\end{equation}
with $\beta_m(nl,nj,k-\{ni\},k)=n^2(l^2+j^2)+(mk+k-\{ni\})^2$.
Then as before, it is easy to see that $\zeta_m\equiv 0$.

(ii).
For $ 1\leq n\leq n_0$,
it follows from  Lemma 2.1 and the definition of the integer
 $n_0$  that  $(nl,nj,\{ni\},k)$ satisfies
 (\ref{4-tuple}).
Thus Assertion (ii) is a direct consequence of  Theorem 2.1.
\end{proof}

\setcounter{equation}{0}
\section{Invariant Subspaces and Properties of the Steady-State Solutions}
In the previous section, we found that the
eigenspaces of  the  spectral problem
$\Delta u -\lambda Au=\rho u$ in $\dbH^2_\sigma$
are  even-dimensional.
Thus it will be difficult  to obtain steady-state
bifurcations in the whole
space $\dbH^2_\sigma$.
Fortunately,  (\ref{ns1}-\ref{ns4}) admit
many flow invariant subspaces, in which it will be convenient to
examine  the bifurcation phenomena.

First, by the divergence free condition $\mbox{dvi} u =0$,
we observe that the bilinear term $B(u,u)=P[u\cdot \nabla u]$ can
be alternatively  given by
\begin{eqnarray*}
B(u,u)&=&u\cdot\nabla u-
\nabla
\Delta^{-1}(
\partial_x(u\cdot\nabla u_1)
+\partial_y(u\cdot\nabla u_2)
+\partial_z(u\cdot\nabla u_1) ) \\
&=& \nabla\cdot (u\otimes u)-\nabla
\Delta^{-1}\nabla^2\cdot (u\otimes u)
\end{eqnarray*}
with
$$\nabla\cdot (u\otimes u) =\left(
\sum_{n=1}^3\partial_n(u_n u_1),    \sum_{n=1}^3\partial_n(u_n u_2),
\sum_{n=1}^3\partial_n(u_n u_3)\right)$$
and
$$\nabla^2\cdot (u\otimes u) =
\sum_{n=1}^3\sum_{m=1}^3\partial_n\partial_m(u_n u_m)$$
for $(\partial_1,\partial_2,\partial_3)= (\partial_x,\partial_y,\partial_z)$.



Let the integers $l\geq 0$, $j\in {\mathbb Z}$ and $0\leq i<k$.
We introduce  the following
subspace of the scalar function space  $H^2({\mathbb T}^3)$:
\begin{eqnarray*}
H^2_{l,j,i,k}({\mathbb T}^3)&=&\left\{v \in
H^2({\mathbb T}^3)\,|
\,v=  \sum_{n\in {\mathbb N}}\xi_n\sin nkz\right.\\
&&+\left.\sum_{m\in {\mathbb N}}\sum_{n\in {\mathbb Z}}
\eta_{m,n}\sin (mlx+mjy+\{mi\}z+knz)\right\},
\\
\tilde{H}^2_{l,j,i,k}({\mathbb T}^3)&=&\left\{v
\in H^2({\mathbb T}^3)\,| \,v= \sum_{n\in {\mathbb N}}\xi_n\sin nkz\right.\\
&&+\left.\sum_{m\in {\mathbb N}}\sum_{n\in {\mathbb Z}}
\eta_{m,n}\cos((2m-1)(lx+jy)+\{2mi-i\}z+ nkz)\right. \\
&&+\left.\sum_{m\in {\mathbb N}}\sum_{n\in {\mathbb Z}}
\zeta_{m,n}\sin(2m(lx+jy)+\{2mi\}z+ nkz)\right\}.\\
\end{eqnarray*}
Here
$H^2_{l,j,i,k}({\mathbb T}^3)$ is a Hilbert space with the norm
\begin{eqnarray*}
\|v\|_{H^2}&=&\|\Delta v\|_{L^2} \\
&=&\left(\sum_{n\in {\mathbb N}}n^4k^4|\xi_n|^2
+\sum_{m\in {\mathbb N}}\sum_{n\in {\mathbb Z}}(m^2(l^2+j^2)+(nk+\{mi\})^2)^2|\eta_{m,n}|^2
\right)^{1/2},
\end{eqnarray*}
while
$\tilde{H}^2_{l,j,i,k}({\mathbb T}^3)$ is a Hilbert space
with the norm
\begin{eqnarray*}
\lefteqn{ \|v\|_{H^2} = \|\Delta v\|_{L^2} }\\
&=&\left(\sum_{n\in {\mathbb N}}n^4k^4|\xi_n|^2
+\sum_{m\in {\mathbb N}}\sum_{n\in {\mathbb Z}}((2m-1)^2(l^2+j^2)+(nk+\{2mi-i\}
)^2)^2|\eta_{m,n}|^2\right.\\
&&\left.+\sum_{m\in {\mathbb N}}\sum_{n\in {\mathbb Z}}((2m)^2(l^2+j^2)
+(nk+\{2mi\})^2)^2|\zeta_{m,n}|^2 \right)^{1/2}.
\end{eqnarray*}

Then we are ready to define the invariant subspaces for
the Navier-Stokes equations by
\begin{align}
& \dbH^2_{l,j,i,k}=\left\{u=(u_1,u_2,u_3)  \in
\dbH^2_\sigma
|\,\,u_1,\,u_2,\,u_3\in H^2_{l,j,i,k}({\mathbb T}^3)\right\},
                 \label{3.2} \\
& \tilde{\dbH}^2_{l,j,i,k}=\left\{u=(u_1,u_2,u_3)
\in \dbH^2_\sigma
|\,\, u_1,\,u_2,\,u_3\in
\tilde{H}^2_{l,j,i,k}({\mathbb T}^3)\right\}.  \label{3.3}
\end{align}

\begin{Lemma}
Let ${\mathbb X}$ be either $ \dbH^2_{l,j,i,k}$ or $\tilde{\dbH}^2_{l,j,i,k}$
with integers  $0\leq l$, $j\in {\mathbb Z}$  and $0\leq i<k$.     Then
for any $u\in {\mathbb X}$, $\Delta^{-1} B(u,u) \in {\mathbb X}$. Consequently, ${\mathbb X}$ is
invariant under the steady-state Navier-Stokes equations (\ref{2.1}).
\end{Lemma}

\begin{proof}
For $u\in {\mathbb X}$, it follows from
H\"older's inequality and Sobolev's inequality  that
$$\|B(u,u)\|_{L^2}\leq  c\|u\|_{L^3}\|\nabla u\|_{L^6}
\leq  c\|\Delta u\|_{L^2}^2.$$
This together with the divergence free condition
$\nabla\cdot \Delta^{-1}B(u,u)=0$
gives \\     $\Delta^{-1}B(u,u) \in \dbH^2_\sigma$.

If   $u=(u_1,u_2,u_3)\in \dbH^2_{l,j,i,k}$,
we see that $u$ is in the form
\begin{eqnarray*}
(u_1,u_2,u_3)&=&\sum_{n\in {\mathbb N}}(\xi_{n,1},\xi_{n,2},\xi_{n,3})\sin nkz\\
&&+\sum_{m\in {\mathbb N}}\sum_{n\in {\mathbb Z}}
(\eta_{m,n,1},\eta_{m,n,2},\eta_{m,n,3})\sin (mlx+mjy+\{mi\}z+knz).
\end{eqnarray*}
This gives, for $i',\,j'=1,\,2,\,3$,
\begin{eqnarray*}
u_{i'}u_{j'}&=&
\big(\sum_{n\in {\mathbb N}}\xi_{n,i'}\sin nkz
+\sum_{m\in {\mathbb N}}\sum_{n\in {\mathbb Z}}
\eta_{m,n,i'}\sin (mlx+mjy+\{mi\}z+knz)\big) \\
&&\times
\big(\sum_{n\in {\mathbb N}}\xi_{n,j'}\sin nkz
+\sum_{m\in {\mathbb N}}\sum_{n\in {\mathbb Z}}
\eta_{m,n,j'}\sin (mlx+mjy+\{mi\}z+knz)\big)\\
&=&c_{i',j'}+\sum_{n\in {\mathbb N}}\xi_{n,i',j'}\cos nkz \\
&&+\sum_{m\in {\mathbb N}}\sum_{n\in {\mathbb Z}}
\eta_{m,n,i',j'}\cos(mlx+mjy+\{mi\}z+knz)
\end{eqnarray*}
for some constants $c_{i',j'}$, $\xi_{n,i',j'}$ and
$\eta_{m,n,i',j'}$.
This yields
\begin{eqnarray*}
B(u,u)&=&\sum_{n\in {\mathbb N}}(\xi'_{n,1},\xi'_{n,2},\xi'_{n,3})\sin nkz\\
&&+\sum_{m\in {\mathbb N}}\sum_{n\in {\mathbb Z}}
(\eta'_{m,n,1},\eta'_{m,n,2},\eta'_{m,n,3})\sin (mlx+mjy+miz+knz).
\end{eqnarray*}
for some constants   $\xi_{n,i'}$ and $\eta'_{m,n,i'}$.
We thus have $\Delta^{-1}B(u,u)\in \dbH^2_{l,j,i,k}$.

Likewise, if  $u=(u_1,u_2,u_3)\in \tilde\dbH^2_{l,j,i,k}$ with
\begin{eqnarray*}
\lefteqn{ (u_1,u_2,u_3) }\\
&=&\sum_{n\in {\mathbb N}}(\xi_{n,1},\xi_{n,2},\xi_{n,3})\sin nkz\\
&&+\sum_{m\in {\mathbb N}}\sum_{n\in {\mathbb Z}}
(\eta_{m,n,1},\eta_{m,n,2},\eta_{m,n,3})
 \cos ((2m-1)(lx+jy)+\{2mi-i\}z+knz)
\\
&&+\sum_{m\in {\mathbb N}}\sum_{n\in {\mathbb Z}}
(\zeta_{m,n,1},\zeta_{m,n,2},\zeta_{m,n,3})\sin
(2m(lx+jy)+\{2mi\}z+knz),
\end{eqnarray*}
we have, for $i',\,j'=1,\,2,\,3$,
\begin{eqnarray*}
u_{i'}u_{j'}&=&\sum_{n\in {\mathbb N}}\xi_{n,i',j'}\cos nkz   \\
&&
+\sum_{m\in {\mathbb N}}\sum_{n\in {\mathbb Z}}
\eta_{m,n,i',j'}\sin ((2m-1)(lx+jy)+\{2mi-i\}z+knz)
\\
&&+c_{i',j'}+\sum_{m\in {\mathbb N}}\sum_{n\in {\mathbb Z}}
\zeta_{m,n,i',j'}\cos (2m(lx+jy)+\{2mi\}z+knz)
\end{eqnarray*}
for some constants $c_{i',j'}$,  $\xi_{n,i',j'}$, $\eta_{m,n,i',j'}$ and
$\zeta_{m,n,i',j'}$.
It thus readily seen that $\Delta^{-1}B(u,u)\in \tilde\dbH^2_{l,j,i,k}$.
\end{proof}

\begin{Remark}
{\rm Both $\dbH^2_{l,j,i,k}$ and $\tilde{\dbH}^2_{l,j,i,k}$
are also invariant
with respect to the time-dependent Navier-Stokes equations
(\ref{ns1}-\ref{ns4}). Hence, (\ref{ns1}-\ref{ns4}) can be restrict to
either of the invariant subspaces.
}
\end{Remark}

\begin{Lemma}
\begin{enumerate}

\item For any $\lambda \in (0, \infty)$, there is at least one steady-state
solution for (\ref{2.1}).

\item There exists an absolute constant $c_1$ independent of $\lambda $ and
$k$ such that
$u_0=(\sin kz,0, 0)$ is the unique solution of (\ref{2.1})
for any $0< \lambda \le c_1/k^2$.

\item
For any solution   $ u\in \dbH^2_\sigma$ of (\ref{2.1}),
there exist two absolute constants independent of $\lambda $, $k$
and $ u$ such that
\begin{equation} \label{3.4}
\gathered
 \| u\|_{H^1}\leq c_2k^2, \\
 \| u\|_{H^2}\leq c_3 k^2(\lambda^2k^4+1).
\endgathered
\end{equation}
\end{enumerate}

\end{Lemma}

\begin{proof}
The proof of the first two parts of the lemma is standard; we omit the details.
The first inequality (\ref{3.4}) is  a direct consequence of
the standard energy estimates, while the second inequality follows from
the following computation:
\begin{eqnarray*}
\|\Delta u\|_{L^2}
&\leq & \lambda \|B(u,u)\|_{L^2}+k^2\| u_0\|_{L^2}\\
&\leq &c\lambda\| u\cdot\nabla u\|_{L^2}+ k^2\| u_0\|_{L^2}\\
&\leq &c\lambda \| u\|_{L^6}\|\nabla u\|_{L^3}+k^2\| u_0\|_{L^2}\\
&\leq &c\lambda \| \nabla u\|_{L^2}^{3/2}\|\Delta u\|_{L^2}^{1/2}+k^2\|
u_0\|_{L^2}\\
&\leq &c\lambda^2 \| \nabla u\|_{L^2}^3+
\frac{1}{2}\|\Delta u\|_{L^2}+k^2\| u_0\|_{L^2}\\
&\leq &c\lambda^2 k^6 \| u_0\|_{L^2}^3+
\frac{1}{2}\|\Delta u\|_{L^2}+k^2\| u_0\|_{L^2}.
\end{eqnarray*}
\end{proof}

\begin{Lemma}
Let ${\mathbb X}$ be either $ \dbH^2_{l,j,i,k}$ or $\tilde{\dbH}^2_{l,j,i,k}$
with integers  $0\leq l$, $j\in {\mathbb Z}$  and $0\leq i<k$.
Then  $\Delta^{-1}A$ and $\Delta^{-1}B$ are compact and
continuous operators  from ${\mathbb X}$ into itself.
\end{Lemma}
\begin{proof}
It suffices to notice that for  $  u\in \dbH^2_\sigma$,
$$\gathered
 \|A u\|_{L^2}\leq c\|\nabla u\|_{L^2}, \\
 \|B( u)\|_{L^2} \leq  c\|u\|_{L^{12}}\|\nabla  u\|_{L^{12/5}}
 \leq c\|\nabla  u\|^2_{L^{12/5}}.
\endgathered
$$
\end{proof}


\setcounter{equation}{0}
\section{Global Branches of Steady-State Solutions}

In this section, we consider (\ref{2.1}) in the flow invariant
subspace $\dbH^2_{n_0l,n_0j,\{n_0i\},k}$  to show  the existence
two global branches of steady-state solutions bifurcating from
the basic flow $u_0$ at a critical
value $\lambda_{l,j,i,k}$.

\begin{Theorem} Let $(l,j,i,k)$ be an integer vector
satisfying (\ref{4-tuple}) and $n_0$ be the integer given in (\ref{n-zero}).
Then there exists a
critical Reynolds number
$\lambda_{l,j,i,k}>0$ such that
problem   (\ref{2.1}) admits two  different global branches of
steady-state solutions
$ u^1_{l,j,i,k,\lambda},\,\, u^2_{l,j,i,k,\lambda}\in
\dbH^2_{n_0l,n_0j,\{n_0i\},k}\subset \dbH^2_{l,j,i,k}\subset \dbH^2_\sigma $
branching  off the bifurcation
point $(\lambda_{l,j,i,k}, u_0)$ continuously such that
for $m=1,\,2$,
 \begin{enumerate}
\item each branch $u^m_{l,j,i,k,\lambda}$ extends to $\lambda=\infty$ in
$\mathcal C_1 \cap  \mathcal C_2$ as shown in  Figures 1.1 and 1.2;

\item
each branch $u^m_{l,j,i,k,\lambda}$  starting from the bifurcation point
$(\lambda_{l,j,i,k},u_0)$
never touches the $\lambda$--axis for $\lambda \neq \lambda_{l,j,i,k}$.
\end{enumerate}    \label{Th4.1}
\end{Theorem}

\begin{proof}
For the integer $n_0$ given in
 (\ref{n-zero}),   we see that
 $\dbH^2_{l,j,i,k}\supset \dbH^2_{n_0l,n_0j,\{n_0i\},k}$,
 which is in the following form
 \begin{eqnarray*}
&&\left\{u\in \dbH^2_\sigma\,|\,u=
\sum_{n\in {\mathbb N}}(\xi_{n,1},\xi_{n,2},\xi_{n,3})\sin nkz \right.\\
&&+\sum_{m\geq 2}\sum_{n\in {\mathbb Z}}
(\eta_{m,n,1},\eta_{m,n,2},\eta_{m,n,3})
\sin (mn_0lx+mn_0jy+\{mn_0i\}z+knz)\\
&&+\left.\sum_{n\in {\mathbb Z}}
(\zeta_{n,1},\zeta_{n,2},\zeta_{n,3})
\sin (n_0lx+n_0jy+\{n_0i\}z+knz)\right\}.
\end{eqnarray*}
It is obvious that the spectral problem
\begin{equation}
\label{4.2}
\Delta u -\lambda A u=\rho u\,\quad \mbox{ for  }\,
\rho >-\beta_0     \mbox{ and } u\in \dbH^2_{n_0l,n_0j,\{n_0i\},k}
\end{equation}
has no eigenfunction in the form
$$\sum_{n\in {\mathbb N}}(\xi_{n,1},\xi_{n,2},\xi_{n,3})\sin nkz. $$
Thus it follows from Theorem 2.3 that (\ref{4.2}) has a unique
eigenvalue  $\rho_{n_0l,n_0j,\{n_0i\},k}$
transversal
across the imaginary axis at the origin, and the choice of the integer
$n_0$ ensures the simplicity of this eigenvalue.

Let us denote by $\lambda_{l,j,i,k}$ the critical Reynolds number
such that
\begin{equation}
\label{4.3}
\rho_{n_0l,n_0j,\{n_0i\},k}(\lambda_{l,j,i,k})=0.
\end{equation}


Taking  $  u'+ u_0$ in place of  $
 u$ in the steady-state equations with respect to
 (\ref{2.1}) reduced in $\dbH^2_{n_0l,n_0j,\{n_0i\},k}$
 and omitting the prime,  we have
\begin{equation}\label{4.4}
u =\lambda\Delta^{-1} A u +\lambda\Delta^{-1} B( u),\quad
u\in \dbH^2_{n_0l,n_0j,\{n_0i\},k},
\end{equation}
of which the linear part
\begin{equation}
\label{4.5}
u -\lambda\Delta^{-1} A u=0,\quad
u\in \dbH^2_{n_0l,n_0j,\{n_0i\},k}
\end{equation}
has a simple eigenvalue $1/\lambda_{l,j,i,k}$ or a simple characteristic
value $\lambda_{l,j,i,k}$ in the notation of Rabinowitz \cite{Rab}.
It follows from \cite[Theorem 1.40]{Rab} that   (\ref{4.4}) admit
two continuous branches of solutions
$(\lambda,u_{l,j,i,k,\lambda}^1-u_0)$ and
$(\lambda, u_{l,j,i,k,\lambda}^2-u_0)$
other than the branch $(\lambda, \,0)$ in
the space $\mathbb  R \times   \dbH^2_{n_0l,n_0j,\{n_0i\},k}$ under the norm
$\|(\lambda,u)\|= (\lambda^2+\|\Delta u\|_{L^2})^{1/2}$. Each of those
two branches meets
 $(\lambda_{l,j,i,k},\,0)$ and either

\begin{enumerate}

\item[(i)]  meets   $\infty $
 in    $\mathbb R\times   \dbH^2_{n_0l,n_0j,\{n_0i\},k}$, or

\item[(ii)] meets $(\tilde\lambda,\,0)$,
 where $\tilde\lambda\neq \lambda_{l,j,i,k}$    is a
 characteristic value
 of  (\ref{4.5}).

\end{enumerate}

Observe that
$(\lambda,0)$ is the unique solution of  (\ref{4.4})
when $0 < \lambda < {c}/{k^2}$.
Thus these two branches of solutions
are contained in the space
$(0,\infty)\times   \dbH^2_{n_0l,n_0j,\{n_0i\},k}$.
Meanwhile, the {\it a priori} bounds of the
the steady-state solutions of  (\ref{2.1})
for every $\lambda >0$ show that both branches of the bifurcation
solutions extends to $\lambda=\infty $
 in $\mathcal C_1 \cap \mathcal C_2$; see Figures 1.1 and 1.2.


Furthermore, the above second alternative is excluded
due to the fact that $\lambda_{l,j,i,k}$ is the only
characteristic value  of  (\ref{4.5}) for $\lambda \in (0,\infty)$.
Hence, each branch $u^m_{l,j,i,k,\lambda}$ intersects with
the $\lambda$--axis
only at the line segment $pq$ shown in Figures 1.1  and 1.2.
In particular, each branch $u^m_{l,j,i,k,\lambda}$ starting
from the bifurcation
point $(\lambda_{l,j,i,k},u_0)$
never touches the $\lambda$--axis for $\lambda \neq \lambda_{l,j,i,k}$.

\end{proof}


\setcounter{equation}{0}
\section{Supercritical Pitchfork Bifurcation}

 In order to reach our main results, it is necessary to give the local
 behavior  of the global branches of the steady-state solutions
 $u^1_{l,j,i,k,\lambda}$ and $u^2_{l,j,i,k,\lambda}$ close to the bifurcation
 point $(\lambda_{l,j,i,k},u_0)$. More precisely,
 we have the following supercritical pitchfork bifurcation theorem.

\begin{Theorem} \label{Th5.1}
    Let the  steady-state solutions
    $u^1_{l,j,i,k,\lambda}$, $u^2_{l,j,i,k,\lambda}\in
    \dbH^2_{n_0l,n_0j,\{n_0i\},k}$ be given in Theorem \ref{Th4.1}.
    Then for  $m=1,$ $2$,
\begin{equation} \label{t5.1}
\begin{split}
& u^m_{l,j,i,k,\lambda}= u_0  \qquad  \qquad \hspace{25mm}
     \mbox{ for  } \lambda <\lambda_{l,j,i,k}, \\
& u_0\neq u^1_{l,j,i,k,\lambda}\neq  u^2_{l,j,i,k,\lambda} \neq u_0\qquad
     \mbox{ for }\lambda_{l,j,i,k}<\lambda,
\end{split}
\end{equation}
provided that 
\begin{equation}\label{t5.2}
\|\Delta(u^m_{l,j,i,k\lambda}-u_0)\|_{L^2}+|\lambda-\lambda_{l,j,i,k}|
<\epsilon
\end{equation}
for some  small constant
    $\epsilon >0$.   
\end{Theorem}


For simplicity, we set 
\begin{equation}\label{t5.3}
\lambda_0=\lambda_{l,j,i,k},\,\,\,\,
u^+_\lambda=u^1_{l,j,i,k,\lambda},\,\,\,\,
u^-_\lambda  =u^2_{l,j,i,k,\lambda}.
\end{equation}


The proof of Theorem 5.1 will be achieved by analyzing the local 
asymptotic expansion of  $u^\pm_\lambda$  in terms of $\lambda$. 

First, notice that he bifurcation phenomenon of the Navier-Stokes system 
     (\ref{2.1}) reduced to  the subspace
     $\dbH^2_{n_0l,n_0j,\{n_0i\},k}$ is excited by the external force
     $k^2u_0$ and the unstable subspace of the bifurcation point
     $(\lambda_0,u_0)$ is contained in $  {\mathbb E}_{n_0l,n_0j,\{n_0i\},k} $.

The unstable subspace is given by 
the nontrivial solution $u^* \in {\mathbb E}_{n_0l,n_0j,\{n_0i\},k}$ 
of the spectral problem  $\Delta u^*-\lambda_0Au^*=0$:
\begin{equation}\label{u*}
\gathered
 u^*=\sum_{n\in {\mathbb Z}}
(\xi_n,\eta_n,\zeta_n)\sin(n_0lx+n_0jy+\{n_0i\}z+nkz) 
            \in {\mathbb E}_{n_0l,n_0j,\{n_0i\},k}, \\
 \zeta_0 =1, \\
 \Delta u^*-\lambda_0Au^*=0, 
\endgathered 
\end{equation}
By Theorem 2.3, for any $n$, $\zeta_n \not=0$. Hence 
there exists a unique nontrivial solution  $u^*$  of (\ref{u*}) with 
$\zeta_0 =1.$

Then we examine  the nonlinear interaction of the unstable direction 
given by $u^*$. To this end, we notice that 
$$B(u^*,u^*) \in {\mathbb E}_{0,0,0,k}\oplus {\mathbb E}_{2n_0l,2n_0j,\{2n_0i\},k}.$$
Hence, we can  decompose $B(u^*,u^*)$  as
\begin{equation}\label{interaction}
\begin{split}
&  B(u^*,u^*)= B_0(u^*,u^*)+ B_2(u,u)+B_3(u,u),  \\
&  B_0(u^*,u^*)=\frac{u_0}{4\pi^3} 
      \int_{{\mathbb T}^3} B(u^*,u^*)\cdot u_0\,dx\,dy\,dz ,  \\
& B_2(u^*,u^*) \in  {\mathbb E}_{0,0,0,k}\,\,\mbox{ with }\,\, 
 \int_{{\mathbb T}^3} B_2(u^*,u^*)\cdot u_0\,dx\,dy\,dz =0, \\
& B_3(u^*,u^*) \in  {\mathbb E}_{2n_0l,2n_0j,\{2n_0i\},k}.
\end{split}
\end{equation}

To prove Theorem 5.1, it suffices for us to prove the following 
stronger result, providing a detailed local asymptotic expansion of 
the bifurcation solutions in terms of $\lambda$.

\begin{Theorem} If (\ref{t5.2})  holds true for some small 
number $\epsilon>0$, then 
the two branches of steady-state solutions $(\lambda, u^\pm_\lambda)$ 
close to the bifurcation 
point $(\lambda_0,u_0)$ undergo  the supercritical pitchfork bifurcation 
in the following sense: 
\begin{equation} \label{expansion}
u^\pm_\lambda  = \cases 
 u_0   & \lambda<\lambda_0, \\[2pt]
 u_0 \pm \alpha \frac{\sqrt{\lambda-\lambda_0}}{\lambda} u^*
 +\frac{\lambda-\lambda_0}{\lambda}
 [ -u_0 + \alpha^2 \Delta^{-1} B_2(u^*,u^*) \\
 +\alpha^2 (-\Delta +\lambda_0 A)^{-1}B_3(u^*,u^*) ]+o(\lambda-\lambda_0) 
    & \lambda>\lambda_0\,.\endcases
\end{equation}
Here  the eigenfunction 
$u^*$ is given by (\ref{u*})  and  the constant  
$\alpha$ is  defined 
\begin{equation}\label{alpha}
 \alpha= \frac{2k\pi^{3/2}}{\displaystyle \left( 
\int_{{\mathbb T}^3} B(u^*,u^*)\cdot u_0 \,dx\,dy\,dz \right)^{1/2}}.
\end{equation}
\end{Theorem}

The definition of the constant $\alpha$ in (\ref{alpha}) 
is justified in the following lemma.

\begin{Lemma}
 $$\int_{{\mathbb T}^3}B(u^*,u^*)\cdot u_0 \,dx\,dy\,dz >0.$$
\end{Lemma}

\begin{proof}[Proof of Lemma 5.1]
   By the divergence free condition and integration by parts,
   \begin{eqnarray}
  \label{n5.6}
  \lefteqn{ \int_{{\mathbb T}^3}B(u^*,u^*)\cdot u_0 \,dx\,dy\,dz }\\
   &=&-\int_{{\mathbb T}^3}B(u^*,u_0)\cdot u^*\,dx\,dy\,dz \nonumber\\
&=& -k\sum_{n,\,m\in {\mathbb Z}}\xi_n\zeta_m\int_{{\mathbb T}^3}
\sin(n_0lx+n_0jy+\{n_0i\}z+nkz)\nonumber\\
&&\times \sin(n_0lx+n_0jy+\{n_0i\}z+mkz) \cos kz    \,dx\,dy\,dz\nonumber \\
&=& -2k\pi^3 \sum_{n\in {\mathbb Z}} \xi_n(\zeta_{n-1}+\zeta_{n+1}). \nonumber
\end{eqnarray}

On the other hand, recall from (\ref{2.15})-(\ref{2.17}) 
that $u^*$ satisfies 
the following set of coupled algebraic equations
\begin{align}
& \frac{2\beta_n^2}{\lambda_{l,j,i,k} n_0l}\zeta_n
+(\beta_{n-1}-k^2)\zeta_{n-1}-(\beta_{n+1}-k^2)\zeta_{n+1}
=0,\label{n5.7} \\
&  \label{n5.8}
\frac{2\beta_n}{\lambda_{l,j,i,k} }(j\xi_n-l\eta_n)
+n_0l(j\xi_{n-1}-l\eta_{n-1})-n_0l(j\xi_{n+1}-l\eta_{n+1})
\\
& \qquad =-kj (\zeta_{n-1}+\zeta_{n+1}), \nonumber \\
& -n_0j\eta_n=n_0l\xi_n+(\{n_0i\}+nk)\zeta_n,
\label{n5.9}
\end{align}
for any $n\in {\mathbb Z}$, where $\beta_n = n_0^2l^2 +n_0^2j^2 +(\{n_0i\} +nk)^2$.

 Multiplying (\ref{n5.7}) by $\zeta_n$ and summing up the resultant equations
 yield 
 \begin{eqnarray}
 \label{n5.10}
 \sum_{n\in{\mathbb Z}}
\frac{2\beta_n^2}{\lambda_{l,j,i,k} n_0l}\zeta_n^2
&=&-\sum_{n\in {\mathbb Z}}(\beta_n-\beta_{n+1})\zeta_n\zeta_{n+1}\\
&=&\sum_{n\in{\mathbb Z}}k(2\{n_0i\}+2nk+k)\zeta_n\zeta_{n+1}.\nonumber
\end{eqnarray}
 Moreover, 
 multiplying (\ref{n5.8}) by $j\xi_n-l\eta_n$ 
 and summing up the resultant equations
 yield 
\begin{equation} \label{nn5.10}
\sum_{n\in {\mathbb Z}}\frac{2\beta_n}{\lambda_{l,j,i,k} }(j\xi_n-l\eta_n)^2
=-k\sum_{n\in{\mathbb Z}}(j^2\xi_n-lj\eta_n)(\zeta_{n-1}+\zeta_{n+1}).
\end{equation}
We infer then from  (\ref{n5.6}), (\ref{n5.9}), (\ref{n5.10}) and 
(\ref{nn5.10}) that
\begin{eqnarray*}  \lefteqn{
\sum_{n\in {\mathbb Z}}\frac{2\beta_n}{\lambda_{l,j,i,k} }(j\xi_n-l\eta_n)^2 }\\
&=&
   -k\sum_{n\in{\mathbb Z}}(j^2\xi_n+l^2\xi_n +\frac{l}{n_0}
   (\{n_0i\}+nk)\zeta_n)(\zeta_{n-1}+\zeta_{n+1})\\      
&=&
   -k(l^2+j^2)\sum_{n\in{\mathbb Z}}\xi_n(\zeta_{n-1}+\zeta_{n+1})      
-\frac{kl}{n_0}\sum_{n\in{\mathbb Z}}(\{n_0i\}+nk)\zeta_n
(\zeta_{n-1}+\zeta_{n+1})\\      
&=&\frac{l^2+j^2}{2\pi^3}\int_{{\mathbb T}^3}B(u^*,u^*)\cdot u_0\,dx\,dy\,dz
-\frac{kl}{n_0}\sum_{n\in{\mathbb Z}}(2\{n_0i\}+2nk+k)\zeta_n\zeta_{n+1}\\      
&=&\frac{l^2+j^2}{2\pi^3}\int_{{\mathbb T}^3}B(u^*,u^*)\cdot u_0\,dx\,dy\,dz
- \sum_{n\in{\mathbb Z}}\frac{2\beta_n^2}{\lambda_{l,j,i,k} n_0^2 }\zeta_n^2.   
\end{eqnarray*}
We thus have 
$$
 \int_{{\mathbb T}^3}B(u^*,u^*)\cdot u_0\,dx\,dy\,dz                                  
 =\frac{4\pi^3}{(l^2+j^2)\lambda_{l,j,i,k}n_0^2} 
 \sum_{n\in {\mathbb Z}}(\beta_nn_0^2(j\xi_n-l\eta_n)^2 +\beta_n^2\zeta_n^2) >0. 
 $$

\end{proof}



Now in order to facilitate the understanding of the problem, 
we consider the infinite 
mode truncation model by projecting  (\ref{2.1})
onto the unstable space Span$\{u_0\}\oplus {\mathbb E}_{n_0l,n_0j,\{n_0i\},k}$:
\begin{equation}\label{proj}
\gathered 
 \mbox{ Find } u\in  W = 
\mbox{ Span}\{u_0\}\oplus {\mathbb E}_{n_0l,n_0j,\{n_0i\},k}  \mbox{  such that  } \\
 \int_{{\mathbb T}^3}(-\Delta(u-u_0) +\lambda B(u,u) )
\cdot \phi \,dx\,dy\,dz =0,  \qquad 
\forall \,\phi \in  W.
\endgathered 
\end{equation}

\begin{Lemma} The parameter $\lambda=\lambda_0$  is 
a supercritical bifurcation point of (\ref{proj}) with 
the two bifurcation branches given by 
\begin{equation}\label{asmpt}
u= \cases
 u_0 & \mbox{ when } \lambda\leq \lambda_0,\\ 
 u_0 \pm \alpha  \frac{\sqrt{\lambda-\lambda_0}} {\lambda} u^*
       - \frac{\lambda-\lambda_0}{\lambda}u_0
   &  \mbox{ when } \lambda> \lambda_0,
\endcases 
\end{equation}
for 
$ \|\Delta(u -u_0)\|_{L^2} +|\lambda-\lambda_0|<\epsilon$
for some constant $\epsilon >0$. 
\end{Lemma}

\begin{proof}[Proof of Lemma 5.2]
Decompose $u$  as
$$u=\mu u_0 + v  \,\,\mbox{  with } \mu\in {\mathbb R} ,\,\,
v\in {\mathbb E}_{n_0l,n_0j,\{n_0i\},k}.$$
Then (\ref{proj})  becomes
\begin{gather}
 \label{trun1}
-\Delta v +\lambda \mu A v = 0,\\
 k^2(\mu -1)+\frac{\lambda }{4\pi^3} \int_{{\mathbb T}^3} B(v,v)\cdot u_0
\,dx\,dy\,dz =0, 
\label{trun2}
\end{gather}
It follows from  Theorem 2.3 
 that $ u=\mu u_0+v$ solves
(\ref{trun1})-(\ref{trun2}) if and only if
$$\mu\lambda =\lambda_0 \,\,\mbox{ and }
\,\, v=cu^*  \,\,\mbox{ for some }\,\,c\in {\mathbb R},$$
 where the eigenfunction  $u^*$ is defined  by (\ref{u*}).
Thus (\ref{trun2}) becomes
$$    \frac{4k^2\pi^3(\lambda -\lambda_0)}{\lambda^2 }=  c^2
\int_{{\mathbb T}^3} B(u^*,u^*)\cdot u_0 \,dx\,dy\,dz.$$
This together with Lemma 5.1 implies 
$$c=\cases 0 & \mbox{ if } \lambda <\lambda_0,  \\[2pt]
\pm \alpha \frac{\sqrt{\lambda-\lambda_0}}{\lambda}
& \mbox{ if } \lambda >\lambda_0\,. \endcases
$$
Thus (\ref{asmpt}) follows, and the lemma is proved.
\end{proof}

\begin{proof}[Proof of Theorem 5.2]
Lemma 5.2  implies  that
the bifurcation solutions $u^\pm_\lambda $ close to 
the bifurcation point $(\lambda_0,u_0)$ 
are in the following form:
\begin{equation}
\label{comp}
u^\pm_\lambda = u_0+ \frac{\sqrt{|\lambda-\lambda_0|}}{\lambda}w_1
+ \frac{\lambda-\lambda_0}{\lambda}w
+o(|\lambda-\lambda_0|),
\end{equation}
with $w_1,\, w\in  \dbH^2_{n_0l,n_0j,\{n_0i\},k}$ independent of
$\lambda$. Here $ w$ can be decomposed as 
$$\gathered
 w=-u_0 + w_2+w_3  \\
 w_2\in {\mathbb E}_{0,0,0,k}, \quad \int_{{\mathbb T}^3}w_2\cdot u_0 \,dx\,dy\,dz =0, \\
 \int_{{\mathbb T}^3}w_3\cdot \phi\,dx\,dy\,dz =0
\qquad \forall \phi\in {\mathbb E}_{0,0,0,k} \,.
\endgathered$$
Thus
\begin{eqnarray*}
\lefteqn{\lambda B(u^\pm_\lambda,u^\pm_\lambda)}\\
&=&
\sqrt{|\lambda-\lambda_0|}\left(B(u_0 +\frac{\lambda-\lambda_0}{\lambda}
(w_2-u_0),w_1)+B
(w_1,u_0+\frac{\lambda-\lambda_0}{\lambda}(w_2-u_0)) \right) \\
&&
+(\lambda-\lambda_0)\left(B(u_0 +\frac{\lambda-\lambda_0}{\lambda}(w_2-u_0)
,w_3)+
B(w_3,u_0+\frac{\lambda-\lambda_0}{\lambda}(w_2-u_0))\right) \\
&&+
              \frac{|\lambda-\lambda_0|}{\lambda}
              B(w_1,w_1)+ o(|\lambda-\lambda_0|)\\
     &=&
\frac{\sqrt{|\lambda-\lambda_0|}}{\lambda}\lambda_0(B(u_0,w_1)+B(w_1,u_0))\\
&&+\frac{\sqrt{|\lambda-\lambda_0|}}{\lambda}\lambda_0(B(u_0,w_3)+B(w_3,u_0))
            +  \frac{|\lambda-\lambda_0|}{\lambda}
              B(w_1,w_1)+ o(|\lambda-\lambda_0|).
              \end{eqnarray*}
Hence the stationary Navier-Stokes system
\begin{equation}
\label{st}
-\Delta(u^\pm_\lambda-u_0) +\lambda B(u^\pm_\lambda,u^\pm_\lambda)=0
\,\,\,\, \,\mbox{ in }
         \dbH^2_{n_0l,n_0j,\{n_0i\},k}
         \end{equation}
becomes
\begin{eqnarray*}
0&=&
-\frac{\sqrt{|\lambda-\lambda_0|}}{\lambda}\Delta w_1
   -\frac{\lambda-\lambda_0}{\lambda}\Delta (w_2-u_0)
   -\frac{\lambda-\lambda_0}{\lambda}\Delta w_3 \\
   &&
+\frac{\sqrt{|\lambda-\lambda_0|}}{\lambda}\lambda_0 A w_1
+\frac{(\lambda-\lambda_0)}{\lambda}\lambda_0 Aw_3
+
              \frac{|\lambda-\lambda_0|}{\lambda}
              B(w_1,w_1)+ o(|\lambda-\lambda_0|).
\end{eqnarray*}
This  yields
\begin{eqnarray*}
-\Delta w_1 +\lambda_0 Aw_1&=&0,\\
(\lambda-\lambda_0)(\Delta(w_2-u_0+w_3) -\lambda_0Aw_3) &=&
|\lambda-\lambda_0|B(w_1,w_1).
\end{eqnarray*}
By Theorem 2.3, we see  $w_1=cu^*\in {\mathbb E}_{n_0l,n_0j,\{n_0i\},k}$ 
for some constant $c$ 
with $u^*$ defined by (\ref{u*}). Hence
$$B(w_1,w_1)=c^2B(u^*,u^*) 
\in {\mathbb E}_{0,0,0,k}\oplus {\mathbb E}_{2n_0l,2n_0j,\{2n_0i\},k}, $$
with $B(u^*,u^*)$  decomposed by (\ref{interaction}).

Then (\ref{st}) yields 
\begin{eqnarray}
  \nonumber     0&=&-\Delta u^* +\lambda_0 Au^*,\\
      \label{n5.16}            (\lambda-\lambda_0)    k^2 u_0 &=&  c^2
|\lambda-\lambda_0| B_0(u^*,u^*)\\
    \nonumber      (\lambda-\lambda_0)\Delta w_2 &=&c^2
|\lambda-\lambda_0|B_2(u^*,u^*),\\
(\lambda-\lambda_0)(-\Delta w_3 +\lambda Aw_3) &=&  
-c^2|\lambda-\lambda_0|B_3(u^*,u^*).  \nonumber
\end{eqnarray}
By Lemma 5.1 and (\ref{n5.16}), we obtain that 
for $\lambda<\lambda_0$, $c=0$. For $\lambda>\lambda_0$, we have  
\begin{eqnarray}
       \nonumber 0&=&  -\Delta u^* +\lambda_0 Au^*,\\
\label{CCC}        c &=&\pm \alpha,\\
   \nonumber   w_2 &=&c^2 \Delta^{-1} B_2(u^*,u^*)\,\,\in {\mathbb E}_{0,0,0,k},\\
 \nonumber w_3&=&  c^2(\Delta -\lambda_0 A)^{-1}B_3(u^*,u^*)  \,\,
 \in {\mathbb E}_{2n_0l,2n_0j,\{2n_0i\},k}.
\end{eqnarray}
Here we have used Theorem 2.3,
the definition of $n_0$ and the Riesz-Schauder theory
to obtain the boundedness of the operator
$(-\Delta +\lambda_0 A)^{-1} $ in ${\mathbb E}_{2n_0l,2n_0j,\{2n_0i\},k}$. Thus
the two branches of steady-state solutions $(\lambda, u^\pm_\lambda)$ 
close to the bifurcation 
point $(\lambda_0,u_0)$ undergo  the supercritical pitchfork bifurcation 
in the sense given by (\ref{expansion}). The theorem is proved.
\end{proof}

To copy the results of Theorems 4.1 and 5.1 to another space  
$\tilde{\dbH}^2_{n_0l,n_0j,\{n_0i\},k}$.
we notice that
$\dbH^2_\sigma\supset
\tilde{\dbH}^2_{n_0l,n_0j,\{n_0i\},k},$ which is in the form
 $$\aligned
 \bigg\{&u\in \dbH^2_\sigma\,|\,
 u=\sum_{n\in {\mathbb N}} (\xi_{n,1},\xi_{n,2},\xi_{n,3}) \sin nkz \\
&+\sum_{m\geq 2}\sum_{n\in {\mathbb Z}}
 (\eta_{m,n,1},\eta_{m,n,2},\eta_{m,n,3})
 \cos ((2m-1)(n_0lx+n_0jy) \\
&\hspace{18mm} +\{2mn_0i-n_0i\}z+knz)\\
&+\sum_{n\in {\mathbb Z}}
 (\zeta_{n,1},\zeta_{n,2},\zeta_{n,3})
 \sin (2m(n_0lx+n_0jy)+\{2mn_0i\}z+knz) \\
&+\sum_{n\in {\mathbb Z}}
 (\eta_{1,n,1},\eta_{1,n,2},\eta_{1,n,3})
 \cos (n_0lx+n_0jy+\{n_0i\}z+knz)\bigg\},
\endaligned$$
and
$$\tilde{\dbH}^2_{l,j,i,k}\supset
\tilde\dbH^2_{n_0l,n_0j,\{n_0i\},k}\,\mbox{ for } n_0 \mbox{ odd }.$$
%(THE LINES FROM NEXT ONE TO THE ONE  ABOVE THEOREM 4.2 ARE NEW )
We see that for each $(l,j,i,k)$ satisfying  (\ref{4-tuple}),
$\tilde{\mathbb E}_{n_0l,n_0j,\{n_0i\},k}\subset\tilde
\dbH^2_{n_0l,n_0j,\{n_0i\},k}$.
It follows from Theorem 2.3 that
 $\rho_{n_0l,n_0j,\{n_0i\},k}$ is simple and is  unique
eigenvalue of the spectral problem
$$
\Delta u -\lambda A u=\rho u\,\quad \mbox{ for  }\,
\rho >-\beta_0     \mbox{ and } u\in\tilde \dbH^2_{n_0l,n_0j,\{n_0i\},k}$$
which is transversal across the imaginary  axis at the origin.
Then the  argument as that given in
the above proof of Theorems \ref{Th4.1} and \ref{Th5.1} shows that there exist two
steady-state solutions
$$\tilde u_{l,j,i,k,\lambda}^1,\,\,\tilde u_{l,j,i,k,\lambda}^2
\in \tilde \dbH^2_{n_0l,n_0j,\{n_0i\},k}$$
globally branching off the bifurcation point $(\lambda_{l,j,i,k},u_0)$,
where $\lambda_{l,j,i,k}$ is defined by (\ref{4.3}).  More precisely,
we have the following Theorem.


\begin{Theorem}
For the critical Reynolds number
$\lambda_{l,j,i,k}>0$ obtained by (\ref{4.3}), the problem  (\ref{2.1})
with $0<\lambda <\infty$
admits also two  steady-state solutions
$$\tilde{u}^1_{l,j,i,k,\lambda},
\,\,\tilde{ u}^2_{l,j,i,k,\lambda}\in
\tilde{\dbH}^2_{n_0l,n_0j,\{n_0i\},k}\subset \dbH^2_\sigma $$
branching  off the point $(\lambda_{l,j,i,k}, u_0)$ continuously
such that for $m=1,\,2$,
$$
\begin{array}{ll}
\displaystyle\tilde{ u}^m_{l,j,i,k,\lambda}=  u_0 &
    \mbox{ for  } \lambda\leq \lambda_{l,j,i,k}, \\[1pt]
\displaystyle
u_0\neq  \tilde{ u}^1_{l,j,i,k,\lambda}\neq  
\tilde{ u}^2_{l,j,i,k,\lambda}\neq u_0
& 
\mbox{ for }\lambda>\lambda_{l,j,i,k},
\end{array}
$$
provided that, for some small constant $\epsilon>0$,
$$\|\Delta(\tilde {u}^m_{l,j,i,k,\lambda}-u_0)\|_{L^2} +
|\lambda-\lambda_0|<\epsilon.$$
Furthermore,  both of the  global branches of  bifurcation solutions
$ \tilde u^m_{l,j,i,k,\lambda}$  
satisfying the following properties:

\begin{enumerate}
\item each branch $ \tilde u^m_{l,j,i,k,\lambda}$ extends
to $\lambda=\infty$ in
$\mathcal C_1 \cap  \mathcal C_2$ as shown in  Figures 1.1  and 1.2;

\item each branch $ \tilde u^m_{l,j,i,k,\lambda}$ intersects
with the $\lambda$--axis
only at the line segment $pq$. In particular, $ \tilde u^m_{l,j,i,k,\lambda}$
never touches the $\lambda$--axis for $\lambda > \lambda_{l,j,i,k}$.
\end{enumerate}
\end{Theorem}

\setcounter{equation}{0}
\section{Proof of Theorems 1.1 and 1.2}


For an given integer vector $(l,j,i,k)$ satisfying (\ref{4-tuple}), let
$u^m_{l,\pm j,i,k,\lambda}$, $\tilde u^m_{l,\pm j,i,k,\lambda}$,
$u^m_{l,\pm j,k-i,k,\lambda}$  and
$\tilde{u}^m_{l,\pm j,k-i,k,\lambda}$
be the bifurcation solutions obtained in  Theorems 4.1, 5.1  and 5.3, and let
\begin{equation} \label{5.1}
 u_{l,j,i,k,n,\lambda}=
 u^n_{l,j,i,k,\lambda} \in \dbH^2_{n_0l,n_0j,\{n_0i\},k},\mbox{ for }
n=1,2.
\end{equation}

For $\{n(k-i)\}=k-\{n_0i\}$,
we define other bifurcation solutions $u_{l,j,i,k,n,\lambda}$ as follows:

\begin{enumerate}

\item[a).] {\sl Case $\{n_0i\}=j=0$}:
\begin{equation}\label{5.2}
       u_{l,0,i,k,n,\lambda}=
\tilde{u}^{n-2}_{l,0,i,k,\lambda}
       \in \tilde  \dbH^2_{n_0l,0,0,k}
\mbox{ for } n=3,4;
\end{equation}

\item[b).] {\sl Case $\{n_0i\}=0, j\neq 0$}:
\begin{equation}\label{5.3}
 u_{l,j,i,k,n,\lambda}=\left\{\begin{array}{ll}
 u^{n-2}_{l,-j,i,k,\lambda}
  \in \dbH^2_{n_0l,-n_0 j,0,k}
&\mbox{ for }   n=3,4,\\[2pt]
\tilde{u}^{n-4}_{l,j,i,k,\lambda}
          \in      \tilde  \dbH^2_{n_0l,n_0 j,0,k}
        &\mbox{ for } n=5,6,\\[2pt]
\tilde{u}^{n-6}_{l,-j,i,k,\lambda}
          \in   \tilde  \dbH^2_{n_0l,-n_0 j,0,k}
&\mbox{ for } n=7,8.
\end{array}\right.
\end{equation}

\item[c).] {\sl Case $\{n_0i\}\not=0, j=0$}:
\begin{equation}\label{5.4}
 u_{l,0,i,k,n,\lambda}=\left\{\begin{array}{ll}
 u^{n-2}_{l,0,k-i,k,\lambda}
  \in \dbH^2_{n_0l,0, k-\{n_0 i\},k},
&\mbox{ for }
n=3,4,\\[2pt]
\tilde{u}^{n-4}_{l,0,i,k,\lambda}
          \in      \tilde  \dbH^2_{n_0l,0, \{n_0 i\},k}
&\mbox{ for } n=5,6,
\\[2pt]
\tilde{u}^{n-6}_{l,0,k-i,k,\lambda}
          \in          \tilde  \dbH^2_{n_0l,0, k-\{n_0 i\},k}
&\mbox{ for } n=7,8.
\vspace{3mm}
\end{array}\right.
\end{equation}

\item[d).] {\sl Case $\{n_0i\}\not=0, \,j\neq 0$}:
\begin{equation} \label{5.5}
 u_{l,j,i,k,n,\lambda}=\left\{\begin{array}{ll}
 u^{n-2}_{l,j,k-i,k,\lambda}
  \in \dbH^2_{n_0l,n_0 j,k-\{n_0i\},k},
 &\mbox{ for }  n=3, 4,
 \\[2pt]
 u^{n-4}_{l,-j,i,k,\lambda}
   \in \dbH^2_{n_0l,-n_0 j,\{n_0i\},k},
 &\mbox{ for }   n=5, 6,
 \\[2pt]
 u^{n-6}_{l,-j,k-i,k,\lambda}
   \in \dbH^2_{n_0l,-n_0 j,k-\{n_0i\},k},
 &\mbox{ for } n=7, 8,
 \\[2pt]
\tilde{u}^{n-8}_{l,j,i,k,\lambda}
               \in      \tilde  \dbH^2_{n_0l,n_0 j,\{n_0i\},k}
\\[2pt]
\tilde{u}^{n-10}_{l,j,k-i,k,\lambda}
                  \in   \tilde  \dbH^2_{n_0l,n_0 j,k-\{n_0i\},k}
&\mbox{ for } n=11, 12,
\\[2pt]
\tilde{u}^{n-12}_{l,-j,i,k,\lambda}
                 \in \tilde  \dbH^2_{n_0l,-n_0 j,\{n_0i\},k}
&\mbox{ for } n=13, 14,
\\[2pt]
\tilde{u}^{n-14}_{l,-j,k-i,k,\lambda}
               \in           \tilde  \dbH^2_{n_0l,-n_0 j,k-\{n_0i\},k}
&\mbox{ for } n=15, 16.
\end{array}\right.
\end{equation}
\end{enumerate}

\begin{Lemma}
The above defined  solutions $u_{l,j,i,k,n,\lambda}$ bifurcate from the same
critical Reynolds number $\lambda_{l,j,i,k}$.
\end{Lemma}

\begin{proof}
By Theorem 2.3, the fact $\{n_0(k-i)\}=k-\{n_0i\}$
in case of $\{n_0i\}\neq 0$ and
the proof of Theorem 2.1, it is easy to see that
$\rho_{n_0l,n_0j,0,k}=\rho_{n_0l,-n_0j,0,k}$  when $\{n_0i\}=0$, and
$\rho_{n_0l,n_0j,\{n_0i\},k}=\rho_{n_0l,n_0j,k-\{n_0i\},k}
=\rho_{n_0l,-n_0j,\{n_0i\},k}=\rho_{n_0l,-n_0j,k-\{n_0i\},k}$
when $\{n_0i\}\neq 0.$
Therefore  we infer from (\ref{4.3}) that
$\lambda_{l,j,i,k}=\lambda_{l,-j,i,k}$  when  $\{n_0i\}=0$, and
$\lambda_{l,j,i,k}=\lambda_{l,j,k-i,k}       =
\lambda_{l,-j,i,k}=\lambda_{l,-j,k-i,k}$ when $\{n_0i\}\neq 0.$
\end{proof}

\begin{proof}[Proof of Theorem 1.1]
It suffices to prove (\ref{1.3}).
To this end, set
$$(l',j',i')=(n_0l,n_0j,\{n_0i\}).$$
Then we consider the
following function  spaces:

{\sl Case $i'=j=0$}:
\begin{equation*}
{\mathbb E}_{l',j',i',k}\subset \dbH^2_{l',j',i',k},\qquad
\tilde{{\mathbb E}}_{l',j',i',k}\subset \tilde \dbH^2_{l',j',i',k};
\end{equation*}


{\sl Case $j\neq  0,\,i'=0$}:

\begin{equation*}
\begin{split}
& {\mathbb E}_{l',j',i',k}\subset  \dbH^2_{l',j',i',k},
\qquad  \tilde{{\mathbb E}}_{l',j',i',k}\subset \tilde \dbH^2_{l',j',i',k},\\
& {\mathbb E}_{l',-j',i',k}\subset \dbH^2_{l',-j',i',k}, \qquad
\tilde{{\mathbb E}}_{l',-j',i',k}\subset \tilde \dbH^2_{l',-j',i',k};
\end{split}
\end{equation*}

{\sl Case $j=0,\,i'>0$}:
\begin{equation*}
\begin{split}
& {\mathbb E}_{l',j',i',k}\subset  \dbH^2_{l',j',i',k}, \qquad
\tilde{{\mathbb E}}_{l',j',i',k}\subset \tilde \dbH^2_{l',j',i',k}, \\
& {\mathbb E}_{l',j',k-i',k}\subset  \dbH^2_{l',j',k-i',k}, \qquad
\tilde{{\mathbb E}}_{l',j',k-i',k}\subset \tilde \dbH^2_{l',j',k-i',k};
\end{split}
\end{equation*}

{\sl Case $ji'\neq 0$}:
\begin{equation*}
\begin{split}
& {\mathbb E}_{l',j',i',k}\subset  \dbH^2_{l',j',i',k}, \qquad
\tilde{{\mathbb E}}_{l',j',i',k}\subset \tilde \dbH^2_{l',j',i',k}, \\
& {\mathbb E}_{l',-j',i',k}\subset  \dbH^2_{l',-j',i',k}, \qquad
\tilde{{\mathbb E}}_{l',-j',i',k}\subset \tilde \dbH^2_{l',-j',i',k}\\
& {\mathbb E}_{l',j',k-i',k}\subset  \dbH^2_{l',j',k-i',k}, \qquad
\tilde{{\mathbb E}}_{l',j',k-i',k}\subset \tilde \dbH^2_{l',j',k-i',k}, \\
& {\mathbb E}_{l',-j',k-i',k}\subset  \dbH^2_{l',-j',k-i',k}, \qquad
\tilde{{\mathbb E}}_{l',-j',k-i',k}\subset \tilde \dbH^2_{l',-j',k-i',k}.
\end{split}
\end{equation*}
It is easy to see that for any
integer vector $(l,j,i,k)$ satisfying (\ref{4-tuple}) and for  each case,
\begin{enumerate}
\item each of these subspaces  contains an
unstable subspace of (\ref{2.1}), generating the bifurcation solutions; and
\item these   subspaces
are orthogonal to each other in $\dbH^2_\sigma$.
\end{enumerate}
Hence (\ref{1.3}), consequently Theorem 1.1,   follows.
\end{proof}

\begin{proof}[Proof of Theorem 1.2]
1. By  Theorems 4.1, 5.1 and 5.2
and (\ref{5.1}-\ref{5.5}),
it remains to show  (\ref{1.4}).
For the integer vector $(l,j,i,k)$ given by
(\ref{4-tuple}), we set
$${\mathbb Y}=\bigg\{ u\in \dbH^2_{l,j,i,k}|\, u=
\sum^{\infty}_{n=1}(\xi_n,\eta_n,\zeta_n)\sin knz\bigg\}.$$
It is easy to see that $u_0=(\sin kz,0,0)$  is the unique steady-state
solution of (\ref{2.1})   in ${\mathbb Y}$.
With this observation in mind,   the cases where $\{n_0i\}=0$ and
$\{n_0i\}\neq 0$ for (\ref{1.4})  are  shown  respectively as follows:

2. In the case
 where $\{n_0i\}=0$  and $j \neq 0$, by (\ref{i-zero}), $i_0=4$.
Then we need to check, for $i_1=0,\,4$,
\begin{equation}\label{5.6}
u_{l,j,i,k,n,\lambda}\neq u_{l,j,i,k,n',\lambda}, \mbox{ if }
\lambda_{l,j,i,k}<\lambda<\infty,\,
i_1+1\leq n\leq i_1+2<n'\leq i_1+4.
\end{equation}


Indeed, when  $i_1=0$ it is obvious that
$$ {\mathbb Y} = \dbH^2_{n_0l,n_0j,0,k}\cap \dbH^2_{n_0l,-n_0j,0,k},$$
and from (\ref{5.1}) and (\ref{5.3}),
$$u_{l,j,i,k,n,\lambda}\in \dbH^2_{n_0l,n_0j,0,k},\,\,\,\,
u_{l,j,i,k,n',\lambda}\in \dbH^2_{n_0l,-n_0j,0,k}.$$
By Assertion 2 of Theorem 1.2,
$$u_{l,j,i,k,n,\lambda}\neq u_0\neq
u_{l,j,i,k,n',\lambda} \mbox{ for } \lambda>\lambda_{l,j,i,k}.$$
This gives (\ref{5.6}) with $i_1=0$.

For  $i_1=4$ we see that
$$ {\mathbb Y} = \tilde  \dbH^2_{n_0l,n_0 j,0,k}
           \cap \tilde \dbH^2_{n_0l,-n_0 j,0,k},
          $$
and from (\ref{5.3}),
$$ u_{l,j,i,k,n,\lambda}\in
           \tilde  \dbH^2_{n_0l,n_0 j,0,k},$$
           $$u_{l,j,i,k,n',\lambda}\in
          \tilde \dbH^2_{n_0l,-n_0 j, 0,k} .
          $$
By Assertion 2 of Theorem 1.2,
$$u_{l,j,i,k,n,\lambda}\neq u_0\neq
u_{l,j,i,k,n',\lambda} \mbox{ for } \lambda>\lambda_{l,j,i,k}.$$
This gives (\ref{5.6}) with $i_1=4$.

3. In the case  where $j\{n_0i\}\neq 0$,  we see that $i_0=16$.
When $i_1=0$, we notice that
\begin{eqnarray*}
{\mathbb Y}  &=&\dbH^2_{n_0l,n_0j,\{n_0i\},k}\cap
\left(
\dbH^2_{n_0l,-n_0j,\{n_0i\},k}\cup
        \dbH^2_{n_0l,-n_0j,k-\{n_0i\},k} \right)\\
 &=&\dbH^2_{n_0l,n_0j,k-\{n_0i\},k}\cap
        \left( \dbH^2_{n_0l,-n_0j,\{n_0i\},k}\cup
        \dbH^2_{n_0l,-n_0j,k-\{n_0i\},k} \right).
 \end{eqnarray*}
From (\ref{5.4}) it follows that
$$ u_{l,j,i,k,n,\lambda}\in
\dbH^2_{n_0l,n_0j,\{n_0i\},k}\cup \dbH^2_{n_0l,n_0 j,k-\{n_0i\},k}
\mbox{ for } 1\leq n\leq 4,$$
and
$$ u_{l,j,i,k,n',\lambda}\in
\dbH^2_{n_0l,-n_0j,\{n_0i\},k}\cup \dbH^2_{n_0l,-n_0 j,k-\{n_0i\},k}
\mbox{ for } 5\leq n\leq 8.$$
By Assertion 2 of Theorem 1.2, we thus have
 $$ u_{l,j,i,k,n,\lambda}\neq
 u_{l,j,i,k,n',\lambda} \mbox{ for } 1\leq n\leq 4<n'\leq 8.$$

When $i_1=8$,  we see that
 \begin{eqnarray*}
 {\mathbb Y} &=&   \tilde\dbH^2_{n_0l,n_0j,\{n_0i\},k}
\cap \left(
        \tilde\dbH^2_{n_0l,-n_0j,\{n_0i\},k}\cup
        \tilde\dbH^2_{n_0l,-n_0j,k-\{n_0i\},k} \right)\\
 &=&\tilde\dbH^2_{n_0l,n_0j,k-\{n_0i\},k}
\cap \left(
        \tilde\dbH^2_{n_0l,-n_0j,\{n_0i\},k}\cup
        \tilde\dbH^2_{n_0l,-n_0j,k-\{n_0i\},k} \right),
\end{eqnarray*}
and, by (\ref{5.4}),
    $$ u_{l,j,i,k,n,\lambda}\in
        \tilde\dbH^2_{n_0l,n_0j,\{n_0i\},k}\cup
        \tilde\dbH^2_{n_0l,n_0j,k-\{n_0i\},k},$$
    $$ u_{l,j,i,k,n',\lambda}\in
        \tilde\dbH^2_{n_0l,-n_0j,\{n_0i\},k} \cup
        \tilde\dbH^2_{n_0l,-n_0j,k-\{n_0i\},k}$$
for $9\leq n\leq 12<n'\leq 16$.
From  Assertion 2 of Theorem 1.2  it thus follows that
$$ u_{l,j,i,k,n,\lambda}\neq
     u_{l,j,i,k,n',\lambda} \mbox{ if } \lambda >\lambda_{l,j,i,k}
     \mbox{ and } 9\leq n\leq 12 <n'\leq 16.$$

Hence, we obtain  (\ref{1.4}) and thus Theorem 1.2    is proved.
\end{proof}


\paragraph{\bf Acknowledgments} 
The authors are grateful to the anonymous referee for his/her
pointing out two gaps, related to the   
algebraic simplicity of the eigenvalues and the supercriticality  
of the bifurcation branches in  the early version of this article, 
and for his/her insightful comments which helped to improve this paper. 
This work  was partially supported by NNSF of China, by the
Office of Naval Research and by the National Science Foundation.


\begin{thebibliography}{99}

\bibitem{AM}Arnold, V. I., Meshalkin, L. D., A. N. Kolmogorov's seminar
on selected problems of analysis (1958-1959).
Russ. Math. Surv. {\bf 15} (1), 247-250 (1960)

\bibitem{BV} Babin, A. V., Vishik, M. I.,
Attractors of partial differential
evolution equations and estimates of their
dimension. Russ. Math. Surv.
{\bf 38}, 151-213 (1983).
\bibitem{CF}Constantin, P., Foias, C., Navier-Stokes Equations,
Chicago;  University of Chicago Press, 1988.

\bibitem{Con}Constantin, P., Foias, C.,  Temam, T.,
On the large time Galerkin approximation of the Navier-Stokes
equations. SIAM J. Num. Anal.
{\bf 21}, 615-634 (1984).


\bibitem{Nir}Caffarelli, L., Kohn, R., Nirenberg, L.,
Partial reqularity of
suitable weak solution of the Navier-Stokes equations.
Comm. Pure. Appl.
Math. {\bf 35}, 771-837 (1982).


\bibitem{CP1}Chen, Z. M., Price, W. G., Time-dependent
periodic Navier-Stokes
flow in a two-dimensional torus. Commun. Math. Phys.
{\bf 179}, 577-597 (1996).

\bibitem{CP2}Chen, Z. M., Price, W. G., Long time behaviour of
Navier-Stokes
flow on a two-dimensional
torus excited by a sinusoidal force. J. Stat. Phys.
{\bf 86}, 301-335 (1997).

\bibitem{FSV} Friedlander, S., Strauss, W., Vishik, M., Nonlinear
instability in an ideal fluid. Ann. Inst. H. Poincar\'e Anal. Non
Lin\'eaire.
{\bf 14}, 187--209 (1997).

\bibitem{Green}Green, J.S.A., Two-dimensional turbulence near the
viscous limit. J. Fluid. Mech. {\bf 62}, 273-287 (1987).


\bibitem{Iu} Iudovich, V. I., Example of the generation of a
secondary stationary or periodic flow when there is loss of
stability of the
laminar flow of a viscous incompressible fluid. J. Math. Mech.
{\bf 29}(4), 587-603 (1965).

\bibitem{Kir}Kirchg\"assner, K., Bifurcation in nonlinear
hydrodynamic
stability. SIAM Review {\bf 17}, 652-683 (1975).


\bibitem{Lin}Lin, C.C., The Theorem of Hydrodynamic Stability,
Cambridge;
Cambridge University Press, 1955

\bibitem{L1} Liu, V. X., A sharp lower bound for the Hausdorff
dimension
of the global attractors of the 2D Navier-Stokes equations.
Commun. Math. Phys.
{\bf 158} 327-339 (1993).

\bibitem{JLL69}
Lions, J. L.,  Quelques M\'ethodes de R\'esolution des
Probl\`emes aux
  Limites Non Lin\'eaires, Dunod, Paris, 1969.

 \bibitem{PLL96}
 Lions, P. J., Mathematical Topics in Fluid Mechanics,
 Oxford science
  Publications, 1996.


\bibitem{Mar}Marchioro, C.,
An example of absence of turbulence for any Reynolds number,
Commun. Math. Phys. {\bf 105}, 99-106 (1986).

\bibitem{MS} Meshalkin, L. D., Sinai, Ya. G., Investigation of
the stability of a stationary solution of a system of
equations for the
plane
movement of an incompressible viscous fluid. J. Math. Mech.
{\bf 19}(9), 1700-1705 (1961).

\bibitem{Obu}Obukhov, A.M., Kolmogorov flow and laboratory
simulation of it.
Russian Math. Surv. {\bf 38}, 113-126 (1983).
\bibitem{Oka}Okamoto, H., Nearly singular two-dimensional
kolmogrov flows
for large Reynolds numbers. J. Dyn. Diff. Equs. {\bf 8}, 203-220 (1996).

\bibitem{Sirov} Platt, N., Sirovich, L., Fitzmaurice, N., An
investigation
of chaotic Kolmogorov flows. Phys. Fluids {\bf A 3 }, 681-696 (1991).

\bibitem{Rab}Rabinowitz, P. H., Some global results for
nonlinear eigenvalue problems. J. Functional Anal. {\bf 7},
487-513 (1971).

\bibitem{Rab68} Rabinowitz, P. H., Existence  and nonuniqueness of
rectangular
solutions of the B\'enard problem, Arch. Rat. Mech. Anal. {\bf 29}, 32-57
(1968).

\bibitem{Rab73} Rabinowitz, P. H., A priori bounds for some
bifurcation
 problems in fluid mechanics, Arch. Rat. Mech. Anal. {\bf 49},
 270-285 (1973).

\bibitem{serrin} Serrin, J., On the stability of viscous
fluid motions. Arch.
Rat. Mech. Anal., {\bf 3} (1959).


\bibitem{Tem} Temam, R., Navier-Stokes Equations. Amsterdam;
North-Holland,
1984.

\bibitem{TX}Tian, G., Xin, Z.,   Gradient Estimation on
Navier-Stokes  Equations, Commun. Geom. Anal. {\bf 7}, 221-257 (1999).

\bibitem{velte} Velte, W., Stabilit\"ats und Verzweigung
station\"arer
L\"osungen
der Navier-Stokesschen Gleichungen, Arch. Rat. Mech. Anal.,
{\bf 16}, 97-125
(1964).

\end{thebibliography}
\end{document}