text:
# Nematic director slippage: Role of the angular momentum of light

## I Introduction

Light-induced reorientation effects in liquid crystals are important both for fundamental reasons and because of the possibility of using these effects in applications Khoo1993 . In the majority of cases, the director reorientation occurs in the cell bulk and is imposed by the optical torque. This torque can originate either from: orientational nonlinearity Zel'dovich1982 , or from photorefractive effects Khoo1996 ; Rudenko1996 , or as a result of photoinduced nonlinearity in dye-doped LCs Janossy1990 ; Janossy1991 ; Janossy1992 ; Janossy1994 ; Zolot'ko1994 ; MuensterPRL1997  .

Light-induced reorientation effects have usually been observed in liquid crystal cells with strong anchoring. By strong anchoring we imply that the director reorientation occurs in the bulk of liquid crystal and is negligible at the aligning surface. Development of new aligning materials providing weak anchoring BryanBrown1999 ; Andrienko1998 ; Andrienko1999  has allowed the observation of surface director reorientation.

The first observation of such a reorientation, referred hereafter as a slippage effect, was reported by Marusii et al Marusii1996  in a combined liquid crystal cell. One of the cell surfaces with strong planar anchoring provided homogeneous alignment of the nematic director parallel to the substrates. The other surface had rather small azimuthal anchoring. Liquid crystal was doped with azo-dye in order to increase the orientational nonlinearity coefficient, due to photoinduced trans-cis isomerization of azo-dye molecules MuensterPRL1997 ; JanossyPRE1998 . A linearly polarized laser beam, propagating through the birefringent mixture, induced a bulk torque accompanied by director reorientation in the cell bulk and at the weak anchoring surface. The surface director reorientation was observed as a change in the polarization state of a probe laser beam propagating from the side of the strongly anchoring surface.

Recently Francescangeli et al Francescangeli1999  found that the slippage effect competes with a the light-induced anchoring effect discovered earlier by Voloshchenko et al Voloshchenko1994 . At the same time, was noticed that this effect is sensitive to the polarization of the laser beam AndrienkoMCLC1999 ; AndrienkoUFJ1999 .

The key element of the present paper is to combine the physics of surface director reorientation with the use of light-induced bulk director reorientation. The reorientation occurs in a cell in which the director rotates in an elliptically polarized laser beam, subject to boundary conditions which constrain the surface director to remain in the plane of the cell. Light which propagates through a linear birefringent medium in this geometry changes its polarization state in a periodic way, with wavelength $`\lambda /(n_en_o)`$. In the course of one wavelength, the ellipticity of the light changes periodically between -1 and 1, as the light alternates between linear and circular polarizations. A liquid crystal, although no longer strictly a linear birefringent medium, in that there is a feedback on the director from the local fields, nevertheless conserves this property except for very intense light. As a consequence, the field-induced torque on the liquid crystal director is also spatially modulated and periodic. However, if the non-linear feedback increases, the ellipticity  spatial modulation is accompanied by a director reorientation. The reorientation of the director in the cell bulk immediately manifests itself as a director reorientation on the surface, if the surface anchoring is weak.

Since the director torque compensates the deposition of the angular momentum of light into the liquid crystal, the director reorientation depends significantly on the polarization state of the incident light beam. As a result, surface director reorientation is also sensitive to the polarization of the incident light beam.

The paper is organized as follows. In §II we solve the self-consistent problem of propagation of monochromatic light wave in the planar liquid crystal cell. We use approximations in the spirit of geometric optics to obtain the director distribution in the cell. Then in §III we compare existing experimental results to the theoretical predictions and estimate the photoinduced nonlinearity coefficient. Finally in §IV we present some brief conclusions.

## II The director and ellipticity profiles

Consider a nematic liquid crystal cell of thickness $`L`$ confined between the planes $`z=0`$ and $`z=L`$ of a cell uniform in the $`xy`$ plane. Let the director $`𝐧=(\mathrm{cos}\phi ,\mathrm{sin}\phi ,0)`$ describe the average molecular orientation in the cell. We suppose $`\phi =\phi (z)`$ and only consider the equilibrium orientation of the director, i.e. the angle $`\phi `$ does not depend on time.

Let the elliptically polarized monochromatic light wave, with wave number $`𝐤=k𝐳`$, ellipticity $`e_0`$, and major ellipse axis making an angle $`\alpha `$ with the $`x`$ axis, be normally incident on the cell from the side of the surface with weak anchoring. Note that, because the light is propagating in an inhomogeneous anisotropic medium, the Poynting vector and the light polarization both vary in space throughout the sample.

The light illumination produces a bulk torque acting on the liquid crystal. The director distribution is determined by the balance between optical, elastic, and surface torques. To calculate the director configuration $`𝐧(z)`$, it is necessary to write down an effective free energy density $`f`$, which can then be minimized with respect to the director profile. This is:

$$f=f_{\mathrm{el}}+f_{\mathrm{opt}},$$
(1)
where $`f_{\mathrm{el}}=\frac{1}{2}K_{22}\left(\phi /z\right)^2`$ is the nematic liquid crystal elastic energy and $`K_{22}`$ is the twist elastic constant.

To obtain an expression for the electromagnetic energy density $`f_{\mathrm{opt}}`$ we suppose the nematic liquid crystal to be a slowly varying uniaxial medium, so that $`\phi `$ varies appreciably only over a length much greater than the optical wavelength $`\lambda `$. We also assume that the birefringence of the medium is small, i.e. $`\mathrm{\Delta }n=n_en_o<<1`$. In this case the light polarization also varies slowly through the medium, and we can use the geometrical optics approximation (GOA) to solve Maxwell’s equations for the field in the cell Ong1983 . We shall also suppose the medium to be non-absorbing. In this case the beam intensity $`I`$, defined as the $`z`$ component of the average Poynting vector, remains constant.

The approach we shall follow was developed by Santamato et al. Santamato1988 ; Abbate1991 . We first remind the reader of the elements of this theory.

The total electromagnetic energy density $`f_{\mathrm{opt}}`$ of the light wave in a non-absorbing, nonmagnetic medium can be rewritten in terms of the light intensity $`I`$, ellipticity $`e`$, and major axis angle $`\psi `$. In the limit of low birefringence this takes the form

$$f_{\mathrm{opt}}=\frac{In_o}{c}\frac{I\mathrm{\Delta }n}{2c}[1+\left(1e^2\right)^{1/2}\mathrm{cos}2\left(\psi \phi \right)],$$
(2)
where $`I`$ is the $`z`$ component of Poynting vector, or equivalently, the total intensity of the light wave.

The ellipticity $`e`$ and the major axis angle $`\psi `$ can be defined in terms of the Stokes parameters $`\left\{S_i\right\}`$, where

$`S_0`$ $`=`$ $`|E_x|^2+|E_y|^2,`$ (3)
$`S_1`$ $`=`$ $`|E_x|^2|E_y|^2,`$
$`S_2`$ $`=`$ $`2\mathrm{R}\mathrm{e}\left(E_x^{}E_y\right),`$
$`S_3`$ $`=`$ $`2\mathrm{I}\mathrm{m}\left(E_x^{}E_y\right).`$
Then the ellipticity, defined by $`e=S_3/S_0`$, is the ellipticity of the polarization ellipse of the light in the cell, and $`2\psi =\mathrm{arctan}\left(S_2/S_1\right)`$ is the angle which the major axis of the polarization ellipse forms with the $`x`$ axis.

Equations governing the polarization of light can be derived from the electromagnetic energy density $`f_{\mathrm{opt}}`$ by considering $`f_{\mathrm{opt}}`$ as a Hamiltonian function. Then $`\psi `$ is a generalized coordinate and $`l_z=(I/\omega )e`$ is conjugate momentum, where $`l_z`$ is the average angular momentum carried by the optical beam along the propagation direction.

These equations are:

$$\frac{\psi }{z}=\frac{\omega }{2c}\mathrm{\Delta }n\frac{e}{\left(1e^2\right)^{1/2}}\mathrm{cos}2\left(\psi \phi \right),$$
(4)
$$\frac{e}{z}=\frac{\omega }{c}\mathrm{\Delta }n\left(1e^2\right)^{1/2}\mathrm{sin}2\left(\psi \phi \right).$$
(5)
The equilibrium orientation of the molecular director is given by the minimizer of the functional $`_0^Lf\left\{\phi \left(z\right)\right\}𝑑z`$ subject to fixed field intensity $`I`$ and polarization state $`(e,\psi )`$.

Instead of deriving variational equations it is more useful to use the conserved quantities corresponding to the Lagrange density $`f`$. This does not depend explicitly on the $`z`$ coordinate. Noether’s theorem then implies that there is a conserved integral analogous to the energy in classical dynamics. Thus the first integral for this problem can be obtained directly from the free energy density (1):

$$\frac{K_2}{2}\left(\frac{\phi }{z}\right)^2+\frac{I\mathrm{\Delta }n}{2c}(1e^2)^{1/2}\mathrm{cos}2(\psi \phi )=E.$$
(6)
There is also another conserved integral which can be derived from this free energy. This is the total angular momentum flux along the propagation direction, and is the sum of the elastic and optical angular momenta:

$$K_2\frac{\phi }{z}+\frac{I}{\omega }e=M.$$
(7)
Eliminating the director angle in eq.(6), using eq.(4) and (7), yields the following equation for the polarization ellipse:

$$\left(\frac{e}{z}\right)^2=\left(\frac{\omega \mathrm{\Delta }n}{c}\right)^2\left(1e^2\right)\left[E\frac{1}{2K_2}\left(M\frac{Ie}{\omega }\right)^2\right]^2\left(\frac{2\omega }{I}\right)^2.$$
(8)
This is the fundamental equation governing the distribution of ellipticity in the liquid crystal cell. Santamato et. al Santamato1988 , then used it to discuss the optical Fredericksz transition in a planar cell. They observed, in addition, that this equation can in general be solved in terms of elliptic integrals.

We now turn to the application of the equation (8) to our geometry. We note that the physical structure of the general solutions strongly depends on the constants $`M,`$ $`E`$, which in turn depend on boundary conditions in a rather complicated way. However, for the cell geometry in which the slippage effect has been observed, considerable simplifications are possible.

The boundary conditions involve a fixed director orientation at $`z=L,`$ and a free but planar boundary at $`z=0`$:

$$\phi |_{z=L}=0,\frac{\phi }{z}|_{z=0}=0.$$
(9)
Note, that the director can slip over the surface with weak anchoring – no surface torque prevents it from the reorientation. As a result, the weak anchoring maximizes the response of the director to a bulk field.

The incident light wave is elliptically polarized and initially incident on the surface with zero anchoring:

$$e|_{z=0}=e_0,\psi |_{z=0}=\alpha .$$
(10)
In this geometry the constants $`E`$ and $`M`$ in eq.(8) can be easily determined from the boundary conditions:

$`M`$ $`=`$ $`\left(I/\omega \right)e_0,`$ (11)
$`E`$ $`=`$ $`{\displaystyle \frac{I}{c}}\left[n_o+{\displaystyle \frac{\mathrm{\Delta }n}{2}}\left[1+\left(1e_0^2\right)^{1/2}\mathrm{cos}2\left(\alpha \phi _0\right)\right]\right],`$
where $`\phi _0=\phi (z=0)`$. Then, eq.(8) for the ellipticity and the director (7) simplify to

