text:
# Spreading of Block Copolymer Films and Domain Alignment at Moving Terrace Steps1footnote 11footnote 1The following article has been submitted to the Journal of Chemical Physics. After it is published, it will be found at http://ojps.aip.org/jcpo

## I Introduction

Equilibrium structure of block copolymer systems was in the center of theoretical and experimental research for a number of years now, and is mostly understood BatesFredrickson1996 ; BatesFredrickson1990 ; Krausch ; Fasolka ; HarrisonInScience ; HahmDefects ; Morkved ; Hahm . It is known that in the bulk, at temperatures below order-disorder transition (ODT), molten diblock copolymers phase-separate into periodic structures with periodicity dependent on polymer size. The actual structures formed depend upon the proportion of different blocks inside copolymer, and include lamellar, cylindrical, spherical and various bicontinuous morphologies.

The spreading dynamics of copolymer films, on the other hand, is much less investigated Joanny1995 . While little difference besides significantly higher viscosity is expected in the spreading of macroscopic copolymer drop, large deviations may arise on microscopic length scales, when the film thickness approaches domain period. In particular, thin films inherit domain structure of the bulk samples, but their domains are typically aligned parallel or perpendicular to the substrate. Combined with extreme anisotropy of diffusion relative to the domain orientation, spreading dynamics of thin films is expected to be very different from that of simple fluids.

In the present paper we look at the spreading of phase-separated copolymer films on wetting substrates in the regime of complete wetting. It is known that copolymer films at small thicknesses adopt a terraced topology, with thickness of each terrace defined by the underlying domain structure Coulon . The copolymer domains are then localized inside these terraces (Fig. 1), so that terraces can slide on top of each other.

General spreading of layered drops was first investigated by de Gennes and Cazabat Cazabat ; Voue2001 , and consists of viscous flow along each terrace, combined with disjoining pressure driven flow at terrace steps. Here we investigate spreading dynamics of the layered copolymer films. On one hand, the thickness of copolymer terraces (layers) significantly exceeds atomic sizes, so that effect of disjoining pressure is negligible at all but the bottommost terrace. Hence, the flow between terraces is driven by other factors, presumably by jump in the internal stress. Additionally, the thickness of each copolymer terrace needs not be uniform. This results in the forced diffusive flow of material along the terrace with diffusion constant determined by the underlying domain topology and strength of copolymer block interaction. We find that at large times the spreading rate $`v`$ for the copolymer film takes the diffusion like time dependence $`vt^{1/2}`$, as in general layered drops of Refs. Cazabat,  and Voue2001, , and contrary to the Tanner’s law for simple liquids, which predicts $`vt^{0.9}`$.

The activated flow of chains at terrace steps, with chains hopping between copolymer domains, can also account for the growth of aligned domain patterns in spreading films. Specifically, asymmetric copolymers that phase separate into cylindrical domains in bulk, often form striped pattern in the thin films, similar to those shown on Figure 1 Jaeger ; HarrisonInScience ; Morkved . Easy control over periodicity of these patterns has driven large interest in their formation, but lack of ordering and generally large number of defects limited their applicability HarrisonInScience ; HahmDefects ; Morkved ; Hahm .

However, Hahm and Sibener Hahm  recently reported that nearly perfect alignment is observed on annuli structures of polystyrene-block-polymethylmethacrylate (PS-PMMA). In that experiment, a thin film of PS-PMMA was spin cast onto a SiO2 substrate that was initially covered with a minor polar solvent. Dewetting properties of the solvent resulted in the formation of bare circular regions with thick annular rims (Fig. 2). On subsequent annealing the thick rims spread out and flatten. The expected striped domains showed unexpected radial alignment. This strong alignment along the spreading direction was confirmed in later experiments with a different geometry Deepak . Below we will argue that hopping of copolymer chains between domains at terrace steps, necessary for the spreading, may be responsible for this alignment. In particular, chain hopping should lead to the continuous growth of cylindrical domains, while thickness gradient near terrace steps should generate alignment of the growing domains. As a result, one may expect aligned, defect-free domain patterns along the whole area of spreading, similar to the observations by Hahm and Sibener Hahm .

Improving the ordering of copolymer film domains is regarded as a key to many potential applications Zehner ; Mansky ; Harrison1998b . Accordingly, much attention has been given to potential alignment mechanisms: local alignment near a terrace step ThomasPRL , thermal annealing HahmDefects ; HarrisonInScience , electric fields Morkved , growth from an aligned boundary SegalmanKramer ; Deepak  and imposed shear CoulonShear . The alignment under spreading observed by Hahm and Sibener represents a potentially new mechanism, different from these others. To understand this mechanism is important for the overall goal of controlling copolymer film patterns.

This paper is organized as follows. In section II we review the known facts about homopolymer spreading. In sections III \- V we develop our model for the copolymer film spreading. Then, in section VI our model is compared with the experimental results on PS-PMMA, and a new dynamic mechanism of domain alignment is discussed.

## II Spreading of Homopolymer Films

We will start with a brief review of the spreading dynamics of regular homopolymers. This spreading has been widely investigated and is mostly understood deGennes1985 ; Bruinsma ; Leger ; deGennes1984 . The rate of this spreading, as well as the shape of the advancing liquid front, can be derived from the lubrication approximation. Within this approximation, the two-dimensional current density of material flow $`J`$ is given by

$$J=\frac{h^2}{\eta }\left(b+\frac{1}{3}h\right)p,$$
(1)
where $`\eta `$ is the viscosity, $`h`$ is drop thickness, $`b`$ is the slip length and $`p`$ is pressure inside the drop. There are generally two major contributions to this pressure: one is from Young-Laplace capillary effect $`\gamma ^2h`$ and the other is from Van-der-Waals disjoining pressure $`W(h)=A/(6\pi h^3)`$. The coefficient $`A`$ is a material dependent constant. Thus

$$p=\gamma ^2h+W(h).$$
(2)
In the frame moving with velocity $`U`$ of the contact line, the current becomes $`J^{}=JUh`$. For steady-state spreading, this current must vanish, and the spreading profile $`h(x)`$ is determined by

$$\left[\gamma ^2h+W(h)\right]\frac{h}{\eta }\left(b+\frac{1}{3}h\right)=U.$$
(3)
Based on the thickness $`h`$ of the drop, there are three distinct regions that can be distinguished. First of all, there is a region of macroscopic drop where film thickness is much larger than the slip length ($`hb`$), and contribution from disjoining pressure is negligible. Equation (3) in this region reduces to the well known spreading equation investigated by Tanner deGennes1985

$$\frac{\gamma }{3\eta }h^2^3h=U$$
(4)
which at large $`h`$ predicts zero curvature wedge-like profile $`h|x|\mathrm{ln}^{1/3}|x|`$ and cubic dependence of the spreading rate on the apparent contact angle $`U\frac{\gamma }{\eta }\theta ^3`$.

The second distinct region of (3) is called the precursor film deGennes1985 , i.e. the area of the film where thickness is very small so that plug flow takes place and Van-der-Waals interactions dominate. Then (3) reduces to

$$\frac{bh}{\eta }W(h)=\frac{Ab}{2\pi \eta }h^3h=U,$$
(5)
which predicts $`h|x|^{1/2}`$ thickness profile.

Finally, there’s an intermediate region characteristic only to polymeric drops and denoted as the foot region. Existence of this region is associated with very large slip length of polymeric liquids, so that there is a region of the film where Van-der-Waals interactions are already insignificant, but flow is still of the “plug” type. Then (3) reduces to

$$\frac{b\gamma }{\eta }h^3h=U,$$
(6)
with film thickness obeying deGennes1984  $`h|x|^{3/2}`$.

The described model for the polymer spreading sufficiently well describes both precursor and macroscopic film regions. However, $`x^{3/2}`$ dependence for the foot is not generally observed Leger .

Similar effects may be expected in copolymer drops. In the precursor film, for example, Van-der-Waals interactions with the substrates significantly exceed repulsive interactions between copolymer blocks, so that copolymer should obey the same dynamics (5) as a homopolymer. Similarly, in the macroscopic drop region, domain orientation may be expected to be random and macroscopic viscosity remains isotropic. Hence, the domain structure may be accounted for via a macroscopic viscosity so that general spreading dynamics still follows the Tanner’s law (4).

However, in the foot region, one may expect even stronger deviations from (6) than those in homopolymeric liquids, as periodic domain structure of copolymer gets aligned parallel to the substrate and a quantized thickness profile emerges. Then the continuous description of equation (3) becomes invalid. To describe the flow of material inside this terraced structure we will look at two different regions of material flow: on one hand, there’s continuous flow inside each terrace, which can be described from a hydrodynamical standpoint; and, on the other hand, there’s activated flow at the terrace steps, where chains have to hop between discontinuous domains of different terraces.

## III Flow Inside a Single Terrace

Hydrodynamic Approximation. The terraced structure of copolymer films comes as a result of domain alignment parallel or perpendicular to the substrate. Some of these terraces may slide on top of each other (Fig. 1), with flow of material inside each terrace being strictly one or two dimensional and parallel to the substrate. Then it is convenient to characterize profile of each terrace in terms of area density $`\rho h`$, so that the flow can be described by hydrodynamic equations of motion:

$$\frac{\rho }{t}+(\rho \stackrel{}{u})=0,$$
(7)
and

$$\rho \frac{Du_j}{Dt}=_i\sigma _{ij}ku_j$$
(8)
where $`k`$ is the viscous friction opposing the flow, and $`\sigma _{ij}`$ is the stress tensor. In the discussion that follows, we will consider only slow ($`Du/Dt0`$) spreading in the $`x`$-direction. Additionally, we will neglect $`y`$-dependence of the stress tensor $`\sigma `$. This dependence, which arises from the domain-induced anisotropy of phase-separated copolymer and potentially affects such parameters as direction of spreading, is insignificant for the present discussion and will be briefly considered in section VI. Hence, (8) reduces to the familiar Darcy’s law

$$u_x=\frac{1}{k}\left[\frac{\sigma _{xx}}{x}+\frac{\sigma _{xy}}{y}\right]\frac{1}{k}\frac{\sigma _{xx}}{x}$$
(9)
Non-uniform stress $`\sigma _{xx}`$ exists in the film because it is not at its equilibrium thickness $`h_0`$. The stress can be expressed in terms of the free energy per chain $`F_c`$. This free energy evidently has a minimum when $`h=h_0`$, and for $`h`$ near $`h_0`$ $`F_c`$ is quadratic in $`(hh_0)`$ Semenov ; WittenCylinders :

$$F_cF_0+g_0\gamma _{AB}v_0\frac{(hh_0)^2}{h_0^3}$$
(10)
Here $`v_0`$ is the volume displaced by each chain. We treat the polymers as an incompressible liquid; thus $`v_0`$ is a fixed constant. The $`\gamma _{AB}`$ is the interfacial energy of the interface between the A and B blocks. The quantity $`g_0`$ is a numerical constant that depends on the domain morphology and how it is distorted in the process of changing $`h`$ (e.g. $`g_0=3/2`$ for plane lamellar terraces g ).

