---
abstract: 'We investigate the wormlike polymer chains using self-consistent field theory and take into account the Onsager excluded-volume interaction between polymer segments.  The propagator of polymer chain is one of the essential physical quantities used to study the conformation of polymers, which satisfies the modified diffusion equation (MDE) for wormlike chain.  The propagator of wormlike chain is not only dependent on the spatial variables, but also on the orientation.  We separate the variables of propagator by using spherical-harmonic series and then simplify the MDE to a coupled set of equations only depends on spatial variables in this paper.  We expand the MDE by spherical-harmonic functions in cylindrical coordinates and spherical coordinates, respectively.  We find that there are three ways to set the orientation, no matter in cylindrical coordinates or spherical coordinates.  But for the convenience of calculation, we compare these three forms and choose the simplest one to simplify the MDE.  And we get a coupled set of equations only depends on spatial variables.' 
author:
- Baopi Liu
title: 'Spherical-harmonic Expansion of the Modified Diffusion Equation for Wormlike Chain in Curvilinear Coordinates'
---

Introduction
============

Semiflexible polymer chains have been widely used to describe the structures and dynamics of a large variety of synthetic and biological polymers such as liquid crystalline polymers and DNA molecule.  Kratky and Porod [@1] introduced the wormlike chain model to describe semiflexible polymer chains.  In Saito-Takahashi-Yunoki (STY) [@2] treatment, the configuration of the continuous wormlike chain is specified by a space curve $\textbf{r}(t)$ in which $t \in [0,1]$ is a contour variable that describes the location of a segment along the backbone of the chain, as shown in Figure \[fig:1\].  The vector $\textbf{u} (t) \equiv d \textbf{r} (t)/Ldt$ is a tangent vector to the chain at contour location $t$ and is constrained to be a unit vector [@3; @4], where $L$ is the contour length. 

![ The continuous wormlike chain model describes the configuration of a polymer as a space curve $\textbf{r}(t)$, where $t\in [0,1]$ is a contour parameter.  The chain end positions correspond to $\textbf{r}(0)$ and $\textbf{r}(1)$.[]{data-label="fig:1"}](1.pdf)

In the external field $W(\textbf{r}, \textbf{u})$, the propagator $q(\textbf{r}, \textbf{u}; t)$ of wormlike chains satisfies the modified diffusion equation (MDE) [@4].  Onsager [@5] developed the polymer chains interact with each other through the excluded-volume interaction, as a function of the density distribution.  It is interesting to study the polymers confined in a restricted space, either inside a cylindrical [@9; @10; @12; @13; @15] or spherical pore [@16; @17].  In these cases the usage of curvilinear coordinates for representing $\textbf{r}(t)$ is more convenient than the Cartesian coordinates.  Theoretical understanding of polymer chains in confinement, is a topic that has been widely studied [@6; @7; @8; @11; @14; @18; @19; @20; @21].  Previous works on wormlike chains have been carried out in Cartesian coordinates [@19; @20; @21; @22; @23].  However, in Qin Liang’s work [@24], they demonstrated that the mathematical form of MDE in curvilinear coordinates is quite different from that in Cartesian coordinates, and an additional term in the partial differential equation needs to be accounted in curvilinear coordinates. 

In this paper, we expand the MDE using spherical-harmonic series [@23] in cylindrical and spherical coordinates.  In this process, the most difficulty is how to match the derivative terms so that they can be expressed by the product of three spherical harmonics, and then expressed as the Clebsch-Gordan coefficients.  For convenience to proceed our calculation, it is significant to choose appropriate form of tangent vector $\textbf{u}$.

Model and Theory
================

=-1

