\section{Introduction}
The valuation of American-style options poses a challenge for both academic and industrial professionals. One of the difficulties comes from the fact that such a valuation process relies on the identification of an optimal exercise policy. So far, considerable effort has been put into simple settings where the underlying asset price follows a log-normal process and the interest rate is constant (i.e., the standard model, or the Black-Scholes model). Within this context, Kim (1990) decomposed the American option price into two parts: the corresponding European option price and an Early Exercise Premium (EEP) that captures the gains from exercising the option prior to its maturity. Similar results are provided by Jacka (1991) and Carr \emph{et al.} (1992). The EEP representation of the American option price has proved extremely useful because it provides a recursive integral equation for the optimal exercise boundary. Solving the integral equations is key to the valuation process: it identifies the optimal exercise policy, providing a parametric formula for the option price. Such an approach, based on the integral equation, is straightforward to implement and shows significant advantages over other numerical procedures such as methods based on binomial lattices, Monte Carlo simulation, and Partial Differential Equations (PDE). See Brodie and Detemple (2004) for a survey of methods on the valuation of American options.

While the valuation of American options in the standard model has been resolved, empirical evidence suggests that the log-normality assumption does not hold in reality. For example, the “volatility smile” phenomenon is a well-known pattern in option pricing practice. To allow for the consistency of models with empirical regularities, non-constant, or even non-deterministic model parameters should be considered. Unfortunately, analytical results in the standard model can not be generalized to models with stochastic parameters in a straightforward manner. Efforts have been made to solve diffusion models with nonconstant parameters. For example, Jacka and Lynn (1992) considered general contingent claims written on diffusion processes.
Detemple and Tian (2002) presented an integral equation approach for the valuation of American-style derivatives when the underlying asset price follows a general diffusion process and the interest rate is stochastic.
See Rutkowski (1994), and Gukhal (2001) for the valuation of American options for other non-standard models.

All the above-mentioned methods rely on the fact that the transition density of the underlying asset dynamics admits a closed functional form. Such conditions have limited the scope of stochastic processes that can be considered. Furthermore, even when the transition density exists in closed form, the structure may be quite complex, and in turn, the method may be difficult to implement. To overcome these difficulties, this article presents a systematic treatment of the valuation of American options based on Hermite polynomial expansions and the EEP formula. Therefore, our contributions to the literature are threefold.

The first contribution is that our method does not rely on the existence of analytical solutions to the transition density or the characteristic function of the distribution of the underlying asset price. Moreover, there are no requirements for affine structures. We propose using Hermite polynomials to approximate the transition density for a given jump-diffusion model. The Hermite polynomial approximation is based on Ait-Sahalia (2002, 2008), and Yu (2007). This approach gives an explicit sequence of closed-form solutions to the transition density and is shown to converge to the turn density. See Ait-Sahalia (1999), Egorov, Li, and Xu (2003), Ait-Sahalia and Kimmel (2007, 2010), and Xiu (2014) for studies related to this approach.

The second contribution is that, our method is fast and accurate, while the price and exercise policy can be jointly approximated by our approximation scheme. Owing to the inherent nature of the EEP approach, by solving the integral equations with the Hermite polynomial-based approximation to the transition density, we can generate an approximation of the optimal exercise boundary. In turn, we provide a theorem (Theorem \ref{thm:aoh} in Subsection \ref{subsection:hp}) on the convergence of our proposed approximation. When we increase the order of the Hermite polynomial in the approximation of the transition density, the proposed approximations of the price and exercise boundary of American options further improve. We can control the smoothness and accuracy of the exercise boundary by changing the polynomial order in the expansion of the transition density.

Third, our method can be easily extended to jump-diffusion models and multidimensional cases. The extension is straightforward to implement without additional theoretical/modeling complications. Kou (2002) established the analytical solutions for European option pricing in a jump-diffusion model. However, the American option pricing with jump-diffusion processes remains challenging. Gukhal (2001) derived an EEP formula for the value of American options in a jump-diffusion model. Despite all these efforts, one major drawback is that jump-diffusion models usually come without closed-form transition densities. Even if such a density exists, its functional form may be quite complicated in structure and the implementation requires a significant amount of human and computer power. Our method, on the other hand, can overcome these difficulties. By using a Hermite polynomial expansion, we can control the computational cost of the pricing algorithm by specifying the order of the expansion. Moreover, unlike conventional approaches such as the PDE-based method (finite difference, for example), Hermite polynomial expansions can be applied to a vector of stochastic processes, and the results can be directly applied to multi-dimensional models.

The structure of this paper is as follows. Section \ref{section:ao_2} describes an approach to the American option valuation when the underlying asset prices follow a general diffusion process. Section \ref{section:ao_3} describes the generalization of the method to jump-diffusion processes. Section \ref{section:ao_4} presents a numerical algorithm for implementing the proposed approach. Section \ref{section:ao_5} provides several examples to demonstrate the efficiency of the proposed approach. Finally, Section \ref{section:ao_6} concludes this paper.

\section{Valuation of American Options in Diffusion Models}\label{section:ao_2}

\subsection{American Options}

We consider the stock price $S$ defined on a probability space $\left(\varOmega,\mathcal{F},\mathbb{P}^{*}\right)$
with filtration $\mathbb{F}=\left(\mathcal{F}_{t}\right)_{0\leq t\leq T}$
satisfying the usual conditions and following:
\begin{equation}
dS_{t}=\left(r\left(S_{t};\theta\right)-\delta\left(S_{t};\theta\right)\right)dt+\sigma\left(S_{t};\theta\right)dW_{t}\label{eq:aomodel}.
\end{equation}
We also denote $\mu\left(S_{t};\theta\right)=r\left(S_{t};\theta\right)-\delta\left(S_{t};\theta\right)$.
Let $D_{S}=\left(\underline{s},\bar{s}\right)$ be the domain of the
diffusion $S$.

The arbitrage-free price of an American put option with a finite expiration $T>0$ and a strike price $K$ can be expressed as the expected
value of its discounted payoff:
\begin{equation}
P\left(t,S_{t}\right)=\sup_{t\leq\tau\leq T}\mathbb{E}^{*}\left[e^{-\left(\tau-t\right)r}\left(K-S_{\tau}\right)^{+}|S_{t}\right]\label{eq:aovalue}
\end{equation}
under the risk-neutral probability measure $\mathbb{P}^{*}$. Here
$\tau$ is the stopping time.

\subsection{Early Exercise Boundary}

Let $\mathcal{B}=\left\{ B_{t}:B_{t}\geq0,t\in\left[0,T\right]\right\} $
denote the optimal early exercise boundary of the American put option.
Then the arbitrage-free price of the American put option, $P\left(t,S_{t}\right)$, solves
the following free boundary problem:
\[
\mathcal{L}P=0,
\]
\[
P\left(T,S_{T}\right)=\left(K-S_{T}\right)^{+},
\]
\[
\lim_{S_{t}\uparrow\infty}P\left(t,S_{t}\right)=0,
\]
\[
\lim_{S_{t}\downarrow B_{t}}P\left(t,S_{t}\right)=K-B_{t},
\]
\[
\lim_{S_{t}\downarrow B_{t}}\frac{\partial P\left(t,S_{t}\right)}{\partial S_{t}}=-1,
\]
where $\mathcal{L}f=\frac{1}{2}\sigma^{2}S_{t}^{2}\frac{\partial^{2}f}{\partial S_{t}^{2}}+\left(r-\delta\right)S_{t}\frac{\partial f}{\partial S_{t}}-rf+\frac{\partial f}{\partial t}.$

\begin{theorem}[Exercise Premium Representation]\label{thm:eep}
We assume that $r$, $\delta$,
and $\sigma$ are continuously differentiable, and (\ref{eq:aomodel})
has a unique strong solution. Then, in the continuation region $\mathcal{C}$\footnote{The continuation region is the set of pairs $(S,t)$ at which immediate exercise is sub-optimal.},
the value of the American put option, $P_{0}\equiv P\left(0,S_{0}=s_{0}\right)$,
has the following early exercise premium representation:
\begin{equation}
P_{0}=p_{0}+e_{0}\label{eq:ao}
\end{equation}
where $p_{0}$ represents the price of a European put option, that is,
\begin{equation}
p_{0}\equiv p\left(0,S_{0}=s_{0}\right)=\int_{0}^{K}e^{-rT}\left(K-S_{T}\right)\psi\left(S_{T};S_{0}=s_{0}\right)dS_{T}\label{eq:ao_eo},
\end{equation}
and $e_{0}$ is the early exercise premium given by
\begin{equation}
e_{0}\equiv e\left(0,S_{0}=s_{0},B\left(\cdot\right)\right)=\int_{0}^{T}\int_{0}^{B_{t}}\left(rK-\delta S_{t}\right)e^{-rt}\psi\left(S_{t};S_{0}=s_{0}\right)dS_{t}dt\label{eq:ao_eep},
\end{equation}
and $\psi\left(S_{t};S_{0}=s_{0}\right)$ denotes the risk-neutral
transitional density function of $S_{t}$, given $S_{0}=s_{0}$. The
exercise boundary $B_{t}$ solves the recursive nonlinear integral
equation
\begin{equation}
K-B_{t}=p\left(t,B_{t}\right)+e\left(t,B_{t},B\left(\cdot\right)\right)\quad\forall t\in\left[0,T\right),\label{eq:ao_bdy}
\end{equation}
subject to the boundary condition
\[
B_{T-}\equiv\lim_{t\uparrow T}B_{t}=\min\left\{ K,\frac{r\left(B_{T};\theta\right)}{\delta\left(B_{T};\theta\right)}K\right\}.
\]
At maturity, $B_{T}=K\geq B_{T-}$. The functions $p$ and $e$ in
(\ref{eq:ao_bdy}) are defined as following:
\begin{equation}
p\left(t,B_{t}\right)\equiv\int_{0}^{K}\left(K-S_{T}\right)\psi\left(S_{T};S_{t}=B_{t}\right)dS_{T},\label{eq:ao_eo2}
\end{equation}
\begin{equation}
e\left(t,B_{t},B\left(\cdot\right)\right)\equiv\int_{t}^{T}\int_{0}^{B_{s}}\left(rK-\delta S_{s}\right)e^{-r\left(s-t\right)}\psi\left(S_{s};S_{t}=B_{t}\right)dS_{s}ds.\label{eq:ao_eep2}
\end{equation}
\end{theorem}
The equations (\ref{eq:ao_eo})-(\ref{eq:ao_bdy}) in Theorem \ref{thm:eep} for the valuation
of American options can be simplified if we make further assumptions
on the model. For example, if we assume that the stock price follows a geometric
Brownian Motion (GBM), that is, $r\left(S_{t};\theta\right)=r$, $\sigma\left(S_{t};\theta\right)=\sigma$
for constants $r$ and $\sigma$, and $\delta\left(S_{t};\theta\right)=0$,
then we have a Black-Scholes style formula for the valuation of American
options. We summarize this result in the following lemma.