$$\left(\frac{e}{s}\right)^2=1e^2\left[\left(1e_0^2\right)^{1/2}\mathrm{cos}2\left(\alpha \phi _0\right)r\left(e_0e\right)^2\right]^2,$$
(12)
$$\frac{\phi }{s}=r\left(e_0e\right),$$
(13)
where $`s=\left(\omega \mathrm{\Delta }n/c\right)z=2\pi \mathrm{\Delta }n\left(z/\lambda \right)`$. We have also introduced the dimensionless parameter $`r`$, which is proportional to the light intensity and defined by:

$$r=Ic/\left(\mathrm{\Delta }nK_{22}\omega ^2\right).$$
(14)
We now consider the two distinct cases, of linearly and circularly polarized incident light waves.

### II.1 Linear polarization of the incident light

Zero ellipticity of the linearly polarized light implies that the light wave does not carry angular momentum, since the latter is proportional to the ellipticity.

The initial conditions for the polarization state $`e`$ read:

$$e_0=0,\psi |_{z=0}=\alpha ,$$
(15)
and the equation for the ellipticity (12) simplifies to

$$\left(\frac{e}{s}\right)^2=1e^2\left[\mathrm{cos}2\left(\alpha \phi _0\right)re^2\right]^2.$$
(16)
Substitution $`y=re^2`$ allows (16) the standard integral for the elliptic functions Gradshtein1994 . After integration, we obtain

$$e_{\mathrm{lin}}(s)=e_0sd\left(\kappa s\right)$$
(17)
where $`e_0=m^{}\sqrt{e_p/r},`$ $`\kappa =\sqrt{e_pr}/m`$, $`m^{}=\sqrt{1m^2}`$, $`m^2=e_p/\left(e_pe_n\right)`$, $`e_{p,n}`$ are the positive and negative roots of the equation $`1e^2\left[\mathrm{cos}2\left(\alpha \phi _0\right)re^2\right]^2=0`$ correspondingly, $`sd(x)=sn(x)/dn(x)`$, where $`sn`$, $`cn`$, $`dn`$ are Jacobi elliptic functions of index $`m`$.

The director distribution can be obtained by integrating (7) and has the form

$$\phi _{\mathrm{lin}}(s)=arccot\left(\frac{m}{m^{}}cn\left[\kappa s_L\right]\right)arccot\left(\frac{m}{m^{}}cn\left[\kappa s\right]\right)$$
(18)
where $`s_L=\left(\omega \mathrm{\Delta }n/c\right)L`$.

The estimation of the typical experimental values for $`r`$ will be done later. We now find the solution for $`r<<1`$.

As it is seen from the angular momentum conservation (7), to obtain the director distribution in linear order in the dimensionless light intensity $`r`$ we need to know the ellipticity $`e`$ in zero order in $`r`$. The latter can be found by solving eq.(16) putting $`r=0`$:

$$e_{\mathrm{lin}}(s)=\mathrm{sin}\left(2\alpha \right)\mathrm{sin}\left(s\right).$$
(19)
Then the equation for the director (7) can be solved using elementary methods, yielding:

$$\phi _{\mathrm{lin}}(0)=2r\mathrm{sin}\left(2\alpha \right)\mathrm{sin}^2\left(s_L/2\right).$$
(20)
The ellipticity distribution and the director deviation given in eqs.(19,20) might in principle have been obtained by assuming a priori that the ordinary and extraordinary light waves simply follow the distribution of the director. Equivalently, one can say that in the low $`r`$ limit the light beam propagates through the cell in the adiabatic or Mauguin regime deGennes1993 . In the adiabatic approach one can write the expression for the electric field components explicitly since the ordinary and extraordinary light waves follow the distribution of the director. Then this field is used as an external field which reorients the director. However, such an assumption is not theoretically satisfactory, and naive use of this approximation as an initial hypothesis can lead to internally inconsistent approximations Ong1983a .

Let us now summarize the important properties of the light propagation and the director distribution in the planar cell. Eq.(19) yields that the polarization ellipse $`e(z)`$ is a periodic function with spatial period $`l=\lambda /\mathrm{\Delta }n.`$ In more general case (17) the spatial period depends on the light intensity, $`l=4K\left(m\right)\left(\kappa \omega \mathrm{\Delta }n/c\right)^1`$, where $`K\left(m\right)`$ is the complete elliptic integral of the first kind.

The liquid crystal director $`\phi (z)`$ follows the distribution of the ellipticity with the same spatial period, and is shown in Fig. 1. It is seen that, in spite of the modulation, the average deviation of the director tends to be along the polarization of the incident light and increases with the light intensity.

We now make some comments on the director deviation $`\phi _0`$ at $`z=0`$. First, eq. (20) shows that the surface director is sensitive to the sign of the birefringence $`\mathrm{\Delta }n`$. If $`\mathrm{\Delta }n>0`$, the director tends to be parallel to the polarization direction. It reorients perpendicular to the polarization direction if $`\mathrm{\Delta }n<0`$.

Second, in the limit of small intensities, $`r<<0`$, the liquid crystal response is maximal at $`\alpha =\pi /4`$, and does not occur for $`\alpha =0,\pi /2`$ since we are below the threshold of Freederiksz transition. However, in the general case, the value of $`\alpha `$ which provides maximal director response depends on the light intensity. For bigger values of the parameter $`r`$ the maximal response of the director occurs at an angle $`\alpha `$ larger than $`\pi /4`$ (Fig. 2).

The most important comment is that the director deviation depends strongly on the cell thickness $`L`$ and refractive index difference $`\mathrm{\Delta }n`$. Indeed, $`\omega \mathrm{\Delta }nL/\left(2c\right)=\pi \mathrm{\Delta }nL/\lambda ,`$ and for typical experimental conditions $`L/\lambda 100`$. Thus, a change in the refractive index $`\mathrm{\Delta }n10^2`$ leads to a considerable change in the director deviation $`\phi _0`$. In practice this means that changes in the refractive index due to the laser-induced heating of the liquid crystal Khoo1988  or molecular phototransformations Pinkevich1992  influence the amplitude of the surface director deviation. This also means that the cell thickness should be tuned rather precisely in order for the surface director orientation to be observed. The optimal condition is that $`\mathrm{\Delta }nL/\lambda =N/2`$, where $`N`$ is an integer.

### II.2 Circular polarization of the incident light

Contrary to the case of linear polarization, circularly polarized light wave carries angular momentum along the propagation direction. The initial conditions for the polarization state $`e`$ are now:

$$e_{z=0}=1.$$
(21)
The boundary conditions for the director are the same as for the linearly polarized light (9).

The equation for the ellipticity then simplifies to

$$\left(\frac{e}{s}\right)^2=1e^2r^2\left[1e\right]^4.$$
(22)
Substitution $`y=1e`$ allows integration of this equation in terms of elliptical integrals Gradshtein1994 . After integration we get the distribution of the ellipticity in the cell:

$$e_{\mathrm{cir}}(s)=1e_r\frac{q\left[1cn\left(\kappa s\right)\right]}{p+q+\left(pq\right)\left[1cn\left(\kappa s\right)\right]}$$
(23)
where $`e_r`$ is the only real root of the equation $`2er^2e^3=0,`$ $`p^2=r^2+3e_r^2,`$ $`q^2=r^2+e_r^2,`$ $`\kappa =r\sqrt{pq}`$. The polarization ellipse $`e_{\mathrm{cir}}(s)`$ has the same dependence on the spatial coordinate as we have in the case of linear polarization (17 ), i.e. is a periodic function with the same spatial period.

Now we turn to the director distribution. Integrating Eq. (7) we obtain

$$\phi \left(s\right)=\phi _0+\frac{e_r}{pq}\left[2\mathrm{arctan}\left(\sqrt{\frac{q}{p}}\mathrm{tan}\left(\frac{\kappa s}{2}\right)\right)\sqrt{\frac{q}{p}}\kappa s\right].$$
(24)
To perform the integration we used an approximation $`cn(x)\mathrm{cos}(x)`$, which is valid for small enough elliptic function indices or alternatively small enough values of $`r`$.

The expressions (23,24) in leading order in $`r`$ are:

$$e_{\mathrm{cir}}\left(s\right)=\mathrm{cos}\left(s\right),$$
(25)
$$\phi _{\mathrm{cir}}(0)=r\left[s_L\mathrm{sin}\left(s_L\right)\right].$$
(26)
Let us now analyze the ellipticity and the director profile in the cell. Eq. (25) yields that polarization ellipse is exactly the same as for the linearly polarized light. The only change is that, in order to match the condition $`e_{\mathrm{cir}}\left(0\right)=1`$ instead of $`e_{\mathrm{lin}}\left(0\right)=0`$, we have to have cosine dependence on $`s`$ (cf. exp.(19)).

However, the liquid crystal director is no longer a periodic function of the cell thickness. Even more, now it is monotonically growing function of $`s`$ (Fig. 3). This fact can be easily proved if one considers the angular momentum conservation (7). For the linear polarization it reads: $`\phi /s=re`$. Since $`e\left(s\right)`$ is a periodic function of $`s`$, $`\phi `$ is also a periodic function of $`s`$ with the same spatial period. However, for circular polarization we have $`\phi /s=r(1e)0`$. Therefore, $`\phi `$ grows monotonically with $`s`$. In fact, the director distribution is a superposition of monotonically growing and oscillating functions, since the ellipticity of the polarization ellipse is a periodic function of $`s`$.

It is also seen from eq. (26) that the director deviation $`\phi _0`$ at $`z=0`$ consists of two parts. The first part, $`r\mathrm{sin}\left(s_L\right)`$, originates from the ellipticity modulation. The other part, $`rs_L`$, is due to the constant injection of light angular momentum into the cell bulk. For typical experimental conditions this contribution dominates, since $`s_L=2\pi \mathrm{\Delta }nL/\lambda >>1>\mathrm{sin}\left(s_L\right)`$ . For the same reasons $`\phi _0^{\mathrm{cir}}2\pi \mathrm{\Delta }nL/\lambda \phi _0^{\mathrm{lin}}`$. Therefore, the surface director deviation in circularly polarized light is $`2\pi \mathrm{\Delta }nL/\lambda `$ times bigger than in the linearly polarized light.

We can make even more interesting remarks.

First, the resulting twist angle $`\phi _0`$ is not sensitive to the variation of $`\mathrm{\Delta }n`$ for typical experimental conditions ($`r<<1`$) and liquid crystal cell parameters ($`s_L>>1`$). This is very different from the situation with linearly polarized light wave. Second, $`\phi _0`$ is proportional to the cell thickness $`L`$, i.e., in thick enough cells, its value can exceed $`\pi `$ and one can get supertwisted distribution of the director or even distribution that realizes in chiral nematics.

### II.3 Elliptically polarized light

Now we try to solve the eqs. (12,13) for the case of elliptically polarized light. Solving these equations for arbitrary value of $`r`$ is a complicated task. However, as one can see from the angular momentum conservation (7), to obtain the director distribution in linear order in the dimensionless light intensity $`r`$ we need to know the ellipticity $`e\left(s\right)`$ in zero order in $`r`$. Solving eq.(12) in this limit we obtain that $`e\left(s\right)`$ is a superposition of the solutions $`e_{\mathrm{lin}}\left(s\right)`$ and $`e_{\mathrm{cir}}\left(s\right)`$:

