text:
# Bernoulli potential at a superconductor surface

\[

## Abstract

The electrostatic Bernoulli potential measured at the surface of a superconductor via Kelvin capacitive coupling is shown to be independent of the pairing mechanism. This contrasts with the Bernoulli potential in the bulk where contributions due to pairing dominate close to $`T_c`$.

\] In 1968 Bok and Klein have first observed non-zero Hall voltage in the Meissner state of a type-I superconductor. Instead of Ohmic contacts, they have used the Kelvin capacitive coupling between the surface and the electrodes. Similar measurements have been performed by Brown and Morris and more recently by Chiang and Shevchenko.

Among these experiments, the highest accuracy has been achieved by Morris and Brown. They found a perfect agreement with the prediction of van Vijfeijken and Staas in which the original Bernoulli potential (derived by Bopp and made familiar by London) is extended by the fountain term of the two-fluid model. In particular, this measurement shows no traces of effects predicted later within more sophisticated thermodynamical approaches.

In this paper we show that the agreement of the two-fluid theory with experiment is in some sense fortuitous. This theory has been developed for the electrostatic potential inside the superconductor but it is the surface potential that is measured by the Kelvin method. The surface potential is influenced by the surface dipole which brings a step-like contribution to the electrostatic potential. With the surface dipole included, the measured surface potential agrees with the prediction of the thermodynamic approaches.

The Hall voltage in superconductors is a difference between electrostatic potentials, $`\phi `$, at the contacts. All theories agree that this potential has the form of the Bernoulli potential,

$$e\phi \frac{1}{2}mv^2,$$
(1)
where $`v`$ is the velocity of the condensate and $`m`$ is the mass of electrons. The point of disagreement concerns the amplitude of this potential.

London assumed that the electric force balances the inertial and Lorentz forces. The energy-conserving integral of motion then yields the Bernoulli potential,

$$e\phi =\frac{1}{2}mv^2.$$
(2)
This formula applies only at zero temperature. To cover finite temperatures, van Vijfeijken and Staas took into account the fact that the normal electrons feel the electric field. Since they stay at rest, there must be a balancing force from the interaction with the condensate. The reaction force acts on the condensate and adds to the electric force. This picture leads to the electrostatic potential which is proportional to the fraction of superconducting electrons,

$$e\phi =\frac{n_s}{n}\frac{1}{2}mv^2.$$
(3)
The formula (3) represents the level of theory before the first successful experiment appeared. With respect to experiments it is profitable to express the Bernoulli potential (3) as a function of the magnetic field $`B`$. At the surface, $`v=\frac{e}{m}\lambda B`$, where $`\lambda =\sqrt{\frac{m}{e^2n_s\mu }}`$ is the London penetration depth, therefore the two-fluid potential (3) reads

$$e\phi =\frac{1}{n}\frac{B^2}{2\mu }.$$
(4)
Morris and Brown experimentally confirmed the validity of (4) in a wide temperature range including the vicinity of the critical temperature.

As discussed by Bok and Klein, there is also a simple theoretical argument supporting formula (4). In the absence of dissipation, the only force acting on the lattice is due to the electric field, $`\phi `$. The total force per unit area acting on a superconducting slab,

$$F=\underset{\mathrm{right}}{\overset{\mathrm{left}}{}}𝑑xen\phi =en\phi _{\mathrm{left}}en\phi _{\mathrm{right}},$$
(5)
is proportional to the difference of the surface potentials. Here, $`en`$ stands for the charge density of the ionic lattice. In the magnetic field $`B`$, the slab with the current density $`J`$ \[in A/m\] feels the Lorentz force $`F=BJ`$. Using the Ampere rule, $`B_{\mathrm{right}}B_{\mathrm{left}}=\mu J`$, and the mean magnetic field, $`B=\frac{1}{2}(B_{\mathrm{right}}+B_{\mathrm{left}})`$, the Lorentz force can be expressed as the difference between the left and right magnetic pressures,

$$F=BJ=\frac{B_{\mathrm{right}}^2}{2\mu }\frac{B_{\mathrm{left}}^2}{2\mu }.$$
(6)
One can see that the Lorentz force (6) agrees with (5) for potential (4).

In spite of its internal consistency and the satisfactory agreement with the experimental data of Bok and Klein, the two-fluid theory came under suspicion when studies based on thermodynamic relations and BCS analyses have been accomplished. These approaches result in a potential

$$e\phi =\frac{n_s}{n}\frac{1}{2}mv^2.$$
(7)
For the phenomenological density of condensate,

$$n_s=n\left(1\frac{T^4}{T_c^4}\right),$$
(8)
potential (7) reads

$$e\phi =\frac{n_s}{n}\frac{1}{2}mv^2+4\frac{n_n}{n}\frac{\mathrm{ln}T_c}{\mathrm{ln}n}\frac{1}{2}mv^2,$$
(9)
where $`n_n=nn_s`$ is the density of normal electrons. The first term is identical to the two-fluid potential (3), the second is the difference which we will call the thermodynamic correction.

Formulas (7) and (3), are very different close to $`T_c`$. Within the BCS relations and the parabolic band approximation, one can rearrange (9) as

$$e\phi =\frac{n_s}{n}\frac{1}{2}mv^2\frac{4}{3}\frac{n_n}{n}\mathrm{ln}\frac{T_D}{T_c}\frac{1}{2}mv^2.$$
(10)
Since the Debye temperature $`T_D`$ is much larger than the critical temperature $`T_c`$, the second term dominates at higher temperatures, approximately for $`T>\frac{2}{3}T_c`$.

The thermodynamic correction arises from the pairing mechanism. According to the BCS critical temperature, $`T_cT_D\mathrm{e}^{1/𝒟V}`$, where $`𝒟`$ is the single-spin density of states and $`V`$ is the BCS interaction parameter, the potential (9) can be written as

$$e\phi =\frac{n_s}{n}\frac{1}{2}mv^2\frac{n_n}{n}\frac{1}{𝒟V}\frac{2}{3}mv^2.$$
(11)
As speculated by Rickayzen, the Hall voltage at temperatures close to $`T_c`$ can offer an experimental access to the BCS interaction $`V`$. Rickayzen expressed his surprise that Bok and Klein did not observe any trace of the thermodynamic correction at the temperature $`T=4.2`$ K for lead with $`T_c=7.2`$ K. For $`T_D=105`$ K the factors in (10) are comparable, $`\frac{n_s}{n}=0.88`$ and $`\frac{4}{3}\frac{n_n}{n}\mathrm{ln}\frac{T_D}{T_c}=0.42`$.

While the measurement of Bok and Klein was not sufficiently precise and the thermodynamic correction could be lost in the experimental error, later measurement of Morris and Brown clearly proved that the Hall voltage is independent of temperature. They also used lead and explored the temperature range from $`4.2`$ to $`7.0`$ K. At $`T=7`$ K the thermodynamic correction dominates the potential, since the factors in (10) are $`\frac{n_s}{n}=0.1`$ and $`\frac{4}{3}\frac{n_n}{n}\mathrm{ln}\frac{T_D}{T_c}=3.2`$. If present, the thermodynamic correction would amplify the signal thirty times. This is far beyond the experimental error, estimated to be less than 10%. Briefly, the thermodynamic correction turned out to be absent. It should be mentioned that Morris and Brown expected a non-zero thermodynamic correction indicated by their less precise preliminary measurement at $`T=4.2`$ K.

One might feel that the experimental results imply that the thermodynamic correction is an ill-defined concept. It is, however, not the case. The electrostatic potential given by (7) applies in the bulk. To obtain the potential at the surface, one cannot simply extrapolate this bulk value. At the surface there is a surface dipole which leads to a discontinuity in the potential.

The surface dipole in a metal introduces a step of several eV’s in the electrostatic potential. In the jellium model this step equals to the Fermi energy $`E_F`$. With the exchange and correlation energy included, the surface dipole is reduced by half which shows the importance of correlations. The phase transition into the superconducting state increases the correlation energy by the condensation energy. The surface dipole of the superconductor is thus slightly smaller than the dipole of the same crystal in the normal state. This difference is of the order of $`\mathrm{\Delta }^2/eE_F`$.

When a supercurrent flows close to the surface, the condensation energy is reduced. The surface dipole thus slightly increases, reaching the normal state value at the critical current. This current dependent part of the surface dipole contributes to the potential observed by the capacitive coupling. Now we show that the surface dipole cancels the thermodynamic contribution (second term of (9)) to the Hall voltage.

In our formulation we follow Rickayzen. The increase of the free energy due to the current equals the kinetic energy of the condensate,

$$\delta =n_s\frac{1}{2}mv^2.$$
(12)
The electrochemical potential, $`\nu =E_F+\delta \nu _v+e\phi `$, is constant in the whole system, therefore $`e\phi =\delta \nu _v`$. Since $`\nu =\frac{}{n}`$, the velocity variation of the local chemical potential is $`\delta \nu _v=\frac{}{n}\delta `$. Accordingly, the electrostatic potential induced by the current is given by

$$e\phi =\frac{}{n}\delta .$$
(13)
With $`\delta `$ given by (12), from (13) one recovers (7).

The step of the electrostatic potential at the surface is given by the Budd-Vannimenus theorem,

$$e\phi _{\mathrm{surf}}e\phi =n\frac{}{n}\frac{\delta }{n}.$$
(14)
We note that Budd-Vannimenus theorem has been derived for the surface potential of the normal metal. Formula (14) is a modification restricted to the potential caused by the supercurrent.

From (14) and the bulk potential (13) one finds the electrostatic potential at the surface,

$$e\phi _{\mathrm{surf}}=\frac{n_s}{n}\frac{1}{2}mv^2.$$
(15)
This potential is identical to the experimentally confirmed two-fluid potential (3). Accordingly the surface dipole cancels exactly the thermodynamic corrections.

In conclusion, the electrostatic potential observed at the surface of superconductors is determined by the magnetic pressure. The only material parameter accessible by its measurement is the charge density of the ionic lattice. In the bulk, the electrostatic potential depends also on the pairing mechanism.

Since the thermodynamic corrections cannot be measured on the surface, it is desirable to observe the internal electric field directly in the bulk of a superconductor. A new experiment in this direction has been reported by Kumagai et al who have measured the effect of charge density deviations on the nuclear quadrupole resonance.

It is not obvious how the surface potential works in layered structures of high-$`T_c`$ materials. A measurement of the Hall voltage on these materials, however, involves enormous technical difficulties.

The authors thank Nai Kwong for helpful discussions. This work was supported by GAČR 202000643 and 202990410, GAAV A1010806 and A1010919 grants. The European ESF program VORTEX is also acknowledged.