\begin{lemma}[Exercise Premium Representation under GBM]\label{lemma:gbm}
If the stock price
$S$ follows geometric Brownian Motion, then the value of the American put option,
$P_{0}$, can be written as:
\begin{equation}
P_{0}=Ke^{-rT}N\left(k_{2}\left(S_{0},K,T\right)\right)-S_{0}N\left(k_{1}\left(S_{0},K,T\right)\right)+rk\int_{0}^{T}e^{-rt}N\left(b_{2}\left(S_{0},B_{t},t\right)\right)dt,\label{eq:ao_gbm}
\end{equation}
where
\[
k_{1}\left(S_{0},K,T\right)\equiv\frac{\log\left(K/S_{0}\right)-\rho_{1}T}{\sigma\sqrt{T}},
\]
\[
k_{2}\left(S_{0},K,T\right)\equiv\frac{\log\left(K/S_{0}\right)-\rho_{2}T}{\sigma\sqrt{T}},
\]
\[
\rho_{1}\equiv\rho_{2}+\sigma^{2}=r+\frac{\sigma^{2}}{2},
\]
\[
b_{2}\left(S_{0},B_{t},t\right)\equiv\frac{\log\left(B_{t}/S_{0}\right)-\rho_{2}t}{\sigma\sqrt{t}},
\]
and $B_{t}$ solves the following integral equation:
\begin{multline}
K-B_{t}=Ke^{-r\left(T-t\right)}N\left(k_{2}\left(B_{t},K,T-t\right)\right)-B_{t}N\left(k_{1}\left(B_{t},K,T-t\right)\right)\\
+rk\int_{t}^{T}e^{-r\left(s-t\right)}N\left(b_{2}\left(B_{t},B_{s},s-t\right)\right)ds.\label{eq:ao_gbm_bdy}
\end{multline}
\end{lemma}

\subsection{Hermite Polynomial-based Approximation}\label{subsection:hp}

Theorem \ref{thm:eep} provides an intuitive approach to the valuation of American
options in diffusion models; however, we still have two difficulties.
First, most of the diffusion models do not admit a closed-form solution
for the transition density. Second, the exercise boundary $\mathcal{B}$
is unknown in (\ref{eq:ao_eep}), and we need to solve the integral
equation (\ref{eq:ao_bdy}) recursively to compute $\mathcal{B}$.

In this study, we propose the use of the Hermite polynomials to approximate
the transition density. Our approach is based on the work of Ait-Sahalia (2002, 2006) and
Yu (2007). The Hermite polynomial approach by Ait-Sahalia (2002) provided
an explicit sequence of closed-form functions to approximate the unknown
transition density. Ait-Sahalia (2006) and Yu (2007) extended the approach
to multivariate case and jump-diffusion models.

To approximate the transition density of the stock price $S$, we
first transform $S$ into a new random variable $Y$ by defining $Y\equiv\gamma\left(S\right)=\int^{S}du/\sigma\left(u\right)$.
We know that $Y$ has a unit diffusion, that is,
\[
dY_{t}=\mu_{Y}\left(Y_{t};\theta\right)dt+dW_{t},
\]
where
\begin{equation}
\mu_{Y}\left(y;\theta\right)=\frac{\mu\left(\gamma^{-1}\left(y;\theta\right);\theta\right)}{\sigma\left(\gamma^{-1}\left(y;\theta\right);\theta\right)}-\frac{1}{2}\frac{\partial\sigma}{\partial S}\left(\gamma^{-1}\left(y;\theta\right);\theta\right).\label{eq:herm_miuy}
\end{equation}
We denote the domain of $Y$ as $D_{Y}=\left(\underline{y},\bar{y}\right)$.
According to Ait-Sahalia (2002), the transition density of $Y$ can
be approximated using Hermite polynomials, and the transition
density of $S$ can then be derived from that of $Y.$ More specifically,
the transition density of $S$ with time interval $\varDelta$ can be
approximated up to order $m$ as following:
\begin{multline}
\tilde{\psi}^{\left(m\right)}\left(S_{t+\varDelta}=S^{\prime};S_{t}=S\right)=\sigma^{-1}\left(S^{\prime};\theta\right)\varDelta^{-\frac{1}{2}}\phi\left(\frac{\gamma\left(S^{\prime};\theta\right)-\gamma\left(S;\theta\right)}{\varDelta^{\frac{1}{2}}}\right)\times\\
\exp\left(\int_{\gamma\left(S;\theta\right)}^{\gamma\left(S^{\prime};\theta\right)}\mu\left(w;\theta\right)dw\right)\times\sum_{k=0}^{m}c_{k}\left(\gamma\left(S^{\prime};\theta\right)|\gamma\left(S;\theta\right);\theta\right)\frac{\varDelta^{k}}{k!},\label{eq:herm_density}
\end{multline}
where $\phi\left(z\right)\equiv\exp\left(-z^{2}/2\right)/\sqrt{2\pi}$
denotes the density function of standard normal distribution, and
for all $j\geq1$,
\begin{multline}
c_{j}\left(\gamma\left(S^{\prime};\theta\right)|\gamma\left(S;\theta\right);\theta\right)=j\left(S^{\prime}-S\right)^{-j}\int_{\gamma\left(S;\theta\right)}^{\gamma\left(S^{\prime};\theta\right)}\left(w-\gamma\left(S;\theta\right)\right)^{j-1}\\
\times\left\{ \lambda\left(w;\theta\right)c_{j-1}\left(w|\gamma\left(S;\theta\right);\theta\right)+\left(\partial^{2}c_{j-1}\left(w|\gamma\left(S;\theta\right);\theta\right)/\partial w^{2}\right)/2\right\} dw\label{eq:herm_c}
\end{multline}
where $\lambda\left(x;\theta\right)\equiv-\left(\mu_{Y}^{2}\left(x;\theta\right)+\partial\mu_{Y}\left(x;\theta\right)/\partial x\right)/2$
with $\mu_{Y}$ defined in (\ref{eq:herm_miuy}), and $c_{0}=1.$

Once we obtain the approximation of the transition density of $S$ in
(\ref{eq:herm_density}), we can plug $\tilde{\psi}^{\left(m\right)}$
into Theorem \ref{thm:eep} and obtain the approximation of the valuation of American
options. More specifically, we have the following approximated early
exercise premium representation for the value of the American put option up to
order $m$:
\begin{equation}
\tilde{P}_{0}^{\left(m\right)}=\tilde{p}_{0}^{\left(m\right)}+\tilde{e}_{0}^{\left(m\right)}\label{eq:ao_appr}
\end{equation}
where $\tilde{p}_{0}^{\left(m\right)}\equiv\tilde{p}^{\left(m\right)}\left(0,S_{0}=s_{0}\right)=\int_{0}^{K}e^{-rT}\left(K-S_{T}\right)\tilde{\psi}^{\left(m\right)}\left(S_{T};S_{0}=s_{0}\right)dS_{T}$
represents the approximated price of a European put option and $\tilde{e}_{0}^{\left(m\right)}$
is the approximated early exercise premium given by
\begin{equation}
\tilde{e}_{0}^{\left(m\right)}\equiv\tilde{e}^{\left(m\right)}\left(0,S_{0}=s_{0},\tilde{B}^{\left(m\right)}\left(\cdot\right)\right)=\int_{0}^{T}\int_{0}^{\tilde{B}_{t}^{\left(m\right)}}\left(rK-\delta S_{t}\right)e^{-rt}\tilde{\psi}^{\left(m\right)}\left(S_{t};S_{0}=s_{0}\right)dS_{t}dt.\label{eq:ao_eep_appr}
\end{equation}
The approximated exercise boundary up to order $m$, $\tilde{B}_{t}^{\left(m\right)}$,
solves the following recursive nonlinear integral equation:
\begin{equation}
K-\tilde{B}_{t}^{\left(m\right)}=\tilde{p}^{\left(m\right)}\left(t,\tilde{B}_{t}^{\left(m\right)}\right)+\tilde{e}^{\left(m\right)}\left(t,\tilde{B}_{t}^{\left(m\right)},\tilde{B}^{\left(m\right)}\left(\cdot\right)\right)\quad\forall t\in\left[0,T\right).\label{eq:ao_bdy_appr}
\end{equation}
Similarly to (\ref{eq:ao_eo2})--(\ref{eq:ao_eep2}), $\tilde{p}^{\left(m\right)}$
and $\tilde{e}^{\left(m\right)}$ are defined as:
\begin{equation}
\tilde{p}^{\left(m\right)}\left(t,\tilde{B}_{t}^{\left(m\right)}\right)\equiv\int_{0}^{K}\left(K-S_{T}\right)\tilde{\psi}^{\left(m\right)}\left(S_{T};S_{t}=\tilde{B}_{t}^{\left(m\right)}\right)dS_{T},\label{eq:ao_eo2_appr}
\end{equation}
\begin{equation}
\tilde{e}^{\left(m\right)}\left(t,\tilde{B}_{t}^{\left(m\right)},\tilde{B}^{\left(m\right)}\left(\cdot\right)\right)\equiv\int_{t}^{T}\int_{0}^{\tilde{B}_{s}^{\left(m\right)}}\left(rK-\delta S_{s}\right)e^{-r\left(s-t\right)}\tilde{\psi}^{\left(m\right)}\left(S_{s};S_{t}=\tilde{B}_{t}^{\left(m\right)}\right)dS_{s}ds,\label{eq:ao_eep2_appr}
\end{equation}
subject to the boundary condition
\[
\tilde{B}_{T-}^{\left(m\right)}\equiv\lim_{t\uparrow T}\tilde{B}_{t}^{\left(m\right)}=\min\left\{ K,\frac{r\left(B_{T};\theta\right)}{\delta\left(B_{T};\theta\right)}K\right\} ,
\]
and $\tilde{B}_{T}^{\left(m\right)}=B_{T}=K\geq\tilde{B}_{T-}^{\left(m\right)}$.

The following theorem guarantees that the proposed approach in (\ref{eq:ao_appr})--(\ref{eq:ao_eep2_appr})
is a well-behaved approximation of the value of American options.
\newpage
\begin{theorem}\label{thm:aoh}
Under Assumptions \ref{ass_ao_1}--\ref{ass_ao_3} given in Appendix A, as $m\rightarrow\infty,$
we have
\begin{enumerate}
\item $\tilde{p}_{0}^{\left(m\right)}\rightarrow p_{0}$,
\item $\tilde{B}_{t}^{\left(m\right)}\rightarrow B_{t}$ for any $t\in\left[0,T\right]$,
\item $\tilde{e}_{0}^{\left(m\right)}\rightarrow e_{0}$,
\item $\tilde{P}_{0}^{\left(m\right)}\rightarrow P_{0}$.
\end{enumerate}
\end{theorem}

\begin{proof}
\textit{In Appendix B}
\end{proof}

\section{Valuation of American Options in Jump-Diffusion Models}\label{section:ao_3}

\subsection{Valuation of American Options}