$$e\left(s\right)=\sqrt{1e_0^2}\mathrm{sin}\left(2\alpha \right)\mathrm{sin}s+e_0\mathrm{cos}s=\sqrt{1e_0^2}e_{\mathrm{lin}}\left(s\right)+e_0e_{\mathrm{cir}}\left(s\right)$$
(27)
Correspondingly, in the leading order in $`r`$, the director deviation at $`z=0`$ is also a superposition of $`\phi _0^{\mathrm{lin}}`$ and $`\phi _0^{\mathrm{cir}}`$ and has the form:

$$\phi _0=r\left[\sqrt{1e_0^2}\mathrm{sin}\left(2\alpha \right)\left(1\mathrm{cos}s_L\right)+e_0\left(\mathrm{sin}s_Ls_L\right)\right]=\sqrt{1e_0^2}\phi _0^{\mathrm{lin}}+e_0\phi _0^{\mathrm{cir}}$$
(28)
It is seen from Eq. (28) that, for circularly polarized light wave, the director at the surface rotates according to the sign of the ellipticity: clockwise polarization leads to the clockwise reorientation of the director and vice versa.

## III Discussion

We now turn to the estimation of the parameter $`r`$ and the director deviation angles one can expect for typical experimental conditions.

Quantitatively, in typical experiments, we have light intensities $`I<10W/cm^2`$, and liquid crystal constants $`\mathrm{\Delta }n=0.1`$, $`K_{22}=310^7dyn`$. With these fundamental parameter values, the dimensionless intensity $`r<0.001<<1`$, which is fairly small. The maximal deviation of the director which can be observed under these circumstances is far less than $`0.1^0`$ for linearly polarized light and is less than $`1^0`$ for circularly polarized light. Such reorientations can hardly be detected in real experiments.

However, the situation is different for the dye-doped liquid crystals. All previous experimental results Janossy1990 ; Janossy1991  and theoretical predictions Janossy1994 ; Marrucci1997  are consistent with the following expression for the dye contribution to the free energy density:

$$f_{\mathrm{dye}}=\eta f_{\mathrm{opt}},$$
(29)
where $`f_{\mathrm{opt}}`$ is the total electromagnetic energy density of the incident light wave (2) and the parameter $`\eta `$ characterizes the efficiency of the dye-induced director reorientation. The parameter $`\eta `$ is proportional to the dye concentration, but also depends on the molecular structures of both dye and liquid crystal.

Eq.(29) implies that the effective torque imposed by the light wave increases proportionally to $`\eta `$. In typical experiments, we obtain a dye assisted nonlinearity coefficient $`\left|\eta \right|200`$ Janossy1992 . This gives $`r0.1`$ and typical director deviations about $`10^0`$, which can be measured experimentally without any difficulties.

Our theoretical framework also permits the evaluation of the photo-induced orientational nonlinearity coefficient $`\eta `$ from the existing experimental data. $`\eta `$ is the crucial phenomenological parameter which governs the interaction between light and the liquid crystal/azo-dye mixture. This quantity is an input parameter in our theory, and must either be measured or calculated using a microscopic theory Janossy1994 ; Marrucci1997 .

Trying to estimate $`\eta `$ we are fully aware that the cell is absorbing, since the effect of the photoinduced nonlinearity is due to the absorption of light by the dye molecules. For the absorbing media, we are no longer able to apply the Lagrange approach, which involves energy and angular momentum conservation. Therefore, estimating coefficient $`\eta `$ we are neglecting effects related to the light absorption.

For the linearly polarized incident light, the value of $`\eta `$ can be determined applying eq.(20), which connects $`\phi _0`$ to the dimensionless intensity $`r`$, to the experimental results presented in Marusii1996 .

Typical experiments were performed in the cell with one isotropic surface and one surface consisting of a substrate covered with a rubbed polyimide film. The cell was filled with the nematic liquid crystal 4-n-pentyl-4-cyanobiphenyl (5CB) doped with the azo dye methyl red (MR) at weight concentration $`0.1`$. The mixture was exposed from the side of the isotropic surface by a linearly polarized beam of a cw Ar<sup>+</sup> laser, irradiation wavelength of which $`\lambda =488\mathrm{n}\mathrm{m}`$ is near the maximum of the absorption band of MR dissolved in 5CB. Typical intensity of the Ar<sup>+</sup> laser was in the range $`0.11\mathrm{W}/\mathrm{cm}^2`$. Its polarization vector was set at $`45^o`$ with respect to the initial director orientation. The cell response was observed with a probe He-Ne laser beam with wavelength $`\lambda =638\mathrm{n}\mathrm{m}`$ which is in the region of MR transparency.

Irradiation with a beam of a Ar<sup>+</sup> laser led to the director reorientation in the liquid crystal bulk and its slippage on the isotropic surface. The latter was detected as a rotation of the polarization of the probe beam. The experimental dependence of the rotation angle on the intensity of the incident light is presented in Fig.4.

Eq. (20) predicts that the surface director reorientation angle $`\phi _0`$ will be a function of the laser intensity $`I`$ and refractive index anisotropy $`\mathrm{\Delta }n=n_en_o`$. The birefringence of the dye-liquid crystal mixture can in principle also depend on the exciting light intensity indirectly both as a result of laser-induced heating of the liquid crystal Khoo1988 , and because molecular phototransformation can change the dipole moment of the phototransformed molecules Pinkevich1992 . The intensity $`I`$ thus has a direct and an indirect effect on the reorientation of the director. In order to interpret the results correctly it is necessary to correct for the indirect effects. To do this we studied the dependence of the dye-liquid crystal birefringence $`\mathrm{\Delta }n`$ on the exciting light intensity independently using a light-induced birefringence technique Pinkevich1992 . We found that in our case:

$$\mathrm{\Delta }n=\mathrm{\Delta }n_0+n_2I,$$
(30)
where $`n_20.0026\left(\mathrm{W}/\mathrm{cm}^2\right)^1`$.

Then, fitting of the data (Fig. 4) yields $`\eta =560`$. The photo-induced orientational nonlinearity is considerably enhanced even with respect to the giant orientational non-linearity of pure liquid crystals Zel'dovich1982  and is consistent with what might be expected as a result of measurements in related systems Zolot'ko1994 ; JanossyPRE1998  using the $`z`$-scan technique.

The theory developed in eq.(20) also requires that the director response to the laser illumination light intensity be linear. In order to check whether this is the case, we measured the dependence $`\phi _0\left(\alpha \right)`$. The results of measurements for $`I=1\mathrm{W}/\mathrm{cm}^2`$ are presented in Fig.5. We find, in agreement with eq. (20), that the liquid crystal response is maximal at $`\alpha =`$$`\pi /4`$ and absent at $`\alpha =0,`$$`\pi /2`$.

As far as we are aware, there have been no direct observations of the slippage effect for circularly or elliptically polarized light. The effect has probably been implicitly observed in azo-dye doped liquid crystal cells AndrienkoMCLC1999 ; AndrienkoUFJ1999 . In these experimetns it was found that when a cell is irradiated by circularly polarized light, an easy orientation axis is induced on an initially isotropic aligning surface. The direction of this axis correlated with the direction of the polarization rotation. The appearance of the easy axis can be explained as follows. First there is a light-induced surface director slippage. The director then freezes as a result of light-induced adsorption of the dye molecules on the aligning surface. By contrast, this reorientation cannot be explained if we use a model which depends only on the selective adsorption of the dye molecules on the aligning surface.

## IV Conclusions

In this paper we have tried to interpret existing experimental facts on the surface director reorientation by assuming the dye-induced torque to be proportional to the optical torque, and treating the nonlinear interaction between the light wave and the mixture in the geometrical optics approximation.

We now turn to a brief discussion of the implications of this work. We believe that the experimental technique we introduce in this paper can be used in further studies of the liquid crystal interface. Relatively straightforward extensions of the technique will, for example, permit measurements of fundamental surface parameters such as in-plane surface viscosities and weak anchoring coefficients.

There are also possible applications of the method involving the writing of the dynamic holographic polarization gratings. Furthermore, the present investigation, together with our results published earlier Francescangeli1999 ; Voloshchenko1994 ; Simoni1997 ; Slussarenko1997  have shown that it is possible to control the light-induced alignment memory effect by changing the concentration of azo-dye in the liquid crystal. There may also be scope for systematic investigations of surface memory effects, the basic mechanisms of which are still little understood.

###### Acknowledgements.

D. A. acknowledges support by EPSRC and support through the Grant No. ORS/99007015 of the Overseas Research Students Awards scheme.


format: text

text:
# Introduction

## Introduction

In  Kauffman defines virtual knot diagrams which are, in some sense, a natural extension of classical knot diagrams. Surprisingly, there can be found quite easily examples of virtual knots with properties that are unknown for classical knots, e.g., a knot with trivial Jones polynomial, see .

In this article, Vassiliev invariants for virtual knots and links are investigated. They arise from the Conway polynomial for virtual links defined in . It is shown that the lowest order coefficient of the Conway polynomial vanishes on virtual knots and, furthermore, that it is a zeroth order Vassiliev invariant for virtual links which cannot detect a change of orientation. In contrast, two examples are given that the corresponding Vassiliev invariant of order one distinguishes a virtual knot from its inverse. For Vassiliev invariants of classical knots it is still unknown whether they are orientation-sensitive or not.

## 1 Virtual Links and Vassiliev Invariants

In the following, standard terminology from classical knot theory will be used and some definitions that can be extended to virtual links in the obvious way will not be given explicitly.

Definition  A virtual link diagram is an oriented 4-valent planar graph embedded in the plane with appropriate orientations of edges and additional crossing information at each vertex (see Fig. 1).

Denote the set of virtual link diagrams by $`𝒱𝒟`$. Two diagrams $`D,D^{}𝒱𝒟`$ are called equivalent if one can be transformed into the other by a finite sequence of extended Reidemeister moves (see Fig. 2)

combined with orientation preserving homeomorphisms of the plane to itself. A virtual link is an equivalence class of virtual link diagrams. A virtual link with one component is called virtual knot.

Though the notion of components of a virtual link, as in the previous definition, comes from the idea to image the virtual link in 3-space as if the virtual crossings were classical ones, a virtual link diagram does not correspond to an object in 3-space in the same way as in classical knot theory, see  for more details. Nevertheless there exist several geometric interpretations for virtual links, see , , , .

In the main part of this article, Vassiliev invariants for virtual links, as introduced in  (and not as in ), will be investigated. For an introduction to classical Vassiliev theory, see  or .

Definition  A singular virtual link diagram is a virtual link diagram that may contain vertices of an additional type called double points, see Fig. 3.

An equivalence relation is defined by adding the rigid vertex moves depicted in Fig. 4

to the extended Reidemeister moves for virtual link diagrams.

In the following, often diagrams are considered which are identical except within a small disk where they differ as depicted in Fig. 5.