To determine the stresses, we consider a small section of the terrace of width $`\mathrm{\Delta }x`$ and length $`\mathrm{\Delta }y`$ with thickness $`hh_0`$. The section is prevented from relaxing to height $`h_0`$ by forces from the adjacent film. If the section were detached from adjacent terrace, it would contract or expand laterally, allowing the height to return to its equilibrium value. To prevent this, the adjacent terrace must exert a tensile stress on the section. In particular, during spreading in $`x`$ direction, $`\mathrm{\Delta }y`$ stays constant and stress $`\sigma _{xx}`$ is

$$\sigma _{xx}=\frac{1}{h\mathrm{\Delta }y}\frac{F_T}{(\mathrm{\Delta }x)}|_{\mathrm{\Delta }y}$$

Using incompressibility of polymer $`h\mathrm{\Delta }x\mathrm{\Delta }y=V`$ and the fact that the terrace’s free energy $`F_T=(V/v_0)F_c`$ we finally get

$$\sigma _{xx}=\frac{h}{v_0}\frac{F_c}{h}|_{\mathrm{\Delta }y}=2g\frac{(hh_0)}{h_0^2}\gamma _{AB}$$
(11)
The numerical constant $`g`$ depends on the $`g_0`$ of Eq. 10, and in the trivial case of lamellar domain $`g=g_0`$. In general, the tensile stress also includes a Laplace-Young contribution $`\gamma ^2h`$ from the curvature of the free surface. This part of the pressure is negligible in the long, flat terraces we are considering. Then the profile of the spreading terrace is given by

$$\frac{h}{t}=\frac{}{x}\left[\frac{2g\gamma _{AB}}{k}\frac{h}{h_0^2}\frac{h}{x}\right].$$
(12)
Away from the terrace steps, the thickness of the terrace is close to $`h_0`$, so that to the first order in $`hh_0`$ the thickness profile is given by

$$\frac{h}{t}=D_h\frac{^2h}{x^2},$$
(13)
where $`D_h=2g\gamma _{AB}/h_0k_{h=h_0}`$.

Dissipative Forces at the Substrate. We now consider the physical origin of the dissipative coefficient $`k`$. The most common type of viscous forces are the frictional forces at the substrate interface. These forces are most easily expressed in terms of the slip length $`b`$, i.e. the distance below the substrate at which velocity of liquid $`u(z)`$ extrapolates to zero: $`b=\left[u/(du/dz)\right]_{z=0}=\eta /\beta `$. Here $`\eta `$ is liquid viscosity, and $`\beta `$ is friction coefficient at the substrate. The magnitude of the friction force does not depend on terrace thickness, so that coefficient $`k`$ entering (9) is $`k=k_s=\beta /h=\eta /bh`$.

While this expression for $`k_s`$ is believed to be true for macroscopic drops, it fails in thin films due to rapid decrease in monomer mobility near a substrate. For instance, near a high-energy substrate monomer mobility may drop $`10^4`$ times compared to its bulk value Bruinsma , leading to a similar drop in $`k`$. Introducing phenomenological parameter $`\xi =\mu _{bulk}/\mu _{substrate}`$, which relates bulk ($`\mu _{bulk}`$) and substrate ($`\mu _{substrate}`$) mobilities, the friction coefficient for the sliding chains becomes $`k=\xi k_s`$.

This decrease in mobility falls off rapidly away from the substrate. Then the chains that are affected most are those in immediate contact with the substrate, while remaining chains can freely diffuse between them. This leads to a two-liquid model of flow, first introduced for polymer spreading by Bruinsma Bruinsma . According to this model, the first liquid consists of chains that come in direct contact with the substrate and slowly slide with velocity $`\stackrel{}{u}_s`$, while the other liquid consists of chains that diffuse between them with velocity $`\stackrel{}{u}_d`$. The net velocity of material flow is then

$$\stackrel{}{u}=\tau \stackrel{}{u}_s+(1\tau )\stackrel{}{u}_d,$$
(14)
where $`\tau `$ is the fraction of the sliding chains, i.e. chains that are in contact with the substrate. If the fraction $`\tau `$ of sliding chains is small, we can temporarily switch to the frame moving with velocity $`U_d`$ where sliding chains become diffusing, with the drag force on them given by Einstein relation $`\stackrel{}{f_0}=\frac{k_BT}{D_c}(\stackrel{}{u}_s\stackrel{}{u}_d)`$. Here $`D_c`$ is the self-diffusion coefficient for chains moving along domain walls. The net force between sliding and diffusing chains is then $`\stackrel{}{f}=\tau \stackrel{}{f_0}/v_0=\tau \frac{k_BT}{v_0D_c}(\stackrel{}{u}_s\stackrel{}{u}_d)`$, and the Darcy’s law for each of these components is

$`(1\tau )\sigma +\tau k_d(\stackrel{}{u}_d\stackrel{}{u}_s)=0`$ (15)
$`\tau \sigma +\xi k_s\stackrel{}{u}_s\tau k_d(\stackrel{}{u}_d\stackrel{}{u}_s)=0`$ (16)
where $`k_d=k_BT/v_0D_c`$, $`k_s=\eta /bh`$ and $`v_0`$ is the volume of one chain. These equations, together with (14) give the final expression for the velocity:

$$\stackrel{}{u}=\left[\frac{(1\tau )^2}{\tau }\frac{1}{k_d}+\frac{1}{\xi k_s}\right]\sigma .$$
(17)
As expected, this expression reduces to (9) in the case when sliding dominates diffusive flow ($`\xi k_sk_d`$). With velocity (17), one-dimensional film spreading is still described by (13), with collective diffusion coefficient $`D_h`$ being

$$D_{slid}=\frac{3\gamma _{AB}}{h_0}\left[\frac{(1\tau )^2}{\tau }\frac{1}{k_d}+\frac{1}{\xi k_s}\right].$$
(18)
Dissipative Forces Between Sliding Terraces. Along with dissipation at the substrate, an alternative dissipation mechanism arises when terraces slide over each other. During this sliding motion, brushes from adjacent terraces interpenetrate and viscous dissipation takes place. The areal density of the dissipated power $`T\mathrm{\Sigma }`$ is given by

$$T\mathrm{\Sigma }=_\zeta ^\zeta \eta ^{}\left(\frac{dv}{dz}\right)^2𝑑z\eta ^{}\frac{v^2}{2\zeta },$$
(19)
where $`v`$ is spreading velocity, $`\eta ^{}`$ is the effective viscosity of one of the copolymer blocks, and $`\zeta `$ is the characteristic depth of interpenetration of two adjacent brushes. In the plane lamellar morphology the brush thickness $`h_AN_A^{2/3}`$ so that the interpenetration depth is WittenInterpenetration :

$$\zeta \left(R_A^4/h_A\right)^{1/3}aN_A^{4/9}.$$
(20)
where $`R_A`$ is radius of gyration of a free chain, and $`a`$ is a microscopic length scale close to the monomer size $`a_0`$:Semenov ; WittenInterpenetration  $`a(kT/\gamma _{AB}a_0^2)^{1/9}a_0`$. The effective viscosity $`\eta ^{}`$ entering (19) is qualitatively different from the homopolymer melt viscosity. This is due to the fact that chains in a segregated copolymer are bound by domain walls and may no longer undergo reptation motion characteristic of linear homopolymers. Instead, the copolymer system becomes increasingly similar to the star polymer melt, where chains are also constrained at their vertex WittenInterpenetration . The characteristic mode of motion of star polymers is repetitive contractions and extensions of one arm of the star, with time scales exponentially increasing with the length of this arm. Correspondingly, the viscosity of star polymers also exhibits exponential increase Fetters  $`\eta ^{}\eta _{arm}e^{N_{arm}/N_e}`$, where $`N_e`$ is the average entanglement length, $`N_{arm}=2\zeta N_A/h_A2N_A^{7/9}`$ is average arm length, and $`\eta _{arm}`$ is the viscosity of the corresponding linear polymer of length $`2N_A^{7/9}`$. Thus the friction coefficient between terraces becomes

$$\beta =\frac{T\mathrm{\Sigma }}{v^2}\frac{\eta _{arm}}{2a}N_A^{4/9}\mathrm{exp}\left(2N_A^{7/9}/N_e\right),$$
(21)
so that final dynamics of the terrace is described by (13) with collective diffusion coefficient $`D_h`$ given by

$$D_{visc}=\frac{2g\gamma _{AB}}{\beta }\frac{6a\gamma _{AB}}{\eta _{arm}}N_A^{4/9}\mathrm{exp}\left(2N_A^{7/9}/N_e\right).$$
(22)
## IV Flow At The Terrace Steps

Besides continuous flow of material inside each of the terraces, flow of material between terraces is needed to sustain spreading. However, topologically different domains of adjacent terraces (Fig. 3) are either disconnected, or connected via some intermediate morphology Harrison1998 . In either case, activated events are needed to create flow at the step: in the former case, single chain hopping is necessary, while, maintaining intermediate morphology in a dynamic system, such as spreading system of Fig. 3 with step moving to the right, requires many chain activation events. As the time scales for many chains events are significantly longer than for single chain hopping, from here on we will assume that single chain hopping dominates flow at terrace steps and intermediate morphologies do not form.

Qualitatively, the decrease in chain self-diffusion due to activation events has the form Barrat ; Kramer

$$DD_0e^{\mathrm{\Delta }E/kT},$$

where $`\mathrm{\Delta }E`$ is the activation free energy for the minority block to hop through the stationary block, and $`D_0`$ is the self-diffusion coefficient of a homopolymer melt with chains assembled from majority block monomers and having the same polymerization as the copolymer chain.

There are two distinct contributions to the above activation energy: the energy cost $`\mathrm{\Delta }E_0`$ for the hopping block to pass through the incompatible stationary block, slightly decreased by the excess energy, or chemical potential, of the chains inside domains. In stationary systems, chemical potentials for all domains are the same and no net flux occurs, while in the spreading systems chemical potentials may differ between terraces, producing net flow of chains.

To estimate this flow, we notice that the average rate at which chains hop between domains with chemical potential $`\mu `$ and spacing $`h_0`$ is LandauVol6

$$(\mathrm{\Delta }t)^1=\frac{2D}{h_0^2}\frac{2D_0}{h_0^2}e^{(\mathrm{\Delta }E_0+\mu )/kT}.$$
(23)
For two parallel lamellar layers separated by $`h_0`$ the current density $`j_0`$ of chains hopping from one layer to the other is given by $`j_0=\rho _0(\mathrm{\Delta }t)^1`$, where $`\rho _0`$ is equilibrium number of chains per unit area. At equilibrium, an equal current density flows in the opposite direction. However, in our spreading terrace, the thin, leading terrace is not in equilibrium with the thicker terrace behind the step. The two regions have different values of $`\mu `$, and a net current is expected.

The relationship between the emitting and receiving regions may take two general forms, as depicted in Fig. 4. In this figure the regions of different chemical potential are labeled 1 and 2. In Fig. 4a the regions lie one above the other; in Fig. 4b, one lies ahead of the step, the other behind it. In view of Eq 23 the current of polymers must be largest in regions with smallest activation free energy $`\mathrm{\Delta }E`$ and largest $`\mathrm{\Delta }\mu `$. In the configuration 4b, the change of $`\mu `$ is concentrated at the step and, hence, so is the current. Similarly, in configuration 4a domain distortion is maximum near the step and the free energy of the transition state may be reduced. For example, the hopping distance for a B block may be reduced, and thus the amount of stretching required for a B block to make the transition should be reduced Helfand . Thus we expect the current to be concentrated at the step in both situations. We shall assume that this is the case in what follows.