In this section, we discuss the approximation of the value of American
options when the underlying asset price follows a jump-diffusion process.
Because of the discontinuous nature of the asset price path, the exercise
premium representation is different from that without jumps. Specifically, we consider the stock price under the risk-neutral measure, and assume that it follows:
\begin{equation}
d\log S_{t}=\left(r\left(S_{t};\theta\right)-\delta\left(S_{t};\theta\right)-\rho j\right)dt+\sigma\left(S_{t};\theta\right) dW_{t}+\left(J-1\right)dq_{t}\label{eq:aomodel_jump}
\end{equation}
where $dq$ is a Poisson process with rate $\rho t$, $J-1$ is the
proportional change in the price due to a jump with density function $\nu$ as a function of jump size with support $D_{J}$, and $j=E\left(J-1\right)$. We assume $r\left(S_{t};\theta\right)$, $\delta\left(S_{t};\theta\right)$, and $\sigma\left(S_{t};\theta\right)$ are smooth functions of $S_{t}$.
Then, based on Gukhal (2001), the value of the American put option, $P_{0}\equiv P\left(0,S_{0}=s_{0}\right)$, has the following
representation:
\begin{equation}
P_{0}=p_{0}+e_{0}+g_{0}\label{eq:ao_jump}
\end{equation}
where
\begin{equation}
p_{0}\equiv p\left(0,S_{0}=s_{0}\right)=\int_{0}^{K}e^{-rT}\left(K-S_{T}\right)\psi\left(S_{T};S_{0}=s_{0}\right)dS_{T}\label{eq:ao_jump_eo}
\end{equation}
\begin{equation}
e_{0}\equiv e\left(0,S_{0}=s_{0},B\left(\cdot\right)\right)=\int_{0}^{T}\int_{0}^{B_{t}}\left(rK-\delta S_{t}\right)e^{-rt}\psi\left(S_{t};S_{0}=s_{0}\right)dS_{t}dt\label{eq:ao_jump_eep}
\end{equation}
and
\begin{multline}
g_{0}\equiv g\left(0,S_{0}=s_{0},B\left(\cdot\right)\right)=\int_{0}^{T}\int_{0}^{B_{t-}}\int_{B_{t}}^{\infty}e^{-rt}\rho\left(P\left(t,J_{t}S_{t-}\right)-\left(K-J_{t}S_{t-}\right)\right)\\
\times\psi\left(S_{t-};S_{0}=s_{0}\right)\psi\left(J_{t}S_{t-};S_{t-}\right)d\left(J_{t}S_{t-}\right)dS_{t-}dt.\label{eq:ao_jump_jump}
\end{multline}
The exercise boundary $B_{t}$ solves the following integral equation
\begin{equation}
K-B_{t}=p\left(t,B_{t}\right)+e\left(t,B_{t},B\left(\cdot\right)\right)-g\left(t,B_{t},B\left(\cdot\right)\right)\label{eq:ao_jump_bdy}
\end{equation}
where
\[
p\left(t,B_{t}\right)=\int_{0}^{K}\left(K-S_{T}\right)\psi\left(S_{T};B_{t}\right)dS_{T}
\]
\[
e\left(t,B_{t},B\left(\cdot\right)\right)=\int_{t}^{T}\int_{0}^{B_{s}}\left(rK-\delta S_{s}\right)e^{-r\left(s-t\right)}\psi\left(S_{s};B_{t}\right)dS_{s}ds
\]
\begin{multline*}
g\left(t,B_{t},B\left(\cdot\right)\right)=\int_{t}^{T}\int_{0}^{B_{s-}}\int_{B_{s}}^{\infty}e^{-r\left(s-t\right)}\rho\left(P\left(s,J_{s}S_{s-}\right)-\left(K-J_{s}S_{s-}\right)\right)\\
\times\psi\left(S_{s-};B_{t}\right)\psi\left(J_{s}S_{s-};S_{s-}\right)d\left(J_{s}S_{s-}\right)dS_{s-}ds.
\end{multline*}

The representation in (\ref{eq:ao_jump}) has a straightforward interpretation.
As in the case without jumps, $p_{0}$ represents the price of a European
put option, $e_{0}$ is the early exercise premium, and $g_{0}$ is the rebalancing
cost due to the jumps of stock prices from the exercise region (the stock price is below
the exercise boundary) into the continuation region (the stock price
is above the exercise boundary).

\subsection{Hermite Polynomial-based Approximation}

Our approach to studying jump-diffusion models is similar to our approach in
Section \ref{section:ao_2}. We first approximate the transition density using Hermite
polynomials. According to Yu (2007), an approximation of the order $m>0$
is obtained as follows:
\begin{multline}
\tilde{\psi}^{\left(m\right)}\left(S_{t+\varDelta}=S^{\prime};S_{t}=S\right)=\varDelta^{-\frac{1}{2}}\exp\left[-\frac{C^{\left(-1\right)}\left(S,S^{\prime}\right)}{\varDelta}\right]\sum_{k=0}^{m}C^{\left(k\right)}\left(S,S^{\prime}\right)\varDelta^{k}\\
+\sum_{k=1}^{m}D^{\left(k\right)}\left(S,S^{\prime}\right)\varDelta^{k}\label{eq:herm_jump_density}
\end{multline}
where
\begin{equation}
C^{\left(-1\right)}\left(S,S^{\prime}\right)=\frac{1}{2}\left[\int_{S}^{S^{\prime}}\sigma\left(s\right)^{-1}ds\right]^{2},\label{eq:herm_jump_c-1}
\end{equation}
\begin{equation}
C^{\left(0\right)}\left(S,S^{\prime}\right)=\frac{1}{\sqrt{2\pi}\sigma\left(S^{\prime}\right)}\exp\left[\int_{S}^{S^{\prime}}\frac{\mu\left(s\right)}{\sigma^{2}\left(s\right)}-\frac{\sigma^{\prime}\left(s\right)}{2\sigma\left(s\right)}ds\right],\label{eq:herm_jump_c0}
\end{equation}
\begin{multline}
C^{\left(k+1\right)}\left(S,S^{\prime}\right)=-\left[\int_{S}^{S^{\prime}}\sigma\left(s\right)^{-1}ds\right]^{-k+1}\int_{S}^{S^{\prime}}\left\{\exp\left[\int_{u}^{S}\frac{\mu\left(u\right)}{\sigma^{2}\left(u\right)}-\frac{\sigma^{\prime}\left(u\right)}{2\sigma\left(u\right)}du\right]\right.\\
\left.\times\sigma\left(s\right)^{-1}\left[\int_{s}^{S^{\prime}}\sigma\left(u\right)^{-1}du\right]^{k}\left[\rho\left(s\right)-\mathcal{L}\right]C^{\left(k\right)}\left(s,S^{\prime}\right)\right\}ds,\quad for\;k\geq0,\label{eq:herm_jump_ck}
\end{multline}
\begin{equation}
D^{\left(1\right)}\left(S,S^{\prime}\right)=\rho\left(S\right)-\upsilon\left(S^{\prime}-S\right),\label{eq:herm_jump_d1}
\end{equation}
\begin{multline}
D^{\left(k+1\right)}\left(S,S^{\prime}\right)=\frac{1}{1+k}\left[\mathfrak{L}D^{\left(k\right)}\left(S,S^{\prime}\right)+\right.\\
\left.\sqrt{2\pi}\rho\left(S\right)\sum_{r=0}^{k}\frac{M_{2r}^{1}}{\left(2r\right)!}\frac{\partial^{2r}}{\partial w^{2r}}m_{k-r}\left(S,S^{\prime},w\right)|_{w=0}\right],\quad for\;k\geq0,\label{eq:herm_jump_dk}
\end{multline}
where
\begin{equation}
m_{k}\left(S,S^{\prime},w\right)\equiv C^{\left(k\right)}\left(w_{B}^{-1}\left(w\right),S^{\prime}\right)\upsilon\left(w_{B}^{-1}\left(w\right)-S\right)\sigma\left(w_{B}^{-1}\left(w\right)\right),\label{eq:herm_jump_m}
\end{equation}
\begin{equation}
M_{2r}^{1}\equiv1/\sqrt{2\pi}\int_{\mathbb{R}}\exp\left(-s^{2}/2\right)s^{2r}ds,\label{eq:herm_jump_M}
\end{equation}
\begin{equation}
w_{B}\left(S,S^{\prime}\right)=\int_{S^{\prime}}^{S}\sigma\left(s\right)^{-1}ds,\label{eq:herm_jump_w}
\end{equation}
\begin{equation}
\mathscr{L}f\left(s,s^{\prime}\right)=\frac{1}{2}\sigma^{2}s^{2}\frac{\partial^{2}f}{\partial s^{2}}\left(s,s^{\prime}\right)+\left(r-\delta-\rho j\right)s\frac{\partial f}{\partial s}\left(s,s^{\prime}\right),\label{eq:herm_jump_A}
\end{equation}
and
\begin{equation}
\mathfrak{L}f\left(s,s^{\prime}\right)=\mathscr{L}f\left(s,s^{\prime}\right)+\rho\int_{D_{J}}\left[f\left(s+c,s^{\prime}\right)-f\left(s,s^{\prime}\right)\right]\upsilon\left(c\right)dc.\label{eq:herm_jump_L}
\end{equation}