If $`D`$ is a singular virtual link diagram with a chosen double point then replacing the double point with a positive or negative crossing yields its positive and negative resolution, $`D_+`$ and $`D_{}`$, respectively. Smoothing the double point yields the corresponding diagram $`D_0`$. Likewise, diagrams $`D_{\epsilon _1\mathrm{}\epsilon _k}`$ with $`\epsilon _i\{+,,0\}`$ can be defined where $`k`$ chosen double points are replaced. The same notation is used when crossings instead of double points are chosen and replaced.

Definition  A Vassiliev invariant $`v`$ is an invariant of singular virtual link diagrams with values in an abelian group such that

$$v(D_{})=v(D_+)v(D_{})\text{(Vassiliev relation)}$$

where $`D_{}`$ denotes a diagram with a chosen double point and $`D_+`$ and $`D_{}`$ its positive and negative resolution, respectively. $`v`$ is said to be of order $`n`$ if it vanishes on diagrams with more than $`n`$ double points and it is said to be of (exact) order $`n`$ if it is of order $`n`$ but not of order $`n1`$.

Remark  Let $`v`$ denote an invariant of virtual links with values in an abelian group. $`v`$ always can be extended to an invariant of singular virtual links by demanding that the Vassiliev relation is fulfilled. $`v`$ is of order zero iff $`v(D_+)=v(D_{})`$ holds for every virtual link diagram $`D`$ and every crossing of $`D`$. $`v`$ is of order $`1`$ iff the equation

$$v(D_{++})v(D_+)v(D_+)+v(D_{})=\mathrm{\hspace{0.33em}0}$$

holds for every virtual link diagram $`D`$ and every pair of crossings of $`D`$.

## 2 Conway Polynomial and its Coefficients

The Conway polynomial for virtual links is derived from the normalized $`Z`$-polynomial defined in  which is an adaption of the Conway polynomial for links in thickened surfaces that has been introduced by Jaeger, Kauffman, and Saleur in . In the following, the construction of the $`Z`$-polynomial is recapitulated.

Let $`D`$ be a virtual link diagram with $`n1`$ classical crossings $`c_1,\mathrm{},c_n`$. Define

$$M_+:=\left(\begin{array}{cc}1x& y\\ xy^1& 0\end{array}\right),M_{}:=\left(\begin{array}{cc}0& x^1y\\ y^1& 1x^1\end{array}\right)=M_+^1$$

and $`M:=diag(M_+^{\epsilon _1},\mathrm{},M_+^{\epsilon _n})`$ where $`\epsilon _i=\pm 1`$ denotes the sign of $`c_i`$. Consider the graph belonging to the virtual link diagram where the virtual crossings are ignored, i.e., the graph consists of $`n`$ vertices $`v_1,\mathrm{},v_n`$ corresponding to the classical crossings and $`2n`$ edges corresponding to the arcs connecting two classical crossings. Subdivide each edge into two half-edges and label the four half-edges belonging to the vertex $`v_i`$ as depicted in Fig. 6.

The assignment

$$(i,a)(j,b)\text{if }i_a^+\text{ and }j_b^{}\text{ belong to the same edge}$$

gives a permutation of the set $`\{1,\mathrm{},n\}\times \{l,r\}`$. Let $`P`$ denote the corresponding $`2n\times 2n`$ permutation matrix where rows and columns are enumerated $`(1,l)`$, $`(1,r)`$, $`(2,l)`$, $`(2,r)`$, …, $`(n,l)`$, $`(n,r)`$, i.e., the $`(i,a)`$-th column of $`P`$ is the $`(j,b)`$-th unit vector.

Define $`Z_D(x,y):=(1)^ndet(MP)`$ (observe that this yields the same polynomial as the definition in ). If $`D`$ has no classical crossings then $`Z_D(x,y):=0`$. See  for example calculations of $`Z_D`$.

###### Theorem 1

$`Z:𝒱𝒟[x^{\pm 1},y^{\pm 1}]`$ is an invariant of virtual links up to multiplication by powers of $`x^{\pm 1}`$. To be more precise, $`Z`$ is invariant with respect to all virtual Reidemeister moves except moves of type I. Creating an additional (classical) crossing by applying a Reidemeister move of type I yields a change of $`Z`$ by a factor as shown in Fig. 7.

Proof:  see

$`\mathrm{}`$

Remark  Let $`T:=diag(\left(\begin{array}{cc}0& 1\\ 1& 0\end{array}\right),\mathrm{},\left(\begin{array}{cc}0& 1\\ 1& 0\end{array}\right))`$. The permutation corresponding to the matrix $`TP`$ can be read off a virtual link diagram by walking along each component of the link and writing down the numbers of the outgoing arcs. This gives a cycle $`\sigma _i`$ for each link component and $`\sigma _1\mathrm{}\sigma _r`$ is the desired permutation. This describes an alternative method to determine the matrix $`P`$.

Define the normalized polynomial $`\stackrel{~}{Z}_D(x,y)`$ as follows. If $`Z_D(x,y)`$ is a non-vanishing polynomial and $`N`$ is the lowest exponent in the variable $`x`$ then define

$$\stackrel{~}{Z}_D(x,y):=x^NZ_D(x,y).$$

Otherwise let $`\stackrel{~}{Z}_D(x,y):=Z_D(x,y)=0`$.

###### Corollary 2

$`\stackrel{~}{Z}:𝒱𝒟[x,y^{\pm 1}]`$ is an invariant of virtual links.

$`\mathrm{}`$

###### Theorem 3

Let $`D,D^{}`$ denote virtual link diagrams and let $`(D_+,D_{},D_0)`$ be a skein triple of virtual link diagrams. Then the following holds.

1. $`Z_D(x,y)=0`$ if $`D`$ has no virtual crossings
2. $`Z_{DD^{}}(x,y)=Z_D(x,y)Z_D^{}(x,y)`$ (disjoint union)
3. $`Z_D(1,y)`$ does not depend on the over-under information of the diagram’s classical crossings
4. $`Z_{D_+}(x,y)xZ_D_{}(x,y)=(1x)Z_{D_0}(x,y)`$ (skein relation of Conway type)

Proof:  see

$`\mathrm{}`$

Remark  Obviously, the parts a), b), c) of Theorem 3 are valid for the normalized polynomial $`\stackrel{~}{Z}`$, too. Part d) is, in general, not true for $`\stackrel{~}{Z}`$ instead of $`Z`$ (simply consider a standard diagram of a Hopf link with one classical crossing replaced by a virtual crossing).

Define the Conway polynomial $`C_D(y,z)[y^{\pm 1},z]`$ by expanding the normalized $`Z`$-polynomial in a Taylor series

$$\stackrel{~}{Z}_D(x,y)=\underset{k}{}c_k(1x)^k$$

and setting $`z:=1x`$. Then $`C_D(y,z)`$ is an invariant of virtual links which satisfies the Conway skein relation

$$C_{D_+}(y,z)C_D_{}(y,z)=zC_{D_0}(y,z)$$

”up to powers of $`x=1z`$”, i.e., there exist integers $`k_+,k_{},k_0`$ such that

$$(1z)^{k_+}C_{D_+}(y,z)(1z)^k_{}C_D_{}(y,z)=z(1z)^{k_0}C_{D_0}(y,z)$$

holds. Extend the Conway polynomial to singular virtual links via the Vassiliev relation.

When the classical Conway polynomial $`_D(z)`$ is extended to (classical) singular link diagrams by

$$_D_{}(z):=_{D_+}(z)_D_{}(z)$$

then it can easily be seen that the coefficient $`c_k`$ of $`_D(z)=c_kz^k`$ is a Vassiliev invariant of order $`k`$ since

$$_D_{}(z)=z_{D_0}(z)$$

by the Conway skein relation. This does not work as easily for the Conway polynomial $`C_D(y,z)`$, see . For the purpose of this article, it is enough to consider the Vassiliev invariants of orders zero and one which are easy to handle.

###### Lemma 4

Let $`P`$ denote the permutation matrix from the definition of $`Z_D`$, and let $`T=diag(\left(\begin{array}{cc}0& 1\\ 1& 0\end{array}\right),\mathrm{},\left(\begin{array}{cc}0& 1\\ 1& 0\end{array}\right))`$. Then

$$c_0=det(diag(y^1,y,\mathrm{},y^1,y)+TP).$$

Proof:  Setting $`x:=1`$ in the definition of $`Z_D`$ gives:

$$\begin{array}{ccc}c_0\hfill & =\hfill & (1)^ndet\left(diag(\left(\begin{array}{cc}0& y\\ y^1& 0\end{array}\right),\mathrm{},\left(\begin{array}{cc}0& y\\ y^1& 0\end{array}\right))P\right)\hfill \\ & =\hfill & (1)^ndet(T)det(diag(y^1,y,\mathrm{},y^1,y)TP)\hfill \\ & =\hfill & det(diag(y^1,y,\mathrm{},y^1,y)+TP)\hfill \end{array}$$

$`\mathrm{}`$

###### Corollary 5

$`c_0`$ is invariant with respect to a change of orientation.

Proof:  A change of orientation in a virtual link diagram corresponds to replacing the permutation matrix $`TP`$ by the transposed matrix $`(TP)^T`$ which does not change the value of $`c_0`$ by Lemma 4.

$`\mathrm{}`$

In  D. Silver and S. Williams show that their polynomial invariant $`\mathrm{\Delta }_0`$ is related to the $`Z`$-polynomial via an appropriate change of variables (Proposition 3.1) and, furthermore, they investigate the behaviour of $`\mathrm{\Delta }_0`$ with respect to a change of orientation (Corollary 5.2). Combining these results with Corollary 5 immediately yields the following relation.

###### Corollary 6

For a virtual link with $`d`$ components,

$$c_0(y)=(1)^dc_0(y^1)$$

where $`c_0(y^1)`$ denotes $`c_0(y)=c_0`$ with the variable $`y`$ replaced by $`y^1`$.

$`\mathrm{}`$

###### Theorem 7

$`c_0=0`$ for virtual knots.

Proof:  The permutation corresponding to the matrix $`TP`$ can be read off a virtual link diagram by walking along each component of the link and writing down the numbers of the outgoing arcs. This gives a cycle for each component. Therefore, if $`D`$ represents a virtual knot, then $`TP`$ is the permutation matrix of a cycle $`\sigma `$ and Lemma 4 gives:

$$\begin{array}{ccc}c_0\hfill & =\hfill & det(diag(y^1,y,\mathrm{},y^1,y)+TP)\hfill \\ & =\hfill & sign(id)(y^1y)\mathrm{}(y^1y)+sign(\sigma )1\mathrm{}1\hfill \\ & =\hfill & 1+(1)^{2n1}\hfill \\ & =\hfill & 0\hfill \end{array}$$

$`\mathrm{}`$

###### Theorem 8

$`c_0`$ is a Vassiliev invariant for virtual links, and $`c_1`$ is a Vassiliev invariant for virtual knots.

Proof:  In the following, let $`c_i^\epsilon `$ with $`\epsilon \{+,,0\}`$ be an abbreviation for $`c_i(D_\epsilon )`$ and likewise $`c_i^{\epsilon _1\epsilon _2}=c_i(D_{\epsilon _1\epsilon _2})`$ for a given virtual link diagram $`D`$ and crossings chosen.

