\documentclass[12pt]{article}
\usepackage{amssymb}
\font\bi=cmbxti10
\topmargin -.7in
\textheight 9in
\textwidth 7.25in
\oddsidemargin -.5in
\evensidemargin -.75in

\makeatletter\@addtoreset{equation}{section}
\def\theequation{\thesection.\arabic{equation}}

\newtheorem{theorem}{Theorem}[section]
\newtheorem{lemma}{Lemma}[section]
\newtheorem{proposition}{Proposition}[section]
\newtheorem{conjecture}{Conjecture}[section]
\newtheorem{example}{Example}[section]
\newtheorem{definition}{Definition}[section]
\newtheorem{corollary}{Corollary}[section]
\newtheorem{remark}{Remark}[section]

\def\proof{\par{\it Proof}. \ignorespaces}
\def\endproof{{\ \vbox{\hrule\hbox{%
   \vrule height1.3ex\hskip0.8ex\vrule}\hrule
  }}\par}

\usepackage{amsmath}

\begin{document}


\begin{center}
{\bf Estimates and Computations for Melting and Solidification Problems}

\vspace{.25in} 

J. M. Greenberg$^{\dagger}$ 
\vspace{.15in} 

Carnegie Mellon University\\
Department of Mathematical Sciences\\
Pittsburgh, PA  15213
\end{center} 




\section{Introduction}

In this note we focus on rigorous estimates for melting and solidification problems.  These estimates guarantee the convergence of solutions to nontrivial equilibrium patterns.  We also present some interesting numerical simulations which demonstrate the equilibrium structures and the approach of the system to equilibrium.  The novel feature of these calculations is the linking of the small parameter in the system, $\delta$, to the grid spacing.

Our governing equations are an energy balance for the medium and an evolution equation for a nonconservative ``order parameter,'' $p$.  In dimensionless variables the energy equation takes the form 

\begin{equation} 
T_t + \left(1-p^2\right)^m\!\!p_t - a^2 \Delta T = 0, \ \ \ (x,y) \in \Omega \label{E1.1}
\end{equation} 

\noindent while the evolution equation for $p$ may be written as 

\begin{equation} 
\alpha \delta p_t - \lambda^2 \delta^2 \Delta p = \mu T \left(1 - p^2 \right)^m + p 
\left(1-p^2\right),\ (x,y) \in \Omega .    \label{E1.2} 
\end{equation}

\noindent Here $\Delta$ is the two-dimensional Laplacian and $T$ is the dimensionless temperature.  The parameters $\alpha, a, \lambda$, and $\mu$ are positive and order $1, 0 < \delta << 1$, and $m$ is a nonnegative integer. The equations with $m=0$ are the most frequently studied, the case where $m=2$ has been investigated by Sekerka and coworkers [1-3] and Almgren and Almgren [4].  For definiteness we take $m=1$ and restrict our attention to the situation where $\Omega$ is the unit square 

\begin{equation} 
\Omega = \left\{ (x,y)\ |\ -1/2 < x < 1/2\ {\rm and}\ -1/2 < y < 1/2 \right\}. \label{E1.3}
\end{equation} 

\noindent On $\partial \Omega$ we assume that the normal derivatives of $T$ and $p$ vanish. Initial conditions for $T$ and $p$ are prescribed.  Noting that at a local maximum (minimum) of $p$ that $\Delta p \leq 0 (\Delta p \geq 0)$ it is easily checked that if the initial data for $p$ satisfies $-1 \leq p (x,y,0^+) \leq 1,\ (x,y) \in \overline{\Omega}$, then the same estimates obtain for all future $t > 0$.  We shall exploit this fact in what follows.

The system (\ref{E1.1}) and (\ref{E1.2}) is customarily referred to as the diffuse-interface description of melting and solidification.  There is a corresponding sharp-interface description which is obtained formally from (\ref{E1.2}) by letting $\delta$ tend to $0^+$.  The derivation of this description may be found in Caginalp [5-6]. 

In the sharp-interface description $p$ only takes on the values $-1$ and $+1$.  Solid regions are those where $p=-1$ and liquid regions those where $p=+1$. In what follows

\vspace{.10in} 

\footnoterule 

$^{\dagger}${\footnotesize This research was partially supported by the Applied Mathematical Sciences Program, U.S. Department of Energy and the Mathematics and Computer Science Division, Army Research Office.}



\begin{equation} 
\Gamma (t) = \left\{(x,y) \ |\ x = \hat{x} (s,t),\ y= \hat{y} (s,t),\ 0 \leq s \leq L \right\} 
\label{E1.4}
\end{equation} 

\noindent will be a curve separating liquid and solid regions and \mbox{\boldmath $n$}$(s,t)$ will be the unit normal to $\Gamma$ which points into the liquid region.  The unit tangent, 
\mbox{\boldmath $t$}, to $\Gamma$ is chosen so that \mbox{\boldmath $n$}$(s,t) \times \mbox{\boldmath $t$} (s,t) =$ \mbox{\boldmath $e$}$_3$.  $K$ is the curvature of $\Gamma$ at $(\hat{x},\hat{y})(s,t)$, and 

\begin{equation} 
c(s,t) \stackrel{def}{=} \left(\hat{x}_t,\ \hat{y}_t\right) \cdot \mbox{\boldmath $n$} (s,t) 
\label{E1.5} 
\end{equation}

\noindent is the speed at which $\Gamma$ moves in the direction \mbox{\boldmath $n$}.  A singular perturbation analysis of (\ref{E1.1}) and (\ref{E1.2}), with $m =1$, yields the following results: 

\begin{enumerate} 

\item $T$ is continuous across $\Gamma$; 

\item On $\Gamma$ the Gibbs-Thomson relation

\begin{equation}
\alpha c + \lambda^2 \delta K = -2^{1/2} \lambda {\mu}T_M \label{E1.6}
\end{equation} 

\noindent holds; and finally across $\Gamma$ 

\begin{equation}
\frac{4}{3} c = a^2 \left(\frac{\partial T^-}{\!\!\!\!\!\!\partial n} - 
\frac{\partial T^+}{\!\!\!\!\!\! \partial n}\right).  \label{E1.7} 
\end{equation} 

\end{enumerate} 

\noindent 

In (\ref{E1.6}), $T_M(s) = T\left(\hat{x}(s,t),\ \hat{y}(s,t)\right)$ is the local melt temperature while in (\ref{E1.7}), ${\displaystyle \frac{\partial T^+}{\!\!\!\!\!\! \partial n}}$ is the normal derivative of $T$ in the liquid region and ${\displaystyle \frac{\partial T^-}{\!\!\!\!\!\! \partial n}}$ is the normal derivative of $T$ in the solid region.  In the sharp-interface description one solves the heat equation 

\begin{equation} 
T_t - a^2 \Delta T = 0 \label{E1.8} 
\end{equation} 

\noindent in the portion of $\Omega$ away from interfaces.  At equilibrium, $c=0$ and (\ref{E1.6}) implies that $T_M = 0(\delta)$.  In what follows, we shall restrict our attention to initial data $T(x,y,0^+) = 0 (\delta^{\frac{1}{2}})$. 

The sharp-interface description becomes unduly complicated in the presence of multiple interfaces since all must be tracked.  Difficulties also arise when interfaces merge.  For these, and a plethora of other reasons, we focus on (\ref{E1.1}) and (\ref{E1.2}).  Our principal results consist of a sequence of a-priori inequalities satisfied by solutions of the system (\ref{E1.1}) and (\ref{E1.2}).  These estimates guarantee regularity of the solutions and that the solutions converge to equilibria as $t$ tends to infinity.  

One point worthy of note is it is possible to obtain solutions of (\ref{E1.1}) and (\ref{E1.2}) which approximate those generated by the sharp interface descriptions.  These are produced by solving (\ref{E1.1}) and (\ref{E1.2}) with initial data $p^{\delta}(x,y,0^+)$ which approximate the jump discontinuous data of the sharp interface theory but transition smoothly from $-1$ to $+1$ over an interval of width of $0(\delta)$ in the direction \mbox{\boldmath $n$}.  For such initial $p$'s (and initial temperature fields satisfying $T(x,y,0^+) = 0(\delta^{\frac{1}{2}})$) we are guaranteed that terms like $\delta {\displaystyle \int \!\int_{\Omega}} | \bigtriangledown p^{\delta}|^2 dxdy$ and ${\displaystyle \frac{1}{\delta}} 
{\displaystyle \int \!\! \int}_{\Omega} \left(1-p^2\right)^2 dxdy$ are bounded independently of $\delta$ and $t$ and it is these estimates which link the solutions of the two theories.

We also note there is an ``intermediate" model where one replaces (\ref{E1.1}) (here with $m=1$) with 

\begin{equation} 
T(x,y,t) \equiv T_0 (t) = h_0 -
{\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}}\left(p-p^3/3\right) (x,y,t) dxdy \label{E1.9} 
\end{equation} 

\noindent and solves 

\begin{equation} 
\alpha \delta p_t - \lambda^2 \delta ^2 \Delta p = \left(\mu T_0 + p \right) \left(1-p^2\right),\ (x,y) \in \Omega \label{E1.10} 
\end{equation} 
 
