\documentclass[11pt]{article}

\usepackage{fullpage}
\usepackage{amssymb}
\usepackage{amsmath}
\usepackage{amsfonts}
\usepackage{epsfig} 
%\usepackage{draftcopy}
\usepackage{comment} 
\usepackage{float}     % with epsfig
\restylefloat{figure}  % with epsfig

\input{/home/roger/TeX/mymacros}

\title{Long}
\date{Draft - \today}
\begin{comment}
\author{Roger J. Thelwell\dag, John D. Carter\ddag, Bernard Deconinck\dag \\ \\
\begin{small}
\begin{tabular}{ll}
  \dag Department of Applied Mathematics& \ddag Mathematics Department \\
  \, University of Washington & \,  Seattle University\\
  \, Seattle, WA 98195-2420 & \, Seattle, WA 98122
\end{tabular}
\end{small}}
\end{comment}

\begin{document}

\maketitle
\tableofcontents

\begin{abstract}
\label{s:Abstract}

The two dimensional cubic nonlinear Schr\"{o}dinger equation
admits a large family of one-dimensional bounded traveling wave
solutions.  All such solutions may be written in terms of an 
amplitude and a phase.  Although solutions with piecewise constant 
phase have been investigated, solutions with non-constant phase have 
been less studied.  These non-constant phase solutions are said to have
nontrivial phase.  We consider the spectral stability of 
nontrivial-phase solutions using both asymptotic and numerical methods.  
These methods allow details of the spectrum to be explored, and we do so.
This paper grew out of an earlier short paper 
which reported that all nontrivial-phase solutions considered are unstable
with respect to both one and two-dimensional perturbations.  The current paper  
presents a rich picture of the mechanisms that lead to instability.
We consider the elliptic and hyperbolic nonlinear Schr\"{o}dinger equation in both the focusing and 
defocusing cases.
\end{abstract}


\section{Introduction}
\label{s:intro}