Since $`c_0=Z_D(1,y)`$ it is clear from Theorem 3 that $`c_0`$ is a Vassiliev invariant of order zero, i.e., $`c_0^+=c_0^{}`$.

For a virtual knot, $`c_0`$ vanishes by Theorem 7 and therefore a comparison of coefficients in the Conway skein relation

$$(c_0^++(c_1^+k^+c_0^+)z+z^2(\mathrm{}))(c_0^{}+(c_1^{}k^{}c_0^{})z+z^2(\mathrm{}))=c_0^0z+z^2(\mathrm{})$$

immediately yields $`c_1^+c_1^{}=c_0^0`$. Thus

$$c_1^{++}c_1^+c_1^++c_1^{}=c_0^{+0}c_0^0=\mathrm{\hspace{0.33em}0}$$

because $`c_0`$ is of order zero. This shows that $`c_1`$ is a Vassiliev invariant of order one.

$`\mathrm{}`$

Remark  Though formulated only for virtual links, Corollaries 5, 6 and Theorems 7, 8 are valid for singular virtual links, too, because of the Vassiliev relation.

## 3 Examples

Example  In general, $`c_1`$ is not a Vassiliev invariant of order one for virtual links with more than one component. A counterexample is depicted in Fig. 8.

$`c_1`$ does not vanish for the singular virtual link diagram $`D`$ shown in Fig. 8 though $`D`$ has two double points:

$$c_1(D)=y+2+y^1$$

Indeed, it is easy to construct examples which show that $`c_1`$ is not a Vassiliev invariant of any order.

Example  In contrast to $`c_0`$, $`c_1`$ is orientation-sensitive. An example is depicted in Fig. 9.

Let $`D`$ be the diagram with the orientation indicated in Fig. 9 and let $`D^{}`$ denote the diagram with the opposite orientation. Then

$$c_1(D)=y^2y+1+y^1\text{and}c_1(D^{})=y+1y^1y^2.$$

Chirality is detected, too. Let $`\overline{D}`$ denote the mirror diagram of $`D`$, i.e., every classical crossing is changed from positive to negative and vice versa. Then

$$c_1(\overline{D})=y1+y^1+y^2\text{and}c_1(\overline{D}^{})=y^2+y1y^1.$$

The Jones polynomial of $`D`$ (see ) is non-trivial and it takes different values on $`D`$ and $`\overline{D}`$. But from the definition via the bracket polynomial it is clear that the Jones polynomial cannot distinguish a virtual link from its inverse.

Example  Let $`D`$ be the diagram with the orientation indicated in Fig. 10

and let $`D^{}`$ denote the diagram with the opposite orientation. $`D`$ is Kauffman’s example of a virtual knot with trivial Jones polynomial and trivial knot group (see ). Again, $`c_1`$ distinguishes $`D`$ and $`D^{}`$. It happens that the values are the same as in the previous example:

$$c_1(D)=y^2y+1+y^1\text{and}c_1(D^{})=y+1y^1y^2.$$

Remark  In the same way as for classical links, quantum invariants can be defined for virtual links, see . It is well-known (and can be shown analogously for virtual links) that quantum invariants cannot detect a change of orientation (see  or , for example) and therefore $`c_1`$ is a Vassiliev invariant of virtual knots which is not a function of quantum invariants. By a result of Graña (), this implies that $`c_1`$ neither is a function of quandle cocycle invariants as defined in .

## Acknowledgement

The author would like to thank Dan Silver for valuable discussions.


text:
# Diffraction-limited CCD imaging with faint reference stars

## 1 Introduction

In recent observations with the Nordic Optical Telescope (NOT) we demonstrated that selection of the best images from a large dataset of short exposures can provide essentially diffraction-limited images at wavelengths shorter than 1 $`\mu `$m from well-figured ground-based telescopes as large as 2.5 m diameter (Baldwin et al (2001)). Images of bright stars derived from the best 1% of exposures of 5 ms duration under good seeing had Strehl ratios as high as 0.3 at 810 nm. The use of the brightest speckle for image selection, shifting and co-adding provides a significant improvement in resolution over recentring and selection techniques based upon the image centroid (see Nieto & Thouvenot (1991) and references therein).

The applicability of this technique to faint astronomical targets depends on the probability that a sufficiently bright reference star lies within the isoplanatic patch prevailing at the times of the selected exposures. This probability is expected to be very much higher for selected exposures than for adaptive optics (AO) at the same wavelength for the following reasons: the atmospheric phase variations at the times of the selected exposures are about 1 rad rms over the whole mirror aperture ($`67r_0`$) rather than $`r_0`$ for non-conjugate AO, so the isoplanatic patch is expected to be correspondingly large; light from the whole aperture ($`67r_0`$) contributes to the short exposure images rather than $`r_0`$ for AO wavefront sensors; and since selection of exposures is a passive activity, exposures can be at least 10 times longer than those necessary for servo-correction in AO, further increasing the limiting magnitude.

If the telescope mirror itself contains phase errors, then the best exposures will be those where the atmosphere has partially corrected these errors. The atmospheric phase variations may then be larger than 1 rad rms over the whole aperture, giving less of an improvement in isoplanatic patch size over AO. A similar situation may arise if a good exposure occurs when phase variations cancel out between two independent atmospheric layers.

A separate consideration is whether the read-out noise of the detector compromises the quality of the final image obtained by summing the selected short exposures. The development of CCDs with negligible readout noise (Mackay et al (2001)) has eliminated this noise penalty for fast frame rate imaging. The use of such CCDs for observations involving image selection should provide yet further improvement in the limiting magnitude of the reference star which can be used and also the limiting sensitivity to faint target objects in the field around the reference star.

The observations described here were designed to test the performance of the exposure selection technique at low light levels using a low noise CCD and to investigate the size of the isoplanatic patch for the technique. We were also interested in determining the quality of relative photometry and astrometry within globular cluster images using a range of reference star magnitudes.

## 2 Observations and Data Reduction

On the nights of 2001 July 5 and 6 we undertook high frame-rate imaging observations at the Cassegrain focus of the NOT using a camera built around a low-noise Marconi CCD. The camera comprised a front-illuminated $`576\times 288`$ frame-transfer CCD65 with $`20\times 30`$ $`\mu `$m pixels cooled in a liquid nitrogen dewar to 120 $`\mathrm{°}`$K to minimise dark current. The CCD was run by an AstroCam 4100 controller modified to provide a variable voltage clock signal for the output gain register of the CCD. We used frame rates between 18 Hz and 140 Hz, with sub-array readout where necessary. The CCD temperature was raised slightly on the second night (July 6) so as to improve the charge transfer efficiency. Short exposure images taken with high output register gain showed time-varying pattern noise with an amplitude equivalent to $``$0.5 detected photons.

The f/11 beam at the focus was converted to f/60 using a single achromat, giving an image scale of $`27\times 40`$ milliarcsec per pixel and a total imaging area of $`11.5\times 15.4`$ arcsec. In order to investigate the image quality achieved at large angular separations from the reference star, the camera optics were designed so that light from two regions of the sky separated by 25 arcsec could be superimposed on the CCD. This allowed science targets to be imaged at separations of up to 30 arcsec from the reference star. All the observations were made at 810 nm with a top-hat filter of 120 nm bandwidth and with no autoguider in operation.

Good exposures were selected on the basis of Strehl ratios calculated from each short exposure image of the reference star, a technique demonstrated previously in Baldwin et al (2001). The spurious bright spots found in images taken with our CCD65 detector (Mackay et al (2001)) did not occur frequently enough to have a significant effect on the reconstructed images, although a small number of frames with bright spots close to the reference star were excluded from our analyses.

The accuracy of measurement of the Strehl ratio and position of the brightest speckle was improved by convolving the individual frames with the diffraction-limited telescope point-spread function before these parameters were calculated. The peak intensity in the convolved images represents the location where correlation with the diffraction-limited point-spread function is maximised, providing a very good estimate for the location and intensity of a diffraction-limited speckle within the original image. After frames have been selected based on the Strehl ratio for the brightest speckle found, the corresponding unprocessed exposures are sinc resampled, shifted and co-added to produce the final image.

## 3 Results

The size of the isoplanatic patch in Strehl-selected exposures was deduced from observations of the binary stars HD 203991, 8 Lac and 61 Cygnii, which have separations of 0.6, 22 and 30 arcsec respectively. These were all bright enough that the selection and image shift correction for exposures was free of the effects of small number statistics in photon detection. One component of each binary was used as the reference to select the best 1% of the exposures, for which both the image full width at half-maximum intensity (FWHM) and Strehl ratio of the other component was also measured. Fig. 1 shows the fractional reduction in Strehl ratio of the second component as a function of angular separation from the reference star. A Gaussian fit to the points gives 50 arcsec for the FWHM of the isoplanatic patch for these selected frames, taken during a period in which the mean seeing was 0.51 arcsec. The image FWHM for objects 30 arcsec from the reference star using the best 1% from a run of 4000 exposures of 61 Cygnii was 130 milliarcsec, only slightly poorer than the FWHM of $``$80 milliarcsec expected for a well sampled diffraction-limited point-spread function. The measured FWHM is increased to 230 milliarcsec when selecting the best 10% of exposures, and 300 milliarcsec when selecting all the exposures for a shift-and-add image.

The effect of reference star magnitude on the quality of the resulting image was studied using observations with the full field of view of the CCD in a densely populated region of M13. The 3.4 MHz pixel rate of the CCD controller limited the frame rate for these observations to 18 Hz, allowing image motion to blur the exposures. This limited the Strehl ratios for reconstructed images to $``$0.16, with image FWHM of $``$100 milliarcsec. Data was taken from 6 runs of 1000 frames each. The selection of the best 1% of the 6000 exposures worked successfully with reference stars as faint as $`I=15.9`$. The effects of “over-resolution” on faint reference stars (Nieto & Thouvenot (1991)) were avoided by using other objects in the field for image quality measurements. Fig. 2 shows the variation in the Strehl ratio and FWHM of nearby stars when a range of different stars are used as the reference for image selection, shifting and adding. I band stellar magnitudes were taken from Cohen et al (1997). The Strehl ratios with reference stars of $`I=13.8`$ and $`I=15.9`$ are 0.13 and 0.065 respectively, a substantial improvement over $``$0.019 for conventional astronomical images generated from each of the individual runs. Fig. 3 shows a comparison of the image quality in these cases. The elongation in Fig. 3c may be due to image drift over the 55 s unguided exposure.

Given a suitable reference star, what combination of exposures gives the best image? The simple method applied so far is to select the best $`x`$%, $`x`$ being chosen to balance good Strehl ratio against adequate sensitivity. If the noise is Gaussian, a better solution would be to weight each exposure proportionally to the square of the signal-to-noise ratio (snr) for the quantity of interest. For example, detection of a faint stellar image with Strehl ratio $`S`$ against a uniform noise background would require weighting $`S^2`$. Since the value of $`S`$ for a faint target is unmeasurable in a single exposure, how well does the Strehl ratio for the reference star, $`S_{\mathrm{ref}}`$, represent that for the target?