Once we obtain the Hermite polynomial approximation of the transition
density as above, we plug the approximation into (\ref{eq:ao_jump_eo})-(\ref{eq:ao_jump_jump}),
and solve the integral equation (\ref{eq:ao_jump_bdy}) recursively.
More specifically, we have the approximated value of the American put option up
to order $m$, $\tilde{P}_{0}^{\left(m\right)}\equiv\tilde{P}^{\left(m\right)}\left(0,S_{0}=s_{0}\right)$:
\begin{equation}
\tilde{P}_{0}^{\left(m\right)}=\tilde{p}_{0}^{\left(m\right)}+\tilde{e}_{0}^{\left(m\right)}+\tilde{g}_{0}^{\left(m\right)}\label{eq:ao_jump_appr}
\end{equation}
where
\begin{equation}
\tilde{p}_{0}^{\left(m\right)}\equiv\tilde{p}^{\left(m\right)}\left(0,S_{0}=s_{0}\right)=\int_{0}^{K}e^{-rT}\left(K-S_{T}\right)\tilde{\psi}^{\left(m\right)}\left(S_{T};S_{0}=s_{0}\right)dS_{T},\label{eq:ao_jump_eo_appr}
\end{equation}
\begin{equation}
\tilde{e}_{0}^{\left(m\right)}\equiv\tilde{e}^{\left(m\right)}\left(0,S_{0}=s_{0},\tilde{B}^{\left(m\right)}\left(\cdot\right)\right)=\int_{0}^{T}\int_{0}^{\tilde{B}_{t}^{\left(m\right)}}\left(rK-\delta S_{t}\right)e^{-rt}\tilde{\psi}^{\left(m\right)}\left(S_{t};S_{0}=s_{0}\right)dS_{t}dt,\label{eq:ao_jump_eep_appr}
\end{equation}
and
\begin{multline}
\tilde{g}_{0}^{\left(m\right)}\equiv\tilde{g}^{\left(m\right)}\left(0,S_{0}=s_{0},\tilde{B}^{\left(m\right)}\left(\cdot\right)\right)=\int_{0}^{T}\int_{0}^{\tilde{B}_{t-}^{\left(m\right)}}\int_{\tilde{B}_{t}^{\left(m\right)}}^{\infty}e^{-rt}\rho\left(\tilde{P}^{\left(m\right)}\left(t,J_{t}S_{t-}\right)-\left(K-J_{t}S_{t-}\right)\right)\\
\times\tilde{\psi}^{\left(m\right)}\left(S_{t-};S_{0}=s_{0}\right)\tilde{\psi}^{\left(m\right)}\left(J_{t}S_{t-};S_{t-}\right)d\left(J_{t}S_{t-}\right)dS_{t-}dt.\label{eq:ao_jump_jump_appr}
\end{multline}
The approximated exercise boundary up to order $m$, $\tilde{B}_{t}^{\left(m\right)}$,
solves the following recursive nonlinear integral equation:
\begin{equation}
K-\tilde{B}_{t}^{\left(m\right)}=\tilde{p}^{\left(m\right)}\left(t,\tilde{B}_{t}^{\left(m\right)}\right)+\tilde{e}^{\left(m\right)}\left(t,\tilde{B}_{t}^{\left(m\right)},\tilde{B}^{\left(m\right)}\left(\cdot\right)\right)-\tilde{g}^{\left(m\right)}\left(t,\tilde{B}_{t}^{\left(m\right)},\tilde{B}^{\left(m\right)}\left(\cdot\right)\right),\label{eq:ao_jump_bdy_appr}
\end{equation}
where
\begin{equation}
\tilde{p}^{\left(m\right)}\left(t,\tilde{B}_{t}^{\left(m\right)}\right)\equiv\int_{0}^{K}\left(K-S_{T}\right)\tilde{\psi}^{\left(m\right)}\left(S_{T};S_{t}=\tilde{B}_{t}^{\left(m\right)}\right)dS_{T},\label{eq:ao_jump_eo2_appr}
\end{equation}
\begin{equation}
\tilde{e}^{\left(m\right)}\left(t,\tilde{B}_{t}^{\left(m\right)},\tilde{B}^{\left(m\right)}\left(\cdot\right)\right)\equiv\int_{t}^{T}\int_{0}^{\tilde{B}_{s}^{\left(m\right)}}\left(rK-\delta S_{s}\right)e^{-r\left(s-t\right)}\tilde{\psi}^{\left(m\right)}\left(S_{s};S_{t}=\tilde{B}_{t}^{\left(m\right)}\right)dS_{s}ds,\label{eq:ao_jump_eep2_appr}
\end{equation}
\begin{multline}
\tilde{g}^{\left(m\right)}\left(t,\tilde{B}_{t}^{\left(m\right)},\tilde{B}^{\left(m\right)}\left(\cdot\right)\right)\equiv\int_{t}^{T}\int_{0}^{\tilde{B}_{s-}^{\left(m\right)}}\int_{\tilde{B}_{s}^{\left(m\right)}}^{\infty}e^{-r\left(s-t\right)}\rho\left(\tilde{P}^{\left(m\right)}\left(s,J_{s}S_{s-}\right)-\left(K-J_{s}S_{s-}\right)\right)\\
\times\tilde{\psi}^{\left(m\right)}\left(S_{s-};\tilde{B}_{t}^{\left(m\right)}\right)\tilde{\psi}^{\left(m\right)}\left(J_{s}S_{s-};S_{s-}\right)d\left(J_{s}S_{s-}\right)dS_{s-}ds.\label{eq:ao_jump_jump2_appr}
\end{multline}

\section{Numerical Method and Algorithm}\label{section:ao_4}

Following Detemple (2006), we divide the period $\left[0,T\right]$
into $N$ equal subintervals and let $\Delta=T/N$. We then use a
step function to compute the exercise boundary recursively. The algorithm
works as follows: suppose that our step function approximation of the exercise
boundary is $\left\{ \tilde{B}_{n\Delta}^{\left(m,N\right)},n=0,\ldots,N\right\} $.
 The terminal condition tells us that
 \[
 \tilde{B}_{N\Delta}^{\left(m,N\right)}=\min\left\{ K,\frac{r\left(B_{T};\theta\right)}{\delta\left(B_{T};\theta\right)}\times K\right\}.
 \]
Suppose that $\tilde{B}_{l\Delta}^{\left(m,N\right)}$ is known for
all $l>n$, then $\left\{ \tilde{B}_{l\Delta}^{\left(m,N\right)},l=0,\ldots,n\right\} $
can be obtained by discretizing the integral in (\ref{eq:ao_bdy_appr})
for a diffusion model, or (\ref{eq:ao_jump_bdy_appr}) for a jump-diffusion
model using the trapezoidal rule. For example, we obtain the following
equation for diffusion models:
\begin{multline}
K-\tilde{B}_{l\Delta}^{\left(m,N\right)}=\tilde{p}^{\left(m\right)}\left(l\Delta,\tilde{B}_{l\Delta}^{\left(m,N\right)}\right)+\sum_{q=l+1}^{N-l}\tilde{\epsilon}^{\left(m\right)}\left(\left(q-l\right)\Delta,\tilde{B}_{l\Delta}^{\left(m,N\right)},\tilde{B}_{q\Delta}^{\left(m,N\right)}\right)\Delta\\
+\left[\tilde{\epsilon}^{\left(m\right)}\left(0,\tilde{B}_{l\Delta}^{\left(m,N\right)},\tilde{B}_{l\Delta}^{\left(m,N\right)}\right)+\tilde{\epsilon}^{\left(m\right)}\left(\left(N-l\right)\Delta,\tilde{B}_{l\Delta}^{\left(m,N\right)},\tilde{B}_{N\Delta}^{\left(m,N\right)}\right)\right]\frac{\Delta}{2}\label{eq:ao_num}
\end{multline}
where
\begin{equation}
\tilde{\epsilon}^{\left(m\right)}\left(s\Delta,\tilde{B}_{t}^{\left(m\right)},\tilde{B}_{t+s\Delta}^{\left(m\right)}\right)\equiv\int_{0}^{\tilde{B}_{t+s\Delta}^{\left(m\right)}}\left(rK-\delta S_{t+s\Delta}\right)e^{-rs\Delta}\tilde{\psi}^{\left(m\right)}\left(S_{t+s\Delta};S_{t}=\tilde{B}_{t}^{\left(m\right)}\right)dS_{t+s\Delta}\label{eq:ao_num_e}
\end{equation}
And for jump-diffusion models, we have
\begin{multline}
K-\tilde{B}_{l\Delta}^{\left(m,N\right)}=\tilde{p}^{\left(m\right)}\left(l\Delta,\tilde{B}_{l\Delta}^{\left(m,N\right)}\right)+\sum_{q=l+1}^{N-l}\tilde{\epsilon}^{\left(m\right)}\left(\left(q-l\right)\Delta,\tilde{B}_{l\Delta}^{\left(m,N\right)},\tilde{B}_{q\Delta}^{\left(m,N\right)}\right)\Delta\\
+\left[\tilde{\epsilon}^{\left(m\right)}\left(0,\tilde{B}_{l\Delta}^{\left(m,N\right)},\tilde{B}_{l\Delta}^{\left(m,N\right)}\right)+\tilde{\epsilon}^{\left(m\right)}\left(\left(N-l\right)\Delta,\tilde{B}_{l\Delta}^{\left(m,N\right)},\tilde{B}_{N\Delta}^{\left(m,N\right)}\right)\right]\frac{\Delta}{2}\\
+\left[\tilde{\eta}^{\left(m\right)}\left(0,\tilde{B}_{l\Delta}^{\left(m,N\right)},\tilde{B}_{l\Delta}^{\left(m,N\right)}\right)+\tilde{\eta}^{\left(m\right)}\left(\left(N-l\right)\Delta,\tilde{B}_{l\Delta}^{\left(m,N\right)},\tilde{B}_{N\Delta}^{\left(m,N\right)}\right)\right]\frac{\Delta}{2}\\
+\sum_{q=l+1}^{N-l}\tilde{\eta}^{\left(m\right)}\left(\left(q-l\right)\Delta,\tilde{B}_{l\Delta}^{\left(m,N\right)},\tilde{B}_{q\Delta}^{\left(m,N\right)}\right)\Delta\label{eq:ao_jump_num}
\end{multline}
where $\tilde{\epsilon}^{\left(m\right)}$ is defined as in (\ref{eq:ao_num_e}),
and $\tilde{\eta}^{\left(m\right)}$ is defined by
\begin{multline}
\tilde{\eta}^{\left(m\right)}\left(s\Delta,\tilde{B}_{t}^{\left(m\right)},\tilde{B}_{t+s\Delta}^{\left(m\right)}\right)\equiv\\
\int_{0}^{\tilde{B}_{t+s\Delta-}^{\left(m\right)}}\int_{\tilde{B}_{t+s\Delta}^{\left(m\right)}}^{\infty}e^{-rs\Delta}\rho\left(\tilde{P}^{\left(m\right)}\left(s\Delta,J_{t+s\Delta}S_{t+s\Delta s-}\right)-\left(K-J_{t+s\Delta}S_{t+s\Delta s-}\right)\right)\\
\times\tilde{\psi}^{\left(m\right)}\left(S_{t+s\Delta-};\tilde{B}_{t}^{\left(m\right)}\right)\tilde{\psi}^{\left(m\right)}\left(J_{t+s\Delta}S_{t+s\Delta-};S_{t+s\Delta-}\right)d\left(J_{t+s\Delta}S_{t+s\Delta-}\right)dS_{t+s\Delta-}\label{eq:ao_num_g}
\end{multline}

We run the above procedure recursively, and obtain the exercise
boundary $\left\{ \tilde{B}_{n\Delta}^{\left(m,N\right)},n=0,\ldots,N\right\} $.
Finally, the value of the American put option can be computed by substituting the exercise
boundary into (\ref{eq:ao_appr}) for diffusion models or (\ref{eq:ao_jump_appr})
for jump-diffusion models.

\section{Applications}\label{section:ao_5}

In this section, we illustrate how to compute the value of the American put option and the corresponding exercise boundary for diffusion models
and jump-diffusion models. We use these examples to illustrate
the accuracy and speed of the proposed approach. We approximate the transition density using $m=2$ for all the examples in this section.

\subsection{Applications to Diffusion Models}

\subsubsection{Geometric Brownian Motion (GBM) Model}

In the geometric Brownian Motion model, the stock price $S$ follows
\[
dS_{t}=\left(r-\delta\right)S_{t}dt+\sigma S_{t}dW_{t},
\]
 where $r$, $\delta$, and $\sigma$ are constants. To test the efficiency of our recursive algorithm in Section \ref{section:ao_2}, we
compare the results from our approach with those
from four widely used methods: the binomial method by Cox, Ross, and
Rubinstein (1979), the accelerated binomial methods by Breen (1991), the
finite difference method, and the analytical approximation by Geske
and Johnson (1984). We use the results from the binomial method with
10,000 time-steps as a benchmark to measure the accuracy. Following
Huang \emph{et al.} (1996) and Geske and Johnson (1984), we
set $S_{0}=40$, $r=4.88\%$, and $\delta=0$.

Table \ref{gbm} reports the valuations of American options from the six approaches. Columns 1 through 3 represent the values of the parameters, $K$ (strike price), $\sigma$ (volatility), and $T$ (maturity), respectively.
Column 4 gives the numerical results from the binomial method with
10,000 time-steps, and we take this approach as a benchmark. Column
5 includes the results in Table I of Geske and Johnson (1984). Columns
6 through 8 report the results from the binomial method with 150 time-steps, the finite difference method with 200 steps, and the accelerated
binomial method with 150 time-steps. Column 9 shows the results
of the proposed approach with 100 time-steps. The accuracy is measured by the root mean squared error, as shown in the last row. It is clear from this table that the proposed approach achieves the best performance in terms of accuracy compared with the other methods.