We consider a monodisperse solution of $n$ wormlike chains in a structureless solvent, occupying a total volume $V$. Each chain is characterized by persistence length $\lambda$ and its cross-section diameter $d$ ($L, \lambda \gg d$).  The propagator satisfies the MDE [@24]: $$\label{MDE2}
\frac{\partial}{\partial t}q(\textbf{r}, \textbf{u}; t)=\left\{-L\textbf{u}\cdot \nabla_{\textbf{r}}+\frac{L}{2 \lambda} \nabla^2_{\textbf{u}}-W(\textbf{r}, \textbf{u})
+L[(\textbf{u} \cdot \nabla_{\textbf{r}})\textbf{u}] \cdot \nabla_{\textbf{u}}\right\} q(\textbf{r}, \textbf{u}; t)$$ $$\label{MDE3}
\frac{\partial}{\partial t}q^{\ast}(\textbf{r}, \textbf{u}; t)=\left\{-L\textbf{u}\cdot \nabla_{\textbf{r}}-\frac{L}{2 \lambda} \nabla^2_{\textbf{u}}+W(\textbf{r}, \textbf{u})
-L[(\textbf{u} \cdot \nabla_{\textbf{r}})\textbf{u}] \cdot \nabla_{-\textbf{u}}\right\} q^{\ast}(\textbf{r}, \textbf{u}; t)$$ Since the two ends of the polymer chains are different, we introduce the function $q^{*}(\textbf{r}, \textbf{u}; t)$, which is defined as the propagator starting from the opposite end of the polymer, and $q^{\ast}(\textbf{r}, \textbf{u}; t)=q(\textbf{r}, -\textbf{u}; 1-t)$ with initial conditions $$q(\textbf{r}, \textbf{u}; 0)=1$$ $$q^{\ast}(\textbf{r}, \textbf{u}; 1)=1$$