To check this, data on the binary stars $`\zeta `$ Boötis (from May 2000, Baldwin et al (2001)), 61 Cygnii and a close pair of stars in M13 were reanalysed. The exposures in each run were ordered by $`S_{\mathrm{ref}}`$, placed in 100 equal bins starting with the best 1%, then the next 1% and so on. The exposures in each bin were shifted and co-added and the resulting Strehl ratio for the companion star was plotted against $`S_{\mathrm{ref}}`$ (Fig. 4a-c). All three plots show a linear relation between the Strehl ratios. The scatter of the points about the line increases from (a) to (c) mainly because of the decreasing number of exposures in each bin, but also in (c) due to photon statistics for the faint reference star. The latter also biases $`S_{\mathrm{ref}}`$, offsetting the mean line from the origin. The interesting conclusion is that $`S_{\mathrm{ref}}`$ is a good measure of $`S`$ for the target objects, even for comparatively low Strehl ratios; there seems no lower limit at which the correlation collapses. In practice, the limiting sensitivity using exposures weighted by $`S_{\mathrm{ref}}^2`$ was determined by pattern noise produced by the camera electronics. Exposure selection was found to give better rejection of this pattern noise, and has been used for the remainder of our analyses.

Fig. 5a shows a high resolution image of a field near the core of M13 generated by selecting the 10% of exposures with the highest Strehl ratios from 6000 taken at 18 Hz frame rate. The $`I=12.7`$ reference star used for exposure selection and recentring has been labelled. Typical stellar FWHM in this image are $``$120 millarcsec. Slight asymmetry in the faint halos around stars may relate to reduced charge transfer efficiency on the CCD at low signal level. Cross-sections through the stars with magnitudes 13.8 and 14.9 labelled X and Y are shown in Fig. 5b.

In order to asses the image quality, the 60 best exposures of M13 were separated into two groups of 30. The exposures in each group were shifted and added together giving two independent images of the field in M13. These images were then convolved with the diffraction-limited telescope point-spread function. Measurements were made of the location and intensity of the brightest pixel for each star, providing crude relative astrometry and photometry within the images. The rms difference in astrometry between the two independent datasets was found to be 6 milliarcsec, and the rms difference in stellar magnitudes was 0.02. The limiting magnitude for a 5$`\sigma `$ detection of a star above the noise in these images was $`I19`$.

## 4 Conclusions

The observations described here show that: (i) the selection of short exposure images can provide high resolution imaging up to 30 arcsec from a reference star; (ii) Strehl ratios and image shift can be accurately calculated with reference stars as faint as $`I=15.9`$; and (iii) assuming a mean density of 400 stars per square degree at the South Galactic Pole with $`I15.9`$ (Reid & Gilmore (1982)), 8% of the South Galactic Pole region is within range of a suitable reference star. If the mean stellar density over the whole sky is 3.3 times higher, then $``$20% of the sky will be within range of suitable reference stars with our current detector. The recent development of thinned low-noise CCDs and improvements to our camera electronics should improve the limiting reference star magnitude, allowing high resolution imaging over $``$25% of the sky. The faintest point source detected at a 5$`\sigma `$ level from 10% of exposures selected in an hour of observing would be $`I21`$ with our old camera or $`I23`$ with the improved instrument.

###### Acknowledgements.

The Nordic Optical Telescope is operated on the island of La Palma jointly by Denmark, Finland, Iceland, Norway and Sweden, in the Spanish Observatorio del Roque de los Muchachos of the Instituto de Astrofisica de Canarias.


text:
# Behavior of tumors under nonstationary theraphy

## I Introduction

Cancer is one of the leading research areas, since this desease is a main cause of death. Surgery and chemotherapy are unsuccessful in many cases. Today the principal efforts are addressed to search new treatment strategies e.g. in immunotherapy (see Refs. and references therein). In this case we refer to the use of cytokynes to stimulate the immune system. This is a protein hormone produced by activated lymphocytes which mediate both natural and specific immunity. The use of cytokines alone to boost the immune system represent one of the methods more commonly used in immunotherapy. The temporal evolution for this treatment comprises different steps: by supplying a starting dose of cytokines the rate of lymphocytes begins to increase due to the immunological reaction  reaching a maximum value; afterwards the lymphocytes begin to decay because of the decrease of cytokines concentration inside the body until it reaches normal values. In the following we shall call this time the activation period of the immune system. This proccess repeats again between any two succesive dose suplies. We consider these injections separated in time by the dosage period. Otherwise, it would provoke an overdose or failure in the treatment for shorter and longer time, respectively.
This work is devoted to understand the temporal tumor behavior when a periodical immunotherapy treatment is provided. Besides, we want to explore the set of parameter values to reproduce with our model the same features presented in (like short tem oscillation in tumor size and long term tumor relapse).
For this reason we reformulate the predator-prey model including new terms in the model which give account of tumor agressiveness, the diffusion of lymphocytes and the effect caused by cytokines on the tumor (Section 2). We analyse the model from the mathematical and biological points of view (Sections 3 and 4). We study in detail tumor evolution for different treatment regimes (Section 5) and, in the last section, we discuss our results and suggest some options for improving the model.

## II The model

An extensive review of all models of tumor-immune system dynamics  as far as we know was done. We agree with the idea that such dynamics is determined by a competence of interacting species resembling the predator-prey model. The main difference between our model and the clasical deterministic model (e.g. Bell 1973) is the inclusion of new terms taking into account, (1) the death of lymphocytes due to the increase of malignant cells population, (2) the flux of lymphocytes towards the place of local interaction and (3) the effect produced by the application of cytokine doses.
Let $`X`$ and $`Y`$ denote respectively the number of malignant and lymphocyte cells. The rate of malignant cells $`\left(\frac{dX}{dt}\right)`$ is given by:

$$\frac{dX}{dt}=aXbXY$$
(1)
We assume a growth rate proportional to $`X`$ and a decrease rate proportional to the frecuency of interaction with lymphocytes. Coefficients are $`a`$ and $`b`$, respectively, where $`a`$ is tissue dependent.

On the other hand, the growth rate of lymphocytes $`\left(\frac{dY}{dt}\right)`$ is described by:

$$\frac{dY}{dt}=dXYfYkX+u$$
(2)
It is proportional to the interaction with malignant cells and also to the flux per unit time of lymphocytes to the place of interaction. These effects are represented by the first and fourth terms in the right-hand side of equation 2. This last term characterizes the difussion process of lymphocytes that takes place in the surroundings of the tumor assuming a constant lymphocytes flux . On the other hand, the decrease rate depends on two factors: natural death and growth of malignant cells related to the effective area of tumor interacting directly with the lymphocytes. These are given by the second and third terms of this same equation where $`f`$ and $`k`$ are their respective coefficients of proportionality.

In order to introduce the effects produced by the treatment with cytokines in the proccess of activation of the immune system, we add a periodical function that mimics the periodical dosage. As a first approximation, we propose the function $`F\mathrm{cos}^2\omega t`$, where $`\omega `$ is the frecuency of the periodical behavior for cytokines inside the body. The modified Eq 2 will be given by the expression

$$\frac{dY}{dt}=dXYfYkX+u+F\mathrm{cos}^2\omega t$$
(3)
Taking into account Eq.1, we get the following system of differential equations.

$`{\displaystyle \frac{dX}{dt}}`$ $`=`$ $`aXbXY`$ (4)
$`{\displaystyle \frac{dY}{dt}}`$ $`=`$ $`dXYfY+Q`$ (5)
$`Q`$ $`=`$ $`KX+u+Fcos^2(\omega t)`$ (6)
where $`Q`$ is a function that picking up all the news contributions respecting the standard predator-prey model.
with $`X(0)=X_0`$ and $`Y(0)=Y_0`$ as initial conditions.

From the analysis for $`Q=0`$ the system of Eq.4, Eq.5and Eq.6 reduces to the predator-prey model, which has as equilibrium points a saddle point at the phase portrait origin and a center in the first quadrant of the phase diagram. This center reflects the oscillating behavior of the competence between both predator and prey, which is characterized by an oscillation frecuency (number of cycles per unit time around the center of the phase diagram). Taking the reverse of this frecuency as a characteristic time $`t_0=\frac{1}{\sqrt{af}}`$ and rescaling the equations of the system of Eq.4, Eq.5and Eq.6 by means of the following scaling parameters:

$`t`$ $`=`$ $`t_0\tau `$ (7)
$`X`$ $`=`$ $`X^{}x`$ (8)
$`Y`$ $`=`$ $`Y^{}y`$ (9)
the system given by the Eq.4, Eq.5and Eq.6 become:

$`{\displaystyle \frac{X^{}}{t_0}}\dot{x}`$ $`=`$ $`axX^{}bxyX^{}Y^{}`$ (10)
$`{\displaystyle \frac{Y^{}}{t_0}}\dot{y}`$ $`=`$ $`(dxX^{}f)yY^{}+Q`$ (11)
$`Q`$ $`=`$ $`F\mathrm{cos}^2(\omega t_0\tau )+uKX^{}x`$ (12)
sustituting the scaling parameters by the following values:

$`t_0`$ $`=`$ $`{\displaystyle \frac{1}{\sqrt{af}}}`$ (13)
$`X^{}`$ $`=`$ $`{\displaystyle \frac{\sqrt{af}}{d}}`$ (14)
$`Y^{}`$ $`=`$ $`{\displaystyle \frac{\sqrt{af}}{b}}`$ (15)
we get the following rescaled equations

$`{\displaystyle \frac{dx}{d\tau }}`$ $`=`$ $`\alpha xxy`$ (16)
$`{\displaystyle \frac{dy}{d\tau }}`$ $`=`$ $`xy{\displaystyle \frac{1}{\alpha }}ykx+\sigma +V\mathrm{cos}^2(\beta \tau )`$ (17)
with $`x(0)=x_0`$ and $`y(0)=y_0`$ as initial conditions.

where $`V=\frac{Fb}{af}`$, $`k=\frac{K}{bd}\sqrt{af}`$,

$`\alpha =\sqrt{\frac{a}{f}}`$, $`\sigma =\frac{ub}{af}`$ and $`\beta =\frac{\omega }{\sqrt{af}}`$
Analysing in detail the prior expressions we can interpret these parameters as follows:

From $`V`$ it is deduced that its value depends on the net value for doses ($`F`$) and on the action of the immune system on malignant cells ($`b`$). Hence it would represent an effective value of all doses employed in the activation of lymphocytes, namely, those necesary for activating the lymphocytes that take action directly against the tumor cells.
On the one hand, $`k`$ is directly proportional to $`K`$ which accounts for the negative effects exerted on the population of lymphocytes due to the size of the tumor and it is inversely proportional to the recognition ($`d`$) and attack ($`b`$) frequencies of the immune systen to malignant cells. So, we infer that the inverse value of this parameter ($`1/k`$) gives, “in some sense”, the control exerted by the immune system over the aggressiveness of the tumor due to its size.
From the similarity of the expressions for $`\sigma `$ and $`V`$, taking into account the remarks done for $`V`$ we interpret that $`\sigma `$ which depends on the action of the immune system on malignant cells ($`b`$) represent of all lymphocytes flux $`u`$ those effective value involved in the attack on tumor cells. Finally, $`\alpha `$ and $`\beta `$ are directly related to the proliferation of malignant cells and the frecuency of treatment, respectively.