\begin{table}[H]
\centering
\caption{Value of American Options Based on Different Numerical Methods}
\renewcommand{\arraystretch}{1.2}
\begin{threeparttable}
\begin{tabular}{|c|c|c|c|c|c|c|c|c|}
\hline
\emph{K} & $\sigma$ & \emph{T }(yr) & Binomial & G\&J & Binomial II & Accelerated & FD & Hermite\tabularnewline
\hline
\hline
35 & 0.2 & 0.0833 & 0.0062 & 0.0062 & 0.0061 & 0.0061 & 0.0278 & 0.0062\tabularnewline
\hline
35 & 0.2 & 0.3333 & 0.2004 & 0.1999 & 0.1995 & 0.1994 & 0.2382 & 0.2004\tabularnewline
\hline
35 & 0.2 & 0.5833 & 0.4328 & 0.4321 & 0.4340 & 0.4331 & 0.4624 & 0.4329\tabularnewline
\hline
40 & 0.2 & 0.0833 & 0.8522 & 0.8528 & 0.8512 & 0.8517 & 0.9874 & 0.8523\tabularnewline
\hline
40 & 0.2 & 0.3333 & 1.5798 & 1.5807 & 1.5783 & 1.5752 & 1.6244 & 1.5800\tabularnewline
\hline
40 & 0.2 & 0.5833 & 1.9904 & 1.9905 & 1.9886 & 1.9856 & 2.0177 & 1.9906\tabularnewline
\hline
45 & 0.2 & 0.0833 & 5.0000 & 4.9985 & 5.0000 & 4.9200 & 5.0052 & 5.0000\tabularnewline
\hline
45 & 0.2 & 0.3333 & 5.0883 & 5.0951 & 5.0886 & 4.9253 & 5.1327 & 5.0886\tabularnewline
\hline
45 & 0.2 & 0.5833 & 5.2670 & 5.2719 & 5.2677 & 5.2844 & 5.2699 & 5.2673\tabularnewline
\hline
35 & 0.3 & 0.0833 & 0.0774 & 0.0744 & 0.0775 & 0.0772 & 0.1216 & 0.0774\tabularnewline
\hline
35 & 0.3 & 0.3333 & 0.6975 & 0.6969 & 0.6993 & 0.6977 & 0.7300 & 0.6976\tabularnewline
\hline
35 & 0.3 & 0.5833 & 1.2198 & 1.2194 & 1.2239 & 1.2218 & 1.2407 & 1.2199\tabularnewline
\hline
40 & 0.3 & 0.0833 & 1.3099 & 1.3100 & 1.3083 & 1.3095 & 1.3860 & 1.3100\tabularnewline
\hline
40 & 0.3 & 0.3333 & 2.4825 & 2.4817 & 2.4799 & 2.4781 & 2.5068 & 2.4828\tabularnewline
\hline
40 & 0.3 & 0.5833 & 3.1696 & 3.1733 & 3.1665 & 3.1622 & 3.1819 & 3.1699\tabularnewline
\hline
45 & 0.3 & 0.0833 & 5.0597 & 5.0599 & 5.0600 & 5.0632 & 5.1016 & 5.0598\tabularnewline
\hline
45 & 0.3 & 0.3333 & 5.7056 & 5.7012 & 5.7065 & 5.6978 & 5.7193 & 5.7059\tabularnewline
\hline
45 & 0.3 & 0.5833 & 6.2436 & 6.2365 & 6.2448 & 6.2395 & 6.2477 & 6.2440\tabularnewline
\hline
35 & 0.4 & 0.0833 & 0.2466 & 0.2466 & 0.2454 & 0.2456 & 0.2949 & 0.2466\tabularnewline
\hline
35 & 0.4 & 0.3333 & 1.3460 & 1.3450 & 1.3505 & 1.3481 & 1.3696 & 1.3461\tabularnewline
\hline
35 & 0.4 & 0.5833 & 2.1549 & 2.1568 & 2.1602 & 2.1569 & 2.1676 & 2.1551\tabularnewline
\hline
40 & 0.4 & 0.0833 & 1.7681 & 1.7679 & 1.7658 & 1.7674 & 1.8198 & 1.7683\tabularnewline
\hline
40 & 0.4 & 0.3333 & 3.3874 & 3.3632 & 3.3835 & 3.3863 & 3.4011 & 3.3877\tabularnewline
\hline
40 & 0.4 & 0.5833 & 4.3526 & 4.3556 & 4.3480 & 4.3426 & 4.3567 & 4.3530\tabularnewline
\hline
45 & 0.4 & 0.0833 & 5.2868 & 5.2855 & 5.2875 & 5.2863 & 5.3289 & 5.2870\tabularnewline
\hline
45 & 0.4 & 0.3333 & 6.5099 & 6.5093 & 6.5103 & 6.5054 & 6.5147 & 6.5101\tabularnewline
\hline
45 & 0.4 & 0.5833 & 7.3830 & 7.3831 & 7.3897 & 7.3785 & 7.3792 & 7.3833\tabularnewline
\hline
\multicolumn{3}{|c|}{RMSE} & 0.0000 & 5.34e-03 & 2.64e-03 & 3.53e-02 & 4.10e-02 & 2.16e-04\tabularnewline
\hline
\end{tabular}

\end{threeparttable}

\label{gbm}
\end{table}

We report the approximated exercise boundary in Figure \ref{boundarygbm1} for several combinations of strike prices and volatility. The parameters values are the same as those in Table \ref{gbm}, and $T$ is 0.5833. This figure shows the marginal effect of strike prices and volatility on the exercise boundary. For example, as the volatility becomes smaller, the responding exercise boundary becomes flatter. The intuition for this result is that when the volatility is small, the return from withholding American options is limited; thus, the American put option will be exercised at a higher boundary instead of a lower one. This result can also be confirmed by checking the partial difference of the exercise boundary with respect to the volatility in (\ref{eq:ao_gbm_bdy}).

The GBM model is one of the limited cases in which we know the true transition density, and to examine the accuracy of our approximated exercise boundary, we plug in the true transition density into the numerical algorithm in Section \ref{section:ao_4}, and compare the results with our approximated exercise boundary. We report this comparison for different strike prices and volatility in Figure \ref{boundarygbm2}.

In Figure \ref{fdvsh}, we compare the results of our proposed approach with the finite difference method for approximating the exercise boundary.

To further investigate the performance of our approach, we report the approximated value of the American put option with respect to strike prices from 10 to 70 in Figure \ref{strike}. In addition, we computed the approximated value of the American put option for different strikes based on various orders of the Hermite polynomial-based approximation of the transition density. More specifically, the first order approximation means $m=1$; the second order means $m=2$; the third order means $m=3$. We use the results from the binomial method as a benchmark for comparison, and report the relative error of the approximation in Figure \ref{strikeerror}.

\subsubsection{Constant Elasticity Volatility (CEV) Model}

The Constant Elasticity Volatility model assumes that the stock price
$S$ follows
\[
dS_{t}=\left(r-\delta\right)S_{t}dt+\sigma S_{t}^{\alpha/2}dW_{t},
\]
where $r$, $\delta$, $\sigma$, and $\alpha$ are constants. Detemple
and Kitapbayev (2018) applied this model to study the pricing of the American
VIX option. Further extension of this model on the valuation of the VIX
option can be found in Goard and Mazur (2013).

We set $K=100$, $r=6/100$, $\delta=r/2$, $\sigma=\sqrt{10}/5$, $S_{0}=40$, and $T=1$. We report in Figure \ref{boundarycev} the approximated exercise boundary of the CEV model for $\alpha=1.9$, and $\alpha=1.7$, respectively. From Figure \ref{boundarycev}, we find that as $\alpha$ decreases, the volatility of the stock price decreases, and thus the optimal exercise boundary becomes higher. This result is the same as that found in the GBM model.

\subsubsection{Nonlinear Mean Reversion (NMR) Model}

The Nonlinear Mean Reversion model assumes that
\[
dS_{t}=\left(\frac{a}{S_{t}}+b+cS_{t}+vS_{t}^{2}\right)dt+\sigma S_{t}^{\gamma}dW_{t},
\]
where $a$, $b$, $c$, $v$, $\sigma$, and $\gamma$ are constants.
This model was discussed in Ait-Sahalia (1996, 1999), and Gallant
and Tauchen (1998) for modeling the interest rates. Eraker and Wang (2012)
proposed a similar model for the VIX option.

In the NMR model, we set $a=500$, $b=5$, $c=0.05$, $v=-0.05$, $\sigma=0.2$, $\gamma=3/2$, $K=20$, $r=5/100$, $\delta=0$, $S_{0}=20$, and $T=0.0833$. We report the approximated exercise boundary shown in Figure \ref{nmr2}.

\subsubsection{Double Mean Reversion (DMR) Model}

The Double Mean Reversion model assumes that
\[
dS_{t}=\beta\left(y_{t}-S_{t}\right)dt+\sigma\sqrt{S_{t}}dW_{t}
\]
\[
dy_{t}=\xi\left(\alpha-y_{t}\right)dt+\kappa\sqrt{y_{t}}dU_{t}
\]
where $W$ and $U$ are two independent Brownian motions, and $\alpha$,
$\beta$, $\xi$, $\kappa$, and $\sigma$ are constants.

Based on the usual square root model, this DMR model includes an additional
stochastic factor for the mean level of the stock price. In this model,
the speed of mean-reversion towards the short-run stochastic mean
level of the stock price is controlled by $\beta$, and the speed of
mean-reversion towards the long-run mean level of the short-run stochastic
mean is controlled by $\xi$. This model was discussed in Amengual (2008),
Mencia and Sentana (2009), and Egloff \emph{et al.} (2010).

In the DMR model, the optimal exercise boundary is a function of time, $t$, and $y$. We set $K=40/100$, $r=4.88/100$, $\delta=0$, $\sigma=0.25$, $\kappa=0.2$, $\beta=2.5$, $\xi=4$, $\alpha=0.25$, $T=0.5$.
In Figure \ref{dmr}, we report the approximated exercise boundary in the DMR model. The boundary is approximated with 20 steps on time, and 100 steps on $y$ for $y$ in $[0,1]$.

\subsection{Applications to Jump-Diffusion Models}