The length of the step region where the current is concentrated has a length in the $`x`$ direction of order $`h_0`$. The current $`i`$ per unit width of the step is thus approximately the net current density times $`h_0`$:

$$ih_0\mathrm{\Delta }j_0\frac{2D_0\rho _0e^{\mathrm{\Delta }E_0/kT}}{h_0}\left(e^{\mu _h/kT}e^{\mu _{\mathrm{}}/kT}\right).$$
(24)
The chemical potential for a terrace of width $`\mathrm{\Delta }x`$ and length $`\mathrm{\Delta }y`$ to the first order in $`hh_0`$ is given by

$$\mu =\frac{F_T}{n}|_{\mathrm{\Delta }x,\mathrm{\Delta }y}\mu _0+2g_0\frac{hh_0}{h_0}\gamma _{AB}^{},$$
(25)
where $`\gamma _{AB}^{}=\gamma _{AB}v_0/h_0`$ is interfacial energy normalized to one chain, $`F_T=h\mathrm{\Delta }x\mathrm{\Delta }yF_c/v_0`$ is terrace free energy (see Eqs. 10 \- 11), and $`n=h\mathrm{\Delta }x\mathrm{\Delta }v/v_0`$. If terraces can slide on top of each other (as in Fig. 4a), then the spreading pressure from the lower terrace does not penetrate into the higher terraces, thus leaving thicknesses and chemical potentials of the latter at equilibrium. Then $`\mathrm{\Delta }\mu =2g\gamma _{AB}^{}(hh_0)/h_0`$.

On the other hand, if two terrace are bound by domain walls (as in Fig. 3 and Fig. 4b), then spreading pressure penetrates past the terrace step. From the balance of forces at the step, the thickness $`h`$ of the terrace on the left (Fig. 4b), and thickness $`h_2`$ of the double terrace on the right are related via $`(hh_0)=2(h_2h_0)`$. Then $`\mathrm{\Delta }\mu =g\gamma _{AB}^{}(hh_0)/h_0`$.

In either case, difference in chemical potentials is proportional to the deformation $`hh_0`$ of the lower terrace, and net current (24) between terraces becomes

$$i2\stackrel{~}{g}D_0\rho _0\frac{hh_0}{h_0^2}\frac{\gamma _{AB}^{}}{kT}e^{(\mathrm{\Delta }E_0+\mu _0)/kT},$$
(26)
where $`\stackrel{~}{g}`$ is between $`g_0`$ and $`2g_0`$ depending upon the interaction between terraces.

## V Spreading of the Copolymer Films

It was shown above that the spreading of copolymer films consists of two parts: flow of material inside each terrace, governed by diffusion law (13), and flow of material between terraces, with the flow rate described by (26). At early stages of the spreading, when terrace lateral extension is small, one may expect the dominant contribution to the total rate to come from the flow at terrace steps, while at late times, the dominant role should shift to the flow inside terraces.

To find this time dependence, we turn to the spreading of the lowest terrace. The advancing front of this terrace ($`x=L(t)`$, Fig. 5) is assumed to be pulled with a constant pressure from the wetting forces. As pressure and thickness of a terrace are proportional (eq. (11)), this is equivalent to stating that advancing front is held at a constant thickness $`h_f`$:

$$h|_{x=L(t)}=h_f.$$
(27)
The boundary condition at the step between terraces ($`x=0`$, Fig. 5) arises from the balance between the current (26) entering from the higher terrace, and the diffusion current $`i=\frac{\rho _0}{h_0}D\frac{h}{x}`$ inside the terrace, where $`D`$ is collective diffusion constant entering (13) and defined by (18) or (22). Hence

$$\frac{h}{x}+A(hh_0)=0,\mathrm{at}x=0,$$
(28)
where

$$A=2\stackrel{~}{g}\frac{D_0}{D}\frac{\gamma _{AB}^{}}{kT}\frac{1}{h_0}e^{(\mathrm{\Delta }E_0+\mu _0)/kT}.$$
(29)
Finally, the lateral extension $`L`$ of the terrace is related to the current $`i_f`$ at the advancing front via $`i_f=\rho _0(h_f/h_0)(dL/dt)`$. Hence

$$\frac{dL}{dt}=\frac{D}{h_f}\frac{h}{x}|_{x=L(t)},$$
(30)
with

$$L|_{t=0}=0.$$
(31)
Thus we come to a system of equations (13, 27, 28, 30, 31). Exact solution of this system is not trivial, so we will look only on large and small time asymptotics.

The large time behavior is best analyzed in new coordinate $`z=x/L(t)`$, which eliminates the moving boundary condition (27). Differential equation (13) with moving boundary condition, similarly to the case with stationary boundaries, has only one solution Tihonov . Then, in the $`z`$ coordinate, the system is likely to asymptotically approach some stable state $`h_1(z)`$. Assuming that $`h(z,t)=h_1(z)h_2(z)t^\beta `$, $`L=at^\alpha `$, we get $`\alpha =\beta =1/2`$, and

$$h_1(z)=h_f+(h_0h_f)\left[1\frac{\mathrm{erf}Rz}{\mathrm{erf}R}\right],$$
(32)
$$h_2(z)=(h_0h_f)\frac{\mathrm{exp}(R^2z^2)}{A\sqrt{\pi D}\mathrm{erf}R}\left[1\frac{\mathrm{erfi}Rz}{\mathrm{erfi}R}\right],$$
(33)
$$L=\sqrt{2D\frac{h_0h_f}{h_f}}t^{1/2},$$
(34)
where $`R=\sqrt{\frac{h_0h_f}{2h_f}}`$. Thus velocity of spreading at late times is given by

$$v_{\mathrm{late}}=\sqrt{\frac{D}{2}\frac{h_0h_f}{h_f}}t^{1/2}.$$
(35)
The early times behavior, on the other hand, is characterized by nearly constant thickness $`hh_f`$ of the growing terrace (Fig. 5). Then all the material entering a terrace (eqs. 26 and 28) should go towards lateral extension of that terrace with spreading velocity given by

$`v_{\mathrm{early}}`$ $`=`$ $`{\displaystyle \frac{dL}{dt}}i{\displaystyle \frac{h_0}{\rho _0h_f}}`$
$``$ $`2\stackrel{~}{g}D_0{\displaystyle \frac{h_0h_f}{h_0h_f}}{\displaystyle \frac{\gamma _{AB}^{}}{kT}}e^{(\mathrm{\Delta }E_0+\mu _0)/kT}.`$
The characteristic cross-over time between the two velocities is found by setting $`v_{\mathrm{early}}v_{\mathrm{late}}`$:

$$t_c\frac{1}{2DA^2}\frac{h_f}{h_0h_f}.$$
(37)
## VI Spreading of PS-PMMA Films on SiO2 and Domain Alignment

Spreading rate of PS-PMMA. The theory above predicts the absolute spreading rates and thickness profile to be expected in the experiments of Ref Hahm, . We may thus test our theory by comparing these expected features with the observed ones. In these experiments, a thin film of 80,000 amu PS-PMMA was spin cast onto a silicon substrate that was initially covered with a minor polar solvent. The presence of the solvent resulted in the formation of annuli structures, likely caused by the dewetting of the substrate, and consisting of a tall rim, bare substrate inside the rim, and uniform film outside (Fig. 2). Later, as the solvent evaporated, the material started flowing downwards onto the bare substrate area inside the rim. This flow is presumed to be driven by wetting forces, so that its rate should be expressed by (35) - (V).

The structure of asymmetric, cylinder-forming PS-PMMA copolymers on SiO2 used in this experiment is well known Jaeger ; Hahm ; HahmDefects , and consists of plane lamellar layer on the substrate, followed by a series of layers (terraces) with cylindrical domains inside (Figure 1). To estimate the spreading rate of this system, we turn to the bottommost terrace pulled by the wetting forces. The collective diffusion constant for material flow inside this terrace is given by (18) and contains two separate contributions. On one hand, there’s a friction on the terrace as a whole (second term in Eq. 18). However, the hydrogen bonding between PMMA and SiO2 is very strong ($`\xi 1`$), and contribution from this term to the diffusion constant is small. The first term, on the other hand, describes self-diffusion of fast chains through the matrix of substrate-bound slow chains. Since there are no walls obstructing this flow in the bottommost terrace of PS-PMMA, contribution from this term dominates and $`D_{slid}3(1\tau )^2v_0D_c\gamma _{AB}/k_BT\tau h_0`$.

The fraction $`\tau `$ of chains that are in contact with the substrate can be estimated based on a strongly-stretched approximation Semenov , which predicts that the density of chain ends in a flat brush is $`g(r)r/(h^2r^2)^{1/2}`$. Since strongly stretched approximation neglects deviations of chain paths from their mean, and all paths are presumed to go unidirectionally towards the free end of the brush, we can assume that all chains ending within a certain small distance from the substrate are bound by the substrate. We estimate this distance as the penetration length $`\zeta `$ discussed in (20). This is the length scale within which the polymers behave like unperturbed random walks rather than like strongly stretched chains. For a $`N_A=200`$ PMMA chains with radius of gyration $`R_Aa_0(N_A/6)^{1/2}`$, we find $`\zeta (R_A/h_A)^{4/3}h_A.5h_A`$. Then $`\tau `$ can be found via

$$\tau \left(_{h_A\zeta }^{h_A}g(r)𝑑r\right)/\left(_0^{h_A}g(r)𝑑r\right)85\%.$$
(38)
In fact, this is an upper estimate for $`\tau `$, as density of ends $`g(r)`$ has a singularity near the open end of the brush produced by the neglect of entropy of chain paths. The lower estimate for $`\tau `$ would come from the assumption of constant $`g(r)`$, so that $`\tau a_0/h_A10\%`$, with actual value being somewhere in between. Fortunately, this variation in $`\tau `$ produces only an insignificant variation in the prefactor of $`D_{slid}`$.

Finally, the self-diffusion constant for the chain flow along domain walls $`D_c`$ is also not readily available. Still, its value can be estimated from the analogy between copolymer flow along domain walls and flow of star polymers, whose self-diffusion is given by $`D_{star}D_{arm}\mathrm{exp}(N_{arm}/N_e)`$. Here $`N_e180`$ is average number of monomers per entanglement, and $`D_{arm}`$ is self-diffusion of linear polymers forming each arm. Then the self-diffusion coefficient for a PS star polymer with 60,000 amu arms is $`D_{\mathrm{PS}}.410^{15}`$ m2/s, and the self-diffusion coefficient for a PMMA star polymer with 20,000 amu arms is $`D_{\mathrm{PMMA}}.410^{15}`$ m2/s. Watanabe ; Milhaupt  Because of similarity of these two coefficients, it is reasonable to assume that the self-diffusion coefficient for the whole block copolymer moving along domain walls is also $`D_c.410^{15}`$ m2/s. Then the forced collective diffusion coefficient for the flow along bottommost terrace (18) is $`D_{slid}.710^{15}`$ m2/s, where known values for other parameters have been used Wu ; Hahm ; HelfandTagami : $`\gamma _{AB}.610^3\mathrm{N}/\mathrm{m}`$, $`h_022`$ nm, and $`v_0Na^31.710^{25}`$ m3.