## III Dynamics and Stability Analysis without treatment

The system of Eq.16 and Eq.17 corresponding to $`V=0`$( no treatment) is an autonomous system given by the following equations:

$`{\displaystyle \frac{dx}{d\tau }}`$ $`=`$ $`\alpha xxy`$ (18)
$`{\displaystyle \frac{dy}{d\tau }}`$ $`=`$ $`xy{\displaystyle \frac{1}{\alpha }}ykx+\sigma `$ (19)
with $`x(0)=x_0`$ and $`y(0)=y_0`$ as initial conditions.
Susbstituting the Eq.18 in the Eq.19 we get the differential equation

$$\frac{d^2x}{d\tau ^2}+(\frac{1}{\alpha }x\frac{1}{x}\frac{dx}{d\tau })\frac{dx}{d\tau }=(k\alpha )x^2+\sigma x$$
(20)
with $`x(0)=x_0`$, $`\dot{x}(0)=v_0`$ as initial conditions.

This equation is similar to that describing the motion of a particle in a force field , whose potential is:

$$U(x)=\frac{1}{3}(k\alpha )x^3\frac{1}{2}\sigma x^2$$
(21)
This potential has two extremes given by

$`x_1=0\text{and}x_2=\frac{\sigma 1}{k\alpha }`$

These extreme points depend on $`\alpha `$, $`k`$ and $`\sigma `$ as:

$`\sigma >1\{\begin{array}{c}x_1=0\text{minimum}\\ \frac{k}{\alpha }>1,x_2>0\text{maximum}\\ \frac{k}{\alpha }<1,x_2<0\text{maximum}\end{array}`$

and for

$`\sigma <1\{\begin{array}{c}x_1=0\text{maximum}\\ \frac{k}{\alpha }>1,x_2>0\text{minimum}\\ \frac{k}{\alpha }<1,x_2<0\text{minimum}\end{array}`$

We only consider motion for $`x>0`$, the suitable potential fields describing the motion of the particle are depicted in Fig. 1a and Fig. 1b. In Fig. 1a , the maximum represents an unstable point for the particle motion, contrary in Fig. 1b the particle oscillates around the minimum.
The analysis of fixed points in the phase space for Eq.18 and Eq.19 shows two steady-states. A fixed point is $`L_0=(0,\alpha \sigma )`$ with associated eigenvalues

$$\lambda _\pm =\frac{\alpha ^2(1\sigma )1}{2\alpha }\pm \left|\frac{\alpha ^2(1\sigma )+1}{2\alpha }\right|$$
(22)
For $`\sigma <1`$, $`L_0`$ is a saddle point while for $`\sigma >1`$ is a stable node.
The other fixed point is $`L_1=(\frac{1\sigma }{\alpha k},\alpha )`$ with associated eigenvalues given by:

$$\lambda _\pm =\frac{k\alpha \sigma }{2\alpha (\alpha k)}\pm \sqrt{\left[\frac{k\alpha \sigma }{2\alpha (\alpha k)}\right]^2(1\sigma )}$$
(23)
The real part of this eigenvalue is zero for $`\frac{k}{\alpha }=\sigma `$ with $`v_0<1`$.
When condition $`(\alpha k)^2\frac{k}{\alpha }>1`$ is fulfilled, we get two values for $`\sigma _c`$ solutions of

$`\alpha ^2\sigma _c^22\alpha [k2\alpha (\alpha k)^2]\sigma _c+k^24\alpha ^2(\alpha k)^2=0`$

given by

$`\sigma _c={\displaystyle \frac{k}{\alpha }}2\alpha ^2\left(1{\displaystyle \frac{k}{\alpha }}\right)^2`$ (24)
$`\pm 2\alpha |1{\displaystyle \frac{k}{\alpha }}|\sqrt{\alpha ^2\left(1{\displaystyle \frac{k}{\alpha }}\right)^2{\displaystyle \frac{k}{\alpha }}+1}`$ (25)
defining the region of complex eigenvalues and focus-like behavior.

The analysis of these eigenvalues provides a rich dynamics. For the case $`\frac{k}{\alpha }<1`$ we have different situations. The states and its stability for the second fixed point are depicted schematically in Fig. 2.
If $`\sigma <\frac{k}{\alpha }`$ ($`Re\lambda _\pm >0`$), we have an unstable focus or node depending on the parameter value $`\sigma _c`$ relative to those given by Eq.25.
On the contrary, stable behavior (focus or node) appears when $`\frac{k}{\alpha }<\sigma <1`$ ($`Re\lambda _\pm <0`$).
Now, if $`1<\sigma `$ the fixed point corresponds to a negative population of malignant cells, with no physical meaning.
For the case $`\frac{k}{\alpha }>1`$ and $`\sigma <1<\frac{k}{\alpha }`$ ($`Re\lambda _\pm <0`$) the critical point moves to the second quadrant of the phase diagram, being discarded as before.
For values of $`\sigma `$ in the ranges $`1<\sigma <\frac{k}{\alpha }`$ and $`1<\frac{k}{\alpha }<\sigma `$ we get a saddle point ($`\lambda _+>0`$ and $`\lambda _{}<0`$), whose separatrix splits the phase portrait into stable and unstable zones as can be seen in Fig.5. In all the cases the dynamics in the phase diagram is represented by a homeomorfism between two fixed points.

## IV Biological Significance

So far, we have presented a detailed analitical study of the linear stability of our model when $`V`$ is set equal to zero. The interpretation of this preliminary results will give us the esential features of the system.
Let us start with the case $`\frac{k}{\alpha }<1`$.
For $`\sigma <\frac{k}{\alpha }`$ the system evolves towards a state of uncontrolable tumor growth(see Fig. 3a) This case can be interpreted as a recurrence like behavior  very similar to q-switching oscillations observed in physical phenomena as, for example, in lasers. On the contrary, when $`\frac{k}{\alpha }<\sigma <1`$, our system evolves towards a controlable mass of malignant cells in a damped oscillating way (Fig. 3b). This state is considered by some authors as a dormant state.
However in both cases, there exist populations of malignant cells that grow towards a state in which immunological activity has been suppressed. In the first case this happens for any initial conditions, while in the second it only happens for an initially weak inmunological response (Fig. 4a and 4b). Let us now analyse the reverse condition $`\frac{k}{\alpha }>1`$.
In this case there are two possible ranges for $`\sigma `$:

$`1<\frac{k}{\alpha }<\sigma `$ and $`1<\sigma <\frac{k}{\alpha }`$

In both situations we are in the presence of a saddle point which means that for populations of cancer cells below the horizontal separatrix the dynamics is irreversible: this curve represents the critical amount of malignant cells for a fixed population of lymphocytes.
This situation is similar to the case analised before for a weak immnune system as an initial condition: the population of malignant cells grows towards a value such that the immunological response is reduced to zero (Fig. 5). This would represent a state where illness is not the cause of death but leaves the body unprotected against other diseases.
However for malignant cells above the horizontal separatrix it is possible to observe regression of tumor as has been reported in clinical experiments(see Ref. and references therein)

## V Stability analysis with treatment and biological implications

The system represented by Eq.16 and Eq.17 with $`V0`$ (cytokines doses amplitude) can be analised as an autonomous system. The procedure consists in substituting the oscillating function $`\mathrm{cos}\beta t`$ of the driven term $`F\mathrm{cos}^2\beta t`$ in the Eq.17 by a new variable $`u`$, which is a solution of the second order differential equation of a linear oscillator

$$\frac{d^2u}{dt^2}+\beta ^2u=0$$
(26)
(where $`u(0)=1`$, $`\dot{u}(0)=0`$) which can be written as two linear coupled differential equations

$`{\displaystyle \frac{dz}{d\tau }}`$ $`=`$ $`\beta ^2u`$ (27)
$`{\displaystyle \frac{du}{d\tau }}`$ $`=`$ $`z`$ (28)
(with $`u(0)=1`$, $`z(0)=0`$)
Then Eq.16 and Eq.17 become:

$`{\displaystyle \frac{dx}{d\tau }}`$ $`=`$ $`\alpha xxy`$ (29)
$`{\displaystyle \frac{dy}{d\tau }}`$ $`=`$ $`xy{\displaystyle \frac{1}{\alpha }}ykx+\sigma +Vu^2`$ (30)
$`{\displaystyle \frac{dz}{d\tau }}`$ $`=`$ $`\beta ^2u`$ (31)
$`{\displaystyle \frac{du}{d\tau }}`$ $`=`$ $`z`$ (32)
with $`x(0)=x_0`$, $`y(0)=y_0`$, $`u(0)=1`$, $`z(0)=0`$ as initial conditions.
This system presents the critical points $`L_0^{}=(0,\alpha \sigma ,0,0)`$ and $`L_1^{}=(\frac{1\sigma }{\alpha k},\alpha ,0,0)`$ whose projection in the $`y`$-$`x`$ plane coincides with those critical points of the unperturbed system (Eq.18 and Eq.19) with the same eigenvalues given by the Eq.22 and Eq.23 plus the new conjugate pair $`\lambda _\pm =\pm i\beta `$.
In this case we are in presence of a center manifold where solutions can be expanding or contracting, i.e., the asymptotic stability analysis carried out before loses its validity, needing more complex developments.
In order to avoid such complexity and gain a better comprehension, we may consider this system like a couple of one linear (Eq.26) and one nonlinear (Eq.20) oscillators, allowing us a more intuitive interpretation of the different regimes. Thus the changes from recurrent to dormant states of tumor cells in the periodical dosage regime can be interpreted as a lock of the unstable oscillations of the nonlinear oscillator imposed by the linear one. We can understand the complex behavior of coupled oscillators by representing its dynamics as a function of control parameters $`V`$ and $`\beta `$. In order to depict it we plot, in the parameter space ($`V\text{vs}\beta `$), those points for which tumor growth is uncontrolable.
In the case $`\frac{k}{\alpha }<1`$ for $`\sigma <\frac{k}{\alpha }`$, for effective value of doses and frecuencies higher than certain threshold, the system can revert from uncontrolable growth to a treatment controled population. Namely, for every set of parameter values $`\alpha `$,$`k`$ and $`\sigma `$, there are threshold values for $`\beta `$ and $`V`$ which split the parameter space (see Fig 6) into two zones corresponding to uncontrolable and controlable growth of malignant cells. The Fig 6 was generated for a fixed set of initial conditions. Although the behavior of the parameter space for different initial conditions is qualitatively the same, the threshold values show strongh dependence of the initial conditions.
From this result and taking into account the meaning of the parameters $`\beta `$ and $`V`$, we can infer that treatment is specific for each patient and kind of tumor since threshold doses values depend on the immunological response of each individual, on the malignant cells population at the begining of the treatment and also on the rate of proliferation of the tissue. Besides, the existence of threshold values reflects the fact that reaching controlable populations of malignant cells is only possible by mantaining a minimal dose above certain threshold given by the continuous line depicted in Fig 6, which can be well fitted by a hiperbolic function.
For effective doses and frecuencies below these threshold values, the system behaves qualitatively the same as without treatment. However, contrary to this statement it is also found that for low frecuencies the growth of malignant cells can be controlled in spite of being below the threshold values (see, for details, Fig 6).
There exists a ”paradoxical” phenomenon observed in experiment and the clinic, consisting in the fact that the enhacement of the inmune system with immunotherapy stimulates tumor growth, which could be explained, “in some way”, by this result, i.e., why, for fixed doses burden , the growth becomes uncontrolable at given frecuencies above those localized in the region of controlable growth of malignant cells (Fig. 6). Such values would represent an optimal treatment as it reduces doses burden and treatment frecuency.
Also, for higher frecuencies with small $`V`$, growth of cancer cells can be controlled as shown in Fig 6. The observance of these optimal values would be important because of the negative effects produced when cytokines concentration reaches above a critical concentration..
Now, increasing the amplitud of the effective dose value for a fixed frecuency above the threshold, a malignant cells population reduction is obtained, in spite of an uncontrolable growth being observed for some higher doses burden (see Fig 6).
On the other hand, setting $`V`$ to some value and varying the frecuency from zero to higher values, different behaviors are reproduced. At zero treatment frecuency, tumor cells population is controlable with the lowest values of doses burden (Fig. 6). However, for frecuencies different from zero, we find zones of recurrent and dormant growth of tumor cells. The population of cancer cells controled under treatment presents an oscillating behavior (see Fig. 7).
For values of the parameters satisfying $`\frac{k}{\alpha }<\sigma <1`$ (that we interpret as a dormant state), malignant cells population can be reduced by increasing effective doses. On the other hand, varying the frecuency for fixed effective doses values, an oscillating behavior for the population of malignant cells is obtained, as in the previous case.
In all these cases, regrowth of malignant cells takes place after treatment interruption . This can be easily understood if we take into account that the population of malignant cells with zero value ($`x=0`$) represents, in the mechanical analogue (Eq.20), an unstable point (a potential maximum, as that shown in Fig. 1b). This means that any variation will lead the system towards a minimal potential position. Therefore for a residual population slightly greater than zero, a regrowth of tumor cells will take place after the treatment.
Analysing the behavior when $`\frac{k}{\alpha }>1`$ we arrive at the following results. The range of values with physical sense for $`\sigma `$, i.e. $`1<\frac{k}{\alpha }<\sigma `$ and $`1<\sigma <\frac{k}{\alpha }`$ allows only two critical points: a stable node and a saddle point. In this case dynamics in the phase portrait is the same as that without treatment. There are no possible changes in the dynamics, so treatment is useless.
Hence, for $`\frac{k}{\alpha }>1`$, only initial conditions determine the final outcome of tumor evolution, irrespective of the applied treatment.

