\documentclass[twoside]{article}
\pagestyle{myheadings} \setcounter{page}{125}
\markboth{\hfil Initial/final-value  problem for Beechem-Haas equations\hfil}%
{\hfil Valentino Anthony Simpao \hfil}
\begin{document}
\title{\vspace{-1in}\parbox{\linewidth}{\footnotesize\noindent
{\sc 15th Annual Conference on Applied Mathematics, Univ. of Central Oklahoma},
\newline
Electronic Journal of Differential Equations, Conference~02, 1999,
 pp. 125--131. 
\newline
ISSN: 1072-6691. URL: http://ejde.math.swt.edu or http://ejde.math.unt.edu
\newline ftp  ejde.math.swt.edu (login: ftp)}
 \vspace{\bigskipamount} \\
  Operational WKB solution to the initial/final-value  problem for 
  Beechem-Haas equations 
\thanks{ {\em 1991 Mathematics Subject Classifications:} 92C45, 34E20, 58J35.
\hfil\break\indent {\em Key words and phrases:} operational 
methods, WKB, Analytic solution. \hfil\break\indent 
\copyright 2000 Southwest Texas 
State University  and University of North Texas. 
\hfil\break\indent Published January 21, 2000.} } 

\date{}
\author{Valentino Anthony Simpao}
\maketitle

\begin{abstract}
Operational and Wantzel-Kramers-Brillouin (WKB) methods are utilized 
to obtain an analytical approximate
solution of the Beechem-Haas equation for  the initial/final-value  problem
 with system parameters  of radial profile. The Beechem-Haas equation herein
is expanded from its traditional form, which determines its solutions in a
spherical coordinate system as functions of the radial  and time coordinates
only, to one which includes the angular coordinates as well, thus providing
a natural modeling framework for polarization phenomena.
\end{abstract}


\section{Introduction}

Consider the Beechem-Haas equation and its  initial/boundary conditions,
as given in \cite{b1},
\begin{eqnarray} \label{I}
&\partial_t\overline{N}(r,t) = \bigg(\sum_{j=1}^n \frac{\alpha_j}{\tau_j}
(1+(R_0/r)^6)\bigg) \overline{N}(r,t)
+ \frac 1{N_0(r)} \partial_r[N_0(r)D(r)\partial_r\overline{N}(r,t)]\,,&
\nonumber \\
&\overline{N}(r,t)\bigg|_{t=0} = 1\,,\quad  \overline{N}_{t=0}(r,t) = 1\,, &\\
&\partial_r\overline{N}(r,t)\bigg|_{r=r_{\mbox{\scriptsize min}}} = 0\,, \quad
 \partial_r\overline{N}(r,t)\bigg|_{r=r_{\mbox{\scriptsize max}}} = 0\,, &\nonumber
\end{eqnarray}
where $0<r_{\mbox{\scriptsize min}}\leq r \leq r_{\mbox{\scriptsize max}}$,
$t\geq 0$, and $\alpha_j,\tau_j,R_0$
are positive real numbers.
Unless otherwise directed, all variables, parameters, and domains are real.

In the present work,  the full 3-dimensional form of the spatial derivatives
in spherical coordinates shall be applied to (\ref{I}), resulting in an
enhanced  form of the Beechem Haas equation,
which will allow for the inclusion of angular components to the solution in
a natural way. Particularly, this allows (\ref{I}) to be extended to  include the
nontrivial effect of angular coupling to the radial component of the density
function (i.e., polarization phenomena). Finally, the solution of the
initial/final-value problem for  the enhanced version of the Beechem-Haas
equation will be addressed herein. This is basically the \lq\lq on/off"
control problem in the time-domain for the system. As a result of solving
this problem, one then determines the system function necessary to produce a
required output state, given a prescribed input state.
Full technical details of the terms in (\ref{I}), other than those supplied 
herein, are unnecessary for the present work. Therefore,
the interested reader is encouraged to examine the original source \cite{b1}.
The primary purpose of the present work  is to simply supply a derivation of
the results named  above; applications and deeper studies shall appear
elsewhere.