The final parameter influencing spreading dynamics (29) - (V) is the activation free energy $`\mathrm{\Delta }E_0`$. Its magnitude strongly depends on the polymer conformation in the transition state. Two distinctly different conformations are possible Helfand . In the first conformation, the hopping block propagates as a compact blob, thus minimizing its stretching energy at the expense of larger interfacial energy. In the second conformation, the hopping block is stretched between original and destination domains, so that only a fraction of the hopping monomers comes in contact with the incompatible stationary block. Using the calculations by Helfand Helfand , one can show that stretching mechanism dominates in the PS-PMMA films under discussion, and activation energy is roughly $`\mathrm{\Delta }E\chi N_{PMMA}kT`$.

Based on these quantities, we can now estimate the spreading rates of the PS-PMMA film. We made these estimates without any assumptions about the spreading pressure. This pressure depends on the molecular structure and interfacial shape on an atomic scale at the advancing front. It is difficult to estimate it without an accurate characterization of the surface and the monomers at this scale. Accordingly, we used the observed thickness at the front, namely $`h_f`$ instead of the pressure. It is observed to be about $`0.85`$ times the equilibrium terrace thickness $`h_0`$. Using this knowledge, we can apply Eqs 34 \- 37 directly.

First of all, the crossover time $`t_c`$ from the step-limited to the diffusion-limited dynamics, using Eq. 37, is found to be about .25 hours. Hence diffusive flow along the terrace should dominate the spreading dynamics during most of the 24 hours of the experiment reported in Ref Hahm, . Second, the typical distance $`L`$ over which PS-PMMA film spreads in 24 hours at 523 K is predicted to be about 5 microns (Eq. 34). This is indeed close to the observations by Hahm and Sibener Hahm , where rim width was observed to increase by 1 – 3 microns. Finally, the spreading velocity after 24 hours is about 100 nm/hr.

Domain alignment. It was already mentioned that hopping of chains between domains at terrace steps may account for alignment of the cylindrical domains in a spreading film. Specifically, in thin films, domain orientation is influenced by the nearby substrate, with domains usually oriented parallel or perpendicular to the substrate. For example, in PS-PMMA films on SiO2 lying-down cylinders are formed on $`nH`$thick terraces (where $`H`$ is domain period), while vertical cylinders are formed on $`(n+1/2)H`$thick terraces. Spreading of the film and lateral growth of the terraces should then be accompanied by changes in the observed domain patterns. Specifically, here we will discuss the growth and alignment of the lying-down cylindrical domains on the growing terraces.

We first look at domain growth at the steps with horizontal domains on both sides of the step (Fig. 3). This is the case, for example, at the ”double steps“ in the aforementioned PS-PMMA systems. Here term “double step” refers to the steps between $`nH`$ and $`(n+1)H`$thick terraces.

Following section IV, topologically different domains of adjacent terraces may be assumed disconnected from each other (Fig. 3) and the flow of material between terraces is governed by single chain hopping events. The direction of this flow is determined by the difference in terrace chemical potentials, and during spreading goes towards stretched lower terrace. Then, upon arrival to the growing terrace, a copolymer chain has a strong energetic bias to merge into existing domains rather than nucleate new ones. Hence the spreading should lead to the extension of the existing domains ensuring growth of continuous, defect free domain patterns.

The thickness gradient in the step region, and to a smaller extent in the bulk of the terrace, can also account for the alignment of the growing surface domains in the thickness gradient direction. To demonstrate this alignment effect, we first imagine that the cylindrical surface domains meet the transitional region at an oblique angle $`\theta _0`$ to the growth direction, as shown in Figure 3. Since surface domains are not linked to the bulk domains of the thicker terrace, they have to end within the transitional region where the thickness is somewhat greater than that of the lower terrace. Then a small increment of motion of this system consists of two separate processes: first, a small amount of new material is added to the ends of the surface domains; second, these domains, being part of the growing lower terrace, move in the direction of the thickness gradient. Film thinning combined with conservation of volume then leads to a change in the direction of these domains, as described below.

The growth increment extends the domains. This growth, we presume, occurs without domain bending so that $`\theta _0`$ remains unchanged, leaving the domain in the lowest energy unbent state. The translation of the step causes a given element of the film to compress in thickness and to lengthen in the direction of motion in order to conserve its volume (Fig. 6). The resulting deformation resembles that of a piece of dough being rolled under a rolling pin. As domains thin, their angle $`\theta _0`$ decreases. To quantify this decrease we notice that change in thickness of the film from $`h_0`$ to $`h_f`$ leads to its lateral extension $`\mathrm{\Delta }x_f/\mathrm{\Delta }x_0=h_0/h_f`$ along the thickness gradient direction $`x`$ (Fig. 6), while lateral extension $`\mathrm{\Delta }y`$ along the step edge remains unchanged. Then the flattened domains emerge at an angle $`\theta _f`$ with the step edge, which is given by $`\mathrm{tan}\theta _f=\mathrm{\Delta }y_f/\mathrm{\Delta }x_f=\mathrm{tan}\theta _0h_0/h_f`$, and the extruded domains are bent towards the normal of the interface.

This mechanism alone causes domains growing in an arbitrary direction to emerge partially aligned to the normal. Further effects should cause this alignment to increase progressively with time. To demonstrate such an effect, we imagine that domains grow in increments of length $`s`$. During each of these increments domains extend by distance $`s`$, and then propagate the same distance along the gradient direction. The change in thickness that accompanies this motion, results in bending of domain tip from angle $`\theta _0`$ to $`\theta _s`$ with latter given by $`\mathrm{tan}\theta _s=\mathrm{tan}\theta _0h_s/h_0=\mathrm{tan}\theta _0\left(1\frac{\left|\mathrm{\Delta }h\right|}{h}s\right)`$. The relaxation of these bent domains into straight ones would require collective motion of large number of chains and exponentially long times. As this is highly unlikely, the next step of the growth will start with domain ends oriented at a smaller angle $`\theta _s`$ to the normal, leading to progressive decrease in $`\theta _s`$ with each growth step:

$$\frac{\mathrm{\Delta }\mathrm{tan}\theta }{\mathrm{\Delta }s}=\frac{\left|h\right|}{h}\mathrm{tan}\theta ,$$

and, hence,

$$\mathrm{tan}\theta =\mathrm{tan}\theta _0e^{\frac{\left|h\right|}{h}s}.$$
(39)
Thus the growing domains progressively align in the normal to the terrace step direction. The characteristic length for this effect is determined by the size of the step region $`h/|h|24h`$ so that characteristic length for alignment is several domain widths.

To complete this section, we turn to an alternative combination of domain morphologies at terrace steps, with lying-down semicylindrical domains on the growing lower terrace and standing cylinders on the shrinking thicker terrace (Fig. 7). In particular, this setup is often observed in PS-PMMA films Hahm ; Jaeger . While the general dynamics of the spreading should remain the same, some aspects of it may differ. First of all, the flow of material between terraces may proceed faster as some standing cylinders of the thicker terrace may continuously deform into short lying-down cylindrical domains (Fig. 7) that later merge with the existing domains of the lower terrace. While this effect should have little effect on the net rate of spreading, which is primarily determined by collective diffusion in the bottommost terrace (not shown in Fig. 7), it may lead to faster disappearance of the terraces with standing cylindrical domains, in agreement with experimental observations Hahm .

Another difference comes from the hexagonal packing of the vertical domains, as it creates preferred directions for easy merging of vertical domains into horizontal domains. Since these directions do not depend on the wetting forces, there may emerge a misalignment between direction of spreading and direction of terrace growth, which was indeed observed in experiments on annuli structures Hahm .

## VII Discussion

We have shown that layered structure of spreading copolymer films should lead to different type of dynamics compared to macroscopic homopolymer or simple liquid drops. First of all, the spreading of copolymer films is determined by combination of two factors: the rate of chains flow inside constant-thickness terraces of the film and rate of chain flow between terraces. The hydrodynamic flow of chains along the terraces leads to the diffusion type law for spreading (13), and its contribution dominates at late times of spreading, with velocity decreasing as $`vt^{1/2}`$. The hopping mechanism of chain flow at terrace steps, on the other hand, controls the early stages of spreading. The crossover time between the two regimes is given by (37). It is notable that both of these regimes are different from the spreading observed in macroscopic drops, where spreading rate decreases as $`vt^{0.9}`$. However, the $`vt^{1/2}`$ behavior is similar to the generic spreading rate of terraced fluids predicted be de Gennes and Cazabat Cazabat .

Both regimes of spreading are predicted to exhibit strong dependence on chain polymerization. At early times (eq. V), the dominant contribution to the spreading rate comes from the activation free energy $`\mathrm{\Delta }E_0N^\alpha `$, so that $`vN^{2.7}e^{N^\alpha /N_e}`$, as opposed to $`v\eta ^1N^3`$ in macroscopic homopolymer drops. Here $`\alpha =1`$ for chains that hop between domains in a single blob, and $`\alpha =2/3`$ for chains that stretch between domains during hopping Helfand .

Dependence on $`N`$ becomes even more complicated at large times, when diffusion along terraces dominates. It was mentioned above that the spreading rate at large times is primarily determined by chain flow in the bottommost terrace with collective diffusion constant $`D_{slid}`$ (18). For most copolymers, the main contribution is expected from the first term in (18), which describes collective diffusion of fast chains through the matrix of substrate-bound slow chains. Then $`D_{slid}(h_0\tau k_d)^1N^{4/3}e^{N/2N_e}`$. However, in systems with large degree of polymerization or small decrease in mobility $`\xi `$ near substrate, the second term becomes dominant with spreading rate independent of $`N`$: $`D_{slid}(h_0k_s)^1N^0`$. Similarly, spreading rate should be independent of $`N`$ if domains in the bottommost terrace obstruct continuous flow of chains, so that two-liquid approximation is inapplicable and first term in (18) disappears.

The flow of material between domains creates a dynamic process of domain growth at the terrace steps. The chains entering the growing terrace preferentially merge the existing domains in that terrace, leading to their continuous growth. In the previous section we have argued that this growth has preferential directions leading to the dynamic alignment of the domain. Before proceeding with further discussion of this model, it seems useful to compare it with other potential mechanisms for domain alignment, such as static, stress-induced alignments.

In the annuli structures described above, stresses may come from the viscous flow, from the terrace steps, or may be caused by the change in circumference as material flows towards the center. In order to alter domain morphology, these stresses must be comparable to the internal stresses in the film created by domain walls. In films with thickness $`h`$ and surface tension between block $`\gamma _{AB}`$, the internal stresses are roughly $`\sigma _{wall}\gamma _{AB}/h.610^5`$ N/m2.

This is significantly larger than the viscous stress produced by the interpenetrating brushes $`\sigma _{viscous}\eta ^{}v/2\varsigma 210^2`$ N/m2 (compare with (19)), making viscous stress an unlikely candidate for the cause of observed alignment.

The other two stresses, on the other hand, come within one order of internal stresses, and could play some role in the alignment. Specifically, stress from a contact line is $`\sigma _{wetting}\gamma (h_0h_f)/h_0^210^4`$ N/m2 (where $`h_h`$ and $`h_f`$ are thicknesses of relaxed and pulled films), and stress from the change in circumference is $`\sigma _c(\gamma /h_0)(\mathrm{\Delta }R/R)10^4`$ N/m2 (where $`R`$ and $`\mathrm{\Delta }R`$ are radius and width of the rim). Despite the similarity of these stresses to the internal stresses in the film, the alignment caused by them would progressively improve towards the center of the rim, which contradicts the observations on annuli structures Hahm .