\noindent and 

\begin{equation} 
\frac{\partial p}{\partial n} = 0 \ \ \ , \ \ \ (x,y) \in \partial \Omega . \label{E1.11} 
\end{equation} 

\noindent The constant $h_0$ in (\ref{E1.9}) is given by 

\begin{equation} 
h_0 = {\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}}\left(T + p -p^3/3\right) (x,y,0^+) dxdy. \label{E1.12} 
\end{equation} 

\noindent This ``intermediate" model is justified by two facts.  The first is that solutions of (\ref{E1.1}) and (\ref{E1.2}) satisfying the boundary conditions ${\displaystyle \frac{\partial T}{\partial n}} = {\displaystyle \frac{\partial p}{\partial n}} = 0,\ (x,y) \in \partial \Omega$ satisfy the overall energy balance: 

\begin{equation} 
h_0 \equiv{\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}}\left(T+ p-p^3/3\right) (x,y,t) dxdy \label{E1.13} 
\end{equation} 

\noindent for all $t \geq 0$ and the second is that the spatially varying portion of $T$ decays to zero as $t$ tends to infinity (for details see the comments following (\ref{E2.17})). 


\section{Entropies and A-Priori Estimates} 

Our goal in this section is a set of inequalities of the form 


\begin{equation}
\frac{\partial \eta}{\partial t} - {\rm div}\ \mbox{\boldmath $q$} \leq G \label{E2.1} 
\end{equation}

\noindent satisfied by solutions of 

\begin{equation}
T_t + \left(1-p^2\right)p_t - a^2 \Delta T = 0\ \ \ , \ \ \ (x,y) \in \Omega, \label{E2.2}
\end{equation}

\begin{equation}
\alpha \delta p_t - \lambda^2 \delta^2 \Delta p = \left(\mu T + p \right) \left(1-p^2\right),\ (x,y) \in \Omega, \label{E2.3}    
\end{equation}
 
\noindent and 

\begin{equation} 
\frac{\partial T}{\partial n} = \frac{\partial p}{\partial n} = 0\ , \ (x,y) \in \partial \Omega . \label{E2.4} 
\end{equation} 

\noindent Once again 

\begin{equation}
\Omega = \left\{(x,y) \ |\ -1/2 < x < 1/2\ {\rm and}\ -1/2 < y < 1/2 ) \right\}. \label{E2.5}
\end{equation} 

\noindent For (\ref{E2.1}) to be of value we require 

\begin{equation} 
0 \leq \eta,\ (x,y) \in \Omega\ {\rm and}\ \mbox{\boldmath $q$} \cdot \mbox{\boldmath $n$} = 0,\ (x,y) \in \partial \Omega \label{E2.6}
\end{equation}

\noindent and that either $G \leq 0$ in $\Omega$ or $G$ is a-priori bounded in $L_1(\Omega \times  (0,\infty))$.  Throughout we exploit the pointwise estimates 

\begin{equation}
-1 \leq p (x,y,t) \leq 1\ , \ (x,y) \in \Omega\ {\rm and}\ t \geq 0 \label{E2.7} 
\end{equation} 

\noindent satisfied by solutions of (2.3).

Once again we record the overall energy balance satisfied by solutions of (\ref{E2.2}) -- (\ref{E2.4}), namely the identity:

$$
\int^{1/2}_{-1/2}  \int^{1/2}_{-1/2} \left(T+p-p^3/3\right) (x,y,t) dxdy \equiv h_0. \eqno{(1.13)} 
$$

\noindent This is obtained by integrating (\ref{E2.2}) over $\Omega$ and exploiting (\ref{E2.4})$_1$. 

Our first inequality is obtained by multiplying (\ref{E2.2}) by $\left(\frac{\mu T}{\delta}\right)$ and (\ref{E2.3}) by $\frac{p_t}{\delta}$ and adding the resulting identities.  This yields (\ref{E2.1}) with 

\begin{equation} 
\eta_1 = \frac{\mu T^2}{2 \delta} + \frac{\lambda^2 \delta}{2} | \bigtriangledown p|^2 + \frac{\left(1-p^2\right)^2}{4\delta}, \label{E2.8} 
\end{equation} 

\begin{equation} 
\mbox{\boldmath $q$}_1 = \frac{\mu a^2}{\delta} T \bigtriangledown T + \lambda^2 \delta p_t \bigtriangledown p, \label{E2.9} 
\end{equation} 

\noindent and

\begin{equation} 
G_1 = - \alpha p^2_t - \frac{\mu a^2 |\bigtriangledown T|^2}{\delta} \label{E2.10} 
\end{equation} 

\noindent and (\ref{E2.8})-(\ref{E2.10}) yield

\begin{equation} 
\begin{array}{llll}
{\displaystyle \frac{d}{dt}} 
{\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} \left({\displaystyle 
\frac{\mu T^2}{2 \delta}} + {\displaystyle \frac{\lambda^2 \delta}{2}} |\bigtriangledown p|^2 + {\displaystyle \frac{(1-p^2)}{4\delta}}^2 \right) (x,y,t) dxdy \\
\\
\leq -{\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}}\left(\alpha p^2_t + {\displaystyle \frac{\mu a^2}{\delta}} | \bigtriangledown T|^2 \right) (x,y,t) dxdy. \end{array} \label{E2.11}
\end{equation}

\noindent Once again we restrict our attention to initial data satisfying 

$${\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} \left({\displaystyle \frac{T^2}{2\delta}} + {\displaystyle \frac{\lambda^2 \delta}{2}} | \bigtriangledown p|^2 + {\displaystyle \frac{(1-p^2)}{4 \delta}}^2\right)(x,y,0^+) dxdy = 0 (1)$$

\noindent independent of $\delta$. 

If we next multiply (\ref{E2.2}) by $\alpha T_t$ we obtain

$$\frac{\alpha a^2}{2}\frac{\partial}{\partial t} | \bigtriangledown T|^2 - \alpha a^2 \ {\rm div}\ (T_t \bigtriangledown T) = - \alpha \left(T^2_t + \left(1-p^2\right) p_t T_t\right)$$

\noindent and if we add this result to (\ref{E2.1}) with $\eta_1,\ \mbox{\boldmath $q$}_1,\ G_1$ we obtain 

$$
\frac{\partial \eta_2}{\partial t} - {\rm div}\ \mbox{\boldmath $q$}_2 \leq G_2  
$$

\noindent where 

\begin{equation} 
\eta_2 = \frac{\mu T^2}{2 \delta} + \frac{\lambda^2 \delta}{2} | \bigtriangledown p|^2 + \frac{(1-p^2)}{4 \delta}^2 + \frac{\alpha a^2}{2} | \bigtriangledown T|^2 
\label{E2.12} 
\end{equation} 

\begin{equation} 
\mbox{\boldmath $q$}_2 = \frac{\mu a^2}{\delta} T \bigtriangledown T + \lambda^2 \delta p_t \bigtriangledown p + \alpha a^2 T_t \bigtriangledown T \label{E2.13} 
\end{equation} 

\noindent and 

\begin{equation} 
G_2 = - \alpha \left(p^2_t + T^2_t + \left(1-p^2\right) T_t p_t\right) - \frac{\mu a^2}{\delta} | \bigtriangledown T|^2 \leq \frac{-\alpha}{2} \left(p^2_t + T^2_t\right) - 
\frac{\mu a^2}{\delta} | \bigtriangledown T|^2 \label{E2.14} 
\end{equation} 

\noindent and (\ref{E2.12})-(\ref{E2.14}) imply that 

\begin{equation}
\begin{array}{lll}
{\displaystyle \frac{d}{dt}} {\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}}\left(
{\displaystyle \frac{\mu T^2}{2 \delta}} + {\displaystyle \frac{\lambda^2 \delta}{2}} | \bigtriangledown p|^2 + {\displaystyle \frac{(1-p^2)}{4 \delta}}^2 + 
{\displaystyle \frac{\alpha a^2}{2}} | \bigtriangledown T|^2 \right) (x,y,t) dxdy \\
\\
\leq - {\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}}\left({\displaystyle \frac{\alpha}{2}} \left(p^2_t + T^2_t \right) + {\displaystyle \frac{\mu a^2}{\delta}} | \bigtriangledown T|^2 \right) (x,y,t) dxdy. \end{array} \label{E2.15} 
\end{equation} 

If we now write

\begin{equation}
T(x,y,t) = T_0 (t) + T_1 (x,y,t)  \ \ , \ \ (x,y) \in \Omega \label{E2.16} 
\end{equation} 

\noindent where 

\begin{equation}
T_0 (t) = \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2} T(x,y,t)dxdy \ {\rm and}\ 
\int^{1/2}_{-1/2}  \int^{1/2}_{-1/2} T_1 (x,y,t) dxdy \equiv 0, \label{E2.17} 
\end{equation} 