\subsection{Classifying the nonlinear Schr\"{o}dinger equation}
The cubic Nonlinear Schr\"{o}dinger (NLS) equation in two spatial dimensions is
given by
\beq
\label{E:NLS}
i\psi_t + \alpha \psi_{xx} + \beta \psi_{yy} + \gamma |\psi|^2 \psi = 0,
\eeq
for a complex wave function $\psi(x,y)$.  Subscripts in~\eqref{E:NLS} 
are used to denote partial derivatives. Among many other situations, 
NLS arises as an
approximate model of the evolution of a nearly monochromatic wave of
small amplitude in pulse propagation along optical
fibers~\cite{opticswaveguides} where
$\alpha\beta>0$, in gravity waves on deep water~\cite{BR}~\cite{DS}
where $\alpha\beta<0$ and in Langmuir waves in a plasma~\cite{Pecseli}
where $\alpha\beta>0$.  As a description of a
superfluid~\cite{Donnelly}, NLS is
known as the Gross-Pitaevskii equation
\cite{Gross}~\cite{Pitaevskii} with
$\alpha\beta>0$.  Derivations of the NLS equation may be found in Benny \& Roseks~\cite{br69}, Davey and Stewertson (1974) and Djordjevic \& Redekopp (1977). Sulem and Sulem~\cite{sulemsulem} examine NLS in detail.

The NLS equation is said be {\em focusing} or
{\em attractive} in the $x$-dimension if $\alpha \gamma >0$.  If $\alpha \gamma < 0$,
NLS is said to be {\em defocusing} or {\em repulsive} in the 
$x$-dimension.  
Similarly, the sign of $\beta \gamma $ leads to focusing or defocusing in the 
$y$-dimension. The quantity $\alpha \beta$ is also 
used to further classify the problem.  
If $\alpha \beta > 0$, then the NLS equation \eqref{E:NLS}
is called an {\em elliptic} cubic NLS equation, and if $\alpha \beta < 0$, 
it is said to be a {\em hyperbolic} cubic NLS equation.

\subsection{Backround}

The first stability calculations of nonzero NLS solutions were linearizations
about plane waves, also known as Stokes' waves (SW). These solutions  
take the form
\begin{equation}
\psi(x,y,t)=\psi_0\mbox{e}^{ikx+ily+i(\gamma|\psi_0|^2-\alpha k^2-\beta
  l^2)t+i\nu},
\label{E:sw}
\end{equation}
where $\psi_0$, $k$, $l$, and $\nu$ are real parameters. Recongnizing the 
Lie group symmetries~\cite{sulemsulem} of the~\eqref{E:NLS} equation allows one to consider the equivalent 
form
\begin{equation}
\psi = \psi_0 e^{\lambda t},
\label{E:SWls}
\end{equation}
where both $\psi_0$  and $\lambda$ are real constants.
In 1974, Davey and Stewertson~\cite{ds} established that SW solutions are 
stable with respect to perturbations only if both $\alpha<0$ and $\beta<0$, or 
if $\alpha \beta =0$. 

Solitary waves, both the bright soliton (BS) and the dark soliton (DS), are solutions to 
the NLS equation.  The BS solution is given by 
\begin{equation}
\psi = \sqrt(2) \gamma \mbox{sech}(x) e^{it},
\label{E:BS}
\end{equation}
while the DS solution is given by
\begin{equation}
\psi = \sqrt(2) \gamma \tanh(x) e^{i2t},
\label{E:DS}
\end{equation}
again recognizing the Lie group symmetries. Also in 1974, Zakarov and 
Rubenchik~\cite{sr} established that all solitary wave to the NLS equation 
are unstable with respect to long-wavelength transverse perturbations. 
Numerous numerical stability results have been published and reviewed in~\cite{kuzrubzak,rr2,check} and more recently by~\cite{kivpel,JCBDsoliton}.   

More generally, bounded one-dimensional traveling wave 
solutions of the NLS \eqref{E:NLS} equation can be written
\beq
\label{E:tp_soln_form}
\psi(x,t) = \phi(x) exp(i \lambda t + \theta),
\eeq
where $phi(x)$ is a spatially-periodic real valued function and both 
$\lambda$ and $theta$ are real constants.
Solutions of this type are known as trivial phase (TP) solutions, of which  
solitary waves solutions are special cases. 
Martin, Yuen and Saffman~\cite{MYS} were the first?? to 
numerically consider the stability of TP solutions with respect to 
(???).  Other TP solution stability results have been reported in~\cite{FilFraYun,kivpel,JCHStp} and recently in~\cite{JCBDtp}. Results here?  Its hard since the finding are contradictory - or should I report your results?
  
The 1981 paper by Infeld and 
Ziemkiewicz~\cite{infeldz}, is notable in that it considers the stability
of both trivial- and 
nontrivial-phase solutions with respect to perturbations. They examine the damped cubic NLS 
equation of the form
\begin{equation}
\label{E:dNLS}
i\psi_t + \alpha \psi_{xx} + \beta \psi_{yy} + \gamma |\psi|^2 \psi  + \rho \psi = 0,
\end{equation}
and find that damping is necessary to establish stability of the NTP solutions under consideration. Our results~/cite{CDTshort} agree with the work, but 
are extended to consider the entire parameter space of the NTP solutions. 

\subsection{Nontrivial wave solutions}

The NLS equation also admits a large class of one-dimensional traveling wave solutions of
the form
\begin{equation}
\psi(x,y,t)=\phi(ax+by-st)\mbox{e}^{ia\bar{s}x+ib\bar{s}y+i\theta(ax+by-st)+i\lambda
  t+i\nu},
\label{completeform}
\end{equation}
where $\phi$ and $\theta$ are real-valued functions,
$\bar{s}=s/(2\alpha a^2+2\beta b^2)$ and $a$, $b$, $s$ and
$\lambda$ are real parameters.  Once again making use of the symmetries of
the NLS equation~\cite{sulemsulem} and by redefining $\lambda$, 
all solutions of this form can be considered by studying the simplified form
\begin{equation}
\psi(x,t)=\phi(x)\mbox{e}^{i\theta(x)+i\lambda t},
\label{nontrivform}
\end{equation}
where $\phi$ and $\theta$ are real-valued functions and $\lambda$ is
a
real constant. 

Or just something like this?

Nontrivial-phase (NTP) solutions, which represent {\em all} bounded one-dimensional traveling-wave solutions of the NLS equation, take a more general form:
\beq
\label{E:ntp_soln_form}
\psi(x,t) = \phi(x) exp(i \theta(x) + i \lambda t).
\eeq
Here $\phi(x)$ is a real valued function, and $\lambda$ is a real constant. The 
phase contribution, $\theta$, is a real valued function of $x$.  If $\theta(x)$ is a piecewise constant, then the NTP solution is considered a TP solution. 
Several restrictions need to be imposed so that $\psi$ is real valued (traveling wave) solution of the NLS equation.
We require 
\begin{subequations}
\label{E:restrictions}
\begin{eqnarray}
\phi^2(x)&=& \disp{\alpha \left(-2 k^2 \sn^2(x,k)+B \right),}\label{E:phi_def} \\
\theta(x)&=& \disp{c\int_0^x\phi^{-2}(\xi)d\xi,}  \label{E:theta_restrict} \\
\lambda&=& \disp{\frac{1}{2}\alpha (3B-2-2k^2),} \label{E:omega_restrict} \\
c^2&=& \disp{-\frac{\alpha^2}{2} B(B-2k^2)(B-2),} \label{E:c_restrict}
\end{eqnarray}
\end{subequations}
where $c$ is a real constant, and B must take its value from a 
constrained interval to ensure that both $c$ and $\phi$ are real valued.
Here $k \in [0,1]$ is the elliptic modulus of the
Jacobi elliptic sine function, sn$(x,k)$.
The function sn$(x,k)$ is periodic if $k\in[0,1)$, with
period given by $L = 4K$, where $K=K(k)$ is defined by a complete 
elliptic ingtegral of the first kind:
\begin{equation}
K(k)=\int_{0}^{\pi/2}\left(1-k^2\sin^2 x\right)^{-1/2} \, dx.
\end{equation}
When $k=0$, $\sn(x,0) = \sin(x)$ with $L = 2 \pi$.  This fact is later used in the numerical experiments section \ref{s:numerics}. As
$k$ approaches 1, $\sn(x,k)$ approaches $\tanh(x)$ and $L$ approaches
infinity~\cite{byrd}. 

Although $\phi(x)$ inherits the periodicity of $\sn(x,k)$, the solution
$\psi(x,t)$ is
typically {\em not} $L$-periodic in the $x$-dimension, because the periods
of $e^{i\theta}$ and $\phi$ are typically non-commensurate.
An additional constraint on $B$ can be imposed to ensure that 
$\psi$ is $2L$-periodic in the $x$-dimension. 
In addition to the restrictions (\ref{E:restrictions}) above, a second
restriction on $c$ is required:
\begin{equation}
c=\frac{2\pi r}{\int_0^{4K}\phi^{-2}(\tau)d\tau}, \label{E:c_restrict_2}
\end{equation}
where $r$ is any real rational number.  If there are two restrictions on
$c(B)$, then one determines $c$ and the other determines $B$.

Without loss of generality, equation (\ref{E:c_restrict}) can be
rewritten as
\begin{equation}
c=\pm\frac{\alpha}{\gamma\sqrt{2}}\sqrt{-B(B-2k^2)(B-2)},
\label{E:c_restriction_3}
\end{equation}
and equation (\ref{E:c_restrict}) can be rewritten as
\begin{equation}
c=\frac{\alpha}{\gamma}\frac{2\pi r}{\int_0^{4K}(-2k^2\mbox{sn}^2(\tau,k)+B)^{-1}d\tau}.
\end{equation}
In order for these two restrictions on $c$ to be compatible,
\begin{equation}
\frac{2\sqrt{2}\hspace{0.1cm}\pi r}{\int_0^{4K}(-2k^2\mbox{sn}^2(\tau,k)+B)^{-1}d\tau}=\sqrt{-B(B-2k^2)(B-2)},
\label{inteqn}
\end{equation}
where $r$ has been redefined in order to absorb the plus or minus sign.  Given
a rational number $r$, equation (\ref{inteqn}) is an integral equation
that defines $B$.  This equation can be rearranged and written as
\begin{equation}
r=r(k,B)=\frac{1}{2\sqrt{2}\hspace{0.1cm}\pi}\Big{(}\int_0^{4K}(-2k^2 \, \mbox{sn}^2(\tau,k)+B)^{-1}d\tau\Big{)}\sqrt{-B(B-2k^2)(B-2)}.
\label{FR}
\end{equation}
This form allows $r$ to be considered a function of both $k$ and $B$.
Note that a choice of $B$ in (\ref{FR}) may lead to an irrational $r$ (and
therefore a quasi-periodic $\psi$).  In order to obtain a $\psi$ that
is periodic in $x$, one must choose a rational $r$ and then
numerically solve equation (\ref{inteqn}) for the corresponding $B$.
Note that $r$ does not depend on any of the NLS parameters, $\alpha$,
$\beta$, or $\gamma$.

\subsection*{Remark:}
In what follows, the NLS parameter $\gamma$ is fixed as 1, in order to  
simplify the discussion. The remaining
NLS parameters $\alpha$ and $\beta$ are choosen to be $\pm1$.  As a consequence,
we will qualify the NLS equation as;  
\bi
\item focusing in the $x$-dimension if $\alpha = +1$,
\item defocusing in the $x$-dimension if $\alpha = -1$,
\item focusing in the $y$-dimension if $\beta= +1$,
\item defocusing in the $y$-dimension if $\beta = -1$,
\item elliptic if $\alpha \beta = 1$,
\item and hyperbolic if $\alpha \beta = 1$.
\ei

\subsection{$\alpha = +1$}
 
\begin{table}[t]
\bc
\begin{tabular}{|c|c|c|}
\hline
Trivial phase solution &  Non-trivial phase ($\alpha <0$) & Non-trivial phase($\alpha >0$) \\
\hline
\hline
Stokes & $k = 0, B \le 0$ & $k = 0, 0 \leq B \leq 2$ \\
\hline
cn     &    N/A                & $0 \le k < 1,  B = 2 k^2$   \\
\hline
dn     & N/A                   & $0 \le k < 1,  B = 2$   \\
\hline
sn     & $0 \le k < 1,  B = 0$  &  N/A \\
\hline
soliton & Dark $k=1, B=0$      &  Bright $k = 1,  B = 2$   \\
& Grey  $k=1, B<0 $     &   \\
\hline
\end{tabular}
\label{T:tp_limit}
\caption{Parameter values which reduce the nontrivial phase problem to
the trivial phase problem.}
\label{T:tp_limit}
\ec
\end{table}

The phase contribution $\theta(x)$ (defined in~\eqref{E:theta_restrict}) 
implicitly depends on $\alpha$ and $B$ in both~\eqref{E:phi_restrict} 
and~\eqref{E:c_restrict}.  Since $\alpha = +1$, then $c$ is 
real valued if $B \in [2 k^2 , 2]$.  The integral entry of~\eqref{E:theta_restrict}$ \ge 0 $ for all admissible choices of $k$, and therefore $\theta \ge 0$. 
Also, since $\alpha = +1$ and $B \ge 0$, $\phi^2(x) \ge 0$, and 
$\phi(x)$ is real valued. 

\begin{figure}
 \bc
  \begin{tabular}{cc}
   \epsfig{file=./pics/focus_full.eps, width = 5.5 cm} \\
  \end{tabular}
  \caption{Admissible $(k,B)$-parameter space 
           resulting in real valued solution to the NLS $x$-focusing equation.
            The shaded region corresponds to NTP solutions}
  \label{F:f_parameter_space}
 \ec
\end{figure}
 
Figure \ref{F:f_parameter_space} represents the $(k,B)$ parameter 
space corresponding to nontrivial phase solutions of NLS.  As the 
phase $\theta$ approaches zero, the nontrivial-phase solutions limit 
to one of four possible trivial phase solutions; 
({\em i}) Stokes wave, ({\em ii}) Jacobi elliptic $\cn$,
({\em iii}) Jacobi elliptic $\dn$, or    
({\em iv}) soliton type solutions.

The boundaries of Fig. \eqref{F:f_parameter_space} correspond to each of 
these trivial phase solutions.  Table (\ref{T:tp_limit}) summarizes the 
parameter values for which the NTP solutions reduce to trivial phase solutions. 

Since numerous results have been published in the
trivial phase setting, it may be instructive to compare the solution
surfaces from a trivial phase solution a nontrivial soltion which
is 'nearby' in parameter space. In the contour plots of Fig. 
\ref{F:f_solution},
the real part of solution $\psi(x,t)$ has been visualized.
In the gravity wave interpretation, for example,
the $\Re(\psi)$ represents surface amplitude of the wave-form.
The addition of nontrivial phase appears to have a large effect in some
parameter ranges.

\begin{figure}
 \bc
  \begin{tabular}{ccc}
   ?? & \epsfig{file=./pics/fm_soln_surf_1.eps, width = 6 cm} &
   \epsfig{file=./pics/fm_soln_surf_2.eps, width = 6 cm} \\ 
   ?? & \epsfig{file=./pics/fm_soln_surf_3.eps, width = 6 cm} &
   \epsfig{file=./pics/fm_soln_surf_4.eps, width = 6 cm} \\
%   stokes & \epsfig{file=./pics/stokes_wave.eps, width = 6 cm} &
%   \epsfig{file=./pics/near_stokes_waves.eps, width = 6 cm} \\
%   cn & \epsfig{file=./pics/fm_soln_surf_1.eps, width = 6 cm} &
%   \epsfig{file=./pics/near_cn_wave.eps, width = 6 cm} \\
%   dn & \epsfig{file=./pics/fm_soln_surf_2.eps, width = 6 cm} &
%   (Trival-phase) &  (Nontrivial-phase)
  \end{tabular}
  \caption{ This is what I want: Focusing solutions.  
            Contour plots in space and time for solutions with trivial-phase 
            and solutions with nontrivial-phase but which are `close' to the 
            trivial-phase solutions. I could also include the phase component }
  \label{F:f_solution}
 \ec
\end{figure}



\subsection{$\alpha = -1$} 

The phase contribution $\theta(x)$ (defined in~\eqref{E:theta_restrict}) 
implicitly depends on $\alpha$ and $B$ in both~\eqref{E:phi_restrict} 
and~\eqref{E:c_restrict}. If $\alpha = -1$, $c \ge 0$ if $B \le 0$ or $B>2k^2$.
However, the fact that $\alpha = -1$ and restriction~\cite{E:phi_restriction} that $2k^2 \mbox{sn}(x,k) + B \le 0$ implies that $psi$ is  
real valued only if $B \le 0$ if $\alpha = -1$. 

\begin{figure}
 \bc
  \begin{tabular}{cc}
   \epsfig{file=./pics/defocus_full.eps, width = 6 cm} &
  \end{tabular}
  \caption{Admissible $(k,B)$-parameter space resulting in real valued 
           solution to the NLS $x$-defocusing equation. The shaded region
           corresponds to NTP solutions}
  \label{F:df_parameter_space}
 \ec
\end{figure}

The shaded region in Fig. \ref{F:df_parameter_space} represents 
the $(k,B)$ parameter 
space corresponding to nontrivial phase solutions of NLS.  As the 
phase $\theta$ approaches zero, the nontrivial phase solutions should limit 
to one of the three possible trivial phase solutions; 
({\em i}) Stokes wave, 
({\em ii}) Jacobi elliptic $\sn$, or   
({\em iii}) a dark soliton.
It is easily seen that the boundaries of the regions in figure 
\eqref{F:parameter_space} correspond to each of these trivial
phase solutions.  Table (\ref{T:tp_limit}) summarizes the 
parameter values for which the NTP solutions reduce to trivial phase solutions. 
Since numerous results have been published in the
trivial phase setting, it may be instructive to compare the solution
surfaces from a trivial phase solution a nontrivial soltion which
is 'nearby' in parameter space. In the contour plots of Fig. 
\ref{F:f_solution},
the real part of solution $\psi(x,t)$ has been visualized.
In the gravity wave interpretation, for example,
the $\Re(\psi)$ represents surface amplitude of the wave-form.
The addition of nontrivial phase appears to have a large effect in some
parameter ranges.

\begin{figure}
 \bc
  \begin{tabular}{ccc}
   ?? & \epsfig{file=./pics/fm_soln_surf_1.eps, width = 6 cm} &
   \epsfig{file=./pics/fm_soln_surf_2.eps, width = 6 cm} \\ 
   ?? & \epsfig{file=./pics/fm_soln_surf_3.eps, width = 6 cm} &
   \epsfig{file=./pics/fm_soln_surf_4.eps, width = 6 cm} \\
%   stokes & \epsfig{file=./pics/stokes_wave.eps, width = 6 cm} &
%   \epsfig{file=./pics/near_stokes_waves.eps, width = 6 cm} \\
%   cn & \epsfig{file=./pics/fm_soln_surf_1.eps, width = 6 cm} &
%   \epsfig{file=./pics/near_cn_wave.eps, width = 6 cm} \\
%   dn & \epsfig{file=./pics/fm_soln_surf_2.eps, width = 6 cm} &
%   \epsfig{file=./pics/near_dn_wave.eps, width = 6 cm} \\
%   sn & \epsfig{file=./pics/sn_wave.eps, width = 6 cm} &
%   \epsfig{file=./pics/near_sn_wave.eps, width = 6 cm} \\
    & Trivial phase &  Nontrivial-phase  
  \end{tabular}
  \caption{This is what I want: Defocusing solutions. 
            Contour plots in space an time for solutions with trivial-phase 
            and solutions with nontrivial-phase but which are `close' to the 
            trivial-phase solutions.  Animations which capture transition from
            between the trivial-phase solutions may be found in the supplimental
            material.}
  \label{F:parameter_space}
 \ec
\end{figure}

\section*{}
This paper explores the stability 
of nontrivial-phase solutions under one- and two dimensional perturbations.
To do so, the required linear problem is quickly developed in this section. 
In later sections (ref{s:perturbation,s:numerics} the stability of this linear problem will be explored
using both perturbative and numerical  methods.  The results
of the numerics may be found in~\cite{CDTntp}. 
While there is some overlap in the preliminary sections, much more detail will be presents here.

\section{The linear stability problem}

In order to study the linear stability of NTP solutions
of the NLS equation, we consider perturbed solutions of the form
\begin{equation}
\psi_{_p}(x,y,t)=(\phi(x)+\epsilon u(x,y,t)+i \epsilon
v(x,y,t) + \mathcal{O}(\eps^2)) \, e^{i\theta(x) + i\lambda t},
\label{E:npert}
\end{equation}
where $u(x,y,t)$ and $v(x,y,t)$ are real-valued functions, $\epsilon$ is a
small real parameter and $\phi(x) \, e^{i\theta(x)+i\lambda t}$ is
an NTP solution of NLS.
Substituting (\ref{E:npert})
in (\ref{E:NLS}), linearizing and separating real and imaginary
parts leads to
\begin{subequations}
\begin{equation}
\lambda u-3\gamma\phi^2u-\beta
u_{yy}+ \alpha c^2 \frac{1}{\phi^4}u-2\alpha c\frac{1}{\phi^{3}}\phi_x v+2\alpha
c\frac{1}{\phi^{2}}v_x-\alpha u_{xx}=-v_t,\label{nreale}
\end{equation}
\begin{equation}
\lambda v-\gamma\phi^2v-\beta v_{yy}+\alpha c^2\frac{1}{\phi^{4}}v+2\alpha
c\frac{1}{\phi^{3}}\phi_x u-2\alpha c\frac{1}{\phi^{2}}u_x-\alpha v_{xx}=u_t.\label{nimage}
\end{equation}
\label{nevaluepde}
\end{subequations}
Since \eqref{nevaluepde} does not
depend on $y$ or $t$ explicitly, we may assume that $u(x,y,t)$ and
$v(x,y,t)$ have the forms
\begin{subequations}
\begin{equation}
u(x,y,t)=U(x,\rho,\Omega) \, e^{i\rho y+\Omega t}+c.c.,
\end{equation}
\begin{equation}
v(x,y,t)=V(x,\rho,\Omega) \, e^{i\rho y+\Omega t}+c.c.,
\end{equation}
\label{uvnontrivial}
\end{subequations}
where $\rho$ is a real constant, $U(x)$ and $V(x)$ are complex-valued
functions,  $\Omega$ is a complex constant and $c.c.$ denotes complex
conjugate.  Notice that $\rho$ is the transverse wavenumber of the
perturbation and $\Omega$ is the exponential growth rate associated
with $\rho$.
If bounded $U,V$ exist such that $\Omega$ has a positive real part,
then the amplitudes of the perturbations grow
exponentially in time and the unperturbed solution is unstable.
If the real part of $\Omega$ is negative, then the
amplitudes of the perturbations decay in time and the unperturbed
solution is said to be 
linearly asymptotically stable.  If $\Omega$ is purely imaginary, then the
perturbations oscillate in time and the unperturbed solution is said
to be linearly (neutrally) stable.

Upon substitution, (\ref{nevaluepde}) yields the spectral problem
\begin{subequations}
\label{E:neval}
\begin{equation}
\lambda U-3\gamma\phi^2U+\beta\rho^2
U+\alpha c^2\frac{1}{\phi^{4}}U-2\alpha c\frac{1}{\phi^{3}}\phi_xV+2\alpha
c\frac{1}{\phi^{2}}V_x-\alpha U_{xx}=-\Omega V,
\end{equation}
\begin{equation}
\lambda V-\gamma\phi^2V+\beta\rho^2 V+\alpha c^2\frac{1}{\phi^{4}}V+2\alpha
c\frac{1}{\phi^{3}}\phi_x U-2\alpha c\frac{1}{\phi^{2}}U_x-\alpha
V_{xx}=\Omega U.
\end{equation}
\end{subequations}

Written in operator form, (\ref{E:neval}) is
\begin{subequations}
\label{E:opneval}
\begin{equation}
(L_1+\beta\rho^2)U+SV=-\Omega V,\label{nL1U}
\end{equation}
\begin{equation}
(L_0+\beta\rho^2)V-SU=\Omega U,\label{nL0V}
\end{equation}
\label{nevalopform}
\end{subequations}
where
\begin{subequations}
\label{E:opsLandS}
\begin{equation}
L_1=\lambda-3\gamma\phi^2+\alpha
c^2\frac{1}{\phi^{4}}-\alpha\partial_{xx},\label{E:nL1}
\end{equation}
\begin{equation}
L_0=\lambda-\gamma\phi^2+\alpha
c^2\frac{1}{\phi^{4}}-\alpha\partial_{xx},\label{E:nL0}
\end{equation}
\begin{equation}
S=-2\alpha c\frac{1}{\phi^{3}}\phi_x+2\alpha c \frac{1}{\phi^{2}}\partial_x.
\end{equation}
\label{nL0L1S}
\end{subequations}

We can write this linear system of equations in terms of a single linear operator,
$\mathcal{L}$, where  
\begin{equation}
\mathcal{L}(k,B)[u] = \Omega u
\label{E:opform}
\end{equation}
and
\begin{equation}
\mathcal{L}(k,B)=  \begin{bmatrix}-S & L_0 \\ L_1 & S\end{bmatrix}
\label{E:opeq}
\end{equation}
and 
\begin{equation}
u =  \begin{bmatrix} U \\ V\end{bmatrix}
\label{E:opeq}
\end{equation}
If $c=0$, then $S=0$ in \eqref{E:opform} and the problem reduces
to the stability problem for TP solutions. 
 
We now consider the stability problem of the NTP solutions to the NLS equation, 
applying perturbative methods to the system~\eqref{E:neval} (or ~\eqref{E:opneval}).

\section{Perturbation approach to stability}
\label{SS:perturbation}

John's section

\subsection{Small $\rho$}

\subsection{Large $\rho$}

Here is the old stuff that I had already.

Generalizing the work of~\cite{JCBDsolition} in the soliton 
setting, we
construct perturbations to explore the stability of the nontrivial phase 
solutions in the large $\rho$ limit.
For large $\rho$, we consider series of the form
\begin{subequations}
\begin{align}
U(x)&\sim u_0(\eta x) + \rho^{-1} u_1( \eta x) + \rho^{-2} u_2(\eta x) + \ldots\\
V(x)&\sim v_0(\eta x) + \rho^{-1} v_1( \eta x) + \rho^{-2} v_2(\eta x) + \ldots \\
\Omega&\sim \omega_{-2}\rho^2 + \omega_{-1}\rho + \ldots \label{E:omega_series}\\
\eta&\sim \eta_{-1} \rho + \eta_0  \label{E:eta_series}
\end{align}
\end{subequations}
where $u_j$ and $v_j$ are complex valued functions, $\omega_j$ is a complex 
constant, and $\eta_j$ is a real constant. 

\begin{figure}
\label{F:rhoim}
\bc
\begin{tabular}{cc}
\epsfig{file=./pics/fp3001109.rhoim.eps, width = 5cm}  &
\epsfig{file=./pics/dp3001109.rhoim.eps, width = 5cm} \\
(a) & (b)
\end{tabular}
\ec
\caption{Projection of spectrum 
onto the $(\rho,\Im(\Omega))$ plane for (a) $\alpha \beta >0$ and (b) $\alpha \beta<0$.} 
\end{figure}

Substituting into the linear problem generates 
\begin{subequations}
\beq
\alpha \eta_{-1}^2 u_0^{''} - \beta u_0 - \omega_{-2} v_0 = 0 
\eeq
\beq
\alpha \eta_{-1}^2 v_0^{''} - \beta v_0 + \omega_{-2} u_0 = 0
\eeq
\end{subequations}
at leading order.  Transforming to a first order system generates
\beq
\label{E:asy_1}
\begin{pmatrix} u_0 \\ u_0' \\ v_0 \\ v_0' \end{pmatrix} ^\prime
= \begin{pmatrix} 
0 & 1 & 0 & 0 \\
\frac{\beta}{\alpha \eta_{-1}^2}& 0 &\frac{\omega_{-2}}{\alpha\eta_1^2} & 0 \\
0 & 0 & 0 & 1 \\
\frac{\omega_{-2}}{\alpha\eta_{-1}^2} & 0 &\frac{\beta}{\alpha\eta_1^2} & 0 
\end{pmatrix}
\begin{pmatrix} u_0 \\ u_0' \\ v_0 \\ v_0' \end{pmatrix}.
\eeq
The eigenvalues $\lambda$ of \eqref{E:asy_1} satisfy
\beq
\left(\lambda^2 - \frac{1}{\alpha\eta_{-1}^2}\left(\beta+i\omega_{-2}\right) \right)
\left(\lambda^2 - \frac{1}{\alpha\eta_{-1}^2}\left(\beta-i\omega_{-2}\right) \right) =0.
\eeq
If the real part of $\omega_{-2}$ is non-zero, then $\Re(\lambda)$ will also be
non-zero. We seek bounded solutions, and so assume that $\omega_{-2}$ is purely
imaginary. Since it is easier to first find all unbounded solutions and then 
exclude them, we do so.  By requiring 
\[
\frac{1}{\alpha} \left(\beta+i\omega_{-2} \right) > 0 
\mbox{ \quad \bf{AND} \quad }
\frac{1}{\alpha} \left(\beta-i\omega_{-2} \right) < 0 
\]
a parameter space corresponding to unbounded solutions may be found.
If $\frac{\beta}{\alpha} > 0$
then there are no bounded eigenfunctions if $-\beta<i \omega_{-2} < \beta$. 
Since $\Omega \sim \omega_{-2} \rho^2$, we conclude that there will be no 
eigenvalues (at first order) of the form \eqref{E:omega_series} with non-zero 
real part in the region between $\rho^2$ and $-\rho^2$. 
If $\frac{\beta}{\alpha} < 0$
then we need to satisfy $\beta > i \omega_{-2}$ and $\beta < -i \omega_{-2}$, 
which is impossible.  As a consequence, the spectrum will will apparently 
occupy the full complex sheet.\\

I stopped here.

or just the summary:\\

\subsection{Summary of results}
\label{ss:perturbation_results}

These assumptions can be shown to imply that, arguing first order contributions
, the cases in which $ \alpha \beta > 0$ should have a spectrum 
whose imaginary part must not lie in the region bounded between 
$\rho^2$ and $-\rho^2$ in the $(\rho,\Im(\Omega))$ plane. Continuing
the analysis to second order leads to the conclusion that $\omega_{-1} = 0$. 
In the case where $\alpha \beta<0$, it is possible to show that 
no region of unbounded growth may be excluded by the assumptions.  This suggests
that the projection onto the spectrum onto the $(\rho,\Im(\Omega))$ plane will
appear full. These claims are supported by the plots of 
figure~\eqref{F:rhoim}, which are qualitatively similar to the plots for 
all admissible choices of $k$ and $B$. 


\section{Numerics}
\label{s:Numerics}

Using the linear system \eqref{E:neval}, we were able to
investigate the stability of all NTP solution numerically, as reported 
in~\cite{CDTshort}, where we established the instbility of all NTP solutions. 
In this section first develope the numerical framework, and then use it to
 examine several fundemental features
that were not apparent from the summary data. Perhaps the most striking feature 
of the instability properties of the NTP solution is the rapid appearance of
exponential growth rates which were 5 to 7 times larger than growth rates of 
nearby sn- and cn-type TP solutions. This in turn lead to the exploration of the 
bifurcations responsible for this rapid growth effect. 
We begin by presenting our numerical framework, which is used to produce
detailed figures of our approximations to the continous spectra 
of the linear stability operator $\mathcal{L}$. 

\subsection{Hill's method}
\label{SS:numerics}

Hill's method is well suited to analyzing the continuous spectra 
of the linear stablity operator $\mathcal{L}$ of~\cite{E:opform}:  it uses the  
periodicity of the coefficient functions in~\cite{E:neval} to greatly 
reduce the computational effort needed to systematically explore the large 
phase space encountered.

To apply Hill's method, Fourier expansions are needed for all coefficient
functions of \eqref{E:neval}.  Using the complex Fourier form, we have
\begin{equation}
\label{E:phihat}
\begin{array}{rr}
\displaystyle{\phi^2(x) = \sum_{n = -\infty}^{\infty} Q_n \, e^{i 2 n \pi x / L},} &  \qquad
\displaystyle{\phi^{-2}(x) = \sum_{n = -\infty}^{\infty} R_n \, e^{i 2 n \pi x / L},} \\
\displaystyle{\phi^{-4}(x) = \sum_{n = -\infty}^{\infty} S_n \, e^{i 2 n \pi x / L},} & \qquad
\displaystyle{\phi^{-3}(x) \phi'(x) = \sum_{n = -\infty}^{\infty} T_n \, e^{i 2
n \pi x / L}},
\end{array}
\end{equation}
where $Q_n,R_n,S_n$ and $T_n$ are Fourier coefficients.
Note that $\phi^2(x)$ has period $L/2$ and that $\phi(x)$ is never zero except in some TP limit cases.

The periodicity of the coefficients containing powers of $\phi$ allows us 
to decompose the eigenfunction components $U$ and $V$ of~\ref{E:neval} 
in the Fourier-Floquet form
%\begin{subequations}
%\label{E:UVhat}
\beq
\label{E:UVhat}
U(x) := e^{i \mu x} \sum_{l=-\infty}^{\infty} U_l e^{-i l \pi x/ PL}   
\mbox{\quad and \quad}
V(x) := e^{i \mu x} \sum_{l=-\infty}^{\infty} V_l e^{-i l \pi x/ PL}.   
\eeq
%\end{subequations}
Here $\mu$ is a Floquet parameter and the factor of $P$ is included 
for mathematical convenience.  The form \eqref{E:UVhat} results 
from an application of Floquet's theorem and then noting that 
we seek eigenfunctions, which are by definition bounded.  Again, 
see~\cite{BDNKhill} for a more complete explanation.

Substitution of \eqref{E:phihat} and \eqref{E:UVhat} into \eqref{E:nevalopform} 
allows us to write equations for $U_n$ and $V_n$ as a coupled 
bi-infinite system of difference equations given by 

\begin{subequations}
\begin{multline}
\label{E:ndL1}
-\left(\lambda + \beta \rho^2 - \alpha\left(i\mu+ \frac{i n \pi}{LP}\right)^2 \right) 
     U_n
+ 3 \gamma 
         \sum_{k = -\infty}^{\infty} Q_{\frac{n-k}{P}} U_k
- \alpha c^2  
         \sum_{k = -\infty}^{\infty} S_{\frac{n-k}{P}} U_k \\
 + 2 \alpha c 
         \sum_{k = -\infty}^{\infty} T_{\frac{n-k}{P}} V_k
- 2 \alpha c \left(i\mu + \frac{i n \pi}{PL}\right) 
          \sum_{k = -\infty}^{\infty} R_{\frac{n-k}{P}} V_k 
    = 
        \Omega V_n 
\end{multline}
\begin{multline}
\label{E:ndL0}
\left(\lambda + \beta \rho^2 - \alpha\left(i\mu+ \frac{i n \pi}{LP}\right)^2 \right) 
        V_n
- \gamma
         \sum_{k = -\infty}^{\infty} Q_{\frac{n-k}{P}} V_k
+ \alpha c^2  
         \sum_{k = -\infty}^{\infty} S_{\frac{n-k}{P}} V_k \\
+ 2 \alpha c 
         \sum_{k = -\infty}^{\infty} T_{\frac{n-k}{P}} U_k
- 2 \alpha c \left(i\mu + \frac{i n \pi}{PL}\right) 
          \sum_{k = -\infty}^{\infty} R_{\frac{n-k}{P}} U_k 
    = 
        \Omega U_n, 
\end{multline}
\end{subequations}
which hold for all $n$. Here  $\mu \in [\frac{-\pi}{PL},\frac{\pi}{PL}]$ and
$Q_\frac{n-k}{P} = 0$ if $\frac{n-k}{P} \not\in \mathbb{Z}$. 
 
In practice, a 
pre-multiplication of the linear system by $\phi^4$ results
in more rapid convergence of the Fourier coefficients. Also,   
an exact cosine series expansion of $\phi^2, \phi^4$ and 
$\phi^6$ can be used, which follows from the Jacobi's~\cite{Jacobi} 
series expansion 
of $\sn^2$. This pre-multiplication transforms the original 
eigenvalue problem into a generalized eigenvalue problem. Golub and Van Loan~\cite{GolubVanLoan} provide a brief discussion of generalized eigenvalue problems. 


Notice that equations [\ref{E:ndL1}-\ref{E:ndL1}] are  
{\em equivalent} to the original system [\ref{E:nL1}-\ref{E:nL0}]. 
Only when truncated to a finite number of modes does this system 
become an approximation.  If, for example, we take $N$ to be maximimal 
mode number in the $U_n$ and $V_n$ approximations, then 
the exact coupled bi-infinite system reduces to a $(4N+2)$ system of equations.

keep this in?

Work on this:
Gershgorin's theorem suggests that the net effect of truncation also behaves 
spectrally. The exponential decay of Fourier coefficients on sub-diagonals 
overwhelms the polynomial growth effect of the Floquet multiplier 
contribution in the columns.  The eigenvalues that result from the $(4N+2)$ 
eigenvalue problem
should converge spectrally in $N$ to eigenvalues of the continuous system.
Once the off-diagonal terms reach machine precision, an increase of $N$ 
results in marginal increase of spectral information at large 
computational cost.  Increasing $N$ linearly increases:
\bi
\item computation time cubically, 
\item the number of eigenvalues linearly,
\item and the accuracy of these values spectrally. 
\ei

\subsubsection{Numerical experiments} 
\label{SSS:Experiments}

$\alpha > 0$ 
Rework this

A large number of parameter values have been explored numerically.  
Approximately 5.2 million generalized eigenvalue problems were 
considered, the size of each determined by the cutoff
mode of the underlying Fourier series.  A truncation to $N$ positive Fourier 
coefficients results in a $(4N+2)\times(4N+2)$ dimensional problem. 
The value of N as a function of $(k,B)$ was chosen by computing sample problems 
over (various $k$'s and $B$'s) in the parameter region for increasingly large 
$N$ until the eigenvalues that resulted stabilized qualitatively.
A simple polynomial was constructed to fit this data.
This polynomial, as well as other parameter ranges used in the 
numerical experiments, is given in Table \ref{T:p_vals},
where $\mathtt{linspace(a,b,m)}$ is a linearly spaced vector from $a$ to 
$b$ of length $m$, and {\tt ceil} is the ceiling function.

\begin{table}{t}
\caption{Parameter values and ranges used in numerical experiments.}
\bc
\label{T:p_vals}
\begin{tabular}{|c|c|l|}
\hline
\label{T:tp_limit}
Parameter & Description & \qquad value \\
\hline
\hline
$k$ & Elliptic Modulus & \, \tt{linspace(0,1,65)} \\ 
\hline
$B$ & Shift & \begin{tabular}{l}
              For $\alpha <0: \quad \mathtt{linspace(-1,0,65)} $ \\ 
              For $\alpha >0: \quad \mathtt{linspace(0,2,65)} \cap (2k^2 , 2)$
              \end{tabular} \\ 
 
\hline
$N$ & Fourier cutoff & \begin{tabular}{l}
                        For $\alpha<0: \quad 15 + \mathtt{ceil(5 k^{5})}$ \\ 
                        For $\alpha>0: \quad 10 + \mathtt{ceil(25 k^{10})}$ 
                        \end{tabular} \\
\hline
$\rho$ & perturbation wavenumber & \, \tt{linspace(0,4,65)} \\
\hline
$\mu$ & Floquet parameter & \, $\mathtt{linspace(-\frac{\pi}{2K},\frac{\pi}{2K},65)}$  \\
\hline
\end{tabular}
\ec
\end{table}

Four cases were considered numerically:
\bi
\item 
 Defocusing in $x$ with focusing perturbation in $y$. 
$(\alpha<0 \mbox{ and } \beta > 0$)
\item 
 Defocusing in $x$ with defocusing perturbation in $y$. 
$(\alpha<0 \mbox{ and } \beta < 0$)
\item 
 Focusing in $x$ with focusing perturbation in $y$.
$(\alpha>0 \mbox{ and } \beta > 0$)
\item 
 Focusing in $x$ with defocusing perturbation in $y$.
 $(\alpha>0 \mbox{ and } \beta < 0$)
\ei

\subsubsection{Numerical results} 
\label{SSS:Num_results}

Using the Fourier-Floquet-Hill method we numerically considered the growth 
instabilities due to transverse
perturbations over the range of parameter values of table \eqref{T:p_vals}.  
Each numerical experiment consisted of computing the 
eigenvalues of a (generalized) eigenvalue problem.  For each parameter triplet 
$(k,B,\rho)$, for $\alpha,\beta = \pm 1$, a sequence of Floquet parameters 
$\mu \in [\frac{-\pi}{PL},\frac{\pi}{PL}]$ the  
corresponding eigenvalues and eigenvectors were computed.  The corresponding 
eigenfunction may be reconstituted from this information.

The entire set of eigenvalues over for each $(k,B,\rho,\mu)$ quadtuplet.
In~\cite{JCBDRTshort}, this large data set was reduced by considering 
the largest real part of all eigenvalues
generated for each $(k,B)$  pair.  Since a single eigenvalue with postive real
part will lead to (linear) instability, the summary information was useful to
dicuss linear stablility.  However, several interesting observations were 
made from this summary information that led to this more detailed study.
In reducing the dimension, a tremendous amount of information and much 
of the physical structure was lost.    

\subsubsection{Cn to Dn series $k = 0.9, B from 2k^2 -> 2$}
The panels of figure \eqref{F:cn_dn} represent the spectrum restricted
such that the imaginary component is less than five in magnitude. This
reduced three dimensional data set, consisting of
$(\Re(\sigma)), \Im(\sigma),\rho)$ was then projected onto the;(a)complex plane
(b)$(\rho,\Re(\sigma))$ plane, and (c)$(\rho,\Im(\sigma))$ plane.  The plots
indicating $\Re(\sigma)$ are relevent for stability discussions, although the
$(\rho,\Im(\sigma))$ appears relevent in the discussion of the bifurcations
that occur over a large range of $B$ values.  The upper panels are close to
the $\cn$ solution, while the lower panel is close to the $\dn$ solution.
In between, however, the spectrum is quite different qualitatively, with the
large number of both growing and decaying solutions that arise as the spectral
bands in the complex plane separate and collide.
                                                                                
\begin{figure}
\bc
\begin{tabular}{ccc}
\epsfig{file=./pics/fp3001073.reim.eps, width = 5cm} &
\epsfig{file=./pics/fp3001073.rhore.eps , width = 5cm} &
\epsfig{file=./pics/fp3001073.rhoim.eps , width = 5cm} \\
\epsfig{file=./pics/fp3001089.reim.eps, width = 5cm} &
\epsfig{file=./pics/fp3001089.rhore.eps , width = 5cm} &
\epsfig{file=./pics/fp3001089.rhoim.eps , width = 5cm} \\
\epsfig{file=./pics/fp3001097.reim.eps, width = 5cm} &
\epsfig{file=./pics/fp3001097.rhore.eps , width = 5cm} &
\epsfig{file=./pics/fp3001097.rhoim.eps , width = 5cm} \\
\epsfig{file=./pics/fp3001105.reim.eps, width = 5cm} &
\epsfig{file=./pics/fp3001105.rhore.eps , width = 5cm} &
\epsfig{file=./pics/fp3001105.rhoim.eps , width = 5cm} \\
\epsfig{file=./pics/fp3001113.reim.eps, width = 5cm} &
\epsfig{file=./pics/fp3001113.rhore.eps , width = 5cm} &
\epsfig{file=./pics/fp3001113.rhoim.eps , width = 5cm} \\
(a) & (b) & (c) \\
\end{tabular}
\caption{ Projections of the spectrum in the case $\alpha = 1$ and $\beta = 1$
(focusing in $x$ and focusing in $y$) for k = 0.9 and  B = [ ?? check] onto
(a) the complex plane $(\Re(\sigma),\Im(\sigma))$, (b)  ($\Re(\sigma),\rho)$ and (c) $(\Im(\sigma),\rho)$ plane.}
\label{F:cn_dn}
\ec
\end{figure}
                        
\subsubsection{Sn series $k = 0.9, B -> 0 to -1$}
                                                                               
The panels of figure \ref{F:sn} represent the spectrum restricted
such that the imaginary component is less than five in magnitude in the
defocusing case.  Recall that a trivial phase solution occurs as $B$ limits to
zero from below.  The value of $B$ is closest to 0 in the top panel, and becomes
increasing negative as the plots continue.  As in the previous section,
the reduced data set is projected onto
the (a) complex plane
(b) $(\rho,\Re(\sigma))$ plane, and (c) $(\rho,\Im(\sigma))$ plane.
Again the stability information is contained in the plots of the
$\Re(\sigma)$, while the $(\rho,\Im(\sigma))$ appears to indicate that
the movement of the purely imaginary eigenvalues play an important role in
determining the bifurcation structure.
\begin{figure}
\bc
\begin{tabular}{ccc}
\epsfig{file=./pics/dp3001073.reim.eps, width = 5cm} &
\epsfig{file=./pics/dp3001073.rhore.eps , width = 5cm} &
\epsfig{file=./pics/dp3001073.rhoim.eps , width = 5cm} \\
\epsfig{file=./pics/dp3001089.reim.eps, width = 5cm} &
\epsfig{file=./pics/dp3001089.rhore.eps , width = 5cm} &
\epsfig{file=./pics/dp3001089.rhoim.eps , width = 5cm} \\
\epsfig{file=./pics/dp3001097.reim.eps, width = 5cm} &
\epsfig{file=./pics/dp3001097.rhore.eps , width = 5cm} &
\epsfig{file=./pics/dp3001097.rhoim.eps , width = 5cm} \\
\epsfig{file=./pics/dp3001105.reim.eps, width = 5cm} &
\epsfig{file=./pics/dp3001105.rhore.eps , width = 5cm} &
\epsfig{file=./pics/dp3001105.rhoim.eps , width = 5cm} \\
\epsfig{file=./pics/dp3001113.reim.eps, width = 5cm} &
\epsfig{file=./pics/dp3001113.rhore.eps , width = 5cm} &
\epsfig{file=./pics/dp3001113.rhoim.eps , width = 5cm} \\
(a) & (b) & (c) \\
\end{tabular}
\caption{Projections of the spectrum in the case $\alpha = -1$ and $\beta = 1$
(defocusing in $x$ and focusing in $y$) with $k = 0.9$ and B = [ ?? check] onto                                                                               
(a) the complex plane $(\Re(\sigma),\Im(\sigma))$, (b) ($\Re(\sigma),\rho)$ and (c) $(\Im(\sigma),\rho)$ plane.}
\ec
\label{F:sn}
\end{figure}

\subsubsection{Grey Soliton $k -> 1, B -> -1$ }

\subsubsection{Stokes limit $k -> 0 to 0.1 B \neq 0$}

We note that Carter and Segur \ref{JCHS} report that stability of the 
trivial phase problem should increase monotonically in the stokes 
limit of $k=0$.  Although difficult to distinguish from the plots of 
figure ~\ref{F:all_growth}, we include figure~\ref{F:logk_growth} 
which supports this claim.  In fact, in all cases, the minimum growth 
rate is found nearest to the stokes limit of $k=0$.  It is interesting 
to note is that the stability factor does not grow monotonically for 
any other choice of $k \in (0,1)$.  Instead, the maximal growth rate 
is achieved at some point in the interior of the admissible parameter space.    
 

\section{Conclusions}

The linear stability of nontrivial phase solutions to the cubic
 Nonlinear Schr\"{o}dinger equation was explored. 
A numerical method was presented that provided eigenvalue information over a 
wide range of parameter values.  Plots were presented that reported
growth factors and wave-numbers resulting in maximal instability. 
A rapid increase in the rate of instability resulted from   
the addition of nontrivial phase $\theta(x)$.  A movement through the shifted $B$ space   
indicated that values which correspond to the largest growth 
 admits a class of Stokes' wave (plane wave) solutions of the form
\begin{equation}
\psi(x,y,t)=\psi_0\mbox{e}^{ikx+ily+i(\gamma|\psi_0|^2-\alpha k^2-\beta
  l^2)t+i\nu},
\label{swsoln}
\end{equation}
where $\psi_0$, $k$, $l$, and $\nu$ are real parameters.  NLS also
admits a large class of one-dimensional traveling wave solutions of
the form
\begin{equation}
\psi(x,y,t)=\phi(ax+by-st)\mbox{e}^{ia\bar{s}x+ib\bar{s}y+i\theta(ax+by-st)+i\lambda
  t+i\nu},
\label{completeform}
\end{equation}
where $\phi$ and $\theta$ are real-valued functions,
$\bar{s}=s/(2\alpha a^2+2\beta b^2)$ and $a$, $b$, $s$ and
$\lambda$ are real parameters.  By making use of the symmetries of
NLS~\cite{sulemsulem} and by redefining $\lambda$, all solutions of this
form can be considered by studying the simplified form
\begin{equation}
\psi(x,t)=\phi(x)\mbox{e}^{i\theta(x)+i\lambda t},
\label{nontrivform}
\end{equation}
where $\phi$ and $\theta$ are real-valued functions and $\lambda$ is a
real constant.  If $\theta^{\prime}(x)=0$, then the solution is said to have
trivial phase (TP) and if $\theta^{\prime}(x)\ne 0$, then the
solution is said to have nontrivial phase (NTP).  Substituting
(\ref{nontrivform}) into (\ref{NLSE}) and separating into real and
imaginary parts gives

generated in the fully 
nontrivial setting.  Additionally, perturbation results are presented 
which compliment the numerical evidence in the case of large $\rho$ 
(long wavelength) perturbation.  
In summary, numerical and perturbative results suggest that bounded, fully 
nontrivial one dimensional solutions to the cubic NLS equation are 
unstable in the presence of transverse perturbation.

What we did:
\bi
\item Linearization
\item Asymptotics
\item Numerics
\item Compared
\ei

\bibliographystyle{plain}
\bibliography{ntppaper}

\end{document}