The dynamic model of domain growth presented above, on the other hand, predicts better alignment on the grown parts of the terraces, that is in the central parts of the terraces and near the thicker terraces, in agreement with the experiment. This suggests that produced alignment is indeed caused by the dynamic effects accompanying spreading.

Our dynamical model of domain growth may provide a new way of producing aligned domain patterns. In particular, any motion of terrace steps in the phase-separated film should lead to the growth of aligned domains on the growing terraces. One such experiment, besides the aforementioned experiment with annuli structures, is reported in Ref Deepak, . In this experiment, copolymer is deposited onto a corrugated grating rather than a flat substrate. After the solvent has evaporated, the copolymer starts climbing up along the groove walls of the grating leaving aligned domain patterns inside the grooves. As this growth proceeds on time scales much longer than the thickness quantization, the observed alignment may be a result of the terrace growth.

An alternative way to stimulate terrace growth may come via control over initial nucleation of islands and holes in the film. Specifically, under normal conditions, this nucleation happens at random places. However, if one applies periodic external catalyst of the island formation, such as periodic temperature gradient, periodic electric field, or, possibly, a standing sound wave, then the formation and growth of islands and holes may proceed in controlled manner leading to the growth of aligned domains.

## VIII Conclusion

We have shown that spreading of the phase separated copolymer films is different from the spreading of homopolymer and other simple liquids. This difference is driven by domain walls inside copolymer and formation of quantized thickness terraces. As a result, regular isotropic reptation mechanism is replaced by anisotropic flow of chains along domains, and slow activated flow between domains. We have shown that at early stages, activated flow between domains at the terrace steps plays a major role, giving constant velocity of spreading. At the late times, hydrodynamic flow inside terraces dominates, resulting in diffusion type velocity decrease $`vt^{1/2}`$, similar to the generic spreading of terraced fluids predicted be DeGennes and Cazabat Cazabat .

The hopping nature and slow rate of chain flow at the terrace steps were shown to produce aligned and continuous domain patterns on the growing lower terraces. Specifically, as new chains enter the growing terrace, they tend to merge with existing domains thus leading to their continuous extension. Furthermore, thickness gradient in the area of domain growth, i.e. near terrace steps, was shown to align growing domains perpendicular to the steps between terraces, thus leading to formation of aligned domain patterns on the growing terrace.

This work shows that wetting flows offer a powerful means of controlling domains in copolymer films. The effect seems to arise not from flow-induced stress, but rather by flow-induced forced hopping of chains at terrace boundaries. The alignment arises at these boundaries, as in the conventional zone refining of crystals. As suggested above, one may imagine many ways to induce wetting flows in order to achieve a wide variety of morphologies. The potential for this new structuring mechanism remains to be explored.

## Acknowledgments

We gratefully acknowledge Prof. Steven Sibener, Deepak Sundrani and Wendy Zhang for their valuable discussions. This work was supported in part by the National Science Foundation under Award Number DMR-9975533 and in part by its MRSEC Program under Award Number DMR-980859.


format: text

text:
# I=2 pion scattering length with the parametrized fixed point action.

## 1 Introduction

The pion scattering length is an important observable characterizing dynamical effects of the strong interaction and its ab initio calculation on the lattice is an important nonperturbative test of QCD. In the literature several attempts to compute the scattering length with Wilson fermions  and staggered fermions  can be found. However, calculations with a chiral Dirac operator which allows for a better control of the small mass region are still missing. In full QCD the scattering length is a quantity which vanishes in the chiral limit, while it is power divergent in the quenched theory  and eventually a study in the full theory is desirable. In this contribution, we present results of the scattering length calculated using the p-FP action. (Calculations with chirally improved fermions is also underway.) The chiral and scaling properties of the light hadron spectrum using p-FP and chirally improved actions were reported by the collaboration in an earlier publication .

## 2 Method and simulation details

The $`\pi ^+\pi ^+`$ scattering length is extracted using the standard finite volume technique of Lüscher . The energy shift of the two pion state, $`\delta E(L)`$, in a finite volume is related to the scattering length, $`a_0`$, in the infinite volume limit through the relation,

$$\delta E=\frac{4\pi a_0}{M_\pi L^3}\{1+c_1\frac{a_0}{L}+c_2\frac{a_0^2}{L^2}\}+\mathrm{}$$

where $`c_1=2.837297`$ and $`c_2=6.375183`$.

The two pion state is created using a Gaussian smeared source for the quark propagator on a single timeslice. A second smeared source which was displaced in time by one lattice unit was used on small lattices to check for the systematics and the often claimed complications due to Fierz mixing. Point sink operators were used on a single timeslice regardless of the position of the sources.

Cutoff effects were studied in a small volume $`(1.2\mathrm{fm})^3`$ using three different lattice spacings ($`0.08,0.10,0.15\mathrm{fm}`$) and in an intermediate volume of $`(1.8\mathrm{fm})^3`$ with $`\mathrm{a}0.10`$ and $`0.15`$ fm. A wide range of quark masses were studied, with the smallest quark mass resulting in $`M_\pi /M_\rho 0.35`$. In this quenched study, we keep those masses with $`M_\pi L4`$.

## 3 Analysis

We extract the energy shift of the two pion state in a finite volume by taking ratios of two pion correlation functions and single pion correlation functions. We follow the procedure and notation of Ref.  and construct the “direct” and “crossed” contractions of

$$<\underset{\stackrel{}{x}_1}{}𝒪(\stackrel{}{x}_1,t)\underset{\stackrel{}{x}_2}{}𝒪(\stackrel{}{x}_2,t)𝒮_2𝒮_1>$$

where $`𝒪`$ annihilates a $`\pi ^+`$ state and $`𝒮_i`$ are the Gaussian sources for the pions and divide by the two single pion propagators. The large time behaviour of the ratio is given by $`Ze^{\delta Et}`$ in the full theory.

### 3.1 Fitting procedure

Periodic boundary conditions were used in the time direction so that we must take into account a backwards propagating pion. The fitting form of the correlation function was

$$Z\frac{e^{(2M_\pi +\delta E)t}+e^{(2M_\pi +\delta E)(Tt)}+ze^{M_\pi T}}{e^{2M_\pi t}+e^{2M_\pi (Tt)+2e^{M_\pi T}}}.$$

The pion mass was fitted from the pion propagator alone within the same bootstrap sample. The momentum spacings of the pions in the largest volume is quite small which resulted in significant excited state contamination (higher momenta states). In this case, we fit to a sum of exponentials (assuming this is the case) to take this into account. Only the amplitude of the excited state was left as an extra fitting parameter and we vary the energy of this state over a wide range to monitor the effect of the choice of this energy. The results are shown in Fig. 2. It can be seen in this example that the ground state energy is nearly independent of the value of the energy of the excited state, however, with much larger errors than with a single exponential fit.

## 4 Results/Discussion

The effect of displacing the second source by one timeslice is shown in Fig. 1. The asymptotic value agrees within errors but there is less excited state contamination when the two pions are not created on the same timeslice. The effect of the second source is being investigated on the largest lattice where this contamination is the source of the large errorbars as was mentioned in the previous section.

The cutoff effects of $`\delta E`$ in the small volume, $`(1.2\mathrm{fm})^3`$, is shown in Fig. 3. We observe no significant lattice artifacts for pseudo-scalars lighter than $`1`$ GeV. A similar analysis in the large volume for $`a_0M_\pi `$ (Fig. 4) confirms that the lattice artifacts are well under control and that the results from $`a0.15\mathrm{fm}`$ and $`L=2.4\mathrm{fm}`$ are very close to the continuum and infinite volume limit.

In Fig. 5, we plot $`a_0/M_\pi `$, against $`M_\pi `$ in physical units. We also include in this plot some recent results from other groups as well as quenched chiral perturbation theory  and chiral perturbation theory predictions . The mass range of the perturbative calculations shown in this figure were chosen for illustration and does not suggest convergence of the series. The effect of the $`\delta ^2/M_\pi ^2`$ term in the quenched theory is still not entirely visible at the present quark masses. The results also suggest that pseudo-scalars with masses less than $`300\mathrm{MeV}`$ may be needed for a reliable extrapolation to the physical pion mass. It remains a challenge to perform a full QCD simulation with such light quarks.

## ACKNOWLEDGEMENTS

K.J. would like to thank G. Colangelo and J. Gasser for providing details of their calculations and insightful discussions. The simulations were done on the Hitachi SR8000 at the Leibniz Rechenzentrum in Munich and at the Swiss Center for Scientific Computing in Manno. This work was partially supported by the European Community’s Human Potential Programme under contract HPRN-CT-2000-00145.


text:
# TeV Blazars — Observations and Models

## 1 Introduction – TeV Blazars

The EGRET (Energetic Gamma Ray Experiment Telescope) detector on board of the Compton Gamma-Ray Observatory discovered strong MeV $`\gamma `$-ray emission from 66 blazars, mainly from Flat Spectrum Radio Quasars and Unidentified Flat Spectrum Radio Sources Hart:99 . Ground-based Cherenkov telescopes discovered TeV $`\gamma `$-ray emission from 6 blazars, 4 of which have not been detected by EGRET. All the established TeV detections belong to the class of BL Lac objects, blazars with relatively low luminosity but with Spectral Energy Distributions (SEDs) that peak at extremely high energies (see Table 1).

The large detection area of Cherenkov telescopes of $``$105 m2 makes it possible to sample the $`\gamma `$-ray lightcurves with a time resolution of several minutes. Large amplitude flux variability on 30 min time scales implies that the TeV emission originates from a small region very near the supermassive black hole Gaid:96 ; Kraw:01a . The relativistic Doppler Factor of the emitting plasma is defined as $`\delta _\mathrm{j}^1=\mathrm{\Gamma }(1\beta \mathrm{cos}(\theta ))`$, where $`\mathrm{\Gamma }`$ is the bulk Lorentz factor of the emitting plasma, $`\beta `$ is its bulk velocity in units of the speed of light, and $`\theta `$ is the angle between jet axis and the line of sight as measured in the observer frame. The requirement that the emitting volume is optically thin for the TeV $`\gamma `$-radiation, together with assumptions about the optical to UV energy spectrum of co-spatially emitted synchrotron emission gives a lower limit of $`\delta _\mathrm{j}>`$ 9 for the emitting plasma Gaid:96 . The Lorentz boost of the observed jet luminosity by the factor $`\delta _\mathrm{j}^{\mathrm{\hspace{0.17em}\hspace{0.17em}4}}`$ makes the rather weak jet emission readily detectable.
The main scientific drivers for broadband observations of TeV blazars are: (i) the study of the matter, energy content, and structure of jets from Active Galactic Nuclei (AGNs), (ii) the time resolved investigation of particle acceleration processes, (iii) the research on the connection between the accreting supermassive black hole system and the jet parameters, and (iv) the measurement of the total infrared/optical luminosity of the universe.
In this paper we give an update on the status of the observations (Sect. 2) and the four major scientific topics (Sects. 36). We conclude with an outlook in Sect. 7. We will focus on observations and modeling of TeV blazars (blazars detected by TeV instruments). Reviews on observations and models of MeV/GeV blazars (detected by EGRET) have been given in Siko:01b ; Copp:99 ; recent overviews of TeV astronomy have been given in Ong:03 ; Buck:02 .