\noindent then (\ref{E2.15}) implies that if $T$ and $p$ are in $H^1(\Omega)$ at $t=0$, then they are in $H^1(\Omega)$ for all $t > 0$ and that $T$ satisfies ${\displaystyle \lim_{t \rightarrow \infty}} {\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} T^2_1 (x,y,t) dxdy = 0$. \footnote{Recall, $H^k(\Omega)$ consists of functions $f$ satisfying 
$$||f||^2_k\ \stackrel{def}{=} \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2} \sum^k_{j=0} 
\left( \begin{array}{l@{\protect \vspace{-.04in}}l} \sum\\
\hspace*{-.10in} \begin{array}{l@{\protect \vspace{-.04in}}ll}
{\scriptstyle \alpha_1 + \alpha_2 =j}\\
{\scriptstyle \alpha_i \geq 0} \end{array}  
 \end{array}  \left(
\frac{\partial^jf}{\partial x^{\alpha_1}\partial y^{\alpha_2}}\right)^2\right) (x,y) dxdy < \infty $$ 

\ \ \ and $H^k (\Omega \times (0,s_*))$ consists of all functions $f$ on $\Omega \times (0,s_*)$ satisfying 
$$|||f|||^2_k\ {\stackrel{def}{=}} \int^{s_*}_0 
\left( \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}  
\sum^k_{j=0} \left( 
\begin{array}{l@{\protect \vspace{-.04in}}l} \sum\\
\hspace*{-.10in} \begin{array}{l@{\protect \vspace{-.04in}}ll}
{\scriptstyle \alpha_1 + \alpha_2 + \alpha_3 =j}\\
{\scriptstyle \alpha_i \geq 0} \end{array} \end{array} 
\left(\frac{\partial^j f}{\partial x^{\alpha_1} \partial y^{\alpha_2}
\partial s^{\alpha_3}}\right)^2\right) (x,y,s) dxdy \right) ds < \infty.$$}

If we make a similar decomposition of $p$, that is write 

\begin{equation} 
p(x,y,t) = p_0 (t) + p_1 (x,y,t) \label{E2.18} 
\end{equation} 

\noindent where 

\begin{equation}
p_0 (t) = \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2} p (x,y,t) dxdy \ {\rm and}\ 
\int^{1/2}_{-1/2}  \int^{1/2}_{-1/2} p_1(x,y,t) dxdy = 0, \label{E2.19} 
\end{equation}

\noindent then (\ref{E2.3}) implies that 

\begin{equation} 
\alpha \delta \frac{dp_0}{dt} = \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2} ( \mu T + p ) (1-p^2) (x,y,t) dxdy \label{E2.20} 
\end{equation}

\noindent and 

\begin{equation}
\alpha \delta \frac{d^2p_0}{dt^2} = \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}
\left(\mu \left(T_t \left(1-p^2\right) - 2 pTp_t \right) + \left(1-3p^2\right) p_t \right) (x,y,t) dxdy \label{E2.21} 
\end{equation}

\noindent and (\ref{E2.20}) and (\ref{E2.21}) along with the fact that (\ref{E2.15}) implies that $p_t$ and $T_t$ are $L_2(\Omega \times(0,\infty))$ combine to yield the result that ${\displaystyle \frac{dp_0}{dt}}$ and ${\displaystyle \frac{d^2p_0}{dt^2}}$ are in 
$L_2(0, \infty)$ and that 

\begin{equation} 
{\displaystyle \lim_{t \rightarrow \infty}} \alpha \delta 
{\displaystyle \frac{dp_0}{dt}} = {\displaystyle \lim_{t \rightarrow \infty}} 
\int^{1/2}_{-1/2}  \int^{1/2}_{-1/2} (\mu T + p) (1-p^2)(x,y,t) dxdy = 0. \label{E2.22} 
\end{equation} 

With the estimate (\ref{E2.15}) we can glean weak information about the long-term behavior of solutions of (\ref{E2.2})-(\ref{E2.4}).  In what follows we let 

\begin{equation} 
 (T^t,p^t)(x,y,s) = (T,p) (x,y,t+s), \ (x,y) \in \Omega \ {\rm and}\ 0 \leq s \leq s_*.   \label{E2.23} 
\end{equation} 

\noindent The inequality (\ref{E2.15}), together with (\ref{E2.7}), implies that for any $s_*$  the functions $T^t$ and $p^t$ are in \newline
$H^1(\Omega \times (0,s_*))$ with bounds which depend only on the $H^1(\Omega)$ norms of $T(x,y,0^+)$ and $p(x,y,0^+)$ and $s_*$ and that 

\begin{equation} 
{\displaystyle \lim_{t \rightarrow \infty}} \int^{s_*}_0 \left(
\int^{1/2}_{-1/2}  \int^{1/2}_{-1/2} \left(\frac{\alpha}{2}\left(\left(T^t_s\right)^2 + \left(p^t_s\right)^2 \right) + \frac{\mu a^2}{\delta} | \bigtriangledown T^t|^2\right) (x,y,s) dxdy\right) ds = 0. \label{E2.24} 
\end{equation} 

\noindent Moreover, (\ref{E2.3}) implies that $\lambda^2 \delta^2 \Delta p^t \in L_2 (\Omega \times (0,s_*))$ with bounds which are independent of $t$. The boundedness of $T^t $ and $p^t$ in $H^1(\Omega \times (0,s_*))$ implies we can find an increasing sequence $\{t_i\}^{\infty}_{i = 1}$ with ${\displaystyle \lim_{i \rightarrow \infty}} t_i = \infty$ and funcitons $T^{\infty}$ and $p^{\infty}$ in $H^1(\Omega \times (0,s_*))$ and $\Delta p^{\infty}$ in $L_2(\Omega \times (0,s_*))$ with the following properties:

\noindent $(T^{\infty},p^{\infty})$ is the strong $L_2(\Omega \times (0,s_*))$ limit of $(T^{t_i},p^{t_i})$ while $T^{\infty}_s,\ p^{\infty}_s, \bigtriangledown T^{\infty},\bigtriangledown p^{\infty}$, and $\Delta p^{\infty}$ are the weak $L_2(\Omega \times (0,s_*))$ limits of the appropriate derivatives of $T^{t_i}$ and $p^{t_i}$.  Additionally, (\ref{E2.24}) implies

\begin{equation} 
T^{\infty}_s = p^{\infty}_s = 0 \ {\rm and}\ \bigtriangledown T^{\infty} = 0\ {\rm in}\ \Omega \times (0,s_*) \label{E2.25} 
\end{equation} 

\noindent and thus that $T^{\infty}$ is a constant and $p^{\infty}$ is independent of $s$.  The energy balance (\ref{E1.13}) yields 

\begin{equation} 
 T^{\infty} = h_0 - \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2} \left(p^{\infty} - \left(p^{\infty}\right) ^3/3\right) (x,y) dxdy  \label{E2.26}
\end{equation} 

\noindent while (\ref{E2.7}) implies that $-1 \leq p^{\infty} \leq 1$.  (\ref{E2.3}) and (\ref{E2.24}) also imply that $p^{\infty}$ satisifes the equilibrium equation 

\begin{equation} 
\lambda^2 \delta^2 \Delta p^{\infty} + \left( \mu T^{\infty} + p^{\infty}\right) \left(1- \left(p^{\infty}\right)^2\right) = 0,\ (x,y) \in \Omega \label{E2.27} 
\end{equation}

\noindent and 

\begin{equation} 
\frac{\partial p^{\infty}}{\partial n} = 0,\ (x,y) \in \partial \Omega \label{E2.28} 
\end{equation}

\noindent where $T^{\infty}$ is given by (\ref{E2.26}) and finally (\ref{E2.27}) and (\ref{E2.28}) imply that $p^{\infty}$ is a critical point of the functional 

\begin{equation} 
N\left(p^{\infty}\right) \stackrel{def}{=} \frac{\mu \left(T^{\infty}\right)^2}{2 \delta} + \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2} \left(\frac{\lambda^2 \delta}{2} |\bigtriangledown p^{\infty}|^2 + \frac{\left(1-\left(p^{\infty}\right)^2 \right)}{4\delta}^2\right) (x,y) dxdy. \label{E2.29} 
\end{equation} 

\noindent Additionally, the functional 

\begin{equation} 
\eta (t,s_*) \stackrel{def}{=} \frac{1}{s_*} \int^{s_*}_0 \! \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2} \left(\frac{\mu (T^t)^2}{2 \delta} + \frac{\lambda^2\delta}{2} | \bigtriangledown p^t|^2 + \frac{\left(1-\left(p^t\right)^2\right)^2}{4 \delta} + \frac{\alpha a^2}{2} | \bigtriangledown T^t|^2 \right) (x,y,s) dxdy ds \label{E2.30} 
\end{equation} 

\noindent is decreasing in $t$ and has a limit $\eta_{\infty}$ satisfying 

\begin{equation} 
N \left(p^{\infty}\right) \leq \eta_{\infty}. \label{E2.31} 
\end{equation} 

To obtain stronger information we examine the time differentiated versions of (\ref{E2.2})-(\ref{E2.4}), namely the equations:

\begin{equation} 
T_{tt}+ (1-p^2)p_{tt} - 2pp^2_t - a^2 \Delta T_t - 0,\ (x,y) \in \Omega, \label{E2.32} 
\end{equation}