## VI Conclusions

In this work we intend to give a new focus to the dynamics of tumor growth in a periodical regime of immunotherapy. We explain such dynamics considering this system as two coupled oscillators, namely, the competence between the immune system and malignant cells analized as a nonlinear oscillator coupled with a linear one that simulates the treatment.
This simple model allowed us to describe all possible states in wich a tumour can be found. It also presented some of the features found in tumor dynamics, outlined by some authors, such as the existence of short term oscillations of tumor size as well as the long term tumor relapse. On the other hand, this model gives the dependence of tumor growth on some parameter values related to the treatment: the frecuency and amount of applied doses. We conclude from this study that the evolution of tumor submitted to immunotherapy has a strong dependence on these parameter values.
In some cases, growth of malignant cells can be reverted with inmunotherapy treatment. Corresponding threshold values are obtained for treatment frequency and dose above which growth is stopped and malignant cells population reduced.
Also, for those inmunological parameters for which a stable population of malignant cells exists, the size of the dormant tumor can be reduced by increasing dose burden.
It was shown, as well, that for certain relation among the parameter values, tumor presents a recurrent behavior with or without treatment.
In all these cases, when a reduction of the tumor is possible, best results are obtained for low constant values of the dose. Nevertheless in all cases previously analyzed, after the interruption of the treatment, tumor regrowth is observed.
This would confirm the fact that treatment with immunotherapy using cytokines alone is not succesful enough in the treatment against cancer. Therefore another kind of therapy would be required.
As a way to improve the model, we propose the introduction of other terms taking into account effects produced by stochastic perturbations due to enviromental conditions. Some authors atribuit to these perturbations the main cause of possible jumps from stable to unstable behavior in tumor growth dynamics. This will be considered in a future work.


text:
# Kaon-nucleon and D-nucleon scattering in the quark model, including spin-orbit interactions

## 1 Introduction

$`KN`$ scattering is ideally suited for studying the origins of the nonresonant “nuclear” force. Conventional s-channel baryon resonance production is excluded because there is no valence annihilation, and one-pion exchange is forbidden because there is no three-pseudoscalar vetex. We may therefore study the nonresonant part of hadron scattering in relative isolation, uncomplicated by one-pion exchange. In this work we calculate $`KN`$ phase shifts at Born order in a quark exhange model with one gluon exchange (OGE) and linear scalar confinement. The dominant hyperfine term has given very good results in S-wave for I=2 $`\pi \pi `$ , I=3/2 $`K\pi `$ , I=1 $`KK`$, $`KN`$ , and short-ranged $`NN`$  scattering. Here we include the spin-orbit and subdominant spin-independent contributions, and give results for higher-L waves. $`KN`$ scattering is an excellent place to test our model since data exists and there is a large spin-orbit effect, especially evident in the P-waves. Understanding the spin-orbit effect here will be important for applications to other hadronic interactions, such as $`DN`$. $`KN`$ scattering is almost entirely elastic below $`K\mathrm{\Delta }`$ threshold, and the Born approximation is expected to be relevant since the interaction is known experimentally to be relatively weak.

## 2 Experiment

There are some basic experimental features of $`KN`$ scattering which are generally agreed upon. The I=1 channel is relatively well determined from $`K^+P`$ scattering. Using the conventional $`L_{I\mathrm{\hspace{0.33em}2}J}`$ notation, the $`S_{11}`$ and $`P_{11}`$ channels are repulsive and the $`P_{13}`$ channel is attractive. I=0 is more in doubt because of inherent difficulties in the experimental analysis. The I=0 $`KN`$ scattering amplitudes have been extracted from $`Kd`$ scattering, using models of the deuteron breakup process and form factors. It is clear, however, that the $`P_{01}`$ channel is strongly attractive (in fact anomalously so). Inelasticities are known to be large above inelastic thresholds, in particular in the $`P_{01}`$ wave. For reference, $`K\mathrm{\Delta }`$ (I=1 $`KN`$ only) opens at $`P_{lab}=0.86`$ GeV, $`K^{}N`$ opens at $`P_{lab}=1.08`$ GeV and $`K^{}\mathrm{\Delta }`$ opens at a much higher $`P_{lab}=1.74`$ GeV.

## 3 Calculation

The fundamental interactions in QCD are between quark and gluon constituents. The quark-quark OGE $`T_{fi}`$ is derived as the nonrelativistic reduction to order $`P^2/m^2`$ of the Feynman amplitude for two quarks to exchange a gluon with phenomenological strength $`\alpha _s`$. The result is the usual Breit-Fermi Hamiltonian, which is well known from atomic physics, and includes the dominant spin-spin hyperfine interaction as well as Coulomb, spin-orbit, tensor, and subdominant spin-independent contributions. The confining interaction is modelled as a linear scalar interaction with strength b (the string tension) between a pair of quarks, and yields linear and spin-orbit contributions, as well as smaller spin-independent terms. The quark-quark interaction is followed by quark interchange, so that the outgoing hadrons emerge as color singlets .

The scattering arises from the interactions between quarks in different hadrons, and the interaction Hamiltonian $`H_{int}`$ is the sum of quark-quark interactions between all such pairs of quarks. To Born order the scattering amplitude is given by the matrix element of $`H_{int}`$ between external hadron states, which ideally should be eigenstates of the respective free Hamiltonians. We approximate the true wavefunctions by Gaussian forms for calculational simplicity, as these lead to analytical results for the scattering amplitudes. This approximation has been tested in previous calculations . The parameters are relatively well determined from quark-model phenomenology; we use a reasonably standard set (below) for our numerical evaluation.

$$\begin{array}{cccc}\hfill \alpha _s& =& 0.6\hfill & \text{strong coupling constant}\hfill \\ \hfill b& =& 0.18\text{GeV}^2\hfill & \text{string tension}\hfill \\ \hfill \beta & =& 0.4\text{GeV}\hfill & \text{meson wavefunction length scale}\hfill \\ \hfill \alpha & =& 0.4\text{GeV}\hfill & \text{baryon wavefunction length scale}\hfill \\ \hfill m_{u,d}& =& 0.330\text{GeV}\hfill & \text{nonstrange constituent quark mass}\hfill \\ \hfill m_s& =& 0.550\text{GeV}\hfill & \text{strange constituent quark mass}\hfill \end{array}$$

The phase shifts are calculated from the scattering amplitudes and phase space factors in the usual way .

## 4 Results

Hashimoto  has performed the most recent comprehensive single-energy phase shift analysis we are aware of, and we take his data for comparison. We also compare our results with those of the recent resonating group method (RGM) calculation of Lemaire, et al. In the I=1 channel, our results are in reasonable agreement with the RGM calculation. Both calculations agree at least qualitatively with the data except in the $`P_{13}`$ wave, for which they both give the wrong sign. In I=0, the absolute values of the calculated phase shifts are in general too small, especially in the remarkably large $`P_{01}`$ wave. The OGE spin-orbit contributions to our quark model $`KN`$ phase shifts are shown in Figure 3. The large $`P_{01}`$ amplitude evident there is perhaps an indication that we have included at least some of the correct physics. It is, however, not large enough by itself to account for the data, and is partially cancelled by the confinement spin-orbit contribution and the negative spin-spin hyperfine term.

Mukhopadhyay and Pirner  have found that the phase shifts calculated from the symmetric part of the OGE spin-orbit interaction agree reasonably well with the experimental phase shifts at low energies except in the $`P_{01}`$ wave, which they seriously underestimate. The antisymmetric part of the OGE spin-orbit interaction is omitted in their calculation because it gives rise to large splittings in the baryon spectroscopy which are not observed. They similarly conclude that the confinement spin-orbit interaction partially cancels the OGE spin-orbit contribution, and the inclusion of the spin-spin interaction actually makes the $`P_{13}`$ phase shift theoretically repulsive.

The opening of inelastic channels is known to be a very large effect in the $`P_{01}`$ wave, and may be responsible for the large experimental phase shift. Coupling to inelastic channels will be treated in a future calculation.

## 5 $`DN`$ scattering

Knowledge of vacuum $`DN`$ scattering amplitudes and observables is important for understanding open-charm hadronic interactions in hot, dense media, as occur in heavy ion collisions. We applied our $`KN`$ scattering model to obtain $`DN`$ elastic phase shifts, which are shown in Figure 4. (The only new parameter was the charmed quark mass, taken to be $`m_c=1.550`$ GeV.) Inelasticities may well be large in $`DN`$, as in $`KN`$ scattering. In I=1, the linear and spin-spin hyperfine terms are roughly comparable, and dominate scattering in the waves shown. The spin-orbit splitting is small in comparison. It is noteable that in I=0 the spin-independent terms (Coulomb and linear) are identically zero and the spin-orbit interactions dominate the scattering.

## References