## 2 Observations

TeV $`\gamma `$-ray observations have achieved a state of maturity with good agreement in fluxes and reasonable agreement in energy spectra between different experiments. The measurements of the steady flux from the Crab Nebula with the CAT, HEGRA, and Whipple experiments lie within 18.5% from the common mean comparable to the systematic error of about 22% quoted by the groups Piro:01 ; Ahar:00 ; Hill:98 . The spectral indices agree to within 0.2 from the common mean, a difference that is larger than the systematic errors of between 0.04 and 0.06 estimated by the 3 groups. Systematic error on flux and spectral changes from one day to the next are estimated to be about 10% and 0.05, respectively.
High-quality TeV $`\gamma `$-ray energy spectra have now been measured for Mrk 421, Mrk 501, H 1426+428, and 1ES 1959+650. Fig. 1 compares the X-ray and TeV $`\gamma `$-ray SEDs of these 4 sources. Compared to the other 3 blazars, the low-energy (synchrotron) component of Mrk 421 peaks at lower frequencies. Comparison of the high-energy TeV $`\gamma `$-ray energy spectra is not straight forward as it depends on the uncertain extent of extragalactic absorption (see Sect. 6). Spectral variability has been established for the two strongest ones, Mrk 421 Kren:02 ; Ahar:02b  and Mrk 501 Djan:99 ; Ahar:01 . The significance of measurements of spectral changes can not be overestimated. Spectral changes (like flux changes) are free of the still substantially uncertain extent of extragalactic $`\gamma `$-ray absorption. Modeling of simultaneously taken X-ray and $`\gamma `$-ray data with well measured spectral changes in both bands has not yet been performed but has the potential of breaking basically all model degeneracies.
Several extensive multiwavelength campaigns have been performed. The correlation of X-ray and TeV $`\gamma `$-ray fluxes had first been seen in Mrk 421 data Buck:96 ; Taka:96 . In the meantime, correlations with high statistical significances have been observed from Mrk 501 Kraw:00 ; Samb:00  and from Mrk 421 Taka:00 ; Foss:03 . A more recent campaign on the TeV blazar 1ES 1959+650 shows an orphan $`\gamma `$-ray flare without an X-ray or optical counterpart Kraw:03 . The orphan flare shows that the simplest models (1-zone Synchrotron Self-Compton models with a randomly oriented magnetic fields and isotropic electron distributions) are not able to account for the data.

Simultaneous measurements of blazar energy spectra with imaging Cherenkov telescopes and solar array Cherenkov detectors have the potential of producing highly interesting results. First lightcurves from STACEE Boon:02  and a first energy spectrum from CELESTE Piro:03  have been reported recently.

## 3 Emission Models and Implications for the Jet Structure

Synchrotron-Compton Models: In Synchrotron Compton (SC) models of TeV blazars the radio to X-ray emission is produced as synchrotron radiation from a population of non-thermal electrons (and possibly positrons). The same electron population emits $`\gamma `$-rays through Inverse Compton processes by electrons scattering synchrotron photons in so-called Synchrotron Self-Compton (SSC) models or with “external” photons that originate outside the jet in so-called External Compton (EC) models (see Siko:01b  and references therein). In the case of BL Lac objects, the lack of strong emission lines is commonly taken as evidence that ambient photon fields are not important, and that SSC models rather than EC models are more likely to explain the data.
A method to distinguish between SSC and EC models is to measure the time lag $`\mathrm{\Delta }t_{\mathrm{S}\mathrm{IC}}`$ between X-ray and $`\gamma `$-ray flux variability Copp:99b . In SSC models one expects a time-lag of the order of $`R`$ $`c^1`$ $`\delta _\mathrm{j}^1`$ with $`R`$ the radius of the emission volume and $`c`$ the speed of light. The time lag originates from the fact that a population of freshly accelerated electrons immediately radiates synchrotron emission; in contrast, the Inverse Compton component peaks only after the new synchrotron seed photons had time to propagate through the source. No such time lag has yet been observed with high statistical confidence.
TeV blazar data have been modeled intensively with a variety of codes. The more simpler ones assume 1-zone homogeneous emission volumes and use phenomenological electron energy spectra to fit the broadband data (e.g. Inou:96 ; Taka:00 ; Tave:01 ). Modeling with time-dependent self-consistently evolved codes (e.g. Copp:92 ; Mast:97 ; Kraw:02 ; Boet:02 ) assures that the energy spectra are physically realizable from initial acceleration energy spectra. Time dependent modeling has to be used whenever the flare duration is comparable or shorter than the longest time scale of the microprocesses that modify the electron energy spectra (e.g. acceleration energy gains, radiative energy gains and losses).
Several conclusions can be drawn from the SSC and EC modeling performed by various authors: (i) relatively simple 1-component and 2-component SSC models are able to explain most multiwavelength data; (ii) even for the nearest blazars with $`z`$ 0.03, extragalactic absorption has to be taken into account and substantially modifies the regions in the parameter space that correspond to acceptable solutions (e.g. Ahar:99a ; Bedn:99 ); (iii) the data prefer high jet Doppler factor of between 50 and 100 Mast:97 ; Kraw:01a ; (iv) in the jet frame the energy density of non-thermal particles is much higher than the energy density of the magnetic field (the jet plasma is a low $`\sigma `$-plasma) (e.g. Kino:02 ; Kraw:02 ); (v) only a small fraction ($`<10^2`$) of the energy transported by the jet is converted into radiation (e.g. Kino:02 ; Kraw:02 ).
The last two conclusions constrain the structure of the jet at its base. While electromagnetic models of jet formation produce high-$`\sigma `$ jets with large electromagnetic field energy to particle energy densities, SSC models indicate low-$`\sigma `$ jets. The low radiative efficiency indicates that the sub-parsec jet can smoothly evolve into the parsec-scale jet, without the need for plasma re-collimation or re-acceleration that would have to follow a more efficient conversion of bulk plasma energy into random particle energy and radiation (see also Celo:93 ).
Hadronic Models: In hadronic models the continuum emission is explained by hadronic interactions of a highly relativistic baryonic outflow which sweeps up ambient matter Pohl:00 , by interactions of high-energy protons with gas clouds moving across the jet Dar:97 , or by interactions of Ultra High Energy protons with ambient photons Mann:93 , with the jet magnetic field Ahar:00a , or with both Muec:02 . Although attractive since they involve Ultra High Energy Cosmic Ray acceleration, hadronic models have difficulties to explain the good X-ray/TeV $`\gamma `$-ray flux correlations found for Mrk 421 and Mrk 501. Furthermore, little work has yet been done to see if proton models can explain more than single snapshots of data.

## 4 Particle Acceleration

Among the many particle acceleration mechanisms that are discussed, Fermi acceleration at strong Magnetohydrodynamic (MHD) shocks has received most attention. Indeed, the data are largely consistent with an electron acceleration energy spectrum $`dN/dE`$ $`E^p`$ with $`p`$ between 2 and 2.23, expected for non-relativistic and ultra-relativistic shocks, respectively (see Kirk:99a  and references therein). Based on a time-dependent SC code Kirk & Mastichiadis (1999) have shown that Fermi energy gains and radiative energy losses can result in characteristic hardness–intensity correlations during individual flares Kirk:99 . X-ray and recently also TeV $`\gamma `$-ray data do show the predicted signatures: during flares the sources seem to go through clockwise and anti-clockwise loops in the X-ray or $`\gamma `$-ray hardness–intensity planes (e.g. Taka:96 ; Horn:02 ). However, the same source shows a wide range of different signatures (e.g. Taka:00 ), indicating that the relative length of the characteristic time scales of particle acceleration and radiative cooling change from flare to flare. The lack of a prevailing signature has cast doubts on whether “cooling loops” and “particle acceleration loops” have really been observed.
Most SC models of TeV blazars invoke a minimum Lorentz factor $`\gamma _{\mathrm{min}}`$ of accelerated particles on the order of $`10^5`$ to account for the data (e.g. Kraw:02  and references therein). Such high $`\gamma _{\mathrm{min}}`$-values could either result from a yet unknown pre-acceleration mechanism, or, from an upstream bulk Lorentz factor on the order of $`10^5`$. A low-entropy, extremely relativistic upstream plasma could be detectable by its Inverse Compton emission if exposed to a suitable external radiation field.

## 5 Black Hole - Jet Connection

From the theoretical side, the connection between black hole and jet properties is far from being understood. Purely electromagnetic Blan:76 ; Love:76  and GR-electromagnetic models Blan:77  of jet formation face the well known $`\sigma `$-problem, namely that they predict Poynting flux dominated energy transport while SC models indicate strongly particle energy dominated jets. Hydrodynamic, MHD, and GR-MHD scenarios of jet formation (e.g. Bege:84 ; Koid:99 ) face the difficulty that they predict jets with bulk Lorentz factors $`\mathrm{\Gamma }<10`$ while SC models favor $`\mathrm{\Gamma }10`$.
From the observational side, first studies of the black hole–jet connection are just becoming possible based on high-quality black hole mass estimates from galactic stellar velocity dispersion measurements Falo:02 ; Bart:03 . The studies performed so far do not indicate any correlation of the black hole mass and the peak location or peak luminosity of the low-energy or high-energy SEDs, the time scale of flux variability (see Fig. 2), or the flare duty cycle. The correlations may be masked by variation of parameters like jet viewing angle, jet magnetic field, or intensity of the ambient photon field Bart:03 ; Kraw:03 . Alternatively, the black hole mass estimators may not work in the case of blazars, as discussed in Bart:03 .

## 6 Measurement of the Cosmic Infrared Background

Extragalactic absorption of TeV $`\gamma `$-rays by the Cosmic Infrared Background (CIB) and the Cosmic Optical Background (COB) in $`\gamma _{\mathrm{TeV}}+\gamma _{\mathrm{CIB}/\mathrm{COB}}`$ $`e^+e^{}`$ pairproduction processes Niki:62 ; Goul:65 ; Stec:92 ; Prim:01 ; Jage:02  substantially complicates the study of the astrophysics of blazar jets. However, the absorption effect allows us to use TeV $`\gamma `$-ray observations to measure the energy spectrum of the CIB and the COB. These backgrounds are extremely interesting owing to their close relation to the total electromagnetic luminosity of the universe since decoupling of matter and radiation $``$300,000 years after the big bang.
The TeV energy spectra of extragalactic sources carry the imprint of extragalactic absorption in the form of “characteristic high-energy cutoffs” in the 30 GeV to 10 TeV energy range. The measurements of these cutoffs makes it possible to constrain the energy spectrum of the infrared to optical background radiation. These constraints are very valuable, as direct measurements are plagued by large systematic errors, owing to strong foreground emission from our solar system and the Milky Way (e.g. Haus:01 ). While the spectral measurement of the nearby blazars Mrk 421 and Mrk 501 ($`z`$ 0.03) have been used to set upper limits on the CIB intensity in the 1-20 $`\mu `$m wavelength range Bill:98 , the recent measurement of the energy spectrum from the more distant blazar H 1426+428 ($`z=`$ 0.129) has given the first demonstration of the full potential of TeV $`\gamma `$-ray observations to decide between different background models Ahar:03a .
Further progress may best be accomplished with a two-pronged approach: (i) Detailed source modeling makes it possible to disentangle, at least partially, source inherent high-energy cutoffs from those caused by extragalactic extinction and thus to improve on the constraints from each single source Copp:99 ; Kraw:02 ; (ii) since the spectral cutoffs from extragalactic absorption are expected to occur at energies which depend only on source redshift, a statistical analysis of the constraints from many sources can be used to test and strengthen the constraints from individual sources.