\begin{equation} 
\alpha \delta p_{tt} - \lambda^2\delta^2 \Delta p_t = (\mu T_t + p_t)(1-p^2) - 2p(\mu T + p)p_t,\ (x,y) \in \Omega, \label{E2.33} 
\end{equation}


\begin{equation} 
\frac{\partial T_t}{\partial n} = \frac{\partial p_t}{\partial n} = 0\ \ , \ \ (x,y) \in \partial \Omega . \label{E2.34} 
\end{equation} 

\noindent We shall also exploit

\noindent {\bf Lemma 1}.  Suppose $f \in H^1 (\Omega)$.  Then, for every $\epsilon > 0$

$$ \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2} f^4 (x,y) dxdy \leq 4 \left[\frac{(1+\epsilon)}{\epsilon} 
\int^{1/2}_{-1/2}  \int^{1/2}_{-1/2} f^2(x,y) dxdy + \epsilon \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2} | \bigtriangledown f|^2 (x,y) dxdy \right]^2 . $$ 

We multiply (\ref{E2.33}) by $p_t$ and make use of (\ref{E2.7}) to obtain 

\begin{equation} 
\begin{array}{lll}
{\displaystyle \frac{\alpha \delta}{2}} {\displaystyle \frac{\partial}{\partial t}} p^2_t - \lambda^2 \delta^2\ \! {\rm div}\ \! \left(p_t \bigtriangledown p_t\right) \leq - \lambda^2 \delta^2 | \bigtriangledown p_t|^2 + \left({\displaystyle \frac{\mu T^2_t}{2}} +  \left(1+ {\displaystyle \frac{\mu}{2}} \right) p^2_t \right) +  2 \mu |T|p^2_t. \end{array} \label{E2.35}
\end{equation} 

\noindent We note that (\ref{E2.15}) implies the term ${\displaystyle \frac{\mu T_t}{2}}^2 + \left(1+{\displaystyle \frac{\mu}{2}}\right)p^2_t$ is in $L_1 (\Omega \times (0, \infty))$ and thus we confine our attention to $2 \mu |T|p^2_t$.  Schwarz's inequality and Lemma 1 imply 

\begin{equation} 
\begin{array}{lll}
2 \mu {\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} \left(|T|p^2_t\right) (x,y,t) dxdy \leq 2 \mu \left( {\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} T^2 (x,y,t) dxdy \right)^{1/2}\ \left({\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} p^4_t (x,y,t) dxdy \right)^{1/2} \\
\\
\leq \! 2 \mu \left( {\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} T^2 (x,y,t) dxdy \right)^{1/2}\!\! \left(\left({\displaystyle \frac{1\!+ \!\epsilon}{\epsilon}}\right)
{\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} p^2_t (x,y,t) dxdy  \! + \! \epsilon {\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} | \! \bigtriangledown \! p_t|^2 (x,y,t) dxdy \right). \\
\\
\\
\\
\\
\   \end{array} \hspace*{-.25in} \label{E2.36} 
\end{equation} 

\ 
\vspace*{-1in}

\noindent Moreover, (\ref{E2.11}) guarantees that 

\begin{equation} 
\left(\int^{1/2}_{-1/2}  \int^{1/2}_{-1/2} T^2 (x,y,t) dxdy\right)^{1/2} \leq \left(\frac{2 \delta}{\mu} \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2} \eta_1(x,y,0^+) dxdy\right)^{1/2} \label{E2.37} 
\end{equation} 

\noindent and this guarantees that we may choose $\epsilon$ so that

\begin{equation} 
2 \mu \epsilon \left(\int^{1/2}_{-1/2}  \int^{1/2}_{-1/2} T^2 (x,y,t) dxdy \right)^{1/2} \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2} | \bigtriangledown p_t |^2 (x,y,t) \leq \frac{\lambda^2 \delta^2}{2} \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2} | \bigtriangledown p_t|^2 (x,y,t) dxdy. \label{E2.38} 
\end{equation}

\noindent With $\epsilon$ so chosen, (\ref{E2.35}) implies that 

\begin{equation} 
\begin{array}{lll}
{\displaystyle \frac{\alpha \delta}{2} \ \frac{d}{dt} \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} p^2_t (x,y,t) \leq 
{\displaystyle \frac{- \lambda^2 \delta^2}{2} \int^{1/2}_{-1/2}}  \int^{1/2}_{-1/2}| \bigtriangledown p_t|^2 (x,y,t) dxdy + {\displaystyle \frac{\mu}{2}} 
{\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} T^2_t (x,y,t) dxdy \\
\\
+ \left( 1 + {\displaystyle \frac{\mu}{2}} + 2^{3/2} \left( \mu \delta 
{\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} \eta_1 (x,y,0^+) dxdy \right)^{1/2} \left(1 + {\displaystyle \frac{1}{\epsilon}}\right) \right) 
{\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} p^2_t (x,y,t) dxdy  \end{array} \label{E2.39} 
\end{equation} 

\noindent and (\ref{E2.39}) then implies that if 

$$p_t (x,y,0^+) = \frac{1}{\alpha \delta} \left(\lambda^2 \delta^2 \Delta p+ (\mu T + p) (1-p^2)\right) (x,y,0^+) $$

\noindent is in $L_2(\Omega)$, then $p_t$ and thus $\Delta p$ are in $L_2(\Omega)$ independently of $t$ and $\bigtriangledown p_t \in L_2 (\Omega \times (0, \infty))$.  

The preceding estimate points up the advantage of the ``intermediate" model discussed briefly in (\ref{E1.9})-(\ref{E1.12}) over the full system.  For the intermediate model the contankerous term $2 \mu|T|p^2_t$ in (\ref{E2.35}) would have been replaced by $2 \mu|T_0|p^2_t$ and this latter term is easily controlled since we know by our previous estimates that 

$$|T_0| \leq \left(\frac{2 \delta}{\mu} \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2} \eta_1 (x,y,0^+) dxdy \right)^{1/2} $$

\noindent and thus we would not have to invoke Lemma 1 and choose $\epsilon = 0 \left(\delta^{3/2}\right)$.  This latter choice, pumped up the coefficient of 
${\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}}p^2_t (x,y,t) dxdy$ in (\ref{E2.39}) to be $0 \left(1/\delta^{3/2}\right)$ which is far from optimal.  If we chose to add some multiple of (\ref{E2.39}) to (\ref{E2.15}) and wanted the right hand side of the resulting inequality to be negative definite our multiplier would have to be $0(\delta^{3/2})$, an unacceptably small number.  For the ``intermediate" model we could use an $0(1)$ multiplier (for details see (\ref{E2.61})-(\ref{E2.63})). 


On the other hand, if we multiply (\ref{E2.32}) by $(T_t + (1-p^2)p_t)/a^2$ we find that 

\begin{equation}
\begin{array}{llll}
{\displaystyle \frac{1}{2a^2}} {\displaystyle \frac{\partial}{\partial t}} \left(T_t + (1-p^2) p_t \right)^2 = {\displaystyle \frac{a^2}{2}} {\displaystyle \frac{\partial}{\partial t}} (\Delta T)^2\\
\\
 = {\rm div}\ ((T_t + (1-p^2)p_t) \bigtriangledown T_t) - |\bigtriangledown T_t|^2 - (1-p^2) \bigtriangledown T_t \cdot \bigtriangledown p_t + 2 pp_t \bigtriangledown T_t  \cdot \bigtriangledown p \end{array}  \label{E2.40} 
\end{equation} 


\noindent and (\ref{E2.40}) implies that 

\begin{equation} 
\begin{array}{lllll}
{\displaystyle\frac{a^2}{2}} {\displaystyle\frac{d}{dt}} {\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} (\Delta T)^2 (x,y,t) dxdy \leq \\
\\
- {\displaystyle \frac{1}{2}} {\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} | \bigtriangledown T_t |^2 (x,y,t) dxdy + 
{\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} | \bigtriangledown p_t |^2 (x,y,t) dxdy \\
\\
+ 4 {\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} p^2_t | \bigtriangledown p |^2  (x,y,t) dxdy. \end{array}    \label{E2.41} 
\end{equation} 

\noindent The first term on the right hand side of (\ref{E2.41}) is in controllable by virtue of the minus sign and the second term is $L_1 (0, \infty)$ by virtue of (\ref{E2.39}).  Thus we can confine our attention to the third term.  To show this term is $L_1(0,\infty)$ we note that 

\begin{equation} 
{\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} p^2_t |\!\! \bigtriangledown \!\! p|^2 (x,y,t) dxdy \!\! \leq \!\! 
\left({\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} p^4_t (x,y,t) dxdy\right)^{1/2} 
\left({\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} |\! \bigtriangledown \! p|^4 (x,y,t) dxdy\right)^{1/2} \!\!. \label{E2.42} 
\end{equation}

\noindent Lemma 1 then guarantees that 

\begin{equation}
\left( {\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} p^4_t (x,y,t) dxdy \right)^{1/2} 
\leq {\displaystyle \frac{(1+ \sqrt{5})}{2}}
\left( {\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} \left(p^2_t + | \bigtriangledown p_t|^2 \right) (x,y,t) dxdy \right)     \label{E2.43} 
\end{equation} 

\noindent and that