The free energy of the system is [@18] $$\label{eq:freenergy}
\beta F=da^{2}\int d\textbf{r}d\textbf{u}\int d\textbf{u}'\rho(\textbf{r}, \textbf{u})|\textbf{u}\times \textbf{u}'|\rho(\textbf{r}, \textbf{u}')-\int d\textbf{r}d\textbf{u}W(\textbf{r}, \textbf{u})\rho(\textbf{r}, \textbf{u})-\ln(\frac{Q^{n}}{n!})$$ where $a$ represents the Kuhn length and is related to the persistence length $\lambda$ by, $a\equiv 2\lambda$. And $Q=\int d\textbf{r}d\textbf{u}q(\textbf{r}, \textbf{u}; t)q^{\ast}(\textbf{r}, \textbf{u}; t)$ is the single chain partition function.  Taking the saddle point of the free energy given in equation (\[eq:freenergy\]) with respect to $\rho(\textbf{r}, \textbf{u})$ and $W(\textbf{r}, \textbf{u})$, we obtain the self-consistent mean-field equations: $$\label{eq:SCFT1}
W(\textbf{r}, \textbf{u})=2da^{2}\int d\textbf{u}'|\textbf{u}\times \textbf{u}'|\rho(\textbf{r}, \textbf{u})$$ $$\label{eq:SCFT2}
\rho(\textbf{r}, \textbf{u})=\frac{n}{Q}\int_{0}^{1}dtq(\textbf{r}, \textbf{u}; t)q^{\ast}(\textbf{r},\textbf{u}; t)$$

Expanding $\rho(\textbf{r}, \textbf{u})$, $W(\textbf{r}, \textbf{u})$, $q(\textbf{r}, \textbf{u}; t)$ and $q^{\ast}(\textbf{r}, \textbf{u}; t)$ in terms of spherical-harmonic series [@23]: $$\begin{aligned}
\begin{aligned}
      &\rho(\textbf{r}, \textbf{u})=\sum\limits_{l,m}\rho_{l,m}(\textbf{r})Y_{l,m}(\textbf{u}),\qquad
      W(\textbf{r}, \textbf{u})=\sum\limits_{l,m}W_{l,m}(\textbf{r})Y_{l,m}(\textbf{u})\\
      &q(\textbf{r}, \textbf{u}; t)=\sum\limits_{l,m}q_{l,m}(\textbf{r},t)Y_{l,m}(\textbf{u}),\qquad
      q^{\ast}(\textbf{r}, \textbf{u}; t)=\sum\limits_{l,m}q_{l,m}^{\ast}(\textbf{r},t)Y_{l,m}(\textbf{u})
\end{aligned}\end{aligned}$$ Since the research objects are all real functions, the expansion coefficients must obey the following conditions [@23]: $$\rho_{l,m}(\textbf{r})=(-1)^{m}\rho_{l,-m}^{\ast}(\textbf{r})$$ Next, the spherical-harmonic expansion for the kernel $|\textbf{u}\times \textbf{u}{'}|$ by use of the addition theorem [@18] $$|\textbf{u}\times \textbf{u}{'}|=\sum\limits_{l,m}\frac{4\pi}{2l+1}d_{l}\cdot Y_{l,m}(\textbf{u})Y_{l,m}^{\ast}(\textbf{u}{'})$$ with [@18] $$d_{l}=
\begin{cases}
     \displaystyle \frac{\pi}{4}, & l=0,\\
     \displaystyle 0, & l=1,3,5,\dots,\\
     \displaystyle -\frac{\pi(2l+1)l!(l-2)!}{2^{2l+1}(\frac{l}{2}-1)!(\frac{l}{2})!(\frac{l}{2})!(\frac{l}{2}+1)!}, &  l=2,4,6,\dots
\end{cases}$$

The mean-field equations (\[eq:SCFT1\]) and (\[eq:SCFT2\]) can be rewrote as $$W_{l,m}(\textbf{r})=\frac{8\pi}{2l+1}d_{l}\cdot da^{2}\rho_{l,m}(\textbf{r})$$ $$\label{eq:trigonometric}
\begin{split}
\rho_{l,m}(\textbf{r})&=\frac{n}{Q}\int_{0}^{1}dt\sum\limits_{l{'},m{'}}\sum\limits_{l{''},m{''}}q_{l{'},m{'}}^{\ast}(\textbf{r},t)q_{l{''},m{''}}(\textbf{r},t)\int
  d\textbf{u}Y_{l{'},m{'}}(\textbf{u})Y_{l{''},m{''}}(\textbf{u})Y_{l,m}^{\ast}(\textbf{u})\\
  &=\frac{n}{Q}\int_{0}^{1}dt\sum\limits_{l{'},m{'}}\sum\limits_{l{''},m{''}}q_{l{'},m{'}}^{\ast}(\textbf{r},t)q_{l{''},m{''}}(\textbf{r},t)
  \sqrt{\frac{(2l{''}+1)(2l{'}+1)}{4\pi(2l+1)}}\times C_{0,0,0}^{l{''},l{'},l}C_{m{''},m{'},m}^{l{''},l{'},l}
\end{split}$$ where $C_{m{''},m{'},m}^{l{''},l{'},l}$ are Clebsch-Gordan coefficients, it is the result of integral for three spherical harmonics. 

$$q_{l,m}(\textbf{r},t)=\int d\textbf{u}Y_{l,m}^{\ast}(\textbf{u})q(\textbf{r}, \textbf{u}; t)$$

Equation (\[MDE2\]) becomes: $$\label{eq:MDE1}
\begin{split}
\frac{\partial}{\partial t}q_{l,m}(\textbf{r},t)&=\int d\textbf{u}Y_{l,m}^{\ast}(\textbf{u})\frac{\partial}{\partial t}q(\textbf{r}, \textbf{u}; t)\\
&=\int d\textbf{u}Y_{l,m}^{\ast}(\textbf{u})\left\{-L\textbf{u}\cdot \nabla_{\textbf{r}}+\frac{L}{2 \lambda} \nabla^2_{\textbf{u}}-W(\textbf{r}, \textbf{u})
+L[(\textbf{u} \cdot \nabla_{\textbf{r}})\textbf{u}] \cdot \nabla_{\textbf{u}}\right\}\sum\limits_{l{'},m{'}}q_{l{'},m{'}}(\textbf{r},t)Y_{l{'},m{'}}(\textbf{u})
\end{split}$$

Spherical-harmonic expansion
============================

Cylindrical coordinates
-----------------------

In the cylindrical coordinate system, the spatial variable $\textbf{r}$ is represented by variables $\rho$, $\Phi$, $z$, the unit vector are given by the orthonormal unit vectors $\textbf{e}_{\rho}$, $\textbf{e}_{\Phi}$, $\textbf{e}_{z}$, as shown in Figure \[fig:2\]. 

![ The considered polymer is represented by gray line, the spatial variable $\textbf{r}$ is represented by variables $\rho$, $\Phi$, $z$ in cylindrical coordinates.  The local coordinate frame is given by orthonormal unit vectors $\textbf{e}_{\rho}$, $\textbf{e}_{\Phi}$, and $\textbf{e}_{z}$.[]{data-label="fig:2"}](2.pdf)

We establish a spherical coordinate system for $\textbf{u}$ within this local frame [@24], and there are three forms of $\textbf{u}$, as shown in Figure \[fig:3\]. 

![ In the cylindrical coordinate system, the three forms of the orientation vector $\textbf{u}$.[]{data-label="fig:3"}](3.pdf)

For calculation convenience, we set the orientation vector $\textbf{u}$ in the following form (Figure \[fig:3\](1)) $$\textbf{u}=\sin\theta \cos\phi\textbf{e}_{\rho}+\sin\theta \sin\phi\textbf{e}_{\Phi}+\cos\theta\textbf{e}_{z}$$ The propagator $q(\textbf{r}, \textbf{u}; t)$ can be represented by $q(\rho,\Phi,z,\theta,\phi; t)$, where $\theta$ is the polar angle and $\phi$ is the azimuthal angle.  And the derivative terms in equation (\[eq:MDE1\]) become $$\textbf{u}\cdot \nabla_{\textbf{r}}q(\textbf{r}, \textbf{u}; t)=\Big(\sin\theta \cos\phi \frac{\partial}{\partial \rho}+\frac{\sin\theta \sin\phi}{\rho}\frac{\partial}{\partial \Phi}+\cos\theta\frac{\partial}{\partial z}\Big)q(\textbf{r}, \textbf{u}; t)$$ and $$\label{eq:de2}
[(\textbf{u}\cdot \nabla_{\textbf{r}})\textbf{u}]\cdot \nabla_{\textbf{u}} q(\textbf{r}, \textbf{u}; t)=\frac{\sin\theta \sin\phi}{\rho}\frac{\partial}{\partial \phi}q(\textbf{r}, \textbf{u}; t)$$ The last term in the right of equation (\[eq:MDE1\]) is zero in Cartesian coordinates, but as can be seen from equation (\[eq:de2\]), it is non-zero in cylindrical coordinates. 

Expanding equation (\[eq:MDE1\]) using spherical-harmonic series, we obtain the following coupled set of equations in cylindrical coordinates[@25]: $$\label{eq:solution1}
\begin{split}
  \frac{\partial}{\partial t}q_{l,m}(\textbf{r},t)=&-L \sum\limits_{l{'},m{'}}\sqrt{\frac{2l{'}+1}{2l+1}}\left\{C_{0,0,0}^{1,l{'},l} \Big[\frac{1}{\sqrt{2}}(C_{-1,m{'},m}^{1,l{'},l}-C_{1,m{'},m}^{1,l{'},l})\frac{\partial}{\partial \rho}\right. \\ 
  & \left.  +\frac{i}{\sqrt{2}}(C_{-1,m{'},m}^{1,l{'},l}+C_{1,m{'},m}^{1,l{'},l})\frac{1}{\rho}\frac{\partial}{\partial \Phi}+C_{0,m{'},m}^{1,l{'},l}\frac{\partial}{\partial z}\Big] \right\}q_{l{'},m{'}}(\textbf{r},t)\\
  &-\frac{L}{2\lambda}l(l+1)q_{l,m}(\textbf{r},t)\\
  &-\sum\limits_{l{'},m{'}}\sum\limits_{l{''},m{''}}\sqrt{\frac{(2l{''}+1)(2l{'}+1)}{4\pi(2l+1)}}C_{0,0,0}^{l{''},l{'},l} \cdot C_{m{''},m{'},m}^{l{''},l{'},l} \cdot W_{l{'},m{'}}(\textbf{r}) \cdot q_{l{''},m{''}}(\textbf{r},t)\\
  &-\sum\limits_{l{'},m{'}}\frac{L}{\sqrt{2}\rho}m{'}\sqrt{\frac{2l{'}+1}{2l+1}}C_{0,0,0}^{1,l{'},l}(C_{-1,m{'},m}^{1,l{'},l}+C_{1,m{'},m}^{1,l{'},l})q_{l{'},m{'}}(\textbf{r},t)
\end{split}$$ with initial conditions $$\label{eq:condition1}
q_{0,0}(\textbf{r},0)=\sqrt{4\pi}, \qquad q_{l,m}(\textbf{r},0)=0, \qquad \rm otherwise.$$ Similarly, we obtain $$\label{eq:solution2}
\begin{split}
  \frac{\partial}{\partial t}q_{l,m}^{\ast}(\textbf{r},t)=&-L \sum\limits_{l{'},m{'}}\sqrt{\frac{2l{'}+1}{2l+1}}
  \left\{C_{0,0,0}^{1,l{'},l} \Big[\frac{1}{\sqrt{2}}(C_{-1,m{'},m}^{1,l{'},l}-C_{1,m{'},m}^{1,l{'},l})\frac{\partial}{\partial \rho}\right. \\ 
  & \left.  +\frac{i}{\sqrt{2}}(C_{-1,m{'},m}^{1,l{'},l}+C_{1,m{'},m}^{1,l{'},l})\frac{1}{\rho}\frac{\partial}{\partial \Phi}+C_{0,m{'},m}^{1,l{'},l}\frac{\partial}{\partial z}\Big] \right\}q_{l{'},m{'}}^{\ast}(\textbf{r},t)\\
  &+\frac{L}{2\lambda}l(l+1)q_{l,m}^{\ast}(\textbf{r},t)\\
  &+\sum\limits_{l{'},m{'}}\sum\limits_{l{''},m{''}}\sqrt{\frac{(2l{''}+1)(2l{'}+1)}{4\pi(2l+1)}}C_{0,0,0}^{l{''},l{'},l} \cdot C_{m{''},m{'},m}^{l{''},l{'},l} \cdot W_{l{'},m{'}}(\textbf{r}) \cdot q_{l{''},m{''}}^{\ast}(\textbf{r},t)\\
  &-\sum\limits_{l{'},m{'}}\frac{L}{\sqrt{2}\rho}m{'}\sqrt{\frac{2l{'}+1}{2l+1}}C_{0,0,0}^{1,l{'},l}(C_{-1,m{'},m}^{1,l{'},l}+C_{1,m{'},m}^{1,l{'},l})q_{l{'},m{'}}^{\ast}(\textbf{r},t)
\end{split}$$ with initial conditions $$\label{eq:condition2}
q_{0,0}^{\ast}(\textbf{r},1)=\sqrt{4\pi}, \qquad q_{l,m}^{\ast}(\textbf{r},1)=0, \qquad \rm otherwise$$

According to the above expressions, the single chain partition function and the free energy (\[eq:freenergy\]) can be represented by $$Q=\sum\limits_{l,m}(-1)^{m}\int d\textbf{r}q_{l,m}(\textbf{r},t)q_{l,-m}^{\ast}(\textbf{r},t)$$ $$\beta F=-\sum\limits_{l,m}\frac{4\pi}{2l+1}d_{l}\cdot da^{2}\int d\textbf{r}|\rho_{l,m}(\textbf{r})|^{2}-\ln(\frac{Q^{n}}{n!})$$

We now analyze the the simplest case where the densities vary in only one spatial dimension, here is chosen to be the $z$ direction.  Then spatial variable of the propagator independent of $\rho$ and $\Phi$, and the unit vector $\textbf{u}(t)$ only depends on the polar angle $\theta$, not on the azimuth angle $\phi$. Therefor, in equations (\[eq:solution1\]) and (\[eq:solution2\]), $m=m{'}=m{''}=0$, so that we subsequently drop the subscript which makes the Clebsch-Gordan coefficients $C_{-1,0,0}^{1,l{'},l}=C_{1,0,0}^{1,l{'},l}=0$ and equations (\[eq:solution1\]) and (\[eq:solution2\]) can be reduced to: $$\begin{split}
  \frac{\partial}{\partial t}q_{l}(z,t)=
  &-L\sum\limits_{l{'}}\sqrt{\frac{2l{'}+1}{2l+1}}(C_{0,0,0}^{1,l{'},l})^{2}\frac{\partial}{\partial z}q_{l{'}}(z,t)-\frac{L}{2\lambda}l(l+1)q_{l}(z,t)\\
  &-\sum\limits_{l{'},l{''}}\sqrt{\frac{(2l{''}+1)(2l{'}+1)}{4\pi(2l+1)}}(C_{0,0,0}^{1,l{'},l})^{2}\cdot W_{l{'}}(z)\cdot q_{l{''}}(z,t)
\end{split}$$ and $$\begin{split}
  \frac{\partial}{\partial t}q_{l}^{\ast}(z,t)=
  &-L\sum\limits_{l{'}}\sqrt{\frac{2l{'}+1}{2l+1}}(C_{0,0,0}^{1,l{'},l})^{2}\frac{\partial}{\partial z}q_{l{'}}^{\ast}(z,t)-\frac{L}{2\lambda}l(l+1)q_{l}^{\ast}(z,t)\\
  &+\sum\limits_{l{'},l{''}}\sqrt{\frac{(2l{''}+1)(2l{'}+1)}{4\pi(2l+1)}}(C_{0,0,0}^{1,l{'},l})^{2}\cdot W_{l{'}}(z)\cdot q_{l{''}}^{\ast}(z,t)
\end{split}$$ with initial conditions $$q_{0}(z,0)=\sqrt{4\pi}, \qquad q_{l}(z,0)=0, \qquad l>0$$ and $$q_{0}^{\ast}(z,1)=\sqrt{4\pi}, \qquad q_{l}^{\ast}(z,1)=0, \qquad l>0$$

Similarly, according to the above expressions, the one-dimensional single-chain partition function and the free energy in equation (\[eq:freenergy\]) can be rewritten as $$Q=\sum\limits_{l}\int dz q_{l}(z,t)q_{l}^{\ast}(z,t)$$ $$\beta F=-\sum\limits_{l}\frac{4\pi}{2l+1}d_{l}\cdot da^{2}\int dz \rho_{l}^{2}(z)-\ln(\frac{Q^{n}}{n!})$$ This result is consistent with that obtained by Dominik Duchs et al.  [@23]. 

Spherical coordinates
---------------------

In the spherical coordinate system, the spatial variable $\textbf{r}$ is represented by variables $\Theta$, $\Phi$, $r$, the unit vector are given by the orthonormal unit vectors $\textbf{e}_{\Theta}$, $\textbf{e}_{\Phi}$, $\textbf{e}_{r}$, as shown in Figure \[fig:4\]. 

![ The considered polymer is represented by gray line, the spatial variable $\textbf{r}$ is represented by variables $\Theta$, $\Phi$, $r$ in spherical coordinates.  The local coordinate frame is given by orthonormal unit vectors $\textbf{e}_{\Theta}$, $\textbf{e}_{\Phi}$, and $\textbf{e}_{r}$.[]{data-label="fig:4"}](4.pdf)

We establish a spherical coordinate system for $\textbf{u}$ within this local frame [@24], and there are three forms of $\textbf{u}$, as shown in Figure \[fig:5\] . 

![ In the spherical coordinate system, the three forms of the orientation vector $\textbf{u}$.[]{data-label="fig:5"}](5.pdf)

For calculation convenience, we set the orientation vector $\textbf{u}$ in the following form (Figure \[fig:5\](2)) $$\textbf{u}=\sin\theta \cos\phi\textbf{e}_{\Theta}+\sin\theta \sin\phi\textbf{e}_{\Phi}+\cos\theta\textbf{e}_{r}$$ The propagator $q(\textbf{r}, \textbf{u}; t)$ is represented by $q(\Theta, \Phi, r, \theta, \phi; t)$, where $\theta$ is the polar angle and $\phi$ is the azimuthal angle.  And the derivative terms in equation $(15)$ can be rewritten as $$\textbf{u}\cdot \nabla_{\textbf{r}}q(\textbf{r}, \textbf{u}; t)=\Big(\frac {\sin\theta \cos\phi}{r} \frac{\partial}{\partial \Theta}+\frac{\sin\theta \sin\phi}{r \sin\Theta}\frac{\partial}{\partial \Phi}+\cos\theta\frac{\partial}{\partial r}\Big)q(\textbf{r}, \textbf{u}; t)$$ and $$[(\textbf{u}\cdot \nabla_{\textbf{r}})\textbf{u}]\cdot \nabla_{\textbf{u}} q(\textbf{r}, \textbf{u}; t)=\Big(\frac{\sin\theta}{r}\frac{\partial}{\partial \theta}+\frac{\sin \theta \sin\phi \cos\Theta}{r \sin\Theta}\frac{\partial}{\partial \phi}\Big)q(\textbf{r}, \textbf{u}; t)$$

Expanding equation (\[eq:MDE1\]) using spherical-harmonic series, yields the following coupled set of equations in spherical coordinates [@25]: $$\begin{split}
  \frac{\partial}{\partial t}q_{l,m}(\textbf{r},t)=&-L \sum\limits_{l{'},m{'}}\sqrt{\frac{2l{'}+1}{2l+1}}\left\{C_{0,0,0}^{1,l{'},l} \Big[\frac{1}{\sqrt{2}}(C_{-1,m{'},m}^{1,l{'},l}-C_{1,m{'},m}^{1,l{'},l})\frac{1}{r}\frac{\partial}{\partial \Theta}\right. \\ 
  & \left.  +\frac{i}{\sqrt{2}}(C_{-1,m{'},m}^{1,l{'},l}+C_{1,m{'},m}^{1,l{'},l})\frac{1}{r \sin\Theta}\frac{\partial}{\partial \Phi}+C_{0,m{'},m}^{1,l{'},l}\frac{\partial}{\partial r}\Big] \right\}q_{l{'},m{'}}(\textbf{r},t)\\
  &-\frac{L}{2\lambda}l(l+1)q_{l,m}(\textbf{r},t)\\
  &-\sum\limits_{l{'},m{'}}\sum\limits_{l{''},m{''}}\sqrt{\frac{(2l{''}+1)(2l{'}+1)}{4\pi(2l+1)}}C_{0,0,0}^{l{''},l{'},l} \cdot C_{m{''},m{'},m}^{l{''},l{'},l} \cdot W_{l{'},m{'}}(\textbf{r}) \cdot q_{l{''},m{''}}(\textbf{r},t)\\
  &+\sum\limits_{l{'},m{'}}\frac{L}{\sqrt{2} r}\sqrt{\frac{2l{'}+1}{2l+1}}C_{0,0,0}^{1,l{'},l}\Big[\sqrt{(l{'}-m{'})(l{'}+m{'}+1)}C_{-1,m{'}+1,m}^{1,l{'},l}\\
  &-\sqrt{(l{'}+m{'})(l{'}-m{'}+1)}C_{1,m{'}-1,m}^{1,l{'},l}-\frac{\cos \Theta}{\sin \Theta}m{'}(C_{-1,m{'},m}^{1,l{'},l}+C_{1,m{'},m}^{1,l{'},l})\Big]q_{l{'},m{'}}(\textbf{r},t)
\end{split}$$ The initial conditions are $$q_{0,0}(\textbf{r},0)=\sqrt{4\pi}, \qquad q_{l,m}(\textbf{r},0)=0,  \qquad \rm otherwise$$ Similarly, we obtain $$\begin{split}
  \frac{\partial}{\partial t}q_{l,m}^{\ast}(\textbf{r},t)=&L \sum\limits_{l{'},m{'}}\sqrt{\frac{2l{'}+1}{2l+1}}\left\{C_{0,0,0}^{1,l{'},l} \Big[\frac{1}{\sqrt{2}}(C_{-1,m{'},m}^{1,l{'},l}-C_{1,m{'},m}^{1,l{'},l})\frac{1}{r}\frac{\partial}{\partial \Theta}\right. \\ 
  & \left.  +\frac{i}{\sqrt{2}}(C_{-1,m{'},m}^{1,l{'},l}+C_{1,m{'},m}^{1,l{'},l})\frac{1}{r \sin\Theta}\frac{\partial}{\partial \Phi}+C_{0,m{'},m}^{1,l{'},l}\frac{\partial}{\partial r}\Big] \right\}q_{l{'},m{'}}^{\ast}(\textbf{r},t)\\
  &-\frac{L}{2\lambda}l(l+1)q_{l,m}^{\ast}(\textbf{r},t)\\
  &-\sum\limits_{l{'},m{'}}\sum\limits_{l{''},m{''}}\sqrt{\frac{(2l{''}+1)(2l{'}+1)}{4\pi(2l+1)}}C_{0,0,0}^{l{''},l{'},l} \cdot C_{m{''},m{'},m}^{l{''},l{'},l} \cdot W_{l{'},m{'}}(\textbf{r}) \cdot q_{l{''},m{''}}^{\ast}(\textbf{r},t)\\
  &+\sum\limits_{l{'},m{'}}\frac{L}{\sqrt{2} r}\sqrt{\frac{2l{'}+1}{2l+1}}C_{0,0,0}^{1,l{'},l}\Big[\sqrt{(l{'}-m{'})(l{'}+m{'}+1)}C_{-1,m{'}+1,m}^{1,l{'},l}\\
  &-\sqrt{(l{'}+m{'})(l{'}-m{'}+1)}C_{1,m{'}-1,m}^{1,l{'},l}+\frac{\cos \Theta}{\sin \Theta}m{'}(C_{-1,m{'},m}^{1,l{'},l}+C_{1,m{'},m}^{1,l{'},l})\Big]q_{l{'},m{'}}^{\ast}(\textbf{r},t)
\end{split}$$ with initial conditions $$q_{0,0}^{\ast}(\textbf{r},1)=\sqrt{4\pi}, \qquad q_{l,m}^{\ast}(\textbf{r},1)=0, \qquad \rm otherwise$$

According to the above expressions, the single chain partition function and the free energy (\[eq:freenergy\]) can be rewritten as $$Q=\sum\limits_{l,m}(-1)^{m}\int d\textbf{r}q_{l,m}(\textbf{r},t)q_{l,-m}^{\ast}(\textbf{r},t)$$ $$\beta F=-\sum\limits_{l,m}\frac{4\pi}{2l+1}d_{l}\cdot da^{2}\int d\textbf{r}|\rho_{l,m}(\textbf{r})|^{2}-\ln(\frac{Q^{n}}{n!})$$

Results and Discussion
======================

In this paper, using wormlike chain model to describe the semiflexible polymers with the Onsager excluded-volume interaction, we are able to obtain the MDE in an external field.  Taking the saddle point approximation of the free energy function with $\rho(\textbf{r}, \textbf{u})$ and $W(\textbf{r}, \textbf{u})$, we obtained the mean-field equations.  It is more convenient to represent $\textbf{r}$ in curvilinear coordinates when wormlike chains confined in a cylindrical or spherical pore.  The solution of the propagator in an external field is an essential tool to calculate other properties such as the average polymer conformation, the density profile and the free energy of a spatially and orientational inhomogeneous system.  In order to obtain the solution of the propagator, we expanded the MDE using spherical-harmonic series.  And we obtained the expansion of spherical-harmonic series of MDE in cylindrical and spherical coordinates. 

We find that there are three ways to set the orientation vector $\textbf{u}$ in cylindrical coordinates and spherical coordinates, respectively.  Compare the three forms of orientation vector $\textbf{u}$ , we choose the appropriate one to make the expansion form most simplest.  The most difficulty of this investigation is that the gradient which depend on the spatial variable $\textbf{r}$ and the gradient, divergence and Laplace operate which depend on the orientation $\textbf{u}$ are mixed together.  We expressed the terms which depend on the orientational variable $\textbf{u}$ as the product of three spherical harmonics.  Then we simplify the MDE to a coupled set of equations only depend on the spatial variable.  In this paper, we only derive the coupled set of equations but do not do numerical calculation for specific examples.  In the following study, we will investigate specific objects using this coupled set of equations. 

#### Acknowledgments:

This work was supported by the National Natural Science Foundation of China (NSFC) Nos.  21873015, 21434001. 

[99]{} O. Kratky, G. Porod, *Recl.  Trav.  Chim.  Pays-Bas* **1949**, *68*, 1106.  N. Saito, K. Takahashi, Y. Yunoki, *J. Phys.  Soc.  Jpn.* ** 1967**, *22*, 219.  K. F. Freed, *Adv.  Chem.  Phys.* ** 1972**, *22*, 1.  G. H. Fredrickson, *The Equilibrium Theory of Inhomogeneous Polymers*, Clarendon Press, Oxford **2006**.  L. Onsager, *Ann.  NY Acad.  Sci.* ** 1949**, *51*, 627.  T. Odijk, *Macromolecules* **1983**, *16*, 1340.  T. Odijk, *Macromolecules* **1986**, *19*, 2313.  T. W. Burkhardt, *J. Phys.  A: Math.  Gen.* **1995**, *28*, L629.  T. W. Burkhardt, *J. Phys.  A: Math.  Gen.* **1997**, *30*, L167.  D. J. Bicout, T. W. Burkhardt, *J. Phys.  A: Math.  Gen.* **2001**, *34*, 5745.  T. W. Burkhardt, Y. Yang, G. Gompper, *Phys.  Rev. E.* **2010**, *82*, 041801.  Y. Yang, T. W. Burkhardt, G. Gompper, *Phys.  Rev. E.* **2007**, *76*, 011804.  Z. Y. Chen, *Macromolecules* **2013**, *46*, 9837.  G. Morrison, D. Thirumalai, *Phys.  Rev. E.* **2009**, *79*, 011924.  I. Kusner, S. Srebnik, *Chem.  Phys.  Lett.* ** 2006**, *430*, 84-88.  A. Nikoubashman, D. A. Vega, K. Binder, A. Milchev, *Phys.  Rev. Lett.* ** 2017**, *118*, 217803.  A. J. Spakowitz and Z. G. Wang, *Phys.  Rev. Lett.* ** 2003**, *91*, 166102.  Z. Y. Chen, *Prog.  Polym.  Sci.* ** 2016**, *54*, 3.  Z. Y. Chen, D. E. Sullivan, X. Yuan, *Europhys.  Lett.* ** 2005**, *72*, 89.  Z. Y. Chen, D. E. Sullivan, *Macromolecules* **2006**, *39*, 7769.  Z. Y. Chen, D. E. Sullivan, X. Yuan, *Macromolecules* **2007**, *40*, 1187.  G. J. Vroege, T. Odijk, *Macromolecules* **1988**, *21*, 2848.  D. Duchs, D. E. Sullivan, *J. Phys.:  Condens.  Matter* **2002**, *14*, 12189.  Q. Liang, J. Li, P. Zhang and Z. Y. Chen, *J. Chem.  Phys.* ** 2013**, *138*, 244910.  G. B. Arfken, H. J. Weber, L. Ruly, *Mathematical Methods for Physicists*, Academic Press, New York **1985**.