\subsubsection{Merton's Jump-Diffusion Model}

Merton (1976) proposed a jump-diffusion model to incorporate discontinuous
returns, and derived a closed-form vanilla option pricing formula. Merton's jump-diffusion model assumes that:
\[
d\log S_{t}=\left(r-\delta-\lambda j\right)dt+\sigma dW_{t}+\left(J-1\right)dq_{t},
\]
where $dq$ is a Poisson process with rate $\lambda t$, $J$ has a lognormal
distribution with mean $\mu_{J}$ and variance $\sigma_{J}^{2}$,
and $j=E\left[J-1\right]=\exp\left(\mu_{J}+\sigma_{J}^{2}/2\right)-1$.

We set $K=40$, $r=4.88/100$, $\sigma=0.2$, $\mu_{J}=0$, $\sigma_{J}=0.2$, $S_{0}=40$, and $T=0.5$. Additionally, we approximate the exercise boundary for different values of $\lambda$. More specifically, we try $\lambda=1/100, 10/100,$ and $25/100$, and present the results in Figure \ref{merton}. By comparing the exercise boundaries in Figure \ref{merton}, we find that when $\lambda$ is smaller, the exercise boundary is higher. The intuition for this result is that when $\lambda$ is smaller, the jump in the return occurs less frequently, and thus the return becomes less volatile. Similar to the models without jumps in this section, when the stock price or the return is less volatile, the exercise boundary becomes higher.

\subsubsection{Kou's Jump-Diffusion Model}

To incorporate the leptokurtic feature of the return distribution
and ``volatility smile'' phenomenon in option market, Kou (2002)
proposed a double exponential jump-diffusion model. This model assumes
\[
d\log S_{t}=\left(r-\delta\right)dt+\sigma dW_{t}+Jdq_{t}
\]
where $J$ has an asymmetric double exponential distribution with
density:

\[
\upsilon\left(z\right)=p*\eta_{1}e^{-\eta_{1}z}\mathbf{1}_{\left\{ z\geq0\right\} }+q*\eta_{2}e^{-\eta_{2}z}\mathbf{1}_{\left\{ z<0\right\} },
\]
where $\eta_{1}>1$, $\eta_{2}>0$, $p+q=1$, and $0\leq p,q\leq1$.
The mean, variance, and skewness of the jump size in log returns are:
\[
\varphi_{1}=\frac{p}{\eta_{1}}-\frac{q}{\eta_{2}},
\]
\[
\varphi_{2}=pq\left(\frac{1}{\eta_{1}}+\frac{1}{\eta_{2}}\right)^{2}+\frac{p}{\eta_{1}^{2}}+\frac{q}{\eta_{2}^{2}},
\]
\[
\varphi_{3}=\frac{2\left(p^{3}-1\right)\eta_{1}^{3}-2\left(q^{3}-1\right)\eta_{2}^{3}+6pq\eta_{1}\eta_{2}\left(q\eta_{2}-p\eta_{1}\right)}{\left(p\eta_{2}^{2}+q\eta_{1}^{2}+pq\left(\eta_{1}+\eta_{2}\right)^{2}\right)^{\frac{3}{2}}}.
\]
Also, $dq$ is a Poisson process with rate $\lambda t$.

In Figure \ref{kou}, we report the approximated exercise boundary for Kou's jump-diffusion model with $\lambda=1/100, 10/100,$ and $20/100$, respectively. We set $K=40$, $r=4.88/100$, $\delta=0$, $\sigma=0.2$, $p=0.04$, $q=0.96$, $\eta_{1}=3.7$, $\eta_{2}=1.8$, $S_{0}=40$, and $T=0.5$. We approximate the boundary by 50 steps on time. Similarly, we find that the smaller the intensity of the jump is, the higher the exercise boundary becomes.

\section{Conclusion}\label{section:ao_6}

In this study, we develop a new approach to approximate the exercise boundary and the value of the American put option based on Hermite polynomials. We also provide a numerical scheme for implementing the proposed approach. We show theoretically that our approximation will converge to the true exercise boundary and the value of the American put option, and provide evidence for the efficiency of our approach through several numerical examples including diffusion processes and jump-diffusion processes. We only discuss the case of the American put option; however, the value of the American call option can be approximated similarly.

A drawback of our approach is its computational complexity. Although we have a closed-form approximation of the transition density for a given jump-diffusion model, we need to evaluate the integral of the transition density, and the integral usually does not admit a closed-form solution. This incurs a heavy computational burden on the numerical implementation. Other approaches for approximating the transition density, such as finite mixture models, can simplify the integral equation, and thus reduce the computational complexity. We leave this to be explored in future research.

\clearpage

\section*{Abstract (Not appropriate in this style!)}%
   \else \small
   \begin{center}{\bf Abstract\vspace{-.5em}\vspace{\z@}}\end{center}%
   \quotation
   \fi
  }%
 }{%
 }%
\@ifundefined{endabstract}{\def\endabstract
  {\if@twocolumn\else\endquotation\fi}}{}%