## 7 Outlook

Over the last several years, the field of TeV blazar studies has benefited greatly from detections of new sources and from very intensive multiwavelength campaigns. Further progress may derive from extending the time scales of the broadband observations. High-accuracy measurement of Gamma-Ray Burst (GRB) afterglows is boosting the theoretical progress in the case of GRBs. In the same way, detection of increased flux levels in the radio to UV bands following major blazar flaring phases may prove crucial for pinning down the jet properties. For this purpose regular broadband monitoring is required over time scales longer than the typical 10 days of current multiwavelength campaigns. As the quality of the data improves, it now becomes apparent that simple one-zone SC model fit snapshots of the data, but are not able to give us deep insights into the jet structure and the nature of the particle acceleration processes. Most present theoretical work focuses on the emission mechanism and largely neglects the flare origin and the development of the plasma that supports the non-thermal particles during the flare. A notable exception is the recent work concerning the internal shock model Spad:01 ; Tani:03 . In future work, combining 3-dimensional MHD simulations of flares with kinetic simulations of the non-thermal particle population would make it possible to explore more complex observational signatures than possible with present codes.
In view of the next-generation Cherenkov telescopes CANGAROO III, H.E.S.S., MAGIC, and VERITAS presently under construction, the upcoming launches of the X-ray observatories SWIFT and ASTRO-E2, and the $`\gamma `$-ray observatories AGILE and GLAST to be launched in 2005 and 2006, respectively, we can expect a wealth of exciting blazar results in the next 5 years.
Acknowledgements: The author gratefully acknowledges support by NASA through the grant NASA NAG5-12974.


text:
# Influence of Planets on Parent Stars: Angular Momentum

## 1 Introduction

The evolution of planets inside the envelopes of evolved stars was studied long before the detection of extrasolar planets (e.g., Eggleton 1978; Livio 1982; Livio & Soker 1984; Harpaz & Soker 1994). This is also true for the study of the influence of planets and brown dwarfs on their parent stars. One of the main aspects of motivation was to explain the formation of some moderate axisymmetrical, i.e., elliptical, planetary nebulae (PNs; bipolar and extreme elliptical PNs require stellar binary companions; Soker 1997; 1998b; 2002b). The immediate progenitors of PNs, which are asymptotic giant branch (AGB) stars, are expected to rotate extremely slowly. Even an earth-like planet can spin-up its parent star if it enters the envelope as the star is about to leave the AGB; Jupiter-like planets may spin-up their parent upper-AGB stars by several orders of magnitude (section 5). Planets were also suggested to enhance mass loss rate from stars on the red giant branch (RGB), such that when the star turns to the horizontal branch (HB) on the HR diagram, it will have low-mass envelope and be a blue HB star (section 4). As the planet is destructed inside the giant envelope, it may change its composition, e.g., deposits lithium. Other effects due to planets were also suggested. For example, planets outside the envelope may accrete mass and have bursts of emission (Struck, Cohanim, & Willson 2002, 2004), and may influence SiO maser emission (Struck-Marcell 1988; Struck et al. 2004), as well as H2O masers (Rudnitskij 2002). I think it will be very difficult to find planets accreting from AGB stars: in order to accrete at a high rate, or to strongly interact with the AGB wind, a low mass companion must be close to the AGB surface; but if it is too close, tidal interaction will cause it to spiral in during a short period of time (see section 3 below). So only a small fraction of systems will go through this stage, and even then the luminosity of the accreting planet will be very low. Planets spiraling in should accrete at a high rate when close to the envelope, and as was found by Struck et al. (2004), may form an accretion disk (I think Struck et al.  overestimate the accretion rate of angular momentum because they use a 2-dimensional rather than 3-dimensional numerical code.) The plant then may blow two jets. The idea that planets may blow jets was raised in the context of young planetary systems, where the planets accrete from the proto planetary system disk (Quillen & Trilling 1998). In any case, my view is that there is no need to invoke direct influence of planets to explain SiO maser variabilities, or replace pulsation with planets to explain oscillations of AGB stars as suggested by Berlioz-Arthaud (2003); the pulsation and strong convection of AGB stars, with possible weak magnetic activity, can account for these variabilities in AGB stars. (I think the results of some of these papers, which deal with planets around AGB stars are also questionable because they ignore the importance of tidal forces when the planets are close to the AGB envelopes.) I think, though, that planets may have an indirect influence on them (see below).

The main point I would like to make is that planets may have pronounced effects on their parent stars, and by studying them we may solve some puzzles in stellar evolution. The first issue to be considered in this respect is the definition of a single star. From the ‘official’ point of view, a star and a planet are a single-star system. However, for someone simulating stellar evolution, this may not be the case. If the evolutionary code, or other means of calculating the evolution, do not include relevant effects of planets around a specific object, e.g., spinning up the envelope and enhancing mass loss rate, causing mixing at the core-envelope interface, depositing fresh hydrogen-rich material to the nuclear burning shell, etc., wrong results will be obtained. For this specific case, the star is not really a single star, although it has no stellar-companions. In the past I have suggested that a non-single star will be one for which one of its relevant properties, e.g., angular momentum, hydrogen abundance in the core, etc., is determined by a gravitationally bound object.

## 2 Nonlinear Effects of Planets

Enhanced equatorial mass loss rate due to centrifugal forces is not important for envelope spun-up by planets. In addition, the gravitational energy of mass accreted onto the planets are in general not sufficient in order to change the structure of the descendant PN. Therefore, there is a need for non-linear effects. By nonlinear effects I refer to effects that are very sensitive to the tiny effect of planets.

### 2.1 Excitation of Waves in Common Envelopes

One such nonlinear effect is the excitation of p-waves (Soker 1992a, 1993) and g-waves (Soker 1992b) during a common envelope phase. While inside the convective envelope of an AGB or RGB star, a companion will excite p-waves which propagate outward with increasing amplitude, mainly in the equatorial plane. For a planet of mass $`10M_J`$, where $`M_J`$ is the Jupiter mass, the surface pressure amplitude $`P^{}`$ is (Soker 1993, eq. 4.1) $`P^{}/P0.4`$ for a secondary at $`a_20.1R_{}`$, where $`P`$ is the average surface pressure, $`a_2`$ is the location of the companion from the center, and $`R_{}`$ is the AGB stellar radius. The perturbation increases linearly with the companion mass (for low mass companions), and increases somewhat as $`a_2`$ decreases (depending on convective viscosity). The amplitude is much larger in the equatorial plane than in the polar directions. Such excited non-radial oscillation can enhance mass loss rate in the equatorial plane.

### 2.2 Destruction of Planets in the Envelope

A study of the fate of planets in the envelope of AGB stars was conducted by Livio & Soker (1984). They assumed that the planet accretes from the envelope at the Bondi-Hoyle accretion rate. However, it is possible that the planet swells as a result of this accretion and does not accrete much, like low-mass main sequence stars do (Hjellming & Taam 1991). Planets may also be evaporated, in particular when they reach the place in the envelope where the envelope’s temperature exceeds the planet’s virial temperature. For RGB and AGB stars, the orbital separation of a planet from the core where fast evaporation starts is (Soker 1998a)

$$a_2(\mathrm{evaporation})10\left(\frac{M_p}{M_J}\right)^1R_{},$$
(1)
where $`M_p`$ is the planet’s mass. The cool and dense evaporated material is still of low entropy, and fraction of it may spiral-in to the core. More massive planets than Jupiter will survive farther in, until they reach a radius where Roche lobe overflow (RLOF) occurs. For a planet of radius $`R_p=0.1\eta R_{}`$, RLOF occurs when the orbital separation from the core is (see Soker 1998a)

$$a_2(\mathrm{RLOF})1.7\eta \left(\frac{M_p}{M_J}\right)^{1/3}R_{}.$$
(2)
The supplement of the destructed planet (or brown dwarf) material to the core and around it may have several effects. First, the low entropy material can absorb heat, and may reduce for a short period of time the stellar luminosity (Harpaz & Soker 1994). If the material reaches the core, or close to it, the release of gravitational energy and nuclear burning of the fresh hydrogen-rich material may lead to stellar expansion and enhanced mass loss rate (Siess & Livio 1999a,b). Recently, Retter & Marom (2003) proposed that three planets which deposited gravitational, and then nuclear, energy into their parent RGB star, along the calculations of Siess & Livio (1999b), with about a month delay from one planet to the next, can account for the erupting of V838 Mon. Second, the high specific angular momentum of the planet’s (or brown dwarf) material may lead to the formation of an accretion disk around the core; such disk can launch two jets (Soker 1996b).

### 2.3 Spinning-up RGB and AGB Envelopes and Magnetic Activity

The spinning-up of RGB and AGB envelopes was discussed in several papers (see Soker 2001b). Basically, RGB and AGB stars slow down rapidly as they lose mass (see figures 1 and 2 in Soker 2001b). Planets and brown dwarfs can then tidally interact with the expanding star (see next section), enter the envelope, and deposit their orbital angular momentum to spin-up the envelop by a factor up to $`10^4`$, depending on the mass of the planet. However, the envelope will still spin much below the Keplerian speed on the equator (the break-up speed). To influence the mass loss process, a non-linear effect must be incorporated. Such an effect appears in the cool magnetic spots model (Soker 1998c; Soker & Clayton 1999), where it is assumed that a weak magnetic field forms cool stellar spots, which facilitate the formation of dust closer to the stellar surface, hence increasing the mass loss rate. If magnetic spots, due to the dynamo activity, are formed mainly near the equator, an enhanced equatorial mass loss is obtained. The weak magnetic field is assumed to be formed by the strong convection in AGB and RGB stars, together with a very slow rotation, which mainly serves for defining the symmetry axis of the magnetic activity.

## 3 Tidal Interaction

RGB and AGB stars tidally interact with their companion before RLOF or common envelope occur. Most of the orbital angular momentum of the companion is deposited to the envelope before the onset of the common envelope phase. Therefore, the star can lose substantial fraction of its mass before the common envelope phase starts (Soker 2002c; Soker & Harpaz 2003). During the evolution along the RGB or AGB, the star expands, and tidal interaction strength increases steeply with the giant radius. It is mandatory to take into account tidal interaction, with substellar or stellar companions, in studying the interaction of giants with their close companions. The tidal interaction, for example, is crucial in determining the fate of the earth as the Sun becomes an AGB (Rybicki & Denis 2001; see also Rasio et al. 1996). Since I find that too many papers dealing with planets around AGB and RGB stars ignore tidal interaction, hence overestimating the importance of the effects they study, I devote a section to this subject.

In Soker (1996a) I found the maximum orbital separation at which tidal interaction is significant. For planets and brown dwarfs, because they cannot bring the envelope to corotate, this radius is the radius below which they will spiral into the envelope of the giant. This maximum radius is given by