format: text

text:
# Monte Carlo Simulation of Liquid-Crystal Alignment and of Chiral Symmetry-Breaking

## I Introduction

Molecular simulations have emerged as an important tool for investigating liquid crystals. The goal of the simulations is to relate the microscopic structure of molecules to the macroscopic properties of liquid-crystal phases. To do this, researchers begin with some assumptions about molecular structure, simulate a large number of interacting molecules, and determine the large-scale order of the resulting system as a function of temperature, applied fields, and other variables. This technique provides insight into the origin of liquid-crystal phase transitions and textures. It also assists in the development of new liquid-crystal systems for technological applications.

Researchers have used several approaches for molecular simulations of liquid crystals. On the most microscopic scale, some studies begin with the structure of liquid-crystal molecules in atomistic detail, and let the molecules interact through a potential derived from the sum of interatomic interactions.  The advantage of this approach is that it can simulate actual experimental materials and can show how macroscopic properties change with slight variations in the chemical structure. The disadvantage is that this approach is computationally intensive, and hence it can only be done for relatively small systems, which may not be large enough to show all the important features of liquid-crystal ordering. Accuracy of classical potentials can also be an issue.

Because of these limitations, many other studies begin with greatly simplified models of molecular structures and interactions, which require less computer time and hence allow larger system sizes. Some studies represent molecules as cylinders or spherocylinders, which interact through a hard or soft excluded-volume repulsion.  Other studies represent molecules as ellipsoids interacting through the Gay-Berne potential, an anisotropic generalization of the Lennard-Jones potential.  A related approach uses the Corner potential.  To model more complicated molecular structures, further studies represent molecules as composites of two Gay-Berne ellipsoids.  Soft-sphere potentials can also be joined together in a linear chain to represent molecules. Such molecules can be rigid  or semiflexible.

As an alternative approach, we model liquid-crystal molecules as a composite of soft spheres bonded rigidly into a fixed shape, with no intramolecular flexibility. This approach is particularly useful for two reasons. First, the intermolecular interaction can be calculated at low computational cost, thus allowing the simulation of large systems. Second, the spheres may be bonded into any shape to represent molecular structures of any symmetry. These features allow us to explore the complex relationship between molecular shape and phase behavior, without undertaking full-scale atomistic simulations. The approach thus represents a compromise between the simplest and most complex models used in liquid crystal simulation.

In a recent paper,  we used this approach to simulate smectic liquid crystals. Our goal was to model the electroclinic effect, in which an applied electric field induces molecular tilt in the smectic-A phase. For that reason, we considered model molecules with the structure shown in Fig. 1(a). The overall shape is a biaxial zig-zag shape like the letter Z, and a transverse electric dipole moment makes the molecule chiral. This shape was motivated by the fact that many real liquid-crystal molecules exhibiting the electroclinic effect have this general structure. Indeed, the biaxiality of this structure should be important for the electroclinic effect because a tilted smectic phase is biaxial. Our simulation results confirmed that the degree of biaxiality or obliqueness in molecular shape is crucial for spontaneous molecular tilt in the smectic-C phase and for induced molecular tilt in the smectic-A phase.

The purpose of our current study is to explore this simulation method further by applying it to two different problems in the physics of liquid crystals. First, we consider the alignment of a nematic liquid crystal by a ridged surface. Second, we explore chiral symmetry-breaking in the smectic phases of a bent-core (or banana-shaped) liquid crystal. In each case, we construct model mesogens consisting of soft spheres bonded rigidly together with an appropriate shape. Our results for these two problems shows that this simulation method is a versatile approach for modeling different types of liquid-crystal ordering.

The first problem that we consider, the surface alignment of a nematic phase, is a classic problem in liquid-crystal science and is fundamental to all applications of liquid crystals.  Liquid-crystal cells are constructed so that the surfaces favor alignment of the molecules in a particular direction, which may be planar (in a particular direction within the surface plane), homeotropic (normal to the surface), or between these extremes. An applied electric field can cause the molecules in the interior of the cell to tilt away from the alignment favored by the surfaces. Thus, the competition between surface alignment and an applied field allows the molecules in the cell to be switched between different orientations, changing the optical properties of the cell. Key questions in designing new surfaces to optimize liquid-crystal alignment are: What type of surface is needed to align liquid-crystal molecules? How far into a cell does the alignment persist?

To address these questions, we perform Monte Carlo simulations of nematic alignment, using the simulation approach discussed above. For the model molecules, we would like to construct a rigid uniaxial rod-like structure, which is known to give a nematic phase at high enough densities. For that reason, we use seven spheres bonded together in a straight line, as shown in Fig. 1(b). For the aligning surface, we would like to model a rubbed polymer substrate, which has parallel ridges that give planar alignment of the liquid crystal. We use a flat surface decorated with regularly spaced elevated ridges, each composed of spheres fixed in a row, as shown in Fig. 2. Our simulation results show that the nematic director aligns with the ridges for all ridge spacings studied. The nematic order parameter near the substrate increases as the ridge spacing is reduced, indicating that narrowly spaced ridges on the substrate suppress director fluctuations. However, the correlation length associated with the drop-off of nematic order away from the substrate is independent of ridge spacing, and is a characteristic of the bulk nematic phase.

Our second problem, chiral symmetry-breaking in systems of bent-core molecules, concerns a new type of smectic liquid crystal. In recent years, researchers have synthesized liquid-crystal molecules with a bent or “banana-shaped” core.  The individual molecules are not chiral. In the bulk or in freely suspended films, the molecules form smectic phases. In the smectic-A phase, the molecules are on average normal to the smectic layers, and hence the phase is not chiral. However, when the material is cooled into a smectic-C phase, the molecules can tilt with respect to the layers in two inequivalent ways that are mirror images of each other. For that reason, the system undergoes a chiral symmetry-breaking transition, in which it forms right- and left-handed chiral domains.  There are actually several different smectic-C phases, depending on how the domains of alternate chirality and alternate tilt pack together. Key questions for the theory of these systems are: What type of intermolecular interaction drives the chiral symmetry-breaking transition? How do the molecules pack together in the achiral and the chiral phases?

To investigate these issues, we perform Monte Carlo simulations of bent-core liquid crystals, again using the same simulation approach. In this case, we construct model molecules by arranging seven soft spheres in a rigid bent-core shape, as shown in Fig. 1(c). In the simulations, we find that a bulk system of these molecules forms a nonchiral smectic-A phase at high temperatures, and a chiral tilted crystalline phase at low temperatures. The chiral symmetry-breaking transition is driven by a favored packing of the molecules, which can fit together more efficiently in either a right- or a left-handed state than in a symmetric state. Thus, at least in this model, the steric interaction of rigid molecules is sufficient for chiral symmetry-breaking; neither molecular flexibility nor dipole-dipole interactions are necessary. Although this system forms a chiral tilted phase with three-dimensional crystalline order, it does not form a chiral smectic-C phase with only one-dimensional positional order. We believe that the simultaneous onset of chiral order and crystalline order is an artifact of our particular choice of molecular shape, because the same type of chiral order can occur in either a crystalline or a smectic-C phase.

## II Surface Alignment of Nematic Phase

To simulate a nematic phase, we use the rod-like model molecule shown in Fig. 1(b). Each molecule is composed of seven spherical particles, or interaction sites, fixed rigidly in a straight line, with no intramolecular degrees of freedom. The sphere-sphere interaction potential is the truncated Lennard-Jones potential cut off at its minimum so there is no attractive tail, also known as the Weeks-Chandler-Andersen potential,