\section{Main result}

With these notions in place,  the enhanced Beechem-Haas equation equation 
with attendant generalized auxiliary conditions becomes
\begin{eqnarray}
&\partial_t\overline{N}_E(r,\theta,\phi,t) = \bigg(\sum_{j=1}^n \frac{\alpha_j}{\tau_j}
(1+(R_0/r)^6)\bigg) \overline{N}_E(r,\theta,\phi,t) &\nonumber\\
&\hspace{2.5cm}+ \frac {4\pi r^2}{N_0(r)} \nabla\cdot [\frac{N_0(r)}{4\pi r^2}D(r)
\nabla\overline{N}_E(r,\theta,\phi,t)]\,,& \label{Ie}\\
&\overline{N}_E(r,\theta,\phi,t)\bigg|_{t=0} = \overline{N}_{E_{t=0}}(r)
\,,\quad
\overline{N}_E(r,\theta,\phi,t)\bigg|_{t=T} = \overline{N}_{E_{t=T}}(r)&
\nonumber\\
&\partial_r[\overline{N}_E(r,\theta,\phi,t)
\bigg|_{r=r_{\mbox{\scriptsize min}}} = 0\,, \quad
 \partial_r[\overline{N}_E(r,\theta,\phi,t)
\bigg|_{r=r_{\mbox{\scriptsize max}}} = 0\,, & \nonumber
\end{eqnarray}
where $0<r_{\mbox{\scriptsize min}}\leq r \leq r_{\mbox{\scriptsize max}}$,
$0\leq\theta\leq\pi$, $0\leq\phi\leq 2\pi$, $t\geq 0$, and
$\alpha_j,\tau_j,R_0$ are positive real numbers.
 Expanding the above expression,  
\begin{eqnarray} 
\lefteqn{\frac 1{D(r)} \partial_t\overline{N}_E(r,\theta,\phi,t) }\nonumber\\
 &=&
\frac 1{D(r)} \bigg(\sum_{j=1}^n \frac{\alpha_j}{\tau_j}
(1+(R_0/r)^6)\bigg) \overline{N}_E(r,\theta,\phi,t) \label{1} \\
&&+\nabla[\ln\big(\frac{N_0(r)D(r)}{4\pi r^2}\big)] 
\cdot\nabla\overline{N}_E(r,\theta,\phi,t)
+\Delta\overline{N}_E(r,\theta,\phi,t)\,. \nonumber
\end{eqnarray}
Recall that the finite Fourier transform and its inverse \cite{k1} are defined
as
$$\displaylines{
\widehat{f}(n)=F_T[f(t)]=\int_0^Tf(t)e^{-in\omega t}\,dt\,, \cr
f(t)=F_T^{-1}[\widehat{f}(n)]=\frac1T\big(\widehat{f}(0)
+\sum_{n=1}^\infty[\widehat{f}(n)e^{in\omega t}+\widehat{f}(-n)e^{-in\omega t}]\big)\cr
}$$
Applying the finite Fourier transform on the time domain of (\ref{1}),
we obtain 
\begin{eqnarray} 
\lefteqn{ \frac 1{D(r)}[\overline{N}_{E_{t=T}}(r)-\overline{N}_{E_{t=0}}(r)]
}\nonumber\\
&=& \frac 1{D(r)} \bigg(\sum_{j=1}^n \frac{\alpha_j}{\tau_j}
(1+(R_0/r)^6)-in\omega\bigg) \widehat{N}_E(r,\theta,\phi,n) \label{3}\\
&&+\nabla[\ln\big(\frac{N_0(r)D(r)}{4\pi r^2}\big)] 
\cdot\nabla\widehat{N}_E(r,\theta,\phi,n)
+\Delta\widehat{N}_E(r,\theta,\phi,n)\,. \nonumber
\end{eqnarray}
Reducing this to canonical form,
\begin{eqnarray} 
\lefteqn{\frac 1{D(r)}[\overline{N}_{E_{t=T}}(r)-\overline{N}_{E_{t=0}}(r)]
\sqrt{\frac{4\pi r^2}{N_0(r)D(r)}}  }\nonumber \\
&=& \Delta \widetilde{N}_E(r,\theta,\phi,n)
+\bigg(\frac 1{D(r)} \big(\sum_{j=1}^n \frac{\alpha_j}{\tau_j}
(1+(R_0/r)^6)-in\omega\big) \label{4} \\
&&-\frac12 \Delta[\ln\big(\frac{N_0(r)D(r)}{4\pi r^2}\big)] 
-\frac14\bigg|\nabla[\ln\big(\frac{N_0(r)D(r)}{4\pi r^2}\big)]\bigg|^2
\bigg)\widetilde{N}_E(r,\theta,\phi,n)\,,  \nonumber
\end{eqnarray}
where
$$\displaylines{
\widetilde{N}_E(r,\theta,\phi,n) = \widehat{N}_E(r,\theta,\phi,n)
\sqrt{\frac{N_0(r)D(r)}{4\pi r^2}}\,, \cr
\partial_r[\widetilde{N}_E(r,\theta,\phi,n)\sqrt{\frac{4\pi r^2}{N_0(r)D(r)}}]
\Big|_{r=r_{\mbox{\scriptsize min}}} = 0\,, \cr
\partial_r[\widetilde{N}_E(r,\theta,\phi,n)\sqrt{\frac{4\pi r^2}{N_0(r)D(r)}}]
\Big|_{r=r_{\mbox{\scriptsize max}}} = 0\,.\cr
}$$
Now we consider the homogeneous equation corresponding to (\ref{4}),
\begin{eqnarray} \label{5}
0&=& \Delta \widetilde{N}_{E_h}(r,\theta,\phi,n)
+\bigg(\frac 1{D(r)} \big(\sum_{j=1}^n \frac{\alpha_j}{\tau_j}
(1+(R_0/r)^6)-in\omega\big) \\
&&-\frac12 \Delta[\ln\big(\frac{N_0(r)D(r)}{4\pi r^2}\big)] 
-\frac14\bigg|\nabla[\ln\big(\frac{N_0(r)D(r)}{4\pi r^2}\big)]\bigg|^2
\bigg)\widetilde{N}_{E_h}(r,\theta,\phi,n)\,.\nonumber
\end{eqnarray}
Let $\widetilde{N}_{E_h}(r,\theta,\phi,n) = R_c(r,n)\Theta(\theta)\Phi(\phi)$.
Then via separation of variables,
\begin{eqnarray}
\lefteqn{ \frac 1{r^2R_c(r,n)} D_r(r^2D_rR_c(r,n)) }\nonumber\\
\lefteqn{ + \frac 1{r^2\sin(\theta)\Theta(\theta)} D_\theta(\sin(\theta)
D_\theta\Theta(\theta))
+ \frac 1{r^2\sin^2(\theta)\Phi(\theta)} D_\theta^2\Phi(\theta) } \label{7}\\
&=& \frac 1{D(r)} \big(\sum_{j=1}^n \frac{\alpha_j}{\tau_j}
(1+(R_0/r)^6)-in\omega\big) \nonumber \\
&&-\frac12 \Delta[\ln\big(\frac{N_0(r)D(r)}{4\pi r^2}\big)] 
-\frac14\bigg|\nabla[\ln\big(\frac{N_0(r)D(r)}{4\pi r^2}\big)]\bigg|^2
\,.\nonumber
\end{eqnarray}
Hence
\begin{eqnarray} \label{8}
 D_r^2(rR_c(r,n)) 
-\bigg(\frac{\lambda}{r^2}+
\frac 1{D(r)} \big(\sum_{j=1}^n \frac{\alpha_j}{\tau_j}
(1+(R_0/r)^6)-in\omega\big)&& \\
-\frac12 \Delta[\ln\big(\frac{N_0(r)D(r)}{4\pi r^2}\big)] 
-\frac14\bigg|\nabla[\ln\big(\frac{N_0(r)D(r)}{4\pi r^2}\big)]\bigg|^2
\bigg) rR_c(r,n) &=& 0\,,\nonumber
\end{eqnarray}
$$\displaylines{
\frac 1{\Phi(\phi)} D_\phi^2\Phi(\phi) = -m^2\,, \cr 
\frac 1{\sin \theta} D_\theta(\sin(\theta)D_\theta\Theta(\theta))
-\frac{m^2}{\sin^2\theta} \Theta(\theta)+\lambda\Theta(\theta) = 0\,,
}$$
where herein $\lambda, m$  are arbitrary real parameters and $n$ is a
positive integer.

Regarding the solutions of the above angular component equations, they 
are most generally expressed as
\begin{eqnarray} \label{9}
\Theta(\theta)\Phi(\phi)&=& [c_{\phi_1}e^{im\phi}+c_{\phi_2}e^{-im\phi}]
\\
&&\times\bigg(c_{\theta_1}\frac 1{\Gamma(1-m)}
(\frac{1+\cos\theta}{1-\cos\theta})^{m/2}
+c_{\theta_2}\frac{\pi\cot(m\pi)}{2\Gamma(1-m)}
(\frac{1+\cos\theta}{1-\cos\theta})^{m/2} \nonumber \\
&&\qquad-c_{\theta_2}\frac{\pi\Gamma(m+\lambda+1)}{2\sin(m\pi)\Gamma(\lambda-m+1)
\Gamma(1+m)}(\frac{1+\cos\theta}{1-\cos\theta})^{-m/2}\bigg) \nonumber\\
&&\times F_{2,1}(-\lambda,\lambda;1-m;\frac{1-\cos\theta}{2})\,,  \nonumber 
\end{eqnarray}
where $c_{\phi_1}$, $c_{\phi_2}$, $c_{\theta_1}$, $c_{\theta_2}$ are
arbitrary  constants. From (\ref{8}) it is clear that the radial ordinary
differential equation couples to the angular ordinary differential
equation only when
$\lambda$ is non-zero. When $\lambda\neq 0$, (\ref{9}) can take different
forms, depending on the values of $\lambda$ and $m$, which may be
arbitrary real numbers.  Indeed for $\lambda$  and $m$ integers, the angular
dependence may be specified in terms of a finite superposition of spherical
harmonics $Y_\lambda^m(\theta,\phi)$, rather than via the Gauss
hyper-geometric functions $F_{2,1}$ used in in (\ref{9}).

However, even if $m$ is an integer and  $\lambda$  is not so restricted,
the angular solutions will still [in general, by necessity] require the
Gauss hyper-geometric formulation \cite{g1}.
 With an eye toward future results,  the phenomena  modeled by  (\ref{8})
may or may not require $\lambda$ and $m$ to be integer, depending upon the
possible presence of more general non-homogeneous auxiliary conditions.
For economy of notation result (\ref{9}) may be denoted as
$\Lambda(\theta,\lambda:\phi,m)$ in general.
However, since the present  homogeneous auxiliary conditions on the
solutions are just the standard ones for boundedness and single-valuedness,
it is here sufficient to particularly designate the angular components as
spherical harmonics $Y_\lambda^m(\theta,\phi)$  with $\lambda$  and $m$
integers. Hence
\begin{equation}  \label{10}
 \widetilde{N}_{Eh}(r,\theta,\phi,n)=R_c(r,n)\Theta(\theta)\Phi(\phi)
=R_c(r,n)\Lambda(\theta,\lambda:\phi,m) 
=R_c(r,n)Y_\lambda^m(\theta,\phi)\,.
\end{equation}
The bulk of the remaining analysis concerns the radial equation (\ref{8}).
Now applying the WKB admissibility criterion \cite{m1} to (\ref{8}),
\begin{eqnarray}  \label{11}
\Phi(r,n)&=&-\frac{\lambda}{r^2}-
\frac 1{D(r)} \big(\sum_{j=1}^n \frac{\alpha_j}{\tau_j}
(1+(R_0/r)^6)-in\omega\big) \\
&&+\frac12 \Delta[\ln\big(\frac{N_0(r)D(r)}{4\pi r^2}\big)] 
+\frac14\bigg|\nabla[\ln\big(\frac{N_0(r)D(r)}{4\pi r^2}\big)]\bigg|^2
\,, \nonumber\\
&& \qquad \frac{|D_r\Phi(r,n)/2|}{|\Phi(r,n)|^{3/2}} \ll 1\,. \nonumber
\end{eqnarray}
Hence the WKB results are accurate when the above criterion is met by the
coefficient  over the given $(r,n)$-domain.
Recalling (\ref{10}), there obtains the linearly independent WKB analytical 
solutions for the homogeneous radial equation in (\ref{8}):
\begin{eqnarray} 
R_{c_1}(r,n)&\approx& R_{c_1\mbox{\tiny WKB}}(r,n) \nonumber\\
&=& \left\{\begin{array}{ll}
(-\Phi(r,n))^{-1/4}\cos(\int (-\Phi(r,n))^{1/2}\,dr) &\mbox{if $\Phi(r,n)<0$,}\\
(\Phi(r,n))^{-1/4}\exp(\int (\Phi(r,n))^{1/2}\,dr) &\mbox{if $0<\Phi(r,n)$.}
\end{array}\right. \nonumber\\
R_{c_2}(r,n)&\approx& R_{c_2\mbox{\tiny WKB}}(r,n) \label{12}\\
&=& \left\{\begin{array}{ll}
(-\Phi(r,n))^{-1/4}\sin(\int (-\Phi(r,n))^{1/2}\,dr) &\mbox{if $\Phi(r,n)<0$,}\\
(\Phi(r,n))^{-1/4}\exp(-\int (\Phi(r,n))^{1/2}\,dr) &\mbox{if $0<\Phi(r,n)$.}
\end{array}\right. \nonumber
\end{eqnarray}
Hence, the general solution of the radial part of the full non-homogeneous 
equation (\ref{4}) may be supplied via (\ref{12}) as 
\begin{eqnarray}
&&c_1(n)R_{c_1}(r,n)+c_2(n)R_{c_2}(r,n) \label{13}\\
&&+ R_{c_2}(r,n)\int 
\frac {\overline{N}_{E_{t=T}}(r)-\overline{N}_{E_{t=0}}(r)}{D(r)\widetilde{w}(r,n)}
 R_{c_1}(r,n)
 \sqrt{\frac{4\pi r^2}{N_0(r)D(r)}} \,dr \nonumber \\
&& - R_{c_1}(r,n)\int 
\frac {\overline{N}_{E_{t=T}}(r)-\overline{N}_{E_{t=0}}(r)}{D(r)\widetilde{w}(r,n)}
 R_{c_2}(r,n)
 \sqrt{\frac{4\pi r^2}{N_0(r)D(r)}} \,dr \,, \nonumber
\end{eqnarray}
where
$$
\widetilde{w}(r,n)= R_{c_1}(r,n)D_rR_{c_2}(r,n) - R_{c_2}(r,n)D_rR_{c_1}(r,n)
\,.
$$
By ansatz (\ref{10}), (\ref{13}), and (\ref{4}), the time-transformed version 
(\ref{3}) of equation (\ref{12}) has general solution
\begin{eqnarray}
\lefteqn{ \widehat{N}_E(r,\theta,\phi,n)=
\sqrt{\frac{4\pi r^2}{N_0(r)D(r)}} \widetilde{N}_E(r,\theta,\phi,n) }\nonumber\\
&=&\sqrt{\frac{4\pi r^2}{N_0(r)D(r)}} Y_\lambda^m(\theta\phi)  
\bigg(c_1(n)R_{c_1}(r,n)+c_2(n)R_{c_2}(r,n) \label{14}\\
&&+ R_{c_2}(r,n)\int 
\frac {\overline{N}_{E_{t=T}}(r)-\overline{N}_{E_{t=0}}(r)}{D(r)\widetilde{w}(r,n)}
 R_{c_1}(r,n)
 \sqrt{\frac{4\pi r^2}{N_0(r)D(r)}} \,dr \nonumber \\
&& - R_{c_1}(r,n)\int 
\frac {\overline{N}_{E_{t=T}}(r)-\overline{N}_{E_{t=0}}(r)}{D(r)\widetilde{w}(r,n)}
 R_{c_2}(r,n)
 \sqrt{\frac{4\pi r^2}{N_0(r)D(r)}} \,dr \bigg)\,, \nonumber
\end{eqnarray}
Now consider the boundary conditions of (\ref{3}), which are the 
time-transformed boundary conditions of (\ref{I}),
\begin{equation} \label{15}
\partial_r\widehat{N}_E(r,\theta,\phi,n)\bigg|_{r=r_{\mbox{\scriptsize min}}} = 0
\,, \quad
\partial_r\widehat{N}_E(r,\theta,\phi,n)\bigg|_{r=r_{\mbox{\scriptsize max}}} = 0
\,. \end{equation}

Without loss of generality, we may ignore the angular factor 
$Y_\lambda^m(\theta\phi)$ and consider only the radial terms in (\ref{14}). 
Given the particular boundary conditions (\ref{15}), the formula (\ref{14}) 
may be used to calculate the appropriate particular solution matching the 
boundary conditions in the $(r,n)$-domain. In general, there will result a 
2x2 matrix equation for the coefficients $c_1(n)$ and  $c_2(n)$.
\begin{equation} \label{17}
\left( \begin{array}{cc} 
P_{11} & P_{21} \\ P_{21} & P_{22} \end{array}\right)
\left( \begin{array}{c} c_1(n) \\ c_2(n) \end{array}\right)
= \left( \begin{array}{c} Q_1 \\ Q_2 \end{array}\right)\,,
\end{equation}
where
$$ \displaylines{
P_{11}=\partial_t\bigg(\sqrt{\frac{4\pi r^2}{N_0(r)D(r)}} 
       R_{c_1}(r,n)\bigg)\bigg|_{r=r_{\mbox{\scriptsize min}}}\,,\cr
P_{12}=\partial_t\bigg(\sqrt{\frac{4\pi r^2}{N_0(r)D(r)}} 
       R_{c_2}(r,n)\bigg)\bigg|_{r=r_{\mbox{\scriptsize min}}}\,,\cr
P_{21}=\partial_t\bigg(\sqrt{\frac{4\pi r^2}{N_0(r)D(r)}} 
       R_{c_1}(r,n)\bigg)\bigg|_{r=r_{\mbox{\scriptsize max}}}\,,\cr
P_{22}=\partial_t\bigg(\sqrt{\frac{4\pi r^2}{N_0(r)D(r)}} 
       R_{c_2}(r,n)\bigg)\bigg|_{r=r_{\mbox{\scriptsize max}}}\,,\cr
Q_1 =-\partial_r\bigg( R_{c_2}(r,n)\int 
\frac {\overline{N}_{E_{t=T}}(r)-\overline{N}_{E_{t=0}}(r)}{D(r)\widetilde{w}(r,n)}
 R_{c_1}(r,n)
 \sqrt{\frac{4\pi r^2}{N_0(r)D(r)}} \,dr  \cr
 - R_{c_1}(r,n)\int 
\frac {\overline{N}_{E_{t=T}}(r)-\overline{N}_{E_{t=0}}(r)}{D(r)\widetilde{w}(r,n)}
 R_{c_2}(r,n)
 \sqrt{\frac{4\pi r^2}{N_0(r)D(r)}} \,dr \bigg)
\bigg|_{r=r_{\mbox{\scriptsize min}}}\,, \cr
Q_2 =-\partial_r\bigg( R_{c_2}(r,n)\int 
\frac {\overline{N}_{E_{t=T}}(r)-\overline{N}_{E_{t=0}}(r)}{D(r)\widetilde{w}(r,n)}
 R_{c_1}(r,n)
 \sqrt{\frac{4\pi r^2}{N_0(r)D(r)}} \,dr  \cr
 - R_{c_1}(r,n)\int 
\frac {\overline{N}_{E_{t=T}}(r)-\overline{N}_{E_{t=0}}(r)}{D(r)\widetilde{w}(r,n)}
 R_{c_2}(r,n)
 \sqrt{\frac{4\pi r^2}{N_0(r)D(r)}} \,dr \bigg)
\bigg|_{r=r_{\mbox{\scriptsize max}}}\,.\cr
}$$
Then by Cramer's Rule, the functions $c_1(n)$ and $c_1(n)$
have  values
\begin{equation} \label{18}
c_1(n)=\frac{P_{22}Q_1 - P_{21}Q_2}{P_{11}P_{22} - P_{21}P_{12}}\,,\quad
c_2(n)=\frac{P_{11}Q_2 - P_{12}Q_1}{P_{11}P_{22} - P_{21}P_{12}}\,.
\end{equation}

Hence the cefficients in (\ref{18}) are the functions of $n$ necessary for
the given time-transformed boundary conditions of (\ref{15}) to be met in
(\ref{14}). 
Finally, substituting the coefficients in (\ref{18}) into (\ref{14}), 
and applying the inverse finite Fourier transform to $\widehat N_E$ 
in (\ref{14}) yields the operational WKB analytical solution to the 
boundary and initial/final value problem (\ref{Ie}) in the 
$(r,\theta,\phi,t)$-domain.   
\begin{eqnarray*}
\lefteqn{ \overline{N}_E(r,\theta,\phi,t)=
F_T^{-1}[\widehat{N}_E(r,\theta,\phi,n)] }\\
&=&\frac1T\big(\widehat{N}_E(r,\theta,\phi,0)
+\sum_{n=1}^\infty[\widehat{N}_E(r,\theta,\phi,n)e^{in\omega t}
+\widehat{N}_E(r,\theta,\phi,-n)e^{-in\omega t}]\big) \,,\cr
\end{eqnarray*}
with $\widehat{N}_E(r,\theta,\phi,n)$ given by (\ref{14}).
 

\begin{thebibliography}{0}
\bibitem{b1} J. Beechem, E. Haas, {\em Simultaneous determination of
intramolecular distance distributions and conformational dynamics by
global analysis of energy trnasfer measurements}, Biophys J.,
vol. 55, 1989, pp. 1225-1236.

\bibitem{k1} W. Kaplan, {\em Operational Methods for Linear Systems},
Addison-Wesley 1962.

\bibitem{g1} I. S. Gradshteyn, I. M. Ryzhik, {\em Table of Integrals,
Series, and Products}, Academic Press, p. 1013.

\bibitem{m1}  P. Morse, W. Feshbach, {\em Methods of Theoretical Physics},
vol. 2, McGraw-Hill.

\end{thebibliography} \bigskip

\noindent{\sc Valentino Anthony Simpao} \\
Mathematical Consultant Services \\
108 Hopkinsville St.\\
Greenville, Kentucky 42345 USA \\
e-mail: mcs007@muhlon.com \\
Tel.: 502-338-5543
\end{document}