$`a_{\mathrm{max}}=3.9R\left({\displaystyle \frac{\tau _{\mathrm{ev}}}{6\times 10^5\mathrm{yr}}}\right)^{1/8}\left({\displaystyle \frac{L}{2000L_{}}}\right)^{1/24}\left({\displaystyle \frac{R}{200R_{}}}\right)^{1/12}`$
$`\times \left({\displaystyle \frac{M_{\mathrm{env}}}{0.5M_1}}\right)^{1/8}\left({\displaystyle \frac{M_{\mathrm{env}}}{0.5M_{}}}\right)^{1/24}\left({\displaystyle \frac{M_2}{0.01M_1}}\right)^{1/8},`$ (3)
where $`L`$, $`R`$, and $`M_1`$ are the luminosity, radius, and mass of the giant (RGB or AGB star), respectively, $`M_{\mathrm{env}}`$ is the giant’s envelope mass, and $`\tau _{\mathrm{ev}}`$ is the evolution time on the upper AGB or RGB.

## 4 Influencing the Horizontal Branch Morphology

After RGB stars ignite helium in their core they move to the horizontal branch (HB) on the HR diagram. Stars which lose more mass on the RGB become bluer (hotter) HB stars. The HB morphologies, i.e., the distribution of stars on the HB of a stellar system, differs substantially from one globular cluster to another. It has long been known that metallicity is the main factor which determines the location of HB stars on the HR diagram. Metallicity is the first parameter. For more than 30 years, though, it has been clear that another factor is required to explain the variation in HB morphologies among globular clusters with similar metallicity (see reviews by Fusi Pecci & Bellazzini 1997; de Boer 1999). This factor is termed the second parameter of the HB. It seems that stellar companions alone cannot be the second parameter (e.g., Rich et al. 1997), nor any other single factor which has been examined (e.g., Ferraro et al. 1997 and references therein). I think that the presence of low mass stars and of planets (or brown dwarfs) is the main second parameter factor (but probably not the only one), with planets occurring more frequently (Soker 1998a).

In recent years it has become a common view that the second parameter determines the HB morphology by regulating the mass loss on the RGB (e.g., Dorman, Rood, & O’Connell 1993; D’Cruz et al. 1996, 2000; Whitney et al. 1998; Catelan 2000). According to this view, the extreme HB (EHB) stars, for example, lose almost all their envelope while still on the RGB (Dorman et al. 1993; D’Cruz et al. 1996); mass loss on the HB itself can’t account for EHB stars (Vink 2003). It is thought by many people that rotation has a connection with the second parameter through its role in determining the mass loss on the RGB, directly or indirectly. I agree with this assertion, and further claim that the source of angular momentum in many cases is the interaction with a planet (Soker & Harpaz 2000; Livio & Soker 2002). Sweigart & Catelan (1998 Moehler, Sweigart, & Catelan 1999) claim that mass loss on the RGB by itself cannot be the second parameter, and it should be supplied by another process, e.g., rotation, or helium mixing, which requires rotation as well. They term the addition of such a process a “noncanonical scenario”. Behr et al. (2000b) find the second parameter problem to be connected with rotation, and note that single star evolution cannot explain the observed rotation of HB stars, even when fast core rotations are considered. The rich variety of HB morphologies (e.g., Catelan et al. 1998) suggests that there is a rich spectrum in the factor(s) behind the second parameter.

After presenting the idea that planets are the main factor in the second parameter (Soker 1998a), I farther explored this idea in three papers. In Soker & Harpaz (2000) we analyzed the angular momentum evolution from the RGB to the HB and along the HB. Using rotation velocities for stars in the globular cluster M13 (Behr et al. 2000b; similar distribution of rotation is in the globular cluster M15; Behr, Cohen, & McCarthy 2000a), we found that the required angular momentum for the fast rotators is up to $`13`$ orders of magnitude (depending on some assumptions) larger than that of the sun. Planets of masses up to five times Jupiter’s mass and up to an initial orbital separation of $`2`$ AU are sufficient to spin-up the RGB progenitors of most of these fast rotators. Other stars have been spun-up by brown dwarfs or low-mass main sequence stars. Our results show that the fast rotating HB stars have been probably spun-up by planets, brown dwarfs, or low-mass main sequence stars, while they evolved on the RGB.

Support of the planet-second parameter idea comes from sdB binary systems. The field sdB stars are post-RGB stars, which have lost most of their envelope, and are parallel to EHB (very hot) stars in globular clusters (Stark & Wade 2003; Vink 2003). The class of objects named EC14026, which have sdB stars and low mass main sequence companions, was discussed by, e.g., Kilkenny et al. (1997), Koen et al. (1997), and Koen et al. (1998), and their relation to EHB stars by (Bono & Cassisi 1999). PG 1336-018, for example, has a secondary of mass $`0.15M_{}`$ with an orbital period of $`0.1`$ days (Kilkenny et al. 1998). Maxted et al. (2000) find the orbital periods and minimum companion masses of two sdB stars: $`0.63M_{}`$ and $`8.33`$ days for PG 0940+068, and $`0.09M_{}`$ and $`0.599`$ days for PG 1247+554. For others, the companion, if it exists, is limited to a spectral type M0 or later (e.g., PG 1605+072, Koen et al. 1998; PG 1047+003, O’Donoghue et al. 1998). For these systems, I suggest that the companion may be a brown dwarf or a massive planet as well. Allard et al. (1994) estimate that $`60\%`$ of hot B subdwarfs have binary stellar companions. Here again, the stars with no stellar companions may have a substellar companion. Green et al. (1998) argued that their “investigations in open clusters and the field strongly suggest that most metal-rich BHB \[blue HB\] stars are in binaries”. In many cases in the proposed scenario for the formation of blue-HB stars the substellar companion will not survive the evolution (Soker 1998a). Therefore, although most of the envelope was lost as a result of the interaction with the substellar companion, it does not exist anymore.

Soker & Hadar (2001) study the correlations between the the HB morphology and some other properties of globular clusters. Strengthening previous results, we find that a general correlation exists only between HB morphology and metallicity. Correlations with other properties, e.g., central density and core radius, exist only for globular clusters within a narrow metallicity range. We conjecture that the lack of correlations with present properties of globular clusters (besides metallicity), is because the variation of the HB morphologies between globular clusters having similar metallicities is caused by a process, or processes, whose effect was determined at the formation time of globular clusters. This process (or processes) is the second parameter. The ‘planet second parameter’ model fits this conjecture. This is because the processes which determine the presence of planets and their properties occur during the formation epoch of the star and its circumstellar disk.

Even if planets are not the second parameter, still, close Jupiter-like planets must increase the mass loss from RGB stars (Livio & Soker 2002). Therefore, if a large fraction of sun-like stars in a group, like in a globular cluster, possesses close planets, many of the HB stars in this group will be blue (hot). Gilliland et al. (2000) used the HST for 8.3 days and found no planets around main sequence stars in the globular cluster 47 Tucanae (NGC 104). This globular cluster contains only a few blue HB stars (Rich et al. 1997; Moehler, Landsman & Dorman 2000), and therefore I do not expect the stars in this globular cluster to have massive and close planets around them. Therefore, 47 Tuc was a bad choice for planets search. (In my talk at the meeting I clearly made that point. Despite that, in a paper to these proceedings \[posted on astro-ph\] a new \[ground\] search for planets in 47 Tuc is presented, without referring to my point.) In the globular cluster M4, on the other hand, close to half of the HB stars are blue (Harris 1996), and indeed, the presence of a planet in this globular cluster was confirmed recently (Sigurdsson et al. 2003). Ferdman et al. (2003) found no candidates for planet transits in their HST study of M4. I suspect their study was not sensitive enough to eliminate the possibility of planets in M4. I encourage more HST observations of M4 and other globular clusters having large population of blue-HB stars.

## 5 Influencing the Mass Loss Geometry

As mentioned in section 1, in several papers I examined the possible role of planets in influencing the mass loss geometry from AGB stars, with the goal of explaining moderate elliptical PNs. The main process (section 2.3) is the dynamo amplification of magnetic fields in slowly rotating AGB stars. This leads to the formation of cool magnetic spots, which enhance dust formation and mass loss, mainly from the equatorial plane. The most recent papers, where more references can be found, are Soker (2001b) and Livio & Soker (2002).

In Soker (2001b) I examine the implications of the recently found extrasolar planets on the planet-induced axisymmetric mass loss model. Since about half of all planetary nebulae are elliptical, i.e., have low equatorial to polar density contrast, it was predicted that about $`50\%`$ of all solar-like stars have Jupiter-like planets around them, i.e., a mass about equal to that of Jupiter, $`M_J`$, or more massive. In light of the new findings that only $`5\%`$ of sun-like stars have such planets, and the mechanism of dust formation near cool magnetic spots, I revise this prediction. In Soker (2001b) I predict that indeed $`50\%`$ of PNs progenitors do have close planets around them, but the planets can have much lower masses, as low as $`0.01M_J`$, in order to efficiently spin-up the envelopes of AGB stars. To support this claim, I follow the angular momentum evolution of single stars with main-sequence mass in the range of $`1.32.4M_{}`$, as they evolve to the post-AGB phase. I find that single stars rotate much too slowly to possess any significant non-spherical mass loss as they reach the upper AGB. It seems, therefore, that planets, in some cases even Earth-like planets, are sufficient to spin-up the envelope of these AGB stars for them to form elliptical PNs. The prediction that on average several such planets orbit each star, as in the solar system, still holds. In Soker (2001a) I show that wide stellar companions to AGB stars may also accrete mass, form an accretion disk, and blow jets, hence forming elliptical PNs. This reduces the fraction of PN progenitors which are needed to have planetary systems from $`50\%`$ to $`35\%`$, or even less.

Another finding from exoplanets is that metal rich stars are more likely to harbor planetary systems. This implies, in the context of planet-shaping of PNs, that spherical PNs will tend to originate from low metallicity stars. Indeed, when carefully defining spherical PNs, this is the case. In Soker (2002a) I examine the mass loss history and distribution of spherical PNs in the galaxy. I argue there that spherical PNs form a special group among all PNs. The smooth surface brightness of most spherical PNs suggests that their progenitors did not go through a final intensive wind (FIW, also termed superwind) phase. While $`70\%`$ of the PNs of all other PNs groups are closer to the galactic center than the sun is, only $`30\%`$ of spherical PNs are; $`70\%`$ of them are farther away from the galactic center. These, plus the well known high scale height above the galactic plane of spherical PNs, suggest that the progenitors of spherical PNs are low mass stars having low metallicity.

I also examine the possibility of detecting signatures of surviving Uranus, Neptune-like planets inside PNs (Soker 1999). Giant planets that are not too close to the stars (orbital separation larger than $`5`$ AU) are likely to survive the entire evolution of the star. As the star turns into a PN, it has a fast wind and strong ionizing radiation. The interaction of the radiation and wind with a planet may lead to the formation of a compact condensation or tail inside the PN, which emits strongly in H$`\alpha `$, but not in \[OIII\]. The position of the condensation (or tail) will change over a time-scale of $`10`$ yr. Such condensations might be detected with currently existing telescopes. This idea was then repeated for planets around white dwarfs (Chu et al. 2001).

Finally, I note that most of the known stars with extrasolar planets will not form PNs at all. Instead, because they have relatively low mass (most have $`M<1.3M_{}`$), I expect them to lose most of their envelope on the RGB, becoming blue-HB stars, and then fading as WD without an observable nebula.

This research was partly supported by the Israel Science Foundation.