$$U_{mn}^{\mathrm{int}}=\{\begin{array}{cc}4ϵ\left[\left(\frac{\sigma }{r_{mn}}\right)^{12}\left(\frac{\sigma }{r_{mn}}\right)^6\right]+ϵ,\hfill & \text{if }r_{mn}r_c=2^{1/6}\sigma \text{;}\hfill \\ 0,\hfill & \text{otherwise.}\hfill \end{array}$$
(1)
where $`r_{mn}=|𝐫_m𝐫_n|`$ and $`m`$ and $`n`$ are the sphere indices referring to different molecules. In this way, the steric interactions between molecules are represented via a short-range repulsive potential, with no attractive or electrostatic interactions.

The simulation cell consists of 5,070 molecules in a three-dimensional rectangular box of size $`35\times 35\times 41`$, where lengths are measured in units of the Lennard-Jones parameter $`\sigma `$, and each molecule is seven units long and one unit wide. The system has periodic boundary conditions along $`x`$ and $`y`$, and is confined in the $`z`$ direction between parallel substrates.

Each substrate consists of a soft repulsive wall that interacts with the spherical particles in each molecule using the same WCA potential. Each substrate is decorated with regularly spaced elevated ridges composed of spheres fixed in parallel rows along the wall, as shown in Fig. 2. Each of these spheres is identical to the spheres making up the liquid-crystal molecules, and again has the same WCA interactions. Thus, interactions between the liquid crystal and the wall are essentially steric in nature, with no attractive or electrostatic terms.

Each wall has an integer number of ridges, spaced evenly. We make the ridge pattern the same for the top and the bottom substrates, so that the nematic is not twisted, having the same preferred planar orientation on both substrates. However we slightly shift the ridge patterns so that the substrates do not have an exact mirror symmetry.

We perform Monte Carlo (MC) simulations at the fixed temperature $`k_BT=5.0`$, well above the nematic–smectic-A transition point, which is about 3.5, measured in units of the Lennard-Jones parameter $`ϵ`$. The density is held fixed at $`\rho =0.7`$ spheres/unit volume. In preliminary simulations, we either randomize the initial molecular orientations or initialize them in a homeotropic orientation, aligned normal to the substrates. We find that either of these initial conditions leads to a slow ordering process, which eventually gives rise to an aligned planar structure, with the nematic director parallel to the ridges on the substrates. After determining the preferred orientation, we perform further simulations with the molecules initialized in the preferred planar orientation, with the molecules aligned parallel to the ridges and the center-of-mass positions randomized.

Each MC step consists of three attempted translations (along $`x`$, $`y`$, and $`z`$) and two attempted rotations (in $`\theta `$ and $`\varphi `$) per molecule, with appropriate sampling in the polar angle.  We carry out 100,000 MC steps to reach equilibrium and then an additional 200,000 MC steps to collect data.

To measure the spatial variation of liquid-crystal order near the substrates, we calculate the nematic order tensor $`Q_{\alpha \beta }`$, which represents the strength and direction of orientational order,

$$Q_{\alpha \beta }=\frac{1}{N}\underset{j=1}{\overset{N}{}}\left(\frac{3}{2}n_{j\alpha }n_{j\beta }\frac{1}{2}\delta _{\alpha \beta }\right),$$
(2)
where $`𝐧_j`$ is the director along the axis of molecule $`j`$ and $`N=5070`$ is the number of molecules. For any $`z`$, the largest eigenvalue of $`Q_{\alpha \beta }(z)`$ represents the scalar nematic order parameter $`S(z)`$. The corresponding eigenvector identifies the local direction of alignment.

We evaluate $`Q_{\alpha \beta }`$ as a function of distance from the substrate by dividing the system into fifty layers and calculating $`Q_{\alpha \beta }`$ in each layer. Because the number of molecules in each layer is relatively small, around 100, we must average over time to get a statistically meaningful sample. Accordingly, we average over a series of 100 independent configurations separated in time by 2000 MC steps each.

We perform simulations using substrates with three different ridge spacings, but with the same bulk density $`\rho =0.7`$ spheres/unit volume, the same temperature $`k_BT=5.0`$, and the same number of molecules $`N=5070`$ (a total of 35,490 particles, plus those decorating the substrate). The number of ridges on the substrates is chosen as 1, 2, and 7, with spacing between ridges of 35, 17.5, and 5 length units, respectively (taking periodic boundary conditions into account). The box dimensions are set at $`H_x=35`$, $`H_y=35`$, and $`H_z41`$, with $`H_z`$ adjusted slightly to make the bulk densities the same for all three systems, as the total number of particles in the box is slightly increased by the addition of extra ridges.

In each case we find that the nematic director throughout the cell aligns parallel to the ridges. Figure 3 shows the nematic order parameter $`S(z)`$ vs. $`z`$ for all three systems, with one substrate at $`z=0`$ and the other at $`z=H_z`$. The bulk order parameter is about the same in all cases, but the order near the substrate varies with ridge spacing. The narrowest ridge spacing gives the largest enhancement in local order, while the widest ridge spacing gives the least enhancement. Figure 3 also shows that the correlation length associated with nematic order near the interface is independent of the ridge spacing. In other words, no matter how densely spaced the substrate ridges, the thickness of the ordered interface does not change. We fit the data in Fig. 3 to the form

$$S(z)=S_0+A(e^{z/\xi }+e^{(H_zz)/\xi }),$$
(3)
where $`S_0`$ is the bulk order parameter value, $`A`$ is the amplitude of enhancement near the substrates, and $`\xi `$ is the associated correlation length. Using the least-squares algorithm, we fit the data for all three ridge spacings. Figure 3 shows an example of the fit for the ridge spacing of 5 units, and Table I gives the fit parameters.

Our results can be understood through an analytical calculation using a Landau free energy functional, which is similar to the model in Ref. . We assume that the nematic order parameter takes the form,

$$S(z)=S_0+S_1(z)$$
(4)
where $`S_0`$ is the bulk value and $`S_1(z)`$ is the surface-induced contribution. Near a single aligning substrate, the free energy can be written as

$$F=𝑑z\left[\frac{1}{2}K\left(\frac{dS}{dz}\right)^2+\frac{1}{2}a(S(z)S_0)^2\right]hS(z=0).$$
(5)
where $`a`$ and $`K`$ are bulk parameters dependent on temperature and density. The term $`hS(z=0)`$ favors enhanced nematic order at the substrate, and one expects the anchoring strength $`h`$ to vary with substrate structure.

Minimizing the free energy functional gives

$$S(z)=S_0+Ae^{z/\xi },$$
(6)
where $`A=h/\sqrt{Ka}`$ and $`\xi =\sqrt{K/a}`$. These equations show that the amplitude $`A`$ is a surface parameter, which depends on $`h`$ as well as $`K`$ and $`a`$ and thus varies with substrate structure. By contrast, the correlation length $`\xi `$ is a bulk parameter, which depends only on $`K`$ and $`a`$ and thus varies with temperature and density but is independent of substrate structure. This is consistent with the simulation results for $`A`$ and $`\xi `$ in Table I.

Apart from the enhanced nematic order near the substrates, the simulations also show a local tendency to smectic ordering near the substrates, even though the temperature is well above the bulk nematic-smectic transition. Figure 4 shows a snapshot of the molecular configuration near one substrate, for the ridge spacing of 5 units. Clearly the molecular density is not uniform. Instead, there are local regions of incipient molecular layering in the direction along the ridges. This local layering suggests that the system is near a surface nematic-smectic transition. The tendency to smectic ordering in a nematic close to a surface with planar orientation has also been studied theoretically in the case of a bare substrate.  In that system, as in our system, the wall suppresses local director fluctuations, giving rise to a more ordered state, just as an applied electric field can drive a nematic-smectic transition.

In experiments on nematic alignment by buffed surfaces, it has sometimes been observed that a nematic aligns in a direction different from, or even perpendicular to, the buffing axis.  While our simulations did not show that kind of ordering, our results lead us to consider possible explanations. One is that the local smectic order next to the substrate could form with the layers parallel to the ridges, and with the molecules oriented perpendicular to the ridges or at some degree of tilt. It is possible that a wider ridge spacing, a different molecular shape, or a temperature closer to the nematic-smectic transition could produce such effects in a simulation like ours. Another likely explanation is that electrostatic or attractive interactions play an important role in determining the preferred orientation. These effects have not been considered in the present model but could be introduced in future work.

## III Chiral Symmetry-Breaking in Systems of Bent-Core Molecules

As discussed in the Introduction, recent experimental studies have found that liquid crystals composed of certain achiral molecules can undergo a spontaneous symmetry-breaking and form chiral domains.  These “bent-core” or “banana-shaped” molecules form layered structures similar to those observed in ferroelectric smectic liquid crystals, with the molecular axis tilted with respect to the layer normal. Each molecule has a net dipole moment which defines a polar axis perpendicular to the long axis of the molecule. The chirality of the tilted phase is determined by the spatial relationship of the layer normal, the tilt direction, and the polar axis; these can be arranged in either right- or left-handed fashion. The two arrangements are degenerate in energy, leading to the formation of right- and left-handed domains with equal probability.

We perform Monte Carlo simulations to investigate what type of intermolecular interaction can drive the chiral symmetry-breaking transition, and how the molecules pack together in the achiral and chiral phases. In these simulations, we use the same approach discussed above for simulations of nematic alignment, but we change the model molecular structure in order to represent bent-core molecules. The model structure consists of seven spheres arranged in a rigid bent-core shape, like the letter V, as shown in Fig. 1(c). The bend angle between the two straight segments of the molecule is $`140^{}`$. As in the earlier simulations, the molecules interact through the soft repulsive sphere-sphere interaction of Eq. (1). We neglect dipole-dipole interactions in order to see what ordering is induced by purely steric interactions.

To characterize the magnitude and sign of chirality of the system, we introduce a chiral order parameter $`\chi _i`$ for each molecule $`i`$. In a smectic phase, let the unit vector $`\widehat{𝐦}_i`$ represent the local layer normal, $`\widehat{𝐧}_i`$ represent the long axis of the molecule, and $`\widehat{𝐩}_i`$ represent the molecular dipole moment. The chiral order parameter can then be defined as

$$\chi _i=[(\widehat{𝐦}_i\times \widehat{𝐧}_i)\widehat{𝐩}_i][\widehat{𝐦}_i\widehat{𝐧}_i].$$
(7)
Note that this order parameter is invariant under the two symmetry operations $`\widehat{𝐦}_i\widehat{𝐦}_i`$ and $`\widehat{𝐧}_i\widehat{𝐧}_i`$, but not under $`\widehat{𝐩}_i\widehat{𝐩}_i`$, which is not a symmetry of the system. The average chiral order parameter is then

$$\overline{\chi }=\frac{1}{N}\underset{i=1}{\overset{N}{}}\chi _i.$$
(8)
That average order parameter can be compared with the rms chiral fluctuations

$$\chi _{\mathrm{rms}}=\left[\frac{1}{N}\underset{i=1}{\overset{N}{}}(\chi _i\overline{\chi })^2\right]^{1/2}.$$
(9)
We simulate 600 molecules in a rectangular box with periodic boundary conditions. We choose the height of the box so that the molecules can form four smectic-like layers. The density is held fixed at $`\rho =0.7`$ spheres per unit volume. At the temperature $`k_BT=5.0`$, the systems forms a chiral crystalline phase. A snapshot of this phase is shown in Fig. 5. In equilibrium, the average chiral order parameter is $`\overline{\chi }=0.24`$, while the rms average chiral fluctuations are $`\chi _{\mathrm{rms}}=0.06`$. The average chiral order parameter is much greater than the chiral fluctuations, which means the system is in a well-defined chiral state. This spontaneous chiral ordering can be seen in the snapshot, which shows the molecules tilting in a way that breaks reflection symmetry.

The system shown in Fig. 5 has a large spontaneous tilt in each layer, but the direction of the tilt appears to vary randomly from layer to layer. The coupling between adjacent layers due to excluded volume interactions is apparently too weak to order the tilt directions. A similar variation in the tilt from layer to layer was found in the smectic-C phase of similar molecules.  Although the tilt directions vary, the chiral order parameter has the same sign in each layer, indicating that the simulation cell represents a domain of uniform chirality.

The structure of a single layer, shown in Fig. 5(b), shows that the system has periodic order with some disordered regions, so it is a crystalline rather than a smectic phase. When we raise the temperature or lower the density of this system, the tilt vanishes and we find a stable untilted smectic-A phase. We do not observe any tilted chiral smectic-C phase. In this system, the onset of chiral order apparently occurs together with the onset of crystalline order.

The observation of a chiral phase, even though it is crystalline rather than smectic, answers the basic questions raised above. First, the simulations show that the steric interaction of rigid molecules is sufficient to drive the chiral symmetry-breaking transition. Because we neglect dipole-dipole interactions in the simulations, we see that they are not necessary for this transition. Furthermore, because we completely suppress intramolecular degrees of freedom, we see that chiral ordering can occur in achiral materials without the involvement of chiral excited molecular states. Second, the snapshot shows how the molecules pack together in the chiral phase. We see that the molecules tilt so that the spheres of one molecule fit into the spaces between the spheres of adjacent molecules. Thus, a tilted chiral state with either handedness is a more efficient packing of neighboring molecules than an untilted achiral state. This shows an explicit mechanism for chiral symmetry-breaking in a system with purely steric interactions.

In fact, the “bumpy” nature of molecular shape may play a crucial role in favoring the tilted geometry, via the details of molecular packing. In our model, two parallel molecules do not pack together efficiently unless they are tilted enough to allow the close packing of the component spheres. A bent-core molecule composed of two spherocylinders has also been simulated,  but it does not show a chiral symmetry-breaking transition, presumably because close packing of spherocylinders does not require molecular tilt.

We do not think that there is any fundamental connection between the onset of chiral order and the onset of crystalline order. The same type of chiral order described by the order parameter of Eq. (7) can occur in either a crystalline or a smectic-C phase. The fact that we do not observe a chiral smectic-C phase is probably an artifact of the particular molecular structure that we have simulated. To test how sensitively the phase behavior depends on molecular structure, future simulations should make slight variations on the molecular structure. In general, a less sharp bend-angle will favor a smectic phase over a crystalline phase at a given temperature and density. A bump strategically placed at the center of the bent-core molecule will change the packing geometry and might favor twist (of either handedness) over parallel packing. Such a bump can be formed by increasing the interaction distance of the central sphere, or by adding an extra sphere to the molecule. The inclusion of dipole-dipole interactions could also destabilize the crystal and favor the formation of a smectic phase.

In summary, we have extended a Monte Carlo simulation method, originally developed for simulations of the electroclinic effect, to model surface alignment of nematic liquid crystals and chiral symmetry-breaking in bent-core liquid crystals. The results give insight into the local molecular arrangements that lead to both of these macroscopic effects. These two studies show that this simulation method is a versatile technique for investigating a wide range of ordering phenomena in liquid crystals.

###### Acknowledgements.

This research was supported by the Office of Naval Research and the Naval Research Laboratory.


text:
# CMC–Surfaces, 𝜑–Geodesics and Caratheodory’sche Vermutung

## Introduction

The constant mean curvature (CMC–) surfaces in $`E^3`$ are known to admit a continuous family of local, non–trivial, isometric deformations preserving mean curvature of the surface ($`H`$–deformations). In the case when surface is compact Umehara  showed that the converse is also true.

The maximal and minimal curvature lines of CMC–surface form an orthogonal net which is called réseau de Bonnet, cf Cartan . The Bonnet Theorem says that if the CMC–surface is simply connected and umbilic–free, then under $`H`$–deformations the orthogonal net “rotates” through a constant angle which can be taken as a parameter of deformation.

If the CMC–surface is not simply connected or umbilic–free, Cartan seems to be the first to ask about possible scenario of evolution of réseau de Bonnet under the $`H`$–deformations.

In the present note we study 111For the reasons which will be clear later. evolution of the orthogonal nets in the case when CMC–surface is simply connected with a single umbilic or, equivalently, doubly connected and umbilic–free. Namely, if we “pinch” the umbilic, the CMC–surface becomes an annulus whose points undergo $`H`$–deformations accordingly with the Bonnet Theorem. In general, the rotation angle is no longer constant at all the points because annulus cannot be covered by a single chart.

However, the Bonnet Theorem implies that every curve of the orthogonal net is a $`\phi `$–geodesic line whatever $`H`$–deformations are applied to a CMC–surface. Metric $`\phi `$ is given by the linear element $`ds=|\phi ||dz|`$, where $`\phi dz^2`$ is a holomorphic quadratic ”differential” associated to the CMC–surface. Of course, $`\phi (0)=0`$ at the umbilical point.

This observation is crucial, because the $`\phi `$–geodesics near $`n`$–th order zero of a holomorphic quadratic form are well–understood due to the works of Strebel . Roughly speaking, the $`\phi `$–geodesics fill–up the annulus either by “hyperbolas” or “radii”. Therefore, possible configurations of réseau de Bonnet near the umbilic looks like a singularity with the finite number of hyperbolic and parabolic sectors.

Despite independent interest, the orthogonal nets are auxiliary for us. We postulate different fact here: $`H`$–deformations of orthogonal nets give an amazingly simple proof to the Caratheodory’sche Vermutung (Conjecture):

###### Theorem 1

Let $`S^2`$ be a $`C^{\mathrm{}}`$ surface which bounds a convex compact body in the Euclidean space $`E^3`$. Then $`S^2`$ has at least two umbilical points. In other words, the Euler-Poincaré index of isolated umbilical point is at most $`+1`$.

(A short overview of this conjecture can be found in ; see also , , .)

Acknowledgments. I am grateful to the referee for a careful reading of this note.

## 1 $`\phi `$–Geodesics

Until further indications, $`M`$ is a simple domain of the complex parameter $`z`$. Let us consider the holomorphic functions $`\phi (z)`$ vanishing at the unique point of $`M`$ which we identify with $`0`$. An order $`n1`$ is assigned to $`0`$, if there exists a complex constant $`a0`$ such that $`\phi (z)=az^n+O(|z|^{n+1})`$.

Flat metric $`\phi `$ with the cone singularity of angle $`(n+2)\pi `$ is given by the formula

$$|ds|=|\phi ||dz|,$$

provided $`\phi (z)dz^2`$ is a quadratic form on $`M`$. By a $`\phi `$–geodesic line in $`M`$ one understands the line conisting of the shortest arcs relatively metric $`\phi `$. Any two points in $`M`$ (including $`0`$) may be joined by the unique $`\phi `$–geodesic line. Strebel classified the possible types of $`\phi `$–geodesics in the neighborhood of $`n`$–th order zero by proving the following lemma.

###### Lemma 1

() Any two points in a neighborhood $`M`$ of $`n`$–th order zero of holomorphic 2–form $`\phi (z)dz^2`$ can be joined by a unique $`\phi `$–geodesic. Moreover, each $`\phi `$–geodesic is either an arc defined by the equation $`Arg\phi (z)dz^2=Const`$, or is composed of the two radii centered in $`0`$ with the minimal angle $`2\pi /(n+2)`$.

The foliation $``$ on $`M\backslash 0`$ is said to be geodesic if every leaf of $``$ is a $`\phi `$–geodesic line. Before we state the general lemma on the structure of geodesic foliations, let us consider an example when all $``$’s can be obtained by a ”brute force”.

If $`0`$ is a double zero, then the $`\phi `$–metric is given by the linear element $`ds^2=(u^2+v^2)(du^2+dv^2)`$ where $`u+iv`$ is a natural parameter. The metric $`|ds|`$ is Liouville’s and the geodesic lines in this metric are completely integrable. The general integral is known to be of the form

$$\frac{du}{\sqrt{u^2a}}\pm \frac{dv}{\sqrt{v^2+a}}=a^{},$$

where $`a,a^{}`$ are two independent constants. Easy calculations show that no information will be lost if we suppose $`a=0`$. The integral takes the form $`\mathrm{ln}|u|\pm \mathrm{ln}|v|=a^{}`$. The geodesic foliation is described by two ”families of curves”: $`v=Cu`$ and $`v=C/u`$, where C is an arbitrary constant. Thus, $``$ near a double zero is either the ”node” with the geodesics radii tending to $`0`$, or the ”saddle” with four sectors filled–up by the geodesic ”hyperbolas”.

Let $`w`$ be a finite ”word” on the alphabet consisting of two symbols $`h`$ and $`p`$. We introduce the elementary operations on $`w`$:

(i) a cyclic permutation of the symbols in $`w`$, and

(ii) a contraction of the $`p`$–symbol: $`p^2=p`$.

Two words are equivalent $`w_1w_2`$ if and only if $`w_2`$ can be obtained from $`w_1`$ by the elementary operations. The equivalence class of word $`w`$ is denoted by $`[w]`$.

Fix an integer number $`n1`$. To every symbol $`h`$ in $`w`$ we assign a weight $`|h|=2\pi /(n+2)`$. To every symbol $`p`$ we assign the weight $`|p|=\alpha _i`$, where $`\alpha _i`$ is a positive real. The weight of $`w`$ is an additive function equal to the sum of weights of the symbols entering $`w`$. The equivalence class $`[w]`$ is called normalized if $`|w|=2\pi `$ for all $`w[w]`$. (Note that the weight of $`w`$ is one and the same for all $`w[w]`$.)

###### Lemma 2

Let $`h`$ and $`p`$ stay for the hyperbolic and the parabolic sectors of the singularity $`w`$, respectively. We encode the singular point $`w`$ by a sequence of symbols $`h`$ and $`p`$ in the order the $`h`$– and the $`p`$–sectors occur when turning clockwise around the singularity. Then:

(i) each $`\phi `$–geodesic foliation $``$ is topologically equivalent to the singularity $`w`$ of a normalized equivalence class $`[w]`$;

(ii) each normalized equivalence class $`[w]`$ can be realized as a $`\phi `$–geodesic foliation $``$ with the singularity $`w[w]`$ in a neighborhood of $`n`$–th order zero of $`\phi `$ for some $`n1`$.

Proof. Denote by $`M`$ a neighborhood of the $`n`$–th order zero of $`\phi `$. Let us introduce a partial order for the points $`x,yM`$: $`xy`$ if and only if $`ArgxArgy`$. If $`xM`$ is an arbitrary point, then by Lemma 1 the $`\phi `$–geodesic line through $`x`$ is either (i) the hyperbola $`Arg\phi dz^2=Const`$ or (ii) the radius $`Ox`$. Let us consider the first possibility.

(i) The hyperbola $`Arg\phi dz^2=Const`$ must tend to the asymptotic rays $`Oz_1,Oz_2`$ with $`z_1<x<z_2`$, enclosing the angle $`2\pi /(n+2)`$. Clearly, the only possibility to the geodesic foliation $``$ is to form a hyperbolic sector $`z_1Oz_2`$. Of course, along $`Oz_1`$ and $`Oz_2`$ $`Arg\phi dz^2`$ is constant.

(ii) Let $`Ox`$ be the geodesic radius through $`x`$, distinct from the boundary radii of the hyperbolic sector. Then through the nearby points $`|xy|<\epsilon `$ one can draw the geodesic radii $`Oy`$’s. Denote by $`y_1Oy_2`$ the maximal connected parabolic sector filled-up with the geodesic radii. Clearly, $`y_1<x<y_2`$. The angle enclosed between two boundary radii, we denote by $`\alpha `$. In general, $`0\alpha 2\pi `$.

If the hyperbolic sector $`h`$ is followed by another hyperbolic sector $`h`$, we write this as $`hh`$. If $`h`$ is followed by a parabolic sector, we put it as $`hp`$. A parabolic sector $`p`$ followed by the parabolic sector $`p`$, gives a larger parabolic sector $`p=pp`$ and the contraction rule (ii) follows. Of course, the ”weights” of the sectors are equal to the angles swept by the sectors.

Finally, according to the definition of normalized equivalence class, each singularity consists of sequence of parabolic and hyperbolic sectors; every curve in these sectors is a geodesic arc.

The part (ii) of Lemma 2 is proved by the similar argument. $`\mathrm{}`$

## 2 CMC–Surfaces

Every smooth immersion $`f:ME^3`$ of an orientable surface $`M`$ into the Euclidean space $`E^3`$ induces a Riemann structure on $`M`$; let $`z=u+iv`$ be the corresponding local parameter. With respect to $`z`$ the first fundamental form can be written as $`ds^2=e^{2\lambda }|dz|^2`$.

If $`ldu^2+2mdudv+ndv^2`$ is the second fundamental form, we consider a complex quadratic form $`\phi dz^2`$, such that $`\phi (z)=\frac{1}{2}(ln)im`$. The Mainardi–Codazzi equations imply that $`\phi `$ is holomorphic on $`M`$ if and only if $`f(M)`$ is a CMC–surface. Locally, along the lines of minimal and maximal curvature $`Arg\phi dz^2=0`$ and $`\phi (0)=0`$ at the umbilic points.

A continuous deformation $`f_t`$ of the immersion $`f=f_0`$ is the isometry of surface $`M`$ such that $`M\times [0,1]E^3`$ is a continuous mapping. The continuous deformation $`f_t`$ is called an $`H`$–deformation if $`H_t=H`$ for all $`t[0,1]`$, where $`H:MR`$ is the mean curvature function.

The CMC–surfaces are known to admit a non–trivial $`H`$–deformations and in the case of compact surfaces, they are the only ones with such a property. Of course, there are known many examples of compact CMC–surfaces of genus $`g>0`$.

What happens with the lines of principal (i.e. minimal or maximal) curvature of the CMC–surface during an $`H`$–deformation? If $`M`$ is a local CMC–surface without umbilics, the principal curvature lines of $`f_0(M)`$ and $`f_t(M)`$ form two families of the parallel lines intersecting each other with the constant angle proportional to the parameter $`t`$ (the Bonnet Theorem, see e.g. ). Note, that if we fix the $`\phi `$–metric on $`M`$ corresponding to $`f_0(M)`$, then the principal curvature lines of $`f_t(M)`$ coincide with the $`\phi `$–geodesic lines of the inclination $`t`$. If the umbilical points are allowed, then a law is given by the following lemma.

###### Lemma 3

Suppose that $`M_0=f_0(M)`$ is a canonical CMC–surface with the quadratic function $`\phi =z^n,n1`$. Let $`\phi `$ be a metric on $`M`$ corresponding to $`M_0`$. If $`M_t=f_t(M)`$ is an $`H`$–deformation of $`M_0`$, then one of the two principal curvature lines of $`M_t`$ coincide with the $`\phi `$–geodesic lines on $`M`$ for any $`t0`$.

Proof. In the polar coordinate system the coefficients of the second fundamental form of surface $`M_t`$ are given by the equations:

$`l`$ $`=`$ $`He^\lambda +|z|^n\mathrm{cos}(2tnArgz),`$
$`m`$ $`=`$ $`|z|^n\mathrm{sin}(2tnArgz),`$ (1)
$`n`$ $`=`$ $`He^\lambda |z|^n\mathrm{cos}(2tnArgz),`$
where $`t`$ is a parameter of the $`H`$–deformation, cf . The following two cases are possible.

(i) An $`H`$–deformation, such that $`t`$ is constant on $`M`$. It can be immediately seen that in new coordinates $`\stackrel{~}{u}=\mathrm{cos}tu+\mathrm{sin}tv`$, $`\stackrel{~}{v}=\mathrm{sin}tu+\mathrm{cos}tv`$ the first and the second forms of surfaces $`M_0`$ and $`M_t`$ are the same. By the fundamental theorem, surfaces $`M_0`$ and $`M_t`$ may differ only by a rigid motion in $`E^3`$. Thus, the $`H`$–deformation is trivial.

(ii) A non–trivial $`H`$–deformation. By item (i), $`t`$ varies for the points of $`M`$. Thus far, associated to every $`zM\backslash 0`$, there is a chart in which the second fundamental form of surface $`M_t(z)`$ writes as

$$l=He^\lambda +\mathrm{cos}2t,m=\mathrm{sin}2t,n=He^\lambda \mathrm{cos}2t,$$

where $`t`$ is the deformation parameter, cf . A straightforward calculation shows that the principal curvature lines of the surface $`M_t(z)`$ coincide with the $`\phi `$–geodesic lines of the slope $`t`$ on $`M`$. (This fact follows also from the Bonnet Theorem.) Since every regular point $`zM`$ can be endowed with such a chart, Lemma 3 is proved. $`\mathrm{}`$

## 3 Proof of Theorem 1

Take a convex $`C^{\mathrm{}}`$ immersion $`f_0:S^2E^3`$ of the 2–sphere into the Euclidean space $`E^3`$ which is not totally umbilic (i.e. there are no $`US^2`$ such that $`f_0(U)`$ is a part of the round sphere). In other words, umbilics are supposed isolated and their number is finite. Denote by $`ds_0`$ a Riemann metric on $`S^2`$ induced by the immersion $`f_0`$ and by $`H:S^2R`$ the corresponding mean curvature function.

###### Definition 1

By a Hopf spheroid in $`E^3`$ we understand a convex $`C^{\mathrm{}}`$ immersion $`f:S^2E^3`$ such that there exists at least one umbilical point $`p`$ and a small closed disc $`Dp`$ such that $`H(D)=Const`$.

###### Lemma 4

There exist infinitely many Hopf spheroids in $`E^3`$.

Proof. By the results of Wente and Kapouleas any compact orientable surface $`S_g`$ of genus $`g>0`$ admits an immersion into $`E^3`$ which is a CMC-surface with $`H>0`$; cf. , . Fix $`g2`$ and consider the lines of principal curvature of any such immersion. By the index argument, there exists an umbilic $`pS_g`$ and a small closed disc $`Dp`$ which is a convex local surface in $`E^3`$. We separate this local surface from $`S_g`$. To obtain a Hopf spheroid, it remains to complete this piece of CMC-surface to a $`C^{\mathrm{}}`$ immersion $`S^2E^3`$. By Urysohn’s lemma this can be done in an infinite number of ways. $`\mathrm{}`$

###### Lemma 5

For the Hopf spheroids the Caratheodory conjecture is true.

Proof. Without loss of generality we can assume that the umbilic point $`p`$ of Hopf spheroid is unique. (For otherwise, if there are more than one umbilic then we are done.) Since a Hopf spheroid is locally CMC, we apply Lemma 3 to identify the curvature lines in the disc $`Dp`$ with $`\phi `$-geodesic lines in the vicinity of a singularity $`w`$.

Let $`w[w]`$ be a word of the minimal length in the normalized equivalence class $`[w]`$. According to Lemma 2, there exists a singularity of order $`n`$ whose topological type is encoded by the sequence $`w`$ of symbols $`h`$ and $`p`$. Let $`w`$ admit $`h`$ symbols of type $`h`$ and $`p`$ symbols of type $`p`$. By the normalization axiom, $`hn+2`$.

To estimate the Euler–Poincaré index of singularity $`w`$, note that the parabolic sectors make no contribution to the index value and the number $`p`$ can be neglected. To the contrary, if there are no hyperbolic sectors (i.e. $`w=p`$) we necessarily have one parabolic sector. The general formula is true:

$$Indw=\{\begin{array}{cc}1\frac{h}{2}\text{if}wp,\hfill & \\ +1\text{if}w=p.\hfill & \end{array}$$

In either case $`Indw1`$ and by the index argument the conjecture follows. $`\mathrm{}`$

Now we are ready to finish the proof of Theorem 1. But first we wish to outline the main idea. To every convex $`C^{\mathrm{}}`$ immersion $`f_0:S^2E^3`$ one can relate a Hopf spheroid. This spheroid is uniquely defined by $`f_0`$ and is a ‘modification’ of $`f_0`$ which has an interesting ‘mechanical’ interpretation.

Suppose that $`f_0`$ is a convex steel ball filled-up with a gas under a pressure. Let $`p`$ be an isolated umbilic of $`f_0`$. We drill a small hole in $`p`$ and glue-up a soap film $`D`$ into this hole maintaining a pressure 222The absolute value of the pressure depends on how ‘flat’ is the surface at the point $`p`$. Of course, by ‘pressure’ we understand difference of pressures inside and outside the steel ball. inside the ball. We also ‘deform’ slightly the ‘edges’ of the cut in order to keep the modified surface $`f:S^2E^3`$ in the class $`C^{\mathrm{}}`$. We claim that $`f`$ is a Hopf spheroid.

Indeed, $`f(D)`$ is a local CMC-surface with an umbilic point $`pD`$. Moreover, the index of umbilic on the Hopf spheroid is equal to the index of $`p`$ on $`f_0`$. (This is because the foliation by principal curvature lines at the ‘steel part’ of ball remains intact.) In general, if $`_0`$ and $``$ are foliations by the principal curvature lines on $`f_0`$ and $`f`$, respectively, then $``$ is obtained from $`_0`$ by a homotopy of openning of a leaf; cf .

Let $`f_0`$ be as above. If $`p`$ is an isolated umbilic of $`f_0`$ then we take a closed disc $`|D|r`$ centred at the point $`p`$. We are going to define a local CMC-surface $`f(D)`$. Let $`z=u+iv`$ be a local parameter which corresponds to a part of CMC-surface with an umbilic; see the beginning of this section. By the results of Umehara  (see also , ) there exists a family of isometric $`H`$-deformations depending on a real parameter $`t`$:

$$I=e^{2\lambda }|dz|^2,II_t=ldu^2+2mdudv+ndv^2,$$
(2)
with $`l,m`$ and $`n`$ given by equations (2). The Mainardi-Codazzi and Gauss equations for $`I,II_t`$:

$$\frac{\phi }{\overline{}z}=\frac{H}{z},|\phi |^2=e^{4\lambda }(H^2K),$$
(3)
where $`\phi =e^{it}z^n`$ is a complex quadratic form $`\phi dz^2`$, are satisfied for any real $`t`$. (Indeed, the first equation is true since $`H=Const`$ and $`\phi `$ is holomorphic; the second equation follows from $`|\phi e^{it}|=|\phi |`$ and the fact that $`H`$-deformation is an isometry.) Therefore, the fundamental forms (2) are realized by a concrete local CMC-surface for each real number $`t`$.

Let $`f_t(D)`$ be a family of local CMC-surfaces described above. Denote by $`A`$ an annular region which surrounds disc $`D`$:

$$A=\{z=u+iv|r|z|r+\epsilon \}.$$
(4)
To glue-up $`f_t(D)`$ properly, we fix the metric $`\lambda `$ so that $`\lambda |_{A_{r+\epsilon }}=\lambda |_{A_r}`$, where the left part denotes a metric on the exterior boundary of $`A`$ which is induced by metric of the surface $`f_0`$. The boundary condition $`\lambda |_{A_r}`$ gives a unique solution $`f_{t=t^{}}(D)`$ to the Gauss equation, so that a representative in the family $`f_t(D)`$ is fixed.

To obtain a $`C^{\mathrm{}}`$ Hopf spheroid it remains to conjugate $`f_t^{}(D)`$ with the rest of the sphere:

$$f(S^2)=\{\begin{array}{cc}f_t^{}(D),\text{if}zDIntD_{r+\epsilon },\hfill & \\ f_0(S^2),\text{if}zS^2\backslash D_{r+\epsilon }.\hfill & \end{array}$$
(5)
By the Urysohn Lemma, function $`f`$ in formula (5) can be chosen $`C^{\mathrm{}}`$ for an arbitrary small $`\epsilon `$, see formula (4). Moreover, taking $`r`$ sufficiently small we can fix number $`n`$ (see (2)) equal to the order of quadratic form $`\phi `$ at point $`p`$ of the surface $`f_0`$. (Such an order is correctly defined for any $`\phi `$, not necessary holomorphic.)

Thus, the surface $`f`$ given by equation (5) is a Hopf spheroid. By the Lawson-Tribuzy theorem $`f`$ is uniquely defined up to a rigid motion in $`E^3`$; see . To finish the proof of Caratheodory conjecture, it remains to notice that passage from $`f_0`$ to $`f`$ gives us a homotopy $`h(_0)=`$ between foliations induced by curvature lines. In particular, $`Indp_0=Indp`$. By Lemma 5, the Caratheodory conjecture follows. $`\mathrm{}`$


text:
# Market price simulator based on analog electrical circuit

## 1 Introduction

Random multiplicative process (RMP) is attracting researchers as a new mechanism of generating fat tails in the distribution of price changes in open markets recently (Takayasu et al. 1997, Jögi et al. 1998, Sato et al. 2000). One of the mechanisms of the fat tails is a feedback mechanism of both positive and negative in a random manner. The positive feedback plays a role of amplification while the negative one plays a role of damping. In the presence of positive feedback one may consider that the system is unstable, however, this is not always correct. When the damping effects dominates the amplification the system can generate statistically stationary fluctuations.

Recently fat tails in the distribution have been confirmed in the market price fluctuation (Mantegna et al. 1995). Among several approaches to this phenomena (Bak et al. 1997, Lux et al. 1999), the authors demonstrated that RMP can be derived from a microscopic artificial dealer’s market simulation model (Sato et al. 1998). The RMP has also been derived from a macroscopic theoretical analysis of the critical dynamics of the balance of demand and supply (Takayasu et al. 1999). The multiplicative effect of the market comes from the dealer’s forecasting which follows the trend of the latest market price changes, which obviously enhances the fluctuations in price.

In this article we review the analog electrical circuit which generates power law fluctuations based on the RMP theory. The application of this circuit is to generate fluctuation which is statistically similar to market price change in inexpensive way at high speed. The goal of this study is a risk estimator that calculates realizations of market prices directly by parallel processing units and calculates a statistical feature of the fluctuation from these simulated data in the meaning of an ensemble (see. fig. 1.).

The outline of this article is the followings. In sec. 2 we analyze tick data of yen/dollar exchange rates, and show that a probability density function of the market price changes has fat tails of roughly power law. In sec. 3 a random multiplicative process is formalized by a stochastic differential equation with a multiplicative noise. It is indicated that the probability density of a dynamical variable follows a power law distribution. In sec. 4 an example of electrical circuit diagram is shown and an equivalent equation of the proposed circuit is described. In sec. 5 we show the result of observation and we will discuss them. Sec. 6 is devoted to the concluding remarks.

## 2 Statistical properties of foreign currency rate

Mantegna et al. investigated time series of stock market index of S&P500 and reported in their famous paper (Mantegna et al. 1995) that a probability density function of changes obeys a power law distribution. They indicated that the power law exponent is estimated as 1.4. However, it is clarified that the exponent is not universal by other researches. Besides stock market prices it is known that the probability density function of foreign exchange rates also follows a power law distribution (Takayasu et al. 2000). We analyze tick data of yen/dollar exchange rate of 578,508 data points from January to March 1999. Let us denote the rate of the $`s`$th tick as $`r_s`$. The change of the rate is defined as $`\mathrm{\Delta }r_s=r_sr_{s1}`$. Fig. 2 shows typical examples of time series of yen-dollar exchange rates and corresponding changes. The cumulative distribution function (CDF), which is defined by

$$F(|x|)=_{\mathrm{}}^{|x|}p(x^{})𝑑x^{}+_{|x|}^{\mathrm{}}p(x^{})𝑑x^{},$$
(1)
is shown in fig. 3. It is clear that the CDF has a liner slope with a quick decay by 1 yen in log-log plots. From the liner slope in log-log plots we can estimated the exponent of the power law distribution,

$$F(|x|)|x|^\beta .$$
(2)
The best fit value of $`\beta `$ we obtain from the exchange rate is $`\beta =1.8`$.

## 3 Theory

In this section we show a brief outline of a continuous time version of random multiplicative process. The random multiplicative process is given by a stochastic differential equation,

$$\frac{dv}{dt}=\nu (t)v(t)+\xi (t),$$
(3)
where $`v(t)`$ is a dynamic variable and $`\nu (t)`$ and $`\xi (t)`$ represent a multiplicative noise and an additive noise, respectively. Denoting an ensemble average by $`\mathrm{}`$, we assume the following relations for both multiplicative and additive noises based on the Gaussian white noise theory.

$`\nu (t)`$ $`=`$ $`\overline{\nu },`$ (4)
$`[\nu (t_1)\overline{\nu }][\nu (t_2)\overline{\nu }]`$ $`=`$ $`2D_\nu \delta (t_1t_2),`$ (5)
$`\xi (t)`$ $`=`$ $`0,`$ (6)
$`\xi (t_1)\xi (t_2)`$ $`=`$ $`2D_\xi \delta (t_1t_2),`$ (7)
where $`\overline{\nu }`$ represents an average of the multiplicative noise, and $`D_\nu `$ and $`D_\xi `$ represent the strength of the multiplicative noise and that of the additive noise, respectively.

The probability density function of $`v(t)`$ in eq. (3) is denoted as $`p(v,t)`$, is known to follow the generalized Fokker-Planck equation (Deutsch 1994, Venkataramani et al. 1996, Nakao 1998).

$$\frac{}{t}p(v,t)=\frac{}{v}\left[(\overline{\nu }+D_\nu )vp(v,t)+\frac{}{v}(D_\nu v^2+D_\xi )p(v,t)\right].$$
(8)
By solving eq. (8) for a steady state $`\frac{}{t}p(v,t)=0`$ and boundary condition, $`\frac{}{v}p(\pm v_u,t)=0`$, we get the stationary distribution.

$$p(v)(D_\xi +D_\nu v^2)^{\frac{\overline{\nu }}{2D_\nu }\frac{1}{2}},$$
(9)
which has power law tails for large $`v`$,

$$p(v)|v|^{\beta 1},$$
(10)
where $`\beta =\frac{\overline{\nu }}{D_\nu }`$. In other words the power law exponent is given by a simple function of the average and variance of the multiplicative random noise.

## 4 The Circuit

From the assumption for the multiplicative noise in the RMP theory described in sec. 3, it is important that $`\nu (t)`$ takes both positive and negative values to realize the power law tails. With an electrical analog circuit this means that it is necessary for the circuit to include both positive and negative feedbacks in the meaning of probability. We solve this problem by using an analog multiplier. Our block diagram of circuit and an implemented noise generator are shown in fig. 4. $`v_o`$ in the figure represents an output voltage, and $`\mu (t)`$ is an output directly from the noise generator. As it is seen from this figure the circuit contains the noise generator, an analog multiplier (Analog devices, 10MHz, 4-quadrant) and an operational amplifier (National Semiconductor, LF157) for integrator (Sato et al. 2000).

The output of the noise generator plays a role of the multiplicative noise in eq. (3). In the noise generator a shot noise of the zener diode in fig. 4 is amplified by an operational amplifier, and it is output from generator through high pass filter. LF157 has the 20M product of a voltage gain (G) and a bandwidth (B). The bandwidth is given by B=200kHz in the noise generator because we put G=100. Thus we expect a frequency characteristics of the noise generator to be up to 200kHz . Let us explain how to realize both positive and negative feedbacks in the meaning of probability. We realize both positive and negative feedbacks by multiplying $`v_o`$ with the output of the noise generator $`\mu (t)`$ in the analog multiplier and by connecting the product to the negative input of the operational amplifier for integrator. Although the additive noise term, $`\xi (t)`$, is not explicitly added in the circuit, it derives from either thermal noises of the operational amplifier or from an external electro-magnetic noises.

An equation equivalent to the circuit diagram of fig. 4 is given as

$$\frac{dv_o}{dt}=\left(\frac{1}{R_fC}+\frac{k}{R_vC}\mu (t)\right)v_o+\xi (t),$$
(11)
where $`k`$ is a factor of the multiplier and $`k=1/10`$. $`\xi (t)`$ represents the additive noise effect. The strength of multiplicative noise $`\mu (t)`$ depends on the value of a variable resistor $`R_v`$ in fig. 4 because $`R_v`$ is a factor of $`\mu (t)`$ in eq. (11). Therefore, we expect that the power law exponent for the output is a function of $`R_v`$.

## 5 Results of observation and discussion

We measure the output of the circuit $`v_o`$ for $`R_v`$ since we expect that the power law exponent $`\beta `$ is a function of $`R_v`$ from the above discussion. We obtained the output through a 12-bit AD converter (Microscience, ADM-652AT) and processed it as digital data in digital computer. A sampling frequency is 125kHz throughout all the observations.

We show a typical example of time series of $`v_o`$ at $`R_v=50\mathrm{\Omega }`$ and $`100\mathrm{\Omega }`$ in fig. 5. The time series at $`R_v=50\mathrm{\Omega }`$ sometimes exhibits larger fluctuations than $`R_v=100\mathrm{\Omega }`$. We show log-log plots of CDF of $`v(t)`$ at $`R_v=25`$,$`50`$,$`75`$,$`100`$$`\mathrm{\Omega }`$ in fig. 6. Log-log plots of the CDF have liner parts for about a decade. We clearly find that the exponent $`\beta `$ depends on $`R_v`$ as expected qualitatively. Moreover as shown in fig. 7 we find a linear relation between $`\beta `$ and $`R_v`$. The autocorrelation function $`R(\tau )`$ and the volatility autocorrelation function $`R^{(2)}(\tau )`$ are defined by

$`R(\tau )`$ $`=`$ $`v_o(t+\tau )v_o(t)v_o(t+\tau )v_o(t),`$ (12)
$`R^{(2)}(\tau )`$ $`=`$ $`v_o(t+\tau )^2v_o(t)^2v_o(t+\tau )^2v_o^2(t).`$ (13)
We show the autocorrelation function and volatility autocorrelation function in fig. 8. Both autocorrelation function and volatility autocorrelation function have quick decay. It is known that the volatility autocorrelation function of real market price changes has a long time correlation. The disagreement of the volatility autocorrelation occurs because the model equation is too simple to describe time structure of fluctuations. We need to study higher order differential equations with multiplicative noises.

## 6 Conclusions

We proposed an analog electrical circuit as an analog generator of power law fluctuations. We described a theoretical equation representing the electrical circuit and showed that the probability density function of a dynamical variable has power law tails. We measured the output of the circuit and observed its cumulative distribution functions. The cumulative distribution of output voltage has power law tails. The power law exponent can be tuned by controlling the variable resistance $`R_v`$. We expect that the proposed circuit is applicable to generate fluctuations having power law distribution in a much cheaper way than any digital computing methods. Moreover, fluctuations of the circuit may be of use for risk estimation in foreign exchange or stock market in the near future.

Acknowledgment

One of the authors (A.-H. Sato) thanks to Yoshihiro Hayakawa for stimulative discussion. This research is partially supported by Japan Society for the Promotion of Science.


text:
# Evolution of Active Galaxies: Black-hole Mass – Bulge relations for Narrow Line Objects

## 1 Introduction

Last couple of years have seen significant progress in our understanding of the relationship between quasar activity and galaxy formation. In the framework of the hierarchical dark-matter cosmogony (e.g. Haehnelt, Natarajan & Rees 1998, HNR98 hereafter), the formation and evolution of galaxies and their active nuclei is intimately related (Fabian 1999, Granato et al. 1999, Haehnelt & Kauffmann 1999). It is assumed that the supermassive black holes form proportionally to the dark matter halos. The process of formation of a massive BH and the active nucleus is the very process of galaxy formation. The active nucleus and the galaxy evolve together, with BH accreting matter and the galaxy making stars. At one stage the winds from the active nucleus blow away the matter surrounding it and a quasar emerges. This is not only the end of active evolution of the quasar but that of the galaxy as well, which is evacuated of its interstellar medium. The quasar then shines as long as there is fuel in the accretion disk (Fabian 1999). In this scenario, high-redshift quasars represent the early stage of galaxy evolution, BALQSOs at z$`2`$ represent the stage when the gas is being blown away, and z$`1`$ quasars would be the passively evolving population. Massive ellipticals found today might be the dead remnants of what were once luminous quasars.

Observational support to the above scheme of evolution is also accumulating. In a classic paper, Maggorian et al. (1998) found that the mass of a BH in the nucleus of a galaxy is correlated with the mass of its bulge (see also van der Marel et al. 1999 and an excellent review by Richstone et al. 1998). Most recently, Gebhardt et al. (2000a) and Ferrarese & Merritt (2000) found a correlation between the BH mass and the bulge velocity dispersion, a much tighter relation than that between the BH mass and bulge luminosity. The important point to note is that the nuclear BH seems to “know” the galaxy it harbors. Laor (1998) showed that the above correlation for normal galaxies extends to quasar host galaxies as well. The similarity of the two correlations is remarkable given the widely different levels of BH activity in luminous quasars and nearby massive normal galaxies. Franceschini et al. (1999) found an interesting similarity between the evolution rates for the total populations of galaxies and AGNs. By comparing the evolution of star formation rate per comoving volume in galaxies and the volume emissivity by AGNs, they found that the evolution of luminous quasars parallels that of stellar populations in massive spheroidal galaxies.

Much of the above observational evidence relates BH masses of luminous quasars, or what were once luminous quasars, with spheroidal masses of massive galaxies. Wandel (1999) extended Laor’s work to low redshift Seyfert galaxies for which BH masses were available using reverberation mapping. Interestingly, Wandel (1999) found that the average BH- to bulge-mass ratio for Seyferts is systematically lower than that for normal galaxies and luminous quasars. Wandel interpreted this result as a scenario of BH growth in which a BH in a Seyfert grows by accretion to become a more massive BH in a quasar. This, however, is highly unlikely since the Seyferts are known to accrete at a sub-Eddington rate while the luminous quasars in Laor’s sample accrete at close to Eddington rate. So, in the above scenario, the accretion rate has to grow with time, which is difficult to achieve and contrary to theoretical expectations. A possible explanation of the Laor and Wandel results might be that the quasars and Seyferts represent two different populations (Mathur 2000). The quasar phenomenon may be a result of galaxy formation due to primordial density fluctuations or mergers at early cosmic epochs, which were to result in massive ellipticals of today. At low redshifts, when new galaxies are formed due to interactions or mergers, similar evolution may take place. Gebhardt et al. (2000b) have now found that the Seyferts and normal galaxies have the same BH mass–bulge velocity dispersion relationship. Note, however, that it still implies that the Seyferts are different from the normal galaxies in that the bulges of Seyferts have lower mass to light ratio, unless the Wandel result is purely due to observational uncertainty.

While we now have the observational evidence for a relation between an AGN and its host galaxy, we do not know how it got there, or what are the evolutionary steps towards the “final” correlation discussed above. Mathur (2000) has argued that narrow line Seyfert 1 galaxies (NLS1s) may represent an early phase in the evolution of active galaxies at low redshift. In our scenario, the luminous quasars and the low-luminosity Seyferts follow separate but similar evolutionary tracks. In both cases the accretion rate ṁ = Ṁ$`/`$Eddington is large in the early stages of evolution and reduces later on, consistent with the theoretical expectations (e.g. HNR98). To test this hypothesis, we determine the BH- to bulge-mass relation as well as the BH-mass to bulge velocity dispersion relation for NLS1s and compare those with the normal Seyfert relation. In $`\mathrm{\S }2`$ we detail our method and present interpretation in $`\mathrm{\S }3`$.

## 2 BH and Bulge Relations for NLS1s

### 2.1 Determination of BH Mass

Our sample consists of 15 narrow-line AGN (11 NLS1s and 4 NL quasars) for which soft- and hard X-ray spectra are available in literature. The central black hole masses for these objects were calculated by fitting their spectral energy distributions with the accretion disk and corona (ADC) model of Witt, Czerny & Życki (1997) (see Kuraszkiewicz et al. 2000 for details). In this model the corona accretes and generates energy through viscosity, and the division between optically thin accretion flow (corona) and optically thick flow (disk) results from the cooling instability discussed by Krolik, McKee and Tarter (1981). The model is fully described by the following parameters: the mass of the central black hole, accretion rate and viscosity parameter. The calculated black hole masses are quoted in Table 1.

We determine the errors of the black hole mass estimated from spectral fitting as follows. The fit is based on two relations: (1) the asymptotic shape of the disk emission in the optical band is: $`F_\nu \nu ^{1/3}(M\dot{M})^{2/3}cosiH_0^2f^{4/3}`$, where $`M`$ is mass of the black hole, $`\dot{M}`$ \- accretion rate, $`i`$ \- the inclination angle of the disk, $`H_0`$ is the Hubble constant, f - spectral correction given by color temperature to the effective temperature ratio; (2) the bolometric luminosity is: $`F_{bol}\dot{M}\eta H_0^2`$, where $`\eta `$ is the efficiency of accretion. In our fits we adopt the following values: $`H_0=50`$ km/s/Mpc, $`cosi=1`$, $`f=1`$ and $`\eta =1/16`$. If any of these values are different, the mass of the black hole would scale as: $`M_{real}=M_{calc}f^2(cosi)^1(H_0/50)^1`$. The inclination angle of the disk for Seyfert 1s is in the range $`0^o`$ to $`48^o`$ (Schmitt et al. 2000 and assuming the accretion disk is aligned with the dusty torus), the upper limit for f is 1.7 (Różańska, private communication), and if there is strong outflow close to the black hole $`\eta `$ decreases (up to a factor 2 if the outflow is energetically comparable to the outflows in radio-loud objects). Taking all these effects into account we estimate the calculated mass of the black hole to be accurate to $`{}_{0.30}{}^{}{}_{}{}^{+0.64}`$ in logarithm.

For the sake of consistency and comparison with earlier results, we also calculated the BH masses for the 7 Sy1s from Wandel’s sample, for which X-ray spectra were available. The results are presented in Table 2. The calculated black hole masses are compared to the virial masses in Tables 2 and 3. The masses calculated by the ADC model are generally larger than masses estimated by Wandel. We calibrate our BH mass determinations with the virial mass estimates for the objects in Table 2 and 3 and find that: $`\mathrm{log}M_{\mathrm{bh}}^{\mathrm{calib}}=\mathrm{log}M_{\mathrm{bh}}^{\mathrm{calc}}0.70\pm 0.15`$. This relation, together with the calculated and virial BH masses for Seyfert 1s and NLS1s (excluding the outlying point for IC 4329A111The virial mass for IC 4329A is adopted from Wandel 1999 who has adopted it from Wandel, Peterson & Malkan 1999. We find that IC 4329A is an outlier in all the relations discussed in Wandel, Peterson & Malkan 1999 as well as in Peterson et al. 2000. This leads us to suspect the virial BH mass of this object. We ignore IC 4329A in rest of the paper.), is displayed in Figure 1. It can be clearly seen that the above equation is a good description of the relation between BH masses calculated by the ADC model and by the virial method. There is clearly some uncertainty in the ADC mass estimates, as discussed above (Tables 1, 2, 3), and the virial mass estimates are also uncertain by a factor of several because of the uncertainty in the geometry of the broad emission line region. Both of these effects contribute towards the scatter seen in Figure 1. The calibrated black hole masses for NLS1s are given in Table 3. These masses were used to calculate the black hole to bulge mass ratio.

### 2.2 Determination of Bulge Mass

For the NLS1s in our sample, we determined the bulge masses of their host galaxies following the general method used by Laor (1998) and Wandel (1999), which we discuss briefly below. For most of our NLS1s we obtain the bulge absolute blue magnitudes ($`M_\mathrm{B}^{\mathrm{bulge}}`$) from Whittle (1992). A galaxy’s total blue magnitude ($`M_\mathrm{B}^{\mathrm{total}}`$) is related to $`M_\mathrm{B}^{\mathrm{bulge}}`$ by Simien & Vaucouleurs (1986) equation (for $`H_0=50`$ km s-1 Mpc-1, $`q_0=0`$): $`M_\mathrm{B}^{\mathrm{total}}M_\mathrm{B}^{\mathrm{bulge}}=0.324\tau 0.054\tau ^2+0.0047\tau ^3`$ where $`\tau =T+5`$ and T is the Hubble stage (defined in de Vaucouleurs, de Vaucoulers & Corwin 1976).

Mrk 110 and Mrk 335 do not have a Hubble type defined in the Whittle (1992) compilation. However NED222This research has made use of the NASA/IPAC Extragalactic Database (NED) which is operated by the Jet Propulsion Laboratory, California Institute of Technology, under contract with the National Aeronautics and Space Administration. quotes a S0/a Hubble type for Mrk 335, which we use to calculate the bulge magnitude for this object. The resulting $`M_\mathrm{B}^{\mathrm{bulge}}`$ (Table 3) is consistent with that in Ho (1999). For Mrk 110 we adopt a canonical Hubble type of Sa.

For NAB 0205+024, PG 1402+261 and PG 1444+407 we take the host galaxy absolute magnitudes from Bahcall et al. (1997), which have been corrected for nuclear emission. For Mrk 1044 we take the host galaxy magnitude from MacKenty (1990), who has nuclear emission included in the total blue magnitude, hence in Table 3 we quote the bulge mass as an upper limit.

After obtaining the bulge blue magnitudes we use the relation between the bulge B and V magnitude (Worthey 1994) B$``$V=0.95, the standard relation between the bulge luminosity and absolute magnitude: $`\mathrm{log}(L_{\mathrm{bulge}}/L_{})=0.4(M_\mathrm{V}^{\mathrm{bulge}}+4.83)`$ and the relation between bulge mass and bulge luminosity for normal galaxies from Magorrian et al. (1998): $`\mathrm{log}M_{\mathrm{bulge}}/M_{}=1.11+1.18\mathrm{log}L_{\mathrm{bulge}}/L_{}`$. The bulge masses calculated this way are listed in Table 3.

### 2.3 Determination of Bulge Velocity Dispersion

We do not have bulge velocity dispersion ($`\sigma _{}`$) measurements for the objects in our sample. However, Nelson (2000) has shown that the width of the narrow emission line \[OIII\] serves as a good representation of the velocity dispersion with $`\sigma _{}`$ = FWHM \[OIII\]/2.35. Adopting from Nelson (2000), we collected from literature, as well as from our optical data, the values of FWHM \[OIII\] for all the NLS1s in our sample (Osterbrock & Pogge 1985, Wamsteker et al. 1985, Dahari & De Robertis 1988, Busko & Steiner 1988, Nelson & Whittle 1995, Wilkes et al. 1999).

## 3 Results and Discussion

### 3.1 BH- to Bulge-mass relation for NLS1s

In Figure 2 we have plotted the $`\mathrm{log}`$ of BH mass vrs. $`\mathrm{log}`$ of bulge mass obtained using the method discussed in $`\mathrm{\S }2`$. The dashed line represents the relation for normal galaxies and luminous quasars (from Laor 1998) while the dotted line represents the relation for Seyfert galaxies (from Wandel 1999, excluding 3 NLS1s in his sample; we have also re-calculated the bulge masses for these objects, as Wandel uses a slightly different relation between the mass and luminosity of the bulge than we do). The NLS1s are represented by filled squares. It is immediately apparent that the relation for Seyfert galaxies does not hold for NLS1s; NLS1s lie systematically below the line for Seyfert galaxies. In addition, we find that the narrow line quasars also lie systematically below the line for quasars. The mean $`M_{\mathrm{BH}}/M_{\mathrm{bulge}}`$ for NLS1s is 0.00005, lower by a factor of six compared to $`M_{\mathrm{BH}}/M_{\mathrm{Bulge}}=0.0003`$ for Wandel’s Seyferts, excluding the three NLS1s in his sample. Similarly, the BH to bulge mass ratio for NL quasars is 0.0005, almost an order of magnitude lower than 0.006 for normal galaxies and quasars.

We would like to emphasize that, since we have calibrated our BH masses with those from the reverberation mapping, the relative difference in the BH- to bulge-mass relation of NLS1s and normal Sy1s is a robust result, though the absolute values may suffer from various uncertainties. We have not been able to do such a calibration for NL quasars, so we cannot place similar confidence in the difference between NL quasars and the quasars in Laor’s sample.

Comparison of BH mass to bulge mass, using bulge luminosity, is a worthwhile exercise even in the light of the new, stronger, correlation with bulge velocity dispersion. Reasons for the differences in two correlations may be many: e.g. Hubble types of many bright AGNs may have systematic errors, or nuclear starbursts may reduce the mass to light ratio of bulges in AGNs (Gebhardt et al. 2000b). Understanding these differences is important. Moreover, here we are comparing two subclasses of Seyfert galaxies: regular, broad line Seyfert 1s and NLS1s, so any systematic problems in determining the bulge luminosity is likely to be the same in both classes. So, any difference in the two is not likely to be a result of observational uncertainty.

### 3.2 BH-mass to Bulge-velocity dispersion relation for NLS1s

Nelson (2000) has determined BH-mass to bulge velocity dispersion relation for Seyfert galaxies. It is noteworthy that all the three NLS1s in his sample: NGC 4051, Mkn 110, and Mkn 335, lie at the lowest boundary of the scatter. In figure 3 we have plotted the $`\mathrm{log}`$ of BH mass vrs FWHM \[OIII\]/2.35. The open squares are the Seyfert 1s and the filled circles are the NLS1s in our sample. The straight line represents the relation for Seyfert galaxies from Nelson (2000). While the Seyfert 1s follow the Nelson relation, NLS1s lie systematically below the line.

### 3.3 Interpretation

The above results support our conjecture that NLS1s are young AGNs. NLS1s are believed to accrete at close to Eddington rate, as first proposed by Pounds, Done & Osborn (1995). As discussed earlier, in our scenario the accretion rate is highest in the early growing stage of an AGN. The BH would grow by accretion and eventually the accretion rate would slow down. Whether the bulge mass would grow or not would depend upon how young/ rejuvenated the galaxy is, or how complete is the interaction/ merger. So an AGN will follow a vertically upward or diagonally upward (with both $`M_{\mathrm{BH}}`$ and $`M_{\mathrm{Bulge}}`$ increasing) track on Figure 2, as it evolves. Similarly, if the bulge velocity dispersion has settled down to its “final” value, an AGN would follow a vertically upward track on Figure 3.

An alternative interpretation of the BH- to bulge-mass relation might be that the bulge mass to light ratio in NLS1s is lower compared to normal Seyferts (see $`\mathrm{\S }1`$). This might be due to higher bulge luminosity as a result of a nuclear starburst, again supporting our proposal of NLS1s being rejuvenated AGN (Mathur 2000). This, however, cannot explain the result regarding velocity dispersion. Alternatively, the above results might mean that the bulge luminosities and FWHM \[OIII\] are systematically overestimated for NLS1s. While not impossible, it seems unlikely that both the quantities are overestimated.

Assuming that our results are not just due to observational uncertainties, and that our scenario of BH growth is valid, some interesting inferences may be obtained. It can be seen from figure 2 that for NLS1s, the relation between BH mass and bulge mass is not strictly linear. So, the BH mass is not a “universal fraction” of the baryonic mass as suggested by Haiman & Loeb (1998). More likely, the accretion process determines the BH mass (e.g. HNR99) at any given time. The Salpeter time for the growth of a BH, i.e. the e-folding time, is $`t_s=4\times 10^7(L_{Eddington}/L)\eta _{0.1}`$ years where $`\eta _{0.1}`$ is the radiative efficiency normalized to 0.1. $`L/L_{Eddington}`$ calculated for the objects in our sample is listed in Table 1. For a bulge with mass $`10^{11}`$M, we have a BH mass of $`5\times 10^6`$M for NLS1 Mrk 493. (Figure 2, Table 3). Let us assume this to be the initial BH mass. The BH would grow to a “final” Seyfert mass of $`5\times 10^7`$M assuming $`M_{\mathrm{BH}}/M_{\mathrm{Bulge}}=0.00047`$. So to grow by a factor of 10, a BH would need $`2.3t_s=3.3\times 10^8\eta _{0.1}`$ years, for a constant accretion rate. Since the accretion rate decreases with time, the growth time scale t$`\genfrac{}{}{0pt}{}{_>}{^{}}`$$`3.3\times 10^8\eta _{0.1}`$ years.

It is possible that being more luminous, quasars evolve faster than the Seyferts. However, if the evolution of quasars and of Seyferts is indeed similar, then quasar lifetimes of about $`10^8`$ years is inferred. This timescale is much larger than $`10^6`$ yr. proposed by Haiman & Leob (1998). HNR98 predict quasar lifetimes anywhere in the range $`10^6`$ to $`10^8`$ yrs. Our results are consistent with, and close to their upper limit $`t_Q`$$`\genfrac{}{}{0pt}{}{_<}{^{}}`$$`4.5\times 10^8\eta _{0.1}`$ years calculated for $`L/L_{Eddington}=1.0`$. HNR98 have considered three different scenarios for BH growth in a quasar: (i) accretion far above the Eddington rate, (ii) accretion obscured by dust, and (iii) accretion below the critical rate leading to an advection dominated accretion flow lasting for a Hubble time. Our results, together with Laor (1998) make the third option unlikely. For NLS1s, we have evidence for accretion rate close to Eddington rate, but not far above the Eddington rate. There are also several NLS1s which are dusty (e.g. IRAS 13349+2438). So options (i) and (ii) are plausible. Studies of low redshift NLS1s thus provide a unique and relatively easy way to understand the high redshift quasars and their evolution.

We gratefully acknowledge financial support through NASA grant NAG5-8913 (LTSA to SM) and NAG5-3363task22 (JK).

Figure Captions:

Figure 1: $`log`$ of BH masses estimated by reverberation mapping is plotted against that calculated with ADC model. Filled symbols are for Seyfert 1s and the open symbols are for NLS1s. The straight line represents the relation between the two mass estimates (see text).

Figure 2: BH mass plotted against bulge mass. The NLS1s, represented by filled squares lie systematically below the dotted line describing the Seyfert relation.

Figure 3: BH mass plotted against the bulge velocity dispersion. Seyfert 1s are represented by open squares and NLS1s by filled circles. The straingt line describes the AGN relation from Nelson (2000).