\begin{equation}
\begin{array}{lll}
\left({\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} | \bigtriangledown p|^4 (x,y,t) dxdy\right)^{1/2} 
\leq \left({\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} p^4_x (x,y,t) dxdy \right)^{1/2} + 
\left( {\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} p^4_y (x,y,t) dxdy \right)^{1/2} \\
\\
\leq {\displaystyle \frac{\left(1+ \sqrt{5}\right)}{2}} {\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} 
\left(p^2_x + p^2_y + p^2_{xx} + p^2_{yy} + 2 p^2_{xy} \right) (x,y,t) dxdy \\
\\
\leq {\displaystyle \frac{\left(1+ \sqrt{5}\right)}{2}} {\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} 
\left(p^2_x + p^2_y + 4 (\Delta p)^2 \right) (x,y,t) dxdy 
\end{array} \label{E2.44} 
\end{equation} 

\noindent and (\ref{E2.43}) and (\ref{E2.44}) together with (\ref{E2.41}) guarantee that if $T(x,y,0^+)$ and $p(x,y,0^+)$ are in $H^2(\Omega)$ then $T(x,y,t)$ is in $H^2(\Omega)$ for all future times and additionally that $\bigtriangledown T_t \in L_2 (\Omega \times (0,\infty ))$. 

The estimates (\ref{E2.39}) and (\ref{E2.44}) also guarantee that if $T(x,y,0^+)$ and $p(x,y,0^+)$ are in $H^2(\Omega)$, then the resulting orbit is compact in $H^1(\Omega)$ and that the $w$-limit set is nonempty in $H^1(\Omega)$.  The fact that $\bigtriangledown T_t \in L_2 (\Omega \times (0, \infty))$ further guarantees that the function $T_1$ defined in (\ref{E2.16}) and (\ref{E2.17}) converges to zero strongly in $H^1(\Omega)$. 

Our final estimates pertain to solutions generated by initial data $T(x,y,0^+) \in H^2(\Omega)$ and\newline
$p(x,y,0^+) \in H^3(\Omega)$.  These are obtained by multiplying (\ref{E2.33}) by $p_{tt}$.  The resulting inequality is 

\begin{equation} 
\begin{array}{lll}
{\displaystyle \frac{\lambda^2 \delta^2}{2}} {\displaystyle \frac{\partial}{\partial t}} | \bigtriangledown p_t|^2 + \alpha \delta p^2_{tt} - \lambda^2 \delta^2 \ {\rm div}\ \left(p_{tt} \bigtriangledown p_t \right) = \\
\\
\mu \left(T_t (1-p^2) - 2pTp_t \right) p_{tt} + (1-3p^2) p_t p_{tt} \\
\\
\leq {\displaystyle \frac{\alpha \delta}{2}} p^2_{tt} + {\displaystyle \frac{2 \mu^2}{\alpha \delta}} \left(T^2_t + 4T^2p^2_t\right) + {\displaystyle \frac{8}{\alpha \delta}} p^2_t . \end{array} \label{E2.45}
\end{equation}

\noindent Our preceding calculations imply that if $T(x,y,0^+) \in H^2(\Omega)$ and $p(x,y,0^+) \in H^3 (\Omega)$, then 

$$
\frac{2 \mu^2}{\alpha \delta} \left(T^2_t + 4T^2p^2_t\right) + {\displaystyle \frac{8}{\alpha \delta}} p^2_t
$$

\noindent is in $L_1 (\Omega \times (0, \infty))$ and thus (\ref{E2.45}) implies that $\bigtriangledown p_t \in L_2 (\Omega)$ independently of $t$ and that $p_{tt} \in L_2 (\Omega \times (0, \infty))$.  The former conclusion, together with (\ref{E2.3}), implies that $p \in H^3 (\Omega)$ with bounds which are independent of $t$ while the latter conclusion implies that 

$${\displaystyle \lim_{t \rightarrow \infty}} \alpha \delta p_t = {\displaystyle \lim_{t \rightarrow \infty}} \left(\lambda^2 \delta^2 \Delta p + (\mu T+p)(1-p^2) \right) = 0$$

\noindent and that this limit is obtained strongly in $L_2 (\Omega)$.

We conclude this section with some observations about the ``intermediate" model introduced in (\ref{E1.9})-(\ref{E1.12}).  Once again in that model one replaces (\ref{E2.3}) with 

\begin{equation} 
\alpha \delta p_t - \lambda^2 \delta^2 \Delta p = (\mu T_0 + p) (1-p^2), (x,y) \in \Omega \label{E2.46} 
\end{equation}

\noindent and 

\begin{equation}
\frac{\partial p}{\partial n} = 0\ \ \ , \ \ \ (x,y) \in \partial \Omega \label{E2.47} 
\end{equation}

\noindent where

\begin{equation} 
T_0(t) = h_0 \! - \! \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2} \left( p-p^3/3\right) (x,y,t) dxdy\   {\rm and} \ 
h_0 = \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2} \left(T + p - p^3/3\right)(x,y,0^+) dxdy. \label{E2.48} 
\end{equation}

Solutions of this latter system also satisfy the pointwise estimates $-1 \leq p (x,y,t) \leq 1$ and a set of ``entropy" inequalities similar to those satisfied by solutions of the full system.  The most basic of these is

\begin{equation}
\frac{\partial \tilde{\eta}_1}{\partial t} - {\rm div}\ \tilde{\mbox{\boldmath $q$}}_1 = \tilde{G}_1 \label{E2.49}
\end{equation}

\noindent where

\begin{equation}\tilde{\eta}_1 = \frac{\mu T^2_0}{2 \delta} + 
\frac{\lambda^2 \delta}{2} |\bigtriangledown p|^2 + \frac{\left(1-p^2\right)^2}{4 \delta} \label{E2.50}
\end{equation}

\begin{equation}
\tilde{\mbox{\boldmath $q$}}_1 = \lambda^2 \delta p_t \bigtriangledown p \label{E2.51} 
\end{equation}

\noindent and 

\begin{equation}
\tilde{G}_1 = - \alpha p^2_t. \label{E2.52}
\end{equation}

\noindent Moreover, (\ref{E2.49}) - (\ref{E2.52}) imply that 

\begin{equation}
\begin{array}{lll}
{\displaystyle \frac{d}{dt}} {\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} 
\left({\displaystyle \frac{\mu T^2_0}{2\delta}} + {\displaystyle \frac{\lambda^2 \delta}{2}} |\bigtriangledown p|^2 + {\displaystyle \frac{\left(1-p^2\right)}{4\delta}}^2 \right) (x,y,t) dxdy \\
\\
\leq - \alpha {\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} p^2_t (x,y,t) dxdy. \end{array} \label{E2.53} 
\end{equation}



\noindent The time differentiated identities 

\begin{equation}
\alpha \delta p_{tt} - \lambda^2 \delta^2 \Delta p_t = \mu \left(\dot{T}_0(1-p^2) - 2pT_0 p_t\right) + (1-3p^2) p_t, \label{E2.54}
\end{equation}

\begin{equation}
\dot{T}_0 = - \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2} (1-p^2) p_t (x,y,t) dxdy, \label{E2.55}
\end{equation} 

\noindent and 

\begin{equation}
\frac{\partial p_t}{\partial n} = 0 \ \ , \ \ (x,y) \in \partial \Omega . \label{E2.56} 
\end{equation}

\noindent also yields: 

\begin{equation}
\begin{array}{lllll}
{\displaystyle \frac{\alpha \delta}{2}} {\displaystyle \frac{\partial p^2_t}{\partial t}} - \lambda^2 \delta^2 \ {\rm div} \ (p_t \bigtriangledown p_t)= \\
\\
- \lambda^2 \delta^2 | \bigtriangledown p_t|^2 + \mu 
\left(\dot{T}_0 (1-p^2) - 2pT_0 p_t \right) p_t + (1-3p^2) p^2_t \\
\\
\leq -\lambda^2 \delta^2 | \bigtriangledown p_t|^2 + {\displaystyle \frac{\mu}{2}} 
\left({\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} p^2_t (x,y,t) dxdy \right)  \\
\\
+ 2^{3/2} \mu^{1/2} \delta^{1/2} \left({\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} \tilde{\eta}_1 (x,y,0^+) dxdy \right)^{1/2} p^2_t+ (1+{\displaystyle \frac{\mu}{2}}) p^2_t \end{array} \label{E2.57} 
\end{equation} 

\noindent and 

\begin{equation}
\begin{array}{lll}
{\displaystyle \frac{\lambda^2 \delta^2}{2}} 
{\displaystyle \frac{ \partial |\bigtriangledown p_t|^2}{\partial t}} - \lambda^2 \delta^2 \ {\rm div}\ (p_{tt} \bigtriangledown p_t) = \\
\\
- \alpha \delta p^2_{tt} + \mu \left(\dot{T}_0(1-p^2)-2pT_0p_t)p_{tt} + (1-3p^2) p_t p_{tt} \right. \\
\\
\leq - {\displaystyle \frac{\alpha \delta}{2}}p^2_{tt}  + {\displaystyle \frac{1}{\alpha \delta}} 
\left( 2 \mu^2 {\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} p^2_t (x,y,t) dxdy + 
\left(4 + 16 \mu \delta {\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} \tilde{\eta}_1 (x,y,0^+)dxdy  \right) p^2_t \right). \end{array} \label{E2.58} 
\end{equation} 

\noindent The former follows from multiplying (\ref{E2.54}) by $p_t$ and the latter from multiplying (\ref{E2.54}) by $p_{tt}$.\newline
Moreover, if we let 

\begin{equation} 
\beta_2 = \left(1+ \mu + 2^{3/2} \mu^{1/2} \delta^{1/2} {\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} \tilde{\eta}_1 (x,y,0^+) dxdy\right), \label{E2.59} 
\end{equation} 

\begin{equation}
\beta_3 = 4 + 2 \mu^2 + 16 \mu \delta \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2} \tilde{\eta}_1 (x,y,0^+) dxdy, \label{E2.60}
\end{equation}

\begin{equation}
\tilde{\eta}_2 = \tilde{\eta}_1 + \frac{\alpha^2 \delta}{8 \beta_2} p^2_t = \frac{\mu T^2_0}{2 \delta} + \frac{\lambda^2 \delta}{2} | \bigtriangledown p|^2 + \frac{(1-p^2)^2}{4 \delta} + \frac{\alpha^2 \delta}{8 \beta_2} p^2_t, \label{E2.61} 
\end{equation} 

\begin{equation} 
 \tilde{\mbox{\boldmath ${q}$}}_2 =  \tilde{\mbox{\boldmath ${q}$}}_1 + 
\frac{\alpha \lambda^2 \delta^2}{4 \beta_2} p_t \bigtriangledown p_t = \lambda^2 \delta p_t \bigtriangledown p + \frac{\alpha \lambda^2 \delta^2}{4 \beta_2} p_t \bigtriangledown p_t, \label{E2.62} 
\end{equation} 

\begin{equation}
\begin{array}{llll}
\tilde{G}_2 = \tilde{G}_1 + {\displaystyle \frac{\alpha}{4 \beta_2}} 
\left( 1+ {\displaystyle \frac{\mu}{2}} + 2^{3/2} \mu^{1/2} \delta^{1/2} {\displaystyle 
\int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} \tilde{\eta}_1 (x,y,0^+) dxdy \right)p^2_t \\
\\
+ {\displaystyle \frac{\alpha \mu}{8 \beta_2}} {\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} p^2_t (x,y,t) dxdy - {\displaystyle \frac{\alpha \lambda^2 \delta^2}{4 \beta_2}} |\bigtriangledown p_t|^2, \end{array} \label{E2.63} 
\end{equation}

\begin{equation} 
\tilde{\eta}_3 = \tilde{\eta}_2 + \frac{\alpha^2 \lambda^2 \delta^3}{8 \beta_3} | \bigtriangledown p_t |^2 = \frac{\mu T^2_0}{2 \delta} + \frac{\lambda^2 \delta}{2} | \bigtriangledown p|^2 + \frac{(1-p^2)}{4 \delta}+ \frac{\alpha^2 \delta}{8 \beta_2} p^2_t + \frac{\alpha^2 \lambda^2 \delta^3}{8 \beta_3} | \bigtriangledown p^2_t| \label{2.64} 
\end{equation} 

\begin{equation} 
 \tilde{\mbox{\boldmath ${q}$}}_3 =  \tilde{\mbox{\boldmath ${q}$}}_2 + \frac{\alpha^2 \lambda^2 \delta^3}{4 \beta_3} p_{tt} \bigtriangledown p_t = \lambda^2 \delta p_t \bigtriangledown p + \frac{\alpha \lambda^2 \delta^2}{4 \beta_2}p_t \bigtriangledown p_t + \frac{\alpha^2 \lambda^2\delta^3}{4 \beta_3}p_{tt} \bigtriangledown p_t, \label{E2.65} 
\end{equation}

\noindent and 

\begin{equation}
\begin{array}{lllll}
\tilde{G}_3 = \tilde{G}_2 + {\displaystyle \frac{\alpha}{4 \beta_3}} \left( 4 + 16 \mu \delta {\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} \tilde{\eta}_1 (x,y,0^+) dxdy \right) p^2_t \\
\\
+ {\displaystyle \frac{\alpha}{2 \beta_3}} \mu^2 {\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} p^2_t (x,y,t) dxdy - {\displaystyle \frac{\alpha^2 \delta^2}{8 \beta_3}} p^2_{tt} \end{array} 
\label{E2.66}
\end{equation}

\noindent we obtain the inequalities 

\begin{equation}
\frac{\partial \tilde{\eta}_2}{\partial t} - {\rm div}\  \tilde{\mbox{\boldmath ${q}$}}_2 \leq \tilde{G}_2\ {\rm and}\ \frac{\partial \tilde{\eta}_3}{\partial t} = {\rm div}  \tilde{\mbox{\boldmath ${q}$}}_3 \leq \tilde{G}_3 \label{E2.67} 
\end{equation}

\noindent and these imply 

\begin{equation}
\begin{array}{lllll}
{\displaystyle \frac{d}{dt}} {\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} \left( 
{\displaystyle \frac{\mu T^2_0}{2 \delta}} + 
{\displaystyle \frac{\lambda^2 \delta |\bigtriangledown p|^2}{2}} + 
{\displaystyle \frac{\left(1-p^2\right)^2}{4 \delta}} + 
{\displaystyle \frac{\alpha^2 \delta}{8 \beta_2}} p^2_t \right) (x,y,t) dxdy\\
\\
\leq - {\displaystyle \frac{3 \alpha}{4}} {\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} p^2_t (x,y,t) dxdy - {\displaystyle \frac{\alpha \lambda^2 \delta^2}{4 \beta_2}} 
{\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}}| \bigtriangledown p_t|^2 (x,y,t) dxdy, \end{array}   \label{E2.68}
\end{equation}

\noindent and

\begin{equation}
\begin{array}{lll}
{\displaystyle \frac{d}{dt}} {\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} 
\left({\displaystyle \frac{\mu T^2_0}{2 \delta}} + 
{\displaystyle \frac{\lambda^2 \delta |\bigtriangledown p|^2}{2}} + 
{\displaystyle \frac{\left(1-p^2\right)^2}{4 \delta}} + 
{\displaystyle \frac{\alpha \delta}{8 \beta_2}} p^2_t + 
{\displaystyle \frac{\alpha^2 \lambda^2 \delta^3}{8 \beta_3}} | \bigtriangledown p_t|^2\right) (x,y,t) dxdy \\
\\
\leq - {\displaystyle \frac{\alpha}{2}} 
{\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}}p^2_t (x,y,t) - 
{\displaystyle \frac{\alpha \lambda^2 \delta^2}{4 \beta_2}} 
{\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}}| \bigtriangledown p_t|^2 (x,y,t) dxdy \\
\\
- {\displaystyle \frac{\alpha^2 \delta^2}{4 \beta_3}} 
{\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} p^2_{tt} (x,y,t) dxdy. \end{array} \label{E2.69} 
\end{equation}

We leave it to the reader to supply the precise conclusions derivable from the inequalities (\ref{E2.53}), (\ref{E2.63}), and (\ref{E2.69}). The key point of these inequalities is that the parameter $\delta$ appears optimally.

\vspace{.15in}

\section{Computational Experiments}

\vspace{.15in}

In this section we present some numerical simulations of the system (\ref{E2.2}) - (\ref{E2.4}). 
We let $N$ be an integer and set 

\begin{equation}
dx = dy = 1/N. \label{E3.1}
\end{equation}

\noindent Our approximate solutions will be piecewise constant on the cells

\begin{equation} 
\Omega_{i,j} = \left\{(x,y) \left|- \frac{1}{2} +  \frac{(i-1)}{N} < x < -\frac{1}{2} +  \frac{i}{N}\ {\rm and}\ - \frac{1}{2} +  \frac{(j-1)}{N} < y < - \frac{1}{2} + \frac{j}{N} \right\} \right. \label{E3.2} 
\end{equation} 

\noindent where $1 \leq i,j \leq N$ and thus the cell averages $\left(T_{i,j}, p_{i,j}\right)$ may also be thought of as the point values at the cell centers 

\begin{equation} 
\left(x_i,y_j\right) = \left( - \frac{1}{2} + \frac{2i-1}{2 N},\ - \frac{1}{2} + 
 \frac{2j-1}{2N}\right). \label{E3.3} 
\end{equation} 

\noindent We impose the boundary conditions 

\begin{equation} 
\frac{\partial T}{\partial n} = \frac{\partial p}{\partial n} = 0 \ \ \ , \ \ \ (x,y) \in \partial \Omega \label{E3.4}
\end{equation} 

\noindent by insisting that the extended grid functions satisfy 

\begin{equation} 
\left. \begin{array}{lll}
& \left(T_{0,k},p_{0,k}\right) = \left(T_{1,k},p_{1,k}\right), \\
\\
& \left(T_{N+1,k},p_{N+1,k}\right) = \left(T_{N,k}, p_{N,k}\right), \\
\\
& \left(T_{k,0},p_{k,0}\right) = \left(T_{k,1},p_{k,1}\right), \\
\\
{\rm and}\\
\\
& \left(T_{k,N+1},p_{k,N+1}\right) = \left(T_{k,N},p_{k,N}\right) \end{array} \right\} \label{E3.5} 
\end{equation} 

\noindent for $1 \leq k \leq N$. 

For any function $u$ defined on the grid and satisfying (\ref{E3.5}) we let 

\begin{equation} 
{\Delta}_N u_{i,j} = N^2 \left(u_{i+1,j} + u_{i-1,j} + u_{i,j+1} + u_{i,j-1} - 4u_{i,j}\right) \label{E3.6}
\end{equation}

\noindent be the discrete laplacian applied to $u$.  For functions $u$ and $v$ satisfying (\ref{E3.5}) we define 

\begin{equation}
[u,v] = {\displaystyle \sum^N_{i,j=1}} u_{i,j} v_{i,j}/N^2 \label{E3.7}
\end{equation} 

\noindent and note that 

\begin{equation}
-[u,\Delta_Nu] = {\displaystyle \sum^N_{i=1}} {\displaystyle \sum^{N-1}_{j=1}} \left(u_{i,j+1} - u_{i,j}\right)^2 + {\displaystyle \sum^N_{j=1}} {\displaystyle \sum^{N-1}_{i=1}} \left(u_{i+1,j}-u_{i,j}\right)^2 \label{E3.8}. 
\end{equation} 

\noindent Such grid functions also satisfy the discrete version of Lemma 1, namely the inequality 

\begin{equation} 
[u^2,u^2]=[u^4,1] \leq 4 \left(\frac{(1+ \epsilon)}{\epsilon} [u,u] - \epsilon [u, \Delta_N u] \right)^2 \label{E3.9} 
\end{equation} 

\noindent for every $\epsilon >  0$. 

Finally, we replace (\ref{E2.2}) - (\ref{E2.4}) by the system 

\begin{equation} 
\dot{T}_{i,j} + (1-p^2_{i,j})\dot{p}_{i,j} - a^2 \Delta_N T_{i,j} = 0 \label{E3.10} 
\end{equation}

\noindent and 

\begin{equation}
\alpha \delta \dot{p}_{i,j} - \lambda^2 \delta^2 \Delta_N p_{i,j} = (\mu T_{i,j}+ p_{i,j})(1-p^2_{i,j}) \label{E3.11} 
\end{equation} 

\noindent where $T_{i,j}$ and $p_{i,j}$ satisfy the discrete boundary conditions (\ref{E3.5}) and the symbol $\cdot$ denotes differentiation with respect to $t$.  Solutions of (\ref{E3.10}) and (\ref{E3.11}) obey the pointwise estimate that if $-1 \leq p_{i,j}(0^+) \leq 1,\ 1 \leq i,j \leq N$, then $-1 \leq p_{i,j}(t) \leq 1,\ 1 \leq i,j \leq N$ for all $t > 0$.  This system also satisfies estimates which are the exact analogue of those obtained in section (2).  Corresponding to (\ref{E2.8}) we have the identity 

\begin{equation} 
[T+p-p^3/3, 1](t) \equiv h_0 \label{E3.12} 
\end{equation} 

\noindent while the analogues of (\ref{E2.11}) and (\ref{E2.15}) are the inequalities 

\begin{equation}\hspace{-.25in} 
{\displaystyle \frac{d}{dt}} \left( \left[ \left( {\displaystyle \frac{\mu T^2}{2 \delta}} + {\displaystyle \frac{(1-p^2)}{4\delta}}^2 \right),1\right] - 
{\displaystyle \frac{\lambda^2 \delta}{2}} [p,\Delta_Np] \right) (t) = - \alpha [p_t,p_t](t) + {\displaystyle \frac{\mu a^2}{\delta}} [T, \Delta_N T] (t),
\end{equation}

\noindent and

\begin{equation}
\begin{array}{llll}
{\displaystyle \frac{d}{dt}} \left( \left[{\displaystyle \frac{\mu T^2}{2 \delta}} + {\displaystyle \frac{(1-p^2)^2}{4 \delta}},1\right] - {\displaystyle \frac{\lambda^2 \delta}{2}} \left[p,\Delta_Np\right] - {\displaystyle \frac{\alpha a^2}{2}} \left[T, \Delta_N T\right]\right) (t) \\
\\
\leq -{\displaystyle \frac{\alpha}{2}} \left(\left[p_t,p_t\right] + \left[T_t,T_t\right] \right) (t) + {\displaystyle \frac{\mu a^2}{\delta}} \left[T,\Delta_NT\right](t). \end{array} \label{E3.14} 
\end{equation}

We spare the reader the detailed analogues of (\ref{E2.39}), (\ref{E2.41}), and (\ref{E2.45}).  Suffice it to say if we replace terms like ${\displaystyle \int^{1/2}_{-1/2} \int^{1/2}_{-1/2}} f(x,y,t)dxdy$  by $\left[ f,1\right](t)$ and terms like 
${\displaystyle \int^{1/2}_{-1/2}  \int^{1/2}_{-1/2}} \mid \nabla f\mid ^2 (x,y,t)dxdy$ by $- \left[ f, \Delta _N f \right] (t)$, then the resulting inequalities go through intact.

When doing simulations with the discrete system (\ref{E3.10}) and (\ref{E3.11}) we link the small parameter $\delta$ to the grid 
$dx = dy = 1/N$; specifically we set 

\begin{equation}
\delta = 1/N 
\label{E3.15}
\end{equation}
   
\noindent and restrict our attention to initial data satisfying

\begin{equation}
N\left[ \frac{\mu T^2}{2} + \frac{(1-p^2)^2}{4},1 \right] \left(0^+\right) - \frac{\lambda ^2}{2N} \left[ p, \Delta_N p \right]\left(0^+\right) - \frac{\alpha a^2}{2} \left[ T, \Delta_N T\right] \left(0^+\right) = 0(1). \label{E3.16}
\end{equation}

\noindent Interesting data satisfying (\ref{E3.16}) may be achieved in a variety of ways.  Below we outline one such construction.  We let $\Gamma$ be a simple closed curve of finite arclength $L(\Gamma)$ which is wholly contained in $\Omega$, ${\cal S}_c$ be all cells $\Omega_{i,j}$ which are cut by $\Gamma$, ${\cal S}_I$ be all cells $\Omega _{i,j}$ which are interior to the region surrounded by $\Gamma$, and ${\cal S}_E$ be those cells $\Omega_{i,j}$ which are in $\Omega$ and exterior to $\Gamma$.  For definiteness we choose

\begin{equation}
p_{i,j} (0^+) = \left\{ \begin{array}{lll}
1, & {\rm if}\ \Omega_{i,j} \in {\cal S}_E\\
p^*_{i,j} \in \left[ -1,1\right], & {\rm if}\ \Omega_{i,j} \in {\cal S}_c\\
-1, & {\rm if}\ \Omega_{i,j} \in {\cal S}_I 
\end{array} \right. \label{E3.17}
\end{equation}

\noindent and  

\begin{equation}
T_{i,j} \left( 0^+ \right) = h_0 - \left[ \left(p - p^3/3\right) \left( 0^+ \right), 1 \right] = 0 \left(1/N^{1/2}\right). \label{E3.18}
\end{equation}

\noindent Such data models the situation where initially water surrounds ice.  With such data it is easily checked that $N\left[ T^2 \left( 0^+ \right),1 \right] = 0(1)$, while $N\left[ (1-p^2)^2,1\right]$ and $- \frac{1}{N} \left[ p,\Delta _N p \right]$ are $0(L(\Gamma))$ and thus the functional defined in (\ref{E3.16}) is $0(1)$ independently of $N$.  For such data we note that $\delta^3 \left[ \Delta_N p, \Delta _N p \right] = \frac{1}{N^3} \left[ \Delta _N p, \Delta_N p \right]$ is also $0(L(\Gamma))$ independently of $N$.  The last functional enters in the discrete analogue of (\ref{E2.68}).

The simulations we present were all run with the following parameters:

\begin{equation}
\left( \alpha,\lambda,a,\mu \right) = \left(1,1,1,1 \right). \label{E3.19}
\end{equation}

\noindent We chose $N = 200$

\begin{equation}
dx = dy = \delta = 1/200 \label{E3.20}
\end{equation}

\noindent and 

\begin{equation}
dt = dx/10. \label{E3.21}
\end{equation}

\noindent With these parameters

\begin{equation}
m1 \stackrel{def}{=} dt/(dx)^2 = 10\label{E3.22}
\end{equation}

\noindent and

\begin{equation}
m2 \stackrel{def}{=} \frac{\delta dt}{(dx)^2} = \frac{1}{10} \label{E3.23}
\end{equation}

\noindent and thus we were compelled to use an implicit integration scheme for (\ref{E3.10}).  Our choice was ADI working with the fixed tri-diagonal matrix whose non-zero entries are

\begin{equation}
A(i,i) = \left\{ \begin{array}{lll}
1+m1 & {\rm if}\ i = 1\ {\rm and}\ i=200\\
1+2m1 & {\rm if}\ 2 \leq i \leq 199\\
\end{array}\right. \label{E3.24}
\end{equation}

\noindent and 

\begin{equation}
A(i,i+1) = A(i+1,i) = -m1,\ \ 1 \leq i \leq 199. \label{E3.25}
\end{equation}

\noindent We used a first order explicit Euler scheme on (\ref{E3.11}).

Our first simulation was run with the following choice of initial data: 

\begin{equation} 
p_{i,j} (0^+) = \left\{ \begin{array}{llllllll}
1,\ {\rm if} \begin{array}{llllll}
\{ 101 + i \leq j \leq 200 & {\rm and} & 1 \leq i \leq 99\}\ {\rm or} \\
\{302 - i \leq j \leq 200 & {\rm and} & 102 \leq i \leq 200 \} \ {\rm or} \\
\{1 \leq j \leq 100 -i & {\rm and}&  1 \leq i \leq 99\}\ {\rm or}\\
\{1 \leq j \leq -101 + i& {\rm and}& 102 \leq i \leq 200 \},\end{array} \\
\\
0, \ {\rm if}\ \begin{array}{lllll}
\{100 + i=j& {\rm and}& 1 \leq i \leq 100\}\ {\rm or}\\
\{301-i=j& {\rm and}& 101 \leq i \leq 200 \}\ {\rm or}\\
\{101-i=j& {\rm and}& 1 \leq i \leq 100 \}\ {\rm or}\\
\{-100 + i = j&  {\rm and}& 101 \leq i \leq 200 \}, \end{array}\\
\\
-1,\ {\rm if} \begin{array}{lllll}
\{102 -i \leq j \leq 99+i & {\rm and}& 2 \leq i \leq 100\}\ {\rm or} \\
\{-99+i \leq j \leq 300 -i& {\rm and}& 101 \leq i \leq 199\}. \end{array} 
\end{array} \right. \label{E3.26}    
\end{equation}

\noindent This data satisfies $[p,1] = [p^3,1]=0$.   We chose $h_0 = .05$ and $T_{i,j}(0^+) \equiv .05,\ 1 \leq i,j \leq 200$.  

The initially superheated square block of ice converged to a circular block centered at $\left(0,0\right)$ with a radius consistent with the $\delta = 0^+$ asymptotics.  Specifically, if we grant that the ``equilibrium" curve separating the ice from the water is a circle of radius $R$, then the ``Gibbs-Thomson" relation (\ref{E1.6}), with our choice of parameters, implies that 

\begin{equation} 
\frac{.005}{2^{1/2}R} = -T_M \label{E3.27} 
\end{equation}

\noindent while (\ref{E3.12}) implies that 

\begin{equation} 
T_M = .05- [p-p^3/3,1]. \label{E3.28} 
\end{equation} 

\noindent Noting that $p_{i,j} \simeq -1$ inside of the circle and $p_{i,j} \simeq +1$ outside yields 

\begin{equation} 
[p-p^3/3,1] = \frac{2}{3} \left(1-2 \pi R^2\right) \label{E3.29} 
\end{equation} 

\noindent and (\ref{E3.27})-(\ref{E3.29}) imply that $R$ satisfies 

\begin{equation} 
\frac{4 \pi}{3} R^2 + \frac{.005}{2^{1/2}R} = \frac{2}{3} - .05. \label{E3.30} 
\end{equation} 

\noindent ({\ref{E3.30}) has one positive solution close to zero and a second larger solution $R=.3808$ which is the approximate equilibrium radius we obtain in our simulation.  This radius, when substituted into the Gibbs-Thomson relation (\ref{E3.27}) yields a melt temperature consistent with what we observe.  The results of this simulation are shown in Figure 1.

Our final simulation starts with 16 regularly spaced circular blocks of ice surrounded by water.  Specifically 

\begin{equation} 
p_{i,j}(0^+) = \left\{\begin{array}{llllll}
-1,\ {\rm if} &  \sqrt{(i-40m)^2+(j-40n)^2} \leq 15 \\
& {\rm and}\ 1 \leq m, n \leq 4\\
1, & {\rm otherwise}. \end{array} \right. \label{E3.31}
\end{equation}

\noindent For this data 

\begin{equation}
[(p-p^3/3)(0^+),1] = .2885. \label{E3.32}
\end{equation}

\noindent We again chose $h_0$ to be $.05$ and at $t=0$ took the uniform temperature field 

\begin{equation}
T_{i,j}(0^+) = h_0 - [(p-p^3/3)(0^+),1] = -.2385. \label{E3.33}
\end{equation} 

\noindent The solutions using this initial data converged to the same equilibrium as before but the approach was far more interesting.  The solutions at times $t=.05,\ .17,\ .24,\ .4,\ .6,\ .7,\ 2, 3.3,$ and 12.5 are shown in Figures 2-10. 

We also ran both sets of data with the ``intermediate" model and obtained virtually identical results. 

\vspace{.15in}

\centerline{\bf Explanation of Graphics}

\vspace{.15in}

In Figures 1-10 the upper left frame is a contour plot of $p$ and $p$ is 1 in a neighbourhood of the boundaries $x = \pm 1/2$ and $y = \pm 1/2$ and the upper right hand frame is a contour plot of $T$.  In Figure 1 the lower left plot shows $p,T$ and $p-p^3/3 + T$ along  the lines $x=0,\ y=0$, and the 45$^o$ diagonal through the origin.  These plots overlay one another because of the circular symmetry of the final state.  In Figure 1 the bottom right plot again displays the same information.  In Figures 2-10 the bottom left plot shows $p, T$, and $p-p^3/3 + T$ along the $45^o$ diagonal through the origin and the bottom right plot shows the same quantities on the line $x=0$.


\newpage
\begin{figure}[h]
\vspace*{6.7in}
\special{psfile=fig5.ps hoffset=-40 voffset=-200 hscale=100 vscale=100 angle=0}
\caption{}
\end{figure}

\newpage
\begin{figure}[h]
\vspace*{6.7in}
\special{psfile=frame1004.ps hoffset=-40 voffset=-200 hscale=100 vscale=100 angle=0}
\caption{}
\end{figure}

\newpage
\begin{figure}[h]
\vspace*{6.7in}
\special{psfile=frame1016.ps hoffset=-40 voffset=-200 hscale=100 vscale=100 angle=0}
\caption{}
\end{figure}

\newpage
\begin{figure}[h]
\vspace*{6.7in}
\special{psfile=frame1023.ps hoffset=-40 voffset=-200 hscale=100 vscale=100 angle=0}
\caption{}
\end{figure}

\newpage
\begin{figure}[h]
\vspace*{6.7in}
\special{psfile=frame1039.ps hoffset=-40 voffset=-200 hscale=100 vscale=100 angle=0}
\caption{}
\end{figure}

\newpage
\begin{figure}[h]
\vspace*{6.7in}
\special{psfile=frame1059.ps hoffset=-40 voffset=-200 hscale=100 vscale=100 angle=0}
\caption{}
\end{figure}

\newpage
\begin{figure}[h]
\vspace*{6.7in}
\special{psfile=frame1069.ps hoffset=-40 voffset=-200 hscale=100 vscale=100 angle=0}
\caption{}
\end{figure}

\newpage
\begin{figure}[h]
\vspace*{6.7in}
\special{psfile=frame1109.ps hoffset=-40 voffset=-200 hscale=100 vscale=100 angle=0}
\caption{}
\end{figure}

\newpage
\begin{figure}[h]
\vspace*{6.7in}
\special{psfile=frame1122.ps hoffset=-40 voffset=-200 hscale=100 vscale=100 angle=0}
\caption{}
\end{figure}

\newpage
\begin{figure}[h]
\vspace*{6.7in}
\special{psfile=frame1214.ps hoffset=-40 voffset=-200 hscale=100 vscale=100 angle=0}
\caption{}
\end{figure}



\newpage
\noindent {\bf References} 

\begin{enumerate} 

\item S.L. Wang, R.F. Sekerka, A.A. Wheeler, B.T. Murray, S.R. Coriell, R.J. Braun, and G.B. McFadden, ``Thermodynamically-Consistent Phase-Field Models, Physica D 69 (1993) 189-200. 

\item G.B. McFadden, A.A. Wheeler, R.J. Brown, S.R. Coriell and R.F. Sekerka, ``Phase-Field Models for Anisotropic Interfaces," Phys. Rev. E 48 (1993) 2016-2024. 

\item S.L. Wang and R.F. Sekerka, ``Algorithms for Phase Field Computations of the Dendritic Operating State at Large Supercoolings," J. Comp. Physics 127 (1996) 110-117.

\item A.S. Almgren and R.F. Almgren, ``Phase Field Instabilities and Adaptive Mesh Refinement," in Modern Methods for Modeling Microstructure in Materials, TMS/SIAM 1996.

\item G. Caginalp, ``Stefan and Hele-Shaw Type Models as Asymptotic Limits of the Phase Field Equation," Phys. Rev. A 39 (1989) 5887-5896.

\item G. Caginalp, ``Phase Field Models and Sharp Interface Limits:  Some Differences in Subtle Situations," Rocky Mountain J. Math. 21 (1996), 2, 603-616.

\end{enumerate}




\end{document}
\end