\@ifundefined{maketitle}{\def\maketitle#1{}}{}%
\@ifundefined{affiliation}{\def\affiliation#1{}}{}%
\@ifundefined{proof}{\def\proof{\noindent{\bfseries Proof. }}}{}%
\@ifundefined{endproof}{\def\endproof{\mbox{\ \rule{.1in}{.1in}}}}{}%
\@ifundefined{newfield}{\def\newfield#1#2{}}{}%
\@ifundefined{chapter}{\def\chapter#1{\par(Chapter head:)#1\par }%
 \newcount\c@chapter}{}%
\@ifundefined{part}{\def\part#1{\par(Part head:)#1\par }}{}%
\@ifundefined{section}{\def\section#1{\par(Section head:)#1\par }}{}%
\@ifundefined{subsection}{\def\subsection#1%
 {\par(Subsection head:)#1\par }}{}%
\@ifundefined{subsubsection}{\def\subsubsection#1%
 {\par(Subsubsection head:)#1\par }}{}%
\@ifundefined{paragraph}{\def\paragraph#1%
 {\par(Subsubsubsection head:)#1\par }}{}%
\@ifundefined{subparagraph}{\def\subparagraph#1%
 {\par(Subsubsubsubsection head:)#1\par }}{}%
\@ifundefined{therefore}{\def\therefore{}}{}%
\@ifundefined{backepsilon}{\def\backepsilon{}}{}%
\@ifundefined{yen}{\def\yen{\hbox{\rm\rlap=Y}}}{}%
\@ifundefined{registered}{%
   \def\registered{\relax\ifmmode{}\r@gistered
                    \else$\m@th\r@gistered$\fi}%
 \def\r@gistered{^{\ooalign
  {\hfil\raise.07ex\hbox{$\scriptstyle\rm\RIfM@\expandafter\text@\else\expandafter\mbox\fi{R}$}\hfil\crcr
  \mathhexbox20D}}}}{}%
\@ifundefined{Eth}{\def\Eth{}}{}%
\@ifundefined{eth}{\def\eth{}}{}%
\@ifundefined{Thorn}{\def\Thorn{}}{}%
\@ifundefined{thorn}{\def\thorn{}}{}%
\def\TEXTsymbol#1{\mbox{$#1$}}%
\@ifundefined{degree}{\def\degree{{}^{\circ}}}{}%
\newdimen\theight
\def\Column{%
 \vadjust{\setbox\z@=\hbox{\scriptsize\quad\quad tcol}%
  \theight=\ht\z@\advance\theight by \dp\z@\advance\theight by \lineskip
  \kern -\theight \vbox to \theight{%
   \rightline{\rlap{\box\z@}}%
   \vss
   }%
  }%
 }%
\def\qed{%
 \ifhmode\unskip\nobreak\fi\ifmmode\ifinner\else\hskip5\p@\fi\fi
 \hbox{\hskip5\p@\vrule width4\p@ height6\p@ depth1.5\p@\hskip\p@}%
 }%
\def\cents{\hbox{\rm\rlap/c}}%
\def\miss{\hbox{\vrule height2\p@ width 2\p@ depth\z@}}%
\def\vvert{\Vert
\def\tcol#1{{\baselineskip=6\p@ \vcenter{#1}} \Column}  %
\def\dB{\hbox{{}}
\def\mB#1{\hbox{$#1$}
\def\nB#1{\hbox{#1}
\def\note{$^{\dag}}%

\defLaTeX2e{LaTeX2e}
\def\chkcompat{%
   \if@compatibility
   \else
     \usepackage{latexsym}
   \fi
}

\ifx\fmtnameLaTeX2e
  \DeclareOldFontCommand{\rm}{\normalfont\rmfamily}{\mathrm}
  \DeclareOldFontCommand{\sf}{\normalfont\sffamily}{\mathsf}
  \DeclareOldFontCommand{\tt}{\normalfont\ttfamily}{\mathtt}
  \DeclareOldFontCommand{\bf}{\normalfont\bfseries}{\mathbf}
  \DeclareOldFontCommand{\it}{\normalfont\itshape}{\mathit}
  \DeclareOldFontCommand{\sl}{\normalfont\slshape}{\@nomath\sl}
  \DeclareOldFontCommand{\sc}{\normalfont\scshape}{\@nomath\sc}
  \chkcompat
\fi

\def\alpha{\Greekmath 010B }%
\def\beta{\Greekmath 010C }%
\def\gamma{\Greekmath 010D }%
\def\delta{\Greekmath 010E }%
\def\epsilon{\Greekmath 010F }%
\def\zeta{\Greekmath 0110 }%
\def\eta{\Greekmath 0111 }%
\def\theta{\Greekmath 0112 }%
\def\iota{\Greekmath 0113 }%
\def\kappa{\Greekmath 0114 }%
\def\lambda{\Greekmath 0115 }%
\def\mu{\Greekmath 0116 }%
\def\nu{\Greekmath 0117 }%
\def\xi{\Greekmath 0118 }%
\def\pi{\Greekmath 0119 }%
\def\rho{\Greekmath 011A }%
\def\sigma{\Greekmath 011B }%
\def\tau{\Greekmath 011C }%
\def\upsilon{\Greekmath 011D }%
\def\phi{\Greekmath 011E }%
\def\chi{\Greekmath 011F }%
\def\psi{\Greekmath 0120 }%
\def\omega{\Greekmath 0121 }%
\def\varepsilon{\Greekmath 0122 }%
\def\vartheta{\Greekmath 0123 }%
\def\varpi{\Greekmath 0124 }%
\def\varrho{\Greekmath 0125 }%
\def\varsigma{\Greekmath 0126 }%
\def\varphi{\Greekmath 0127 }%

\def\Greekmath 0272 {\Greekmath 0272 }
\def\FindBoldGroup{%
   {\setbox0=\hbox{$\mathbf{x\global\edef\theboldgroup{\the\mathgroup}}$}}%
}

\def\Greekmath#1#2#3#4{%
    \if@compatibility
        \ifnum\mathgroup=\symbold
           \mathchoice{\mbox{\boldmath$\displaystyle\mathchar"#1#2#3#4$}}%
                      {\mbox{\boldmath$\textstyle\mathchar"#1#2#3#4$}}%
                      {\mbox{\boldmath$\scriptstyle\mathchar"#1#2#3#4$}}%
                      {\mbox{\boldmath$\scriptscriptstyle\mathchar"#1#2#3#4$}}%
        \else
           \mathchar"#1#2#3#
        \fi
    \else
        \FindBoldGroup
        \ifnum\mathgroup=\theboldgroup
           \mathchoice{\mbox{\boldmath$\displaystyle\mathchar"#1#2#3#4$}}%
                      {\mbox{\boldmath$\textstyle\mathchar"#1#2#3#4$}}%
                      {\mbox{\boldmath$\scriptstyle\mathchar"#1#2#3#4$}}%
                      {\mbox{\boldmath$\scriptscriptstyle\mathchar"#1#2#3#4$}}%
        \else
           \mathchar"#1#2#3#
        \fi
	  \fi}

\newif\ifGreekBold  \GreekBoldfalse
\let\SAVEPBF=\pbf
\def\pbf{\GreekBoldtrue\SAVEPBF}%

\@ifundefined{theorem}{\newtheorem{theorem}{Theorem}}{}
\@ifundefined{lemma}{\newtheorem{lemma}[theorem]{Lemma}}{}
\@ifundefined{corollary}{\newtheorem{corollary}[theorem]{Corollary}}{}
\@ifundefined{conjecture}{\newtheorem{conjecture}[theorem]{Conjecture}}{}
\@ifundefined{proposition}{\newtheorem{proposition}[theorem]{Proposition}}{}
\@ifundefined{axiom}{\newtheorem{axiom}{Axiom}}{}
\@ifundefined{remark}{\newtheorem{remark}{Remark}}{}
\@ifundefined{example}{\newtheorem{example}{Example}}{}
\@ifundefined{exercise}{\newtheorem{exercise}{Exercise}}{}
\@ifundefined{definition}{\newtheorem{definition}{Definition}}{}

\@ifundefined{mathletters}{%

  \newcounter{equationnumber}
  \def\mathletters{%
     \addtocounter{equation}{1}
     \edef\@currentlabel{\arabic{equation}}%
     \setcounter{equationnumber}{\c@equation}
     \setcounter{equation}{0}%
     \edef\arabic{equation}{\@currentlabel\noexpand\alph{equation}}%
  }
  \def\endmathletters{%
     \setcounter{equation}{\value{equationnumber}}%
  }
}{}

\@ifundefined{BibTeX}{%
    \def\BibTeX{{\rm B\kern-.05em{\sc i\kern-.025em b}\kern-.08em
                 T\kern-.1667em\lower.7ex\hbox{E}\kern-.125emX}}}{}%
\@ifundefined{AmS}%
    {\def\AmS{{\protect\usefont{OMS}{cmsy}{m}{n}%
                A\kern-.1667em\lower.5ex\hbox{M}\kern-.125emS}}}{}%
\@ifundefined{AmSTeX}{\def\AmSTeX{\protect\AmS-\protect\TeX\@}}{}%

\ifx\ds@amstex\relax
   \message{amstex already loaded}\makeatother\endinpu
\else
   \@ifpackageloaded{amstex}%
      {\message{amstex already loaded}\makeatother\endinput}
      {}
   \@ifpackageloaded{amsgen}%
      {\message{amsgen already loaded}\makeatother\endinput}
      {}
\fi
\let\DOTSI\relax
\def\RIfM@{\relax\ifmmode}%
\def\FN@{\futurelet\next}%
\newcount\intno@
\def\iint{\DOTSI\intno@\tw@\FN@\ints@}%
\def\iiint{\DOTSI\intno@\thr@@\FN@\ints@}%
\def\iiiint{\DOTSI\intno@4 \FN@\ints@}%
\def\idotsint{\DOTSI\intno@\z@\FN@\ints@}%
\def\ints@{\findlimits@\ints@@}%
\newif\iflimtoken@
\newif\iflimits@
\def\findlimits@{\limtoken@true\ifx\next\limits\limits@true
 \else\ifx\next\nolimits\limits@false\else
 \limtoken@false\ifx\ilimits@\nolimits\limits@false\else
 \ifinner\limits@false\else\limits@true\fi\fi\fi\fi}%
\def\multint@{\int\ifnum\intno@=\z@\intdots@
 \else\intkern@\fi
 \ifnum\intno@>\tw@\int\intkern@\fi
 \ifnum\intno@>\thr@@\int\intkern@\fi
 \int
\def\multintlimits@{\intop\ifnum\intno@=\z@\intdots@\else\intkern@\fi
 \ifnum\intno@>\tw@\intop\intkern@\fi
 \ifnum\intno@>\thr@@\intop\intkern@\fi\intop}%
\def\intic@{%
    \mathchoice{\hskip.5em}{\hskip.4em}{\hskip.4em}{\hskip.4em}}%
\def\negintic@{\mathchoice
 {\hskip-.5em}{\hskip-.4em}{\hskip-.4em}{\hskip-.4em}}%
\def\ints@@{\iflimtoken@
 \def\ints@@@{\iflimits@\negintic@
   \mathop{\intic@\multintlimits@}\limits
  \else\multint@\nolimits\fi
  \eat@
 \else
 \def\ints@@@{\iflimits@\negintic@
  \mathop{\intic@\multintlimits@}\limits\else
  \multint@\nolimits\fi}\fi\ints@@@}%
\def\intkern@{\mathchoice{\!\!\!}{\!\!}{\!\!}{\!\!}}%
\def\plaincdots@{\mathinner{\cdotp\cdotp\cdotp}}%
\def\intdots@{\mathchoice{\plaincdots@}%
 {{\cdotp}\mkern1.5mu{\cdotp}\mkern1.5mu{\cdotp}}%
 {{\cdotp}\mkern1mu{\cdotp}\mkern1mu{\cdotp}}%
 {{\cdotp}\mkern1mu{\cdotp}\mkern1mu{\cdotp}}}%
\def\RIfM@{\relax\protect\ifmmode}
\def\RIfM@\expandafter\text@\else\expandafter\mbox\fi{\RIfM@\expandafter\RIfM@\expandafter\text@\else\expandafter\mbox\fi@\else\expandafter\mbox\fi}
\let\nfss@text\RIfM@\expandafter\text@\else\expandafter\mbox\fi
\def\RIfM@\expandafter\text@\else\expandafter\mbox\fi@#1{\mathchoice
   {\textdef@\displaystyle\f@size{#1}}%
   {\textdef@\textstyle\tf@size{\firstchoice@false #1}}%
   {\textdef@\textstyle\sf@size{\firstchoice@false #1}}%
   {\textdef@\textstyle \ssf@size{\firstchoice@false #1}}%
   \glb@settings}

\def\textdef@#1#2#3{\hbox{{%
                    \everymath{#1}%
                    \let\f@size#2\selectfont
                    #3}}}
\newif\iffirstchoice@
\firstchoice@true
\def\Let@{\relax\iffalse{\fi\let\\=\cr\iffalse}\fi}%
\def\vspace@{\def\vspace##1{\crcr\noalign{\vskip##1\relax}}}%
\def\multilimits@{\bgroup\vspace@\Let@
 \baselineskip\fontdimen10 \scriptfont\tw@
 \advance\baselineskip\fontdimen12 \scriptfont\tw@
 \lineskip\thr@@\fontdimen8 \scriptfont\thr@@
 \lineskiplimit\lineskip
 \vbox\bgroup\ialign\bgroup\hfil$\m@th\scriptstyle{##}$\hfil\crcr}%
\def\Sb{_\multilimits@}%
\def\endSb{\crcr\egroup\egroup\egroup}%
\def\Sp{^\multilimits@}%
\let\endSp\endSb
\newdimen\ex@
\ex@.2326ex
\def\rightarrowfill@#1{$#1\m@th\mathord-\mkern-6mu\cleaders
 \hbox{$#1\mkern-2mu\mathord-\mkern-2mu$}\hfill
 \mkern-6mu\mathord\rightarrow$}%
\def\leftarrowfill@#1{$#1\m@th\mathord\leftarrow\mkern-6mu\cleaders
 \hbox{$#1\mkern-2mu\mathord-\mkern-2mu$}\hfill\mkern-6mu\mathord-$}%
\def\leftrightarrowfill@#1{$#1\m@th\mathord\leftarrow
\mkern-6mu\cleaders
 \hbox{$#1\mkern-2mu\mathord-\mkern-2mu$}\hfill
 \mkern-6mu\mathord\rightarrow$}%
\def\overrightarrow{\mathpalette\overrightarrow@}%
\def\overrightarrow@#1#2{\vbox{\ialign{##\crcr\rightarrowfill@#1\crcr
 \noalign{\kern-\ex@\nointerlineskip}$\m@th\hfil#1#2\hfil$\crcr}}}%
\let\overarrow\overrightarrow
\def\overleftarrow{\mathpalette\overleftarrow@}%
\def\overleftarrow@#1#2{\vbox{\ialign{##\crcr\leftarrowfill@#1\crcr
 \noalign{\kern-\ex@\nointerlineskip}$\m@th\hfil#1#2\hfil$\crcr}}}%
\def\overleftrightarrow{\mathpalette\overleftrightarrow@}%
\def\overleftrightarrow@#1#2{\vbox{\ialign{##\crcr
   \leftrightarrowfill@#1\crcr
 \noalign{\kern-\ex@\nointerlineskip}$\m@th\hfil#1#2\hfil$\crcr}}}%
\def\underrightarrow{\mathpalette\underrightarrow@}%
\def\underrightarrow@#1#2{\vtop{\ialign{##\crcr$\m@th\hfil#1#2\hfil
  $\crcr\noalign{\nointerlineskip}\rightarrowfill@#1\crcr}}}%
\let\underarrow\underrightarrow
\def\underleftarrow{\mathpalette\underleftarrow@}%
\def\underleftarrow@#1#2{\vtop{\ialign{##\crcr$\m@th\hfil#1#2\hfil
  $\crcr\noalign{\nointerlineskip}\leftarrowfill@#1\crcr}}}%
\def\underleftrightarrow{\mathpalette\underleftrightarrow@}%
\def\underleftrightarrow@#1#2{\vtop{\ialign{##\crcr$\m@th
  \hfil#1#2\hfil$\crcr
 \noalign{\nointerlineskip}\leftrightarrowfill@#1\crcr}}}%

\def\qopnamewl@#1{\mathop{\operator@font#1}\nlimits@}
\let\nlimits@\displaylimits
\def\setboxz@h{\setbox\z@\hbox}

\def\varlim@#1#2{\mathop{\vtop{\ialign{##\crcr
 \hfil$#1\m@th\operator@font lim$\hfil\crcr
 \noalign{\nointerlineskip}#2#1\crcr
 \noalign{\nointerlineskip\kern-\ex@}\crcr}}}}

 \def\rightarrowfill@#1{\m@th\setboxz@h{$#1-$}\ht\z@\z@
  $#1\copy\z@\mkern-6mu\cleaders
  \hbox{$#1\mkern-2mu\box\z@\mkern-2mu$}\hfill
  \mkern-6mu\mathord\rightarrow$}
\def\leftarrowfill@#1{\m@th\setboxz@h{$#1-$}\ht\z@\z@
  $#1\mathord\leftarrow\mkern-6mu\cleaders
  \hbox{$#1\mkern-2mu\copy\z@\mkern-2mu$}\hfill
  \mkern-6mu\box\z@$}

\def\qopnamewl@{proj\,lim}{\qopnamewl@{proj\,lim}}
\def\qopnamewl@{inj\,lim}{\qopnamewl@{inj\,lim}}
\def\mathpalette\varlim@\rightarrowfill@{\mathpalette\varlim@\rightarrowfill@}
\def\mathpalette\varlim@\leftarrowfill@{\mathpalette\varlim@\leftarrowfill@}
\def\mathpalette\varliminf@{}{\mathpalette\mathpalette\varliminf@{}@{}}
\def\mathpalette\varliminf@{}@#1{\mathop{\underline{\vrule\@depth.2\ex@\@width\z@
   \hbox{$#1\m@th\operator@font lim$}}}}
\def\mathpalette\varlimsup@{}{\mathpalette\mathpalette\varlimsup@{}@{}}
\def\mathpalette\varlimsup@{}@#1{\mathop{\overline
  {\hbox{$#1\m@th\operator@font lim$}}}}

\def\tfrac#1#2{{\textstyle {#1 \over #2}}}%
\def\dfrac#1#2{{\displaystyle {#1 \over #2}}}%
\def\binom#1#2{{#1 \choose #2}}%
\def\tbinom#1#2{{\textstyle {#1 \choose #2}}}%
\def\dbinom#1#2{{\displaystyle {#1 \choose #2}}}%
\def\QATOP#1#2{{#1 \atop #2}}%
\def\QTATOP#1#2{{\textstyle {#1 \atop #2}}}%
\def\QDATOP#1#2{{\displaystyle {#1 \atop #2}}}%
\def\QABOVE#1#2#3{{#2 \above#1 #3}}%
\def\QTABOVE#1#2#3{{\textstyle {#2 \above#1 #3}}}%
\def\QDABOVE#1#2#3{{\displaystyle {#2 \above#1 #3}}}%
\def\QOVERD#1#2#3#4{{#3 \overwithdelims#1#2 #4}}%
\def\QTOVERD#1#2#3#4{{\textstyle {#3 \overwithdelims#1#2 #4}}}%
\def\QDOVERD#1#2#3#4{{\displaystyle {#3 \overwithdelims#1#2 #4}}}%
\def\QATOPD#1#2#3#4{{#3 \atopwithdelims#1#2 #4}}%
\def\QTATOPD#1#2#3#4{{\textstyle {#3 \atopwithdelims#1#2 #4}}}%
\def\QDATOPD#1#2#3#4{{\displaystyle {#3 \atopwithdelims#1#2 #4}}}%
\def\QABOVED#1#2#3#4#5{{#4 \abovewithdelims#1#2#3 #5}}%
\def\QTABOVED#1#2#3#4#5{{\textstyle
   {#4 \abovewithdelims#1#2#3 #5}}}%
\def\QDABOVED#1#2#3#4#5{{\displaystyle
   {#4 \abovewithdelims#1#2#3 #5}}}%

\def\tint{\mathop{\textstyle \int}}%
\def\tiint{\mathop{\textstyle \iint }}%
\def\tiiint{\mathop{\textstyle \iiint }}%
\def\tiiiint{\mathop{\textstyle \iiiint }}%
\def\tidotsint{\mathop{\textstyle \idotsint }}%
\def\toint{\mathop{\textstyle \oint}}%
\def\tsum{\mathop{\textstyle \sum }}%
\def\tprod{\mathop{\textstyle \prod }}%
\def\tbigcap{\mathop{\textstyle \bigcap }}%
\def\tbigwedge{\mathop{\textstyle \bigwedge }}%
\def\tbigoplus{\mathop{\textstyle \bigoplus }}%
\def\tbigodot{\mathop{\textstyle \bigodot }}%
\def\tbigsqcup{\mathop{\textstyle \bigsqcup }}%
\def\tcoprod{\mathop{\textstyle \coprod }}%
\def\tbigcup{\mathop{\textstyle \bigcup }}%
\def\tbigvee{\mathop{\textstyle \bigvee }}%
\def\tbigotimes{\mathop{\textstyle \bigotimes }}%
\def\tbiguplus{\mathop{\textstyle \biguplus }}%

\def\dint{\mathop{\displaystyle \int}}%
\def\diint{\mathop{\displaystyle \iint }}%
\def\diiint{\mathop{\displaystyle \iiint }}%
\def\diiiint{\mathop{\displaystyle \iiiint }}%
\def\didotsint{\mathop{\displaystyle \idotsint }}%
\def\doint{\mathop{\displaystyle \oint}}%
\def\dsum{\mathop{\displaystyle \sum }}%
\def\dprod{\mathop{\displaystyle \prod }}%
\def\dbigcap{\mathop{\displaystyle \bigcap }}%
\def\dbigwedge{\mathop{\displaystyle \bigwedge }}%
\def\dbigoplus{\mathop{\displaystyle \bigoplus }}%
\def\dbigodot{\mathop{\displaystyle \bigodot }}%
\def\dbigsqcup{\mathop{\displaystyle \bigsqcup }}%
\def\dcoprod{\mathop{\displaystyle \coprod }}%
\def\dbigcup{\mathop{\displaystyle \bigcup }}%
\def\dbigvee{\mathop{\displaystyle \bigvee }}%
\def\dbigotimes{\mathop{\displaystyle \bigotimes }}%
\def\dbiguplus{\mathop{\displaystyle \biguplus }}%
\def\stackunder#1#2{\mathrel{\mathop{#2}\limits_{#1}}}%
\begingroup \catcode `|=0 \catcode `[= 1
\catcode`]=2 \catcode `\{=12 \catcode `\}=12
\catcode`\\=12
|gdef|@alignverbatim#1\end{align}[#1|end[align]]
|gdef|@salignverbatim#1\end{align*}[#1|end[align*]]

|gdef|@alignatverbatim#1\end{alignat}[#1|end[alignat]]
|gdef|@salignatverbatim#1\end{alignat*}[#1|end[alignat*]]

|gdef|@xalignatverbatim#1\end{xalignat}[#1|end[xalignat]]
|gdef|@sxalignatverbatim#1\end{xalignat*}[#1|end[xalignat*]]

|gdef|@gatherverbatim#1\end{gather}[#1|end[gather]]
|gdef|@sgatherverbatim#1\end{gather*}[#1|end[gather*]]

|gdef|@gatherverbatim#1\end{gather}[#1|end[gather]]
|gdef|@sgatherverbatim#1\end{gather*}[#1|end[gather*]]

|gdef|@multilineverbatim#1\end{multiline}[#1|end[multiline]]
|gdef|@smultilineverbatim#1\end{multiline*}[#1|end[multiline*]]

|gdef|@arraxverbatim#1\end{arrax}[#1|end[arrax]]
|gdef|@sarraxverbatim#1\end{arrax*}[#1|end[arrax*]]

|gdef|@tabulaxverbatim#1\end{tabulax}[#1|end[tabulax]]
|gdef|@stabulaxverbatim#1\end{tabulax*}[#1|end[tabulax*]]

|endgroup

\def\align{\@verbatim \frenchspacing\@vobeyspaces \@alignverbatim
You are using the "align" environment in a style in which it is not defined.}
\let\endalign=\endtrivlist

\@namedef{align*}{\@verbatim\@salignverbatim
You are using the "align*" environment in a style in which it is not defined.}
\expandafter\let\csname endalign*\endcsname =\endtrivlist

\def\alignat{\@verbatim \frenchspacing\@vobeyspaces \@alignatverbatim
You are using the "alignat" environment in a style in which it is not defined.}
\let\endalignat=\endtrivlist

\@namedef{alignat*}{\@verbatim\@salignatverbatim
You are using the "alignat*" environment in a style in which it is not defined.}
\expandafter\let\csname endalignat*\endcsname =\endtrivlist

\def\xalignat{\@verbatim \frenchspacing\@vobeyspaces \@xalignatverbatim
You are using the "xalignat" environment in a style in which it is not defined.}
\let\endxalignat=\endtrivlist

\@namedef{xalignat*}{\@verbatim\@sxalignatverbatim
You are using the "xalignat*" environment in a style in which it is not defined.}
\expandafter\let\csname endxalignat*\endcsname =\endtrivlist

\def\gather{\@verbatim \frenchspacing\@vobeyspaces \@gatherverbatim
You are using the "gather" environment in a style in which it is not defined.}
\let\endgather=\endtrivlist

\@namedef{gather*}{\@verbatim\@sgatherverbatim
You are using the "gather*" environment in a style in which it is not defined.}
\expandafter\let\csname endgather*\endcsname =\endtrivlist

\def\multiline{\@verbatim \frenchspacing\@vobeyspaces \@multilineverbatim
You are using the "multiline" environment in a style in which it is not defined.}
\let\endmultiline=\endtrivlist

\@namedef{multiline*}{\@verbatim\@smultilineverbatim
You are using the "multiline*" environment in a style in which it is not defined.}
\expandafter\let\csname endmultiline*\endcsname =\endtrivlist

\def\arrax{\@verbatim \frenchspacing\@vobeyspaces \@arraxverbatim
You are using a type of "array" construct that is only allowed in AmS-LaTeX.}
\let\endarrax=\endtrivlist

\def\tabulax{\@verbatim \frenchspacing\@vobeyspaces \@tabulaxverbatim
You are using a type of "tabular" construct that is only allowed in AmS-LaTeX.}
\let\endtabulax=\endtrivlist

\@namedef{arrax*}{\@verbatim\@sarraxverbatim
You are using a type of "array*" construct that is only allowed in AmS-LaTeX.}
\expandafter\let\csname endarrax*\endcsname =\endtrivlist

\@namedef{tabulax*}{\@verbatim\@stabulaxverbatim
You are using a type of "tabular*" construct that is only allowed in AmS-LaTeX.}
\expandafter\let\csname endtabulax*\endcsname =\endtrivlist

\def\@@eqncr{\let\@tempa\relax
    \ifcase\@eqcnt \def\@tempa{& & &}\or \def\@tempa{& &}%
      \else \def\@tempa{&}\fi
     \@tempa
     \if@eqnsw
        \iftag@
           \@taggnum
        \else
           \@eqnnum\stepcounter{equation}%
        \fi
     \fi
     \global\@ifnextchar*{\@tagstar}{\@tag}@false
     \global\@eqnswtrue
     \global\@eqcnt\z@\cr}

 \def\endequation{%
     \ifmmode\ifinner
      \iftag@
        \addtocounter{equation}{-1}
        $\hfil
           \displaywidth\linewidth\@taggnum\egroup \endtrivlist
        \global\@ifnextchar*{\@tagstar}{\@tag}@false
        \global\@ignoretrue
      \else
        $\hfil
           \displaywidth\linewidth\@eqnnum\egroup \endtrivlist
        \global\@ifnextchar*{\@tagstar}{\@tag}@false
        \global\@ignoretrue
      \fi
     \else
      \iftag@
        \addtocounter{equation}{-1}
        \eqno \hbox{\@taggnum}
        \global\@ifnextchar*{\@tagstar}{\@tag}@false%
        $$\global\@ignoretrue
      \else
        \eqno \hbox{\@eqnnum
        $$\global\@ignoretrue
      \fi
     \fi\fi
 }

 \newif\iftag@ \@ifnextchar*{\@tagstar}{\@tag}@false

 \def\@ifnextchar*{\@tagstar}{\@tag}{\@ifnextchar*{\@tagstar}{\@tag}}
 \def\@tag#1{%
     \global\@ifnextchar*{\@tagstar}{\@tag}@true
     \global\def\@taggnum{(#1)}}
 \def\@tagstar*#1{%
     \global\@ifnextchar*{\@tagstar}{\@tag}@true
     \global\def\@taggnum{#1
}

\makeatother
\endinput