text:
# Feeding upon negative entropy in a thermal-equilibrium environment

## Abstract

The validity of the Second Law of thermodynamics, indisputable in the macroscopic world, is challenged at the mesoscopic level: a mesoscopic isolated system, possessing spatial dimensions of the order of a few microns, is capable, as shown by a straightforward kinetic analysis, to exhibit a perpetuum mobile behavior associated with large negative variations of the Clausius entropy of the system. This violation of the Second Law is expedient for devising a cyclic process through which an isolated system can extract energy from a surrounding thermal bath.

As first clearly stated by Schroedinger “a living organism tends to approach the dangerous state of maximum entropy, which is death. It can only keep aloof from it, i.e., alive, by continually drawing from its environment negative entropy. What an organism fed upon is negative entropy”.Schroedinger  Is it possible for an organism immersed in a thermal bath to be insulated, and thus avoid thermalization, and still be able to reduce its entropy? Unfortunately, the above scenario contradicts one of the most cherished laws of physics, that is the Second Law of Thermodynamics, which has victoriously resisted all the attempts aimed at finding a particular case where it could be violated. However, one has to consider that the Second Law has been formulated in the context of macroscopic physics and it is in this context that it has been successfully applied and verified. While the Second Law, being inherently statistical in its nature, cannot carry over to microscopic cases where the number of involved particles is too small, its limits of validity are not well understood in the region bridging the macroscopic to the microscopic world. In particular, the intermediate mesoscopic regime, which includes part of the biological realm, is still terra incognita and open to possible surprises. Recent years have actually witnessed a growing interest in the thermodynamics of small-scale non-equilibrium devices, especially in connection with the operation and efficiency of Brownian motors.P1  We believe that we have been able to devise a particular case in which an isolated object, placed in a thermal-equilibrium environment, is nevertheless capable of reducing its entropy or, equivalently, to extract work from the thermal bath surrounding it. We have identified this object with a mesoscopic system, typically of the size of a few microns, consisting of two isolated cavities, filled with a common number of molecules, separated by a movable adiabatic impermeable partition (see Fig.1). Its macroscopic equivalent is the so-called adiabatic piston, sometimes referred to as “enigmatic” because the description of its dynamical evolution toward equilibrium is far from being trivial, the problem being actually undetermined from the point of view of elementary thermodynamics.P3  A careful analysis of the problem, which has been tackled a few years ago by an approach based on gas kinetic theory,P4  shows that thermodynamic fluctuations cannot be neglected at the mesoscopic level and that they actually play a major role in determining the dynamic evolution of this isolated system. As shown in this paper, its behavior implies that the system fails to settle down to a well-defined equilibrium state and undergoes negative entropy fluctuations whose ensemble average over many replicas of the system can be quite large. Thus, on the average, the system is able to feed on negative entropy. It is also possible to show how the system can function as an engine and turn the entropy extracted from the surrounding thermal bath into work. These remarkable results are derived in the frame of a model which does not contain any ad hoc parameters and appears to possess some simple universal feature which makes it a potential candidate for supporting some still unknown ubiquitous biological process.

Our system consists of two isolated cavities A and B, filled with a common number of moles of the same gas, separated by a movable impermeable partition which is assumed to be initially at rest in a position corresponding to equal volumes and temperatures (and thus pressures) in the two cavities. Under these conditions, one would expect the wall to remain at rest, apart from microscopic random oscillations around the equilibrium position.

This is obviously the case in most ordinary situations, where the wall thermal conductivity is large enough to maintain a common temperature in the two cavities, so that any small random displacement of the wall is immediately reversed by the induced pressure change. Conversely, if thermal conductivity is negligible, it will be shown below that the wall can actually undergo sizeable random displacements, during which the pressures on the two sides remain the same. In each realization, the system evolves and the total entropy variation, from the initial state in which both cavities have the same volume $`V_0/2`$ and temperature $`T_0`$ to a generic one $`(V_A,T_A;V_B,T_B)`$, is

$`\mathrm{\Delta }S`$ $`=`$ $`\mathrm{\Delta }S_A+\mathrm{\Delta }S_B=nc_p\mathrm{ln}{\displaystyle \frac{T_A}{T_0}}+nc_p\mathrm{ln}{\displaystyle \frac{T_B}{T_0}}`$ (1)
$`=`$ $`nc_p\mathrm{ln}\left[1\left({\displaystyle \frac{\mathrm{\Delta }V}{V_0/2}}\right)^2\right],`$
where $`\mathrm{\Delta }V=V_AV_0/2`$, $`c_p`$ being the molar heat at constant pressure and $`n`$ the common number of gas moles on each side. The ensemble average $`\mathrm{\Delta }S`$ of this intrinsically negative entropy change can assume large values in specific mesoscopic situations, as it will be proved in this paper. Note that in deriving Eq.(1), we have exploited the relation $`T_A+T_B=2T_0`$, which corresponds to neglecting the kinetic energy of the wall. In fact, under the adiabatic hypothesis, it can be regarded as a body possessing a single (translational) degree of freedom. As a consequence, its energy and entropy are of the order $`kT_0/2`$ and $`k`$, respectively, and can be neglected.

In order to corroborate the above considerations, we investigate the so-called “adiabatic-piston problem”, P3 ; P4  dealing with a system which, in the language of classical thermodynamics, consists of an isolated cylinder divided into two parts by means of a frictionless adiabatic piston, each section containing the same number n of moles of the same perfect gas (see Fig.(1)).

Initially, the piston is held in the fixed position corresponding to equal volumes and pressures in the two sections. Once the piston is released, the system evolves in a labile quasi-equilibrium state, driven by the random motion of the wall induced by the elastic molecular collisions, while the gas pressure on both sides remains equal to the initial one. P4 ; P5  We describe the dynamics of the system in terms of the random variable $`X(t)`$, that is the instantaneous position of the piston, whose ensemble-average value $`X(t)`$ remains equal to $`L/2`$ for obvious symmetry reasons. The normalized mean-square value $`[X(t)L/2]^2/(L/2)^2`$ represents the ensemble-average of the quantity $`(\mathrm{\Delta }V)^2/(V_0/2)^2`$ appearing in Eq.(1) and determines the corresponding entropy decrease undergone by the system.

The system evolution from an initial configuration of macroscopic equilibrium, that is equal pressures and temperatures on both sides and piston at rest, is described by the stochastic equationP5

$`{\displaystyle \frac{d^2X}{dt^2}}`$ $`+`$ $`\sqrt{{\displaystyle \frac{16NkT_0}{\pi \mu ML}}}\left({\displaystyle \frac{1}{\sqrt{X}}}+{\displaystyle \frac{1}{\sqrt{LX}}}\right){\displaystyle \frac{dX}{dt}}`$ (2)
$`+`$ $`{\displaystyle \frac{2}{\mu }}\left[{\displaystyle \frac{XL/2}{X(LX)}}\right]\left({\displaystyle \frac{dX}{dt}}\right)^2=a(t),`$
where $`k`$ is Boltzmann’s constant, $`N`$ the common number of molecules on each side, and $`\mu =M/M_g`$, $`M`$ and $`M_g`$ being respectively the mass of the piston and the common value of the gas mass in each side. Determining the correct expression of $`a(t)`$ is a delicate task, since the standard Langevin approach does not in general carry over to nonlinear dynamical systems P6 , as the one described by Eq.(2). In order to take advantage of the Langevin method, we linearize the above equation: a) by considering small displacements around the starting position, that is $`|(XL/2)|/(L/2)1`$, and b) by approximating the square of the piston velocity $`(dX/dt)^2`$ with its mean-square velocity $`kT_0/M`$ (both hypotheses will be proved consistent a posteriori). Proceeding in this way, a straightforward application of the dissipation-fluctuation theorem yields $`a(t)a(t^{})=[8(2mkT/\pi )^{1/2}PS/M^2]\delta (tt^{})`$, where $`m`$ is the mass of the individual molecule, $`P`$ the common pressure on the two sides of the piston and $`S`$ its area. After introducing the variable $`x=XL/2`$, the linearized form of Eq.(2) reads

$`{\displaystyle \frac{d^2x}{dt^2}}`$ $`+`$ $`8\sqrt{{\displaystyle \frac{2NkT_0}{\pi \mu ML^2}}}{\displaystyle \frac{dx}{dt}}+{\displaystyle \frac{8kT_0}{\mu ML^2}}x=a(t).`$ (3)
We now observe that the above equation is formally identical to the one describing the Brownian motion of a harmonically-bound particle of mass $`M`$, that is

$$\ddot{x}+\beta \dot{x}+\omega ^2x=A(t),$$
(4)
where $`A(t)`$ is the Langevin acceleration, a problem which has been thoroughly described in the literature P7 . By comparing Eqs.(3) and (4), we can obviously apply the results of P7  to our case by identifying $`\beta `$ with $`8(2NkT_0/\pi \mu ML^2)^{1/2}`$ and $`\omega ^2`$ with $`8kT_0/\mu ML^2`$. Whenever $`\beta \omega `$ (“overdamped” case), which in our situation is equivalent to the obviously satisfied relation $`\sqrt{N}1`$, the analysis carried out in P7  naturally highlights the existence of two significant time scales $`t_{th}=1/2\beta `$ and $`t_{as}=\beta /2\omega ^2`$ ($`t_{as}t_{th}`$). They represent the thermalization time $`t_{th}`$, i.e., the time over which the mean-square velocity $`(dx/dt)^2`$ attains its equipartition value $`kT_0/M`$, and the much longer time $`t_{as}`$ over which the mean-square displacement reaches its asymptotic value $`x^2=kT_0/M\omega ^2`$. For our system $`t_{as}=(NL/w)(M/M_g)/\pi ^{1/2}`$ and $`t_{as}/t_{th}=16N/\pi `$, where $`w=(2kT_0/m)^{1/2}`$ is the most probable velocity of the gas Maxwellian distribution function. Therefore, since $`N1`$, $`t_{th}`$ is much smaller than $`t_{as}`$, a circumstance which justifies a posteriori the replacement of $`(dX/dt)^2`$ in Eq.(2) by its average value $`kT_0/M`$. The asymptotic value of the mean-square displacement of the piston from its central position reads

$$<x^2>=\frac{kT_0}{m\omega ^2}=\frac{\mu }{2}(\frac{L}{2})^2\text{,}$$
(5)
so that, in the limit $`\mu 1`$ (that is, small piston mass with respect to gas mass), the above assumptions a) and b) allowing us to linearize Eq.(2) are both satisfied.

We can now evaluate the entropy change $`\mathrm{\Delta }S`$ averaged over many realizations of our system, starting from the initial state corresponding to $`X(t=0)=L/2`$. Since in our case $`(\mathrm{\Delta }V)^2/(V_0/2)^2=x^2/(L/2)^2=\mu /21`$, Eq.(1) approximately yields

$$\frac{\mathrm{\Delta }S}{k}=\frac{c_p}{2R}N\mu ,$$
(6)
which corresponds to a large entropy decrease whenever $`\mu N1`$. The question naturally arises: what are the spatial and temporal scales over which this violation of the second law of thermodynamics can actually occur? Before answering this question, the conceptual meaning of the result implied by Eq.(6) has to be clarified. To this aim, we note that, according to standard thermodynamics, any closed system in an equilibrium state, once an internal constraint is removed, eventually reaches a new equilibrium state characterized by a larger value of the entropy. The case considered in this paper represents a remarkable exception to this statement. More precisely, our closed system is made up of the two gases and the initial internal constraint is provided by the piston held in the central position. After the piston is released, our analysis shows that no final equilibrium state is eventually reached, since the piston keeps wandering, performing random oscillations around X=L/2, with mean-square amplitude $`x^2=(L/2)^2\mu /2`$.

The system behaves as a perpetuum mobile of the second kind. If embedded in a thermal bath at temperature $`T_0`$ and pressure $`P_0`$, one can devise a process through which work can be extracted from this environment. To this end, let us assume our system to be capable of splitting into the two separate cavities $`A`$ and $`B`$ if a large fluctuation $`\mathrm{\Delta }V`$ occurs (see Fig. 2). The two cavities can then undergo a reversible adiabatic process which drives them back to the initial temperature $`T_0`$. If, at this point, the wall of the cavities becomes thermally conductive for a finite amount of time, the two cavities can be brought back to the initial volume $`V_0/2`$ through a reversible isothermal process at the bath temperature $`T_0`$. Let us now assume that the wall recovers its insulating nature and that the two cavities are reunited again: the total work extracted in this cyclic process is precisely $`W=T_0|\mathrm{\Delta }S|kT_0N\mu `$, where $`\mathrm{\Delta }S`$ is the negative entropy variation associated with $`\mathrm{\Delta }V`$.

We now return to the question of the spatial and temporal scales over which our results apply. To this end, let us consider the specific case of a gas under standard conditions of temperature and pressure, for which, expressing hereafter $`L`$ in microns, $`N310^7L^3`$. According to Eq.(6), if we refer as an example to a biatomic gas, we have $`\mathrm{\Delta }S/k510^7L^3\mu `$ and, recalling the expression of $`t_{as}`$, that is $`t_{as}=NL\mu /w\sqrt{\pi }`$, we obtain $`t_{as}510^2\mu L^4sec`$ (having assumed $`w410^8`$ microns/sec, molecular oxygen). This extremely sensitive dependence of $`t_{as}`$ on the linear dimension of the cylinder ($`t_{as}L^4`$) appears to limit the applicability of our model to values of $`L`$ up to a few microns. In fact, beyond this mesoscopic scale, $`t_{as}`$ becomes so large as to render unrealistic the adiabatic piston assumption over this time interval. As an example, for $`\mu =10^2`$ and $`L=1`$ cm, we obtain $`t_{as}510^{12}`$ sec, that is about 1000 centuries! Conversely, by taking $`L=1`$ micron, we get the reasonable value $`t_{as}510^4`$ sec and $`\mathrm{\Delta }S/k510^5`$. This corresponds a violation of the second law in the mesoscopic realm.

The above approach has allowed us to deal with the situation $`\mu 1`$. In order to have an insight into the behavior of our process in the more general case $`\mu 1`$, we assume Langevin’s approach to be approximately valid also in this moderately nonlinear regime, and use the nonlinear Eq.(2) with the same stochastic acceleration $`a(t)`$ worked out in the linear case. After introducing the dimensionless units $`\xi =X/L`$ and $`\tau =t/t_o`$, where $`t_o=4\sqrt{2}t_{th}`$, Eq.(2) reads

$$\ddot{\xi }+\left(\frac{1}{\sqrt{\xi }}+\frac{1}{\sqrt{1\xi }}\right)\dot{\xi }\frac{1}{\mu }\left(\frac{1}{\xi }\frac{1}{1\xi }\right)\dot{\xi }^2=\sigma \alpha (\tau ),$$
(7)
where the dot stands for derivative with respect to $`\tau `$, $`a(\tau )`$ is a unitary-power white noise process and $`\sigma ^2=(\pi /2\sqrt{2})(\mu /N)`$. This last equation can be numerically integrated by adopting a second-order leap-frog algorithm as the one developed in P9 . In particular, we are able to evaluate the ensemble average (over 1000 realizations of the piston) of the time evolution of the mean square-root deviation $`x^2`$ of the piston position. This analysis is reported in Fig.(1) for $`\mu =0.5`$ and $`N=310^4`$. Its inspection clearly shows that the piston undergoes random fluctuations around the central position $`X=L/2`$, which increase with time up to an asymptotic value of the order of $`L/10`$. The ensemble average of the total entropy change $`\mathrm{\Delta }S/k=N(c_p/R)\mathrm{ln}(14x^2/L^2)`$ of our isolated system can also be numerically evaluated. In Fig.(1), $`\mathrm{\Delta }S/k`$ is reported as a function of the normalized time $`\tau `$.

Further numerical analysis shows (see Fig.(3)) that, keeping $`N`$ fixed and varying $`\mu `$, the asymptotic value of $`\mathrm{\Delta }S/k`$ undergoes an approximately linear decrease with $`\mu `$, up to $`\mu =0.6`$, as predicted by Eq.(6). Beyond this value the relative entropy decrease slows down, but we are now in a range where Langevin’s model may not be valid. Conversely, in the limit $`\mu 0`$ our approach furnishes $`(X<X>)^2/(L/2)^2=\mu /20`$ and $`t_{as}\mathrm{}`$, results which are consistent with those of Chernov02 .

We wish to note that the celebrated Boltzmann’s H-theorem P10 , which can be considered as a “proof” of the second law, does not apply in our case since one of its main hypotheses, that is the complete molecular-chaos assumption (Stosszahlansatz), is not valid in our dynamics: in fact, the correlation induced by the random motion of the piston favors a common sign of $`v_x`$ and $`v_x^{}`$ in the two-particle correlation function of the gas near the piston, while the single-particle correlation functions are independent from the sign of $`v_x`$ and $`v_x^{}`$. This prevents the two-particle correlation distribution to factorize into the product of the one-particle distribution functions and, thus, the molecular-chaos assumption does not apply.

Finally, the validity of our conclusions appears to be corroborated by suitable molecular dynamic simulations of the evolution of our system. These typically involve a considerable number of point particles, which model the gas inside the cylinder, separated by a frictionless piston against which they undergo perfect elastic collisions. In particular, a microscopic model consisting of $`N=500`$ hard disks furnishes a relaxation time $`\tau `$ (corresponding to our $`t_{as}`$) in good qualitative and quantitative agreement with out results.Kestemont00

Besides, numerical investigations of a system consisting of a number of particles of the order of $`10^3`$, indicates that the difference between the temperatures on the two sides undergoes relevant oscillations, so that the system does not reach equilibrium, as predicted by our model.P11

Further recent molecular dynamic simulations describing the time evolution of the piston position around $`XL/2`$ are in fairly good agreement with our results, both qualitatively and quantitatively. White02

We wish to thank Noel Corngold for many fruitful discussions.


format: text

text:
# Polyelectrolyte Multilayering on a Charged Planar Surface

## I Introduction

PE multilayer structures are often obtained in a so-called layer-by-layer method by alternating exposure of a charged substrate to solutions of polycations (PCs) and polyanions (PAs). This widely used technique was first introduced by Decher and coworkers Decher\_1992 ; Decher\_1997 , and its simplicity and versatility trigger a large interest in the engineering community. As examples of technological applications, one can mention: biosensing Caruso\_Langmuir\_1998 , catalysis Onda\_1999 , non-linear optical devices Wu\_JACS\_1999 , nanoparticles coating Caruso\_Science\_1998 ; Caruso\_JPCB\_2001 , etc.

On the theoretical side, there exit a few analytical works about PE multilayers on charged planar surfaces based on different levels of approximation Solis\_JCP\_1999 ; Netz\_Macromol\_1999b ; Castelnovo\_2000 . Solis and Olvera de la Cruz considered the conditions under which the spontaneous formation of polyelectrolyte layered structures can be induced by a charged wall Solis\_JCP\_1999 . Based on Debye-Hückel approximations for the electrostatic interactions and including some lateral correlations by the consideration of given adsorbed PE structures, Netz and JoannyNetz\_Macromol\_1999b  found a remarkable stability of the (semi-flexible) PE multilayers supported by scaling laws. For weakly charged flexible polyelectrolytes at high ionic strength qualitative agreements between theory Castelnovo\_2000 , also based on scaling laws, and experimental observations Ladam\_Langmuir\_2000  (such as the predicted thickness and net charge of the PE multilayer) have been provided.

The important driving force for all these PE multilayering processes is of electrostatic origin. More precisely, it is based on an overcharging mechanism, where the first layer overcharges the substrate and, along the PE multilayering process, the top layer overcharges the adsorbed PE layers underneath. Nevertheless, the strong correlations existing between oppositely charged polyions, especially for highly charged PEs, provide a formidable challenge for the understanding of the PE multilayer microstructures. In this respect, numerical simulations are of great help. It is only recently, that MC simulations were carried out to study such PE structures built up on spherical charged surfaces Messina\_Langmuir\_2003 .

In this paper, we provide a detailed study of the PE multilayer structure adsorbed on a charged planar surface and discuss the basic mechanisms that are involved there by means of MC simulations. Our paper is organized as follows: Sec. II is devoted to the description of our MC simulation technique. The measured quantities are specified in Sec. III. The PE monolayering is studied in Sec. IV, and the PE bilayering in Sec. V. Then the PE multilayering process is addressed in Sec. VI. Finally, Sec. VII contains some brief concluding remarks.

## II Simulation details

The setup of the system under consideration is similar to that recently investigated with a spherical substrate Messina\_Langmuir\_2003 . Within the framework of the primitive model we consider a PE solution near a charged hard wall with an implicit solvent (water) of relative dielectric permittivity $`ϵ_r80`$. This charged substrate located at $`z=0`$ is characterized by a negative surface bare charge density $`\sigma _0e`$, where $`e`$ is the (positive) elementary charge and $`\sigma _0>0`$ is the number of charges per unit area. Electroneutrality is always ensured by the presence of explicit monovalent ($`Z_c=1`$) plate’s counterions of diameter $`a`$. PE chains ($`N_+`$ PCs and $`N_{}`$ PAs) are made up of $`N_m`$ monovalent monomers ($`Z_m=1`$) of diameter $`a`$. Hence, all microions are monovalent: $`Z=Z_c=Z_m=1`$. For the sake of simplicity, we only consider here symmetrical complexes where PC and PA chains have the same length and carry the same charge in absolute value.

All these particles making up the system are confined in a $`L\times L\times \tau `$ box. Periodic conditions are applied in the $`(x,y)`$ directions, whereas hard walls are present at $`z=0`$ (location of the charged plate) and $`z=\tau `$ (location of an uncharged wall). To avoid the appearance of image charges Torrie\_JCP\_1982 ; Messina\_image\_2002 , we assume that on both parts of the charged plate (at $`z=0`$) the dielectric constants are the same.

The total energy of interaction of the system can be written as

$`U_{tot}`$ $`=`$ $`{\displaystyle \underset{i}{}}\left[U_{hs}^{(plate)}(z_i)+U_{coul}^{(plate)}(z_i)+U_{vdw}^{(plate)}(z_i)\right]+`$
$`{\displaystyle \underset{i,i<j}{}}\left[U_{hs}(r_{ij})+U_{coul}(r_{ij})+U_{FENE}(r_{ij})+U_{LJ}(r_{ij})\right],`$
where the first (single) sum stems from the interaction between an ion $`i`$ (located at $`z=z_i`$) and the charged plate, and the second (double) sum stems from the pair interaction between ions $`i`$ and $`j`$ with $`r_{ij}=|𝐫_i𝐫_j|`$. All these contributions to $`U_{tot}`$ in Eq. (II) are described in detail below.

Excluded volume interactions are modeled via a hardcore potential defined as follows

$$U_{hs}(r_{ij})=\{\begin{array}{cc}\mathrm{},\hfill & \mathrm{for}r_{ij}<a\hfill \\ 0,\hfill & \mathrm{for}r_{ij}a\hfill \end{array}$$
(2)
for the microion-microion one, except for the monomer-monomer one note\_HS , and

$$U_{hs}^{(plate)}(z_i)=\{\begin{array}{cc}\mathrm{},\hfill & \mathrm{for}z_i<a/2\hfill \\ \mathrm{},\hfill & \mathrm{for}z_i>\tau a/2\hfill \\ 0,\hfill & \mathrm{for}a/2z_i\tau a/2\hfill \end{array}$$
(3)
for the plate-microion one.

The electrostatic energy of interaction between two ions $`i`$ and $`j`$ reads

$$\frac{U_{coul}(r_{ij})}{k_BT}=\pm \frac{l_B}{r_{ij}},$$
(4)
where +(-) applies to charges of the same (opposite) sign, and $`l_B=e^2/4\pi ϵ_0ϵ_rk_BT`$ is the Bjerrum length corresponding to the distance at which two monovalent ions interact with $`k_BT`$. The electrostatic energy of interaction between an ion $`i`$ and the (uniformly) charged plate reads

$$\frac{U_{coul}^{(plate)}(z_i)}{k_BT}=\pm 2\pi l_B\sigma _0z_i$$
(5)
where +(-) applies to positively (negatively) charged ions. An appropriate and efficient modified Lekner sum was utilized to compute the electrostatic interactions with periodicity in two directions Brodka\_MolPhys\_2002 ; note\_Lekner . To link our simulation parameters to experimental units and room temperature ($`T=298`$K) we choose $`a=4.25`$ Å leading to the Bjerrum length of water $`l_B=1.68a=7.14`$ Å. The surface charge density of the planar macroion was chosen as $`\sigma _0e0.165\mathrm{C}/\mathrm{m}^2`$.

The polyelectrolyte chain connectivity is modeled by employing a standard FENE potential in good solvent (see, e.g., Kremer\_FENE\_1993 ), which reads

$$U_{FENE}(r)=\{\begin{array}{cc}\frac{1}{2}\kappa R_0^2\mathrm{ln}\left[1\frac{r^2}{R_0^2}\right],\hfill & \text{for}r<R_0\hfill \\ & \\ \mathrm{},\hfill & \text{for}rR_0\hfill \end{array}$$
(6)
with $`\kappa =27k_BT/a^2`$ and $`R_0=1.5a`$. The excluded volume interaction between chain monomers is taken into account via a purely repulsive Lennard-Jones (LJ) potential given by

$$U_{LJ}(r)=\{\begin{array}{cc}4ϵ\left[\left(\frac{a}{r}\right)^{12}\left(\frac{a}{r}\right)^6\right]+ϵ,\hfill & \text{for}r2^{1/6}a\hfill \\ & \\ 0,\hfill & \text{for}r>2^{1/6}a\hfill \end{array}$$
(7)
where $`ϵ=k_BT`$. These parameter values lead to an equilibrium bond length $`l=0.98a`$.

An important interaction in PE multilayering is the non-electrostatic short ranged attraction, $`U_{vdw}^{(plate)}`$, between the planar macroion and the PC chain. To include this kind of interaction, we choose without loss of generality a (microscopic) van der Waals (VDW) potential of interaction between the planar macroion and a PC monomer that is given by

$$U_{vdw}^{(plate)}(z)=ϵ\chi _{vdw}\left(\frac{a}{z+a/2}\right)^6\text{for}za/2,$$
(8)
where $`\chi _{vdw}`$ is a positive dimensionless parameter describing the strength of this attraction. Thereby, at contact (i.e., $`z=a/2`$), the magnitude of the attraction is $`\chi _{vdw}ϵ=\chi _{vdw}k_BT`$ which is, in fact, the relevant characteristic of this potential. Since it is not straightforward to directly link this strength of adsorption to experimental values, we choose $`\chi _{vdw}=5`$ (also considered among other values in the case of a spherical macroion Messina\_Langmuir\_2003 ), so as to mimic good “anchoring” properties to the planar substrate.

All the simulation parameters are gathered in Table 1. The set of simulated systems can be found in Table 2. The equilibrium properties of our model system were obtained by using standard canonical MC simulations following the Metropolis scheme Metropolis\_JCP\_1953 ; Allen\_book\_1987 . Single-particle moves were considered with an acceptance ratio of $`30\%`$ for the monomers and $`50\%`$ for the counterions. Typically, about $`5\times 10^4`$ to $`\times 10^6`$ MC steps per particle were required for equilibration, and about $`5\times 10^510^6`$ subsequent MC steps were used to perform measurements. To improve the computational efficiency, we omitted the presence of PE counterions when $`N_+=N_{}`$ so that the system is still globally electroneutral. We have systematically checked for $`N_+=N_{}=20`$ (system $`C`$) that the (average) PE configurations (especially the monomer distribution) are indistinguishable, within the statistical uncertainty, from those where PE counterions are explicitly taken into account, as it should be.

## III Measured quantities

We briefly describe the different observables that are going to be measured. In order to characterize the PE adsorption, we compute the monomer density $`n_\pm (z)`$ that is normalized as follows

$$_{a/2}^{\tau a/2}n_\pm (z)L^2𝑑z=N_\pm N_m$$
(9)
where $`+()`$ applies to PCs (PAs). This quantity is of special interest to characterize the degree of ordering in the vicinity of the planar macroion surface.

The total number of accumulated monomers $`\overline{N}_\pm (r)`$ within a distance $`z`$ from the planar macroion is given by

$$\overline{N}_\pm (z)=_{a/2}^zn_\pm (z^{})L^2𝑑z^{}$$
(10)
where $`+()`$ applies to PCs (PAs). This observable will be addressed in the study of PE monolayer (Sec. IV) and PE bilayer (Sec. V).

Another relevant quantity is the global net fluid charge $`\sigma (z)`$ which is defined as follows

$$\sigma (z)=_{a/2}^z\left[\stackrel{~}{n}_+(z^{})\stackrel{~}{n}_{}(z^{})\right]𝑑z^{},$$
(11)
where $`\stackrel{~}{n}_+`$ ($`\stackrel{~}{n}_{}`$) stand for the density of all the positive (negative) microions (i.e., monomers and counterions). Thus, $`\sigma (z)`$ corresponds to the net fluid charge per unit area (omitting the bare macroion surface-charge $`\sigma _0`$) within a distance $`z`$ from the charged wall. At the uncharged wall, electroneutrality imposes $`\sigma (z=\tau a/2)=\sigma _0`$. By simple application of the Gauss’ law, $`\left[\sigma (z)\sigma _0\right]`$ is directly proportional to the mean electric field at $`z`$. Therefore $`\sigma (z)`$ can measure the screening strength of the macroion-plate charge by the neighboring solute charged species.

## IV Monolayer

In this part, we study the adsorption of PC chains for two chain lengths $`N_m=10`$ (system $`A`$) and $`N_m=20`$ (system $`B`$), and for two different couplings $`\chi _{vdw}=0`$ and $`\chi _{vdw}=5`$. Experimentally, this would correspond to the formation of the first polyelectrolyte layer. This is a decisive step to elucidate the even more complex PE multilayer structures where additionally PAs are also present.

Here, where $`N_{}=0`$ (i.e., no polyanions), global electroneutrality is ensured by the presence of explicit PC’s counterions (i.e., monovalent anions) and the macroion-plate’s counterions (i.e., monovalent cations). Also, we recall that the total number of monomers, $`N_+N_m=200`$, is identical for both systems $`A`$ and $`B`$ under consideration (see Table 2). Hence, the total monomer charge is the same for systems $`A`$ and $`B`$.

The profiles of the monomer density $`n_+(z)`$ are depicted in Fig. 1. At $`\chi _{vdw}=0`$, the density $`n_+(z)`$ near contact ($`za/2`$) is basically independent on the chain size $`N_m`$. But away from the surface, the density of monomers is slightly higher for larger $`N_m`$. This is a combined effect of (i) entropy and (ii) electrostatic correlations. These underlying mechanisms at $`\chi _{vdw}=0`$ can be explained with simple ideas as follows:

* At fixed number of total monomers, entropic effects are larger the shorter the chains, and in the limiting case of $`N_m=1`$ (i.e., the electrolyte limit) entropy effects are maximal leading to the highest monomer “release”. It is to say that the chain connectivity lowers the entropy of the system.
* In parallel, electrostatic correlations Shklowskii\_PRE\_1999b ; Messina\_PRE\_2001 ; Levin\_RepProgPhys\_2002  are also higher the higher the valence of the adsorbed particles. In our case $`N_m`$ plays the role of the polyion valence.

The density $`n_+(z)`$ near contact increases considerably with $`\chi _{vdw}`$ (here about one order of magnitude) as expected. It turns out that with $`\chi _{vdw}=5`$, the $`n_+(z)`$-profiles are basically identical for $`N_m=10`$ and $`N_m=20`$. This is due to the sufficiently strong non-electrostatic attractive force that can overcompensate the antagonistic entropic effects that were more efficient at $`\chi _{vdw}=0`$.

The fraction $`\overline{N}_+(z)/(N_+N_m)`$ of adsorbed monomers can be found in Fig. 2. At a $`z`$-distance of $`1.5a`$ from the planar macroion surface (corresponding to a width of two monomers), about 90% of the monomers are adsorbed for $`\chi _{vdw}=5`$ against only $`50\%`$ for $`\chi _{vdw}=0`$. Again, at given $`\chi _{vdw}`$, $`\overline{N}_+(z)/(N_+N_m)`$ is larger for longer chains due to the same coupled effects of entropy and electrostatic correlations explained above.

The (global) net fluid charge $`\sigma (z)`$ is reported in Fig. 3. In all cases we observe a macroion-surface charge reversal (i.e., $`\sigma (z)/\sigma _0>1`$). The position $`z=z^{}`$ at which $`\sigma (z^{})`$ gets its maximal value decreases with $`\chi _{vdw}`$, due to the $`\chi _{vdw}`$-enhanced adsorption of the PCs. Concomitantly, this overcharging increases with $`\chi _{vdw}`$, since the (extra) gain in energy by macroion-monomer VDW interactions can better overcome (the higher $`\chi _{vdw}`$) the cost of the self-energy stemming from the adsorbed excess charge Messina\_Langmuir\_2003 . More quantitatively, we have $`\sigma (z^{})/\sigma _01.7`$ at $`\chi _{vdw}=5`$ against only $`\sigma (z^{})/\sigma _01.25`$ at $`\chi _{vdw}=0`$. Note that the maximal value of charge reversal of $`(20090)/90=122\%`$ (i.e., $`\sigma (z^{})/\sigma _0=2.22`$) allowed by the total charge of PCs can not be reached due to a slight accumulation of microanions. In agreement with the profiles of $`n_+(z)`$ and $`\overline{N}_+(z)`$ (see Fig. 1 and Fig. 2), at given $`\chi _{vdw}`$, the overcharging gets higher the higher the chain length. Those (locally) overcharged states should be the driving force for the building of subsequent PE bilayers when PA chains are added.

Typical equilibrium configurations can be found in Fig. 4. The qualitative difference between $`\chi _{vdw}=0`$ \[Figs. 4(a) and (b)\] and $`\chi _{vdw}=5`$ \[Figs. 4(c) and (d)\] is rather spectacular. Without additional VDW attraction ($`\chi _{vdw}=0`$) the adsorption is much weaker than at $`\chi _{vdw}=5`$, where in the latter situation the $`z`$-fluctuation is very weak within the adsorbed layer. Basically the first layer is glued at $`\chi _{vdw}=5`$, and the excess PC chains float in the solution. It is typically this type of configurations for the first layer that is wanted in experimental situations.

The next section (Sec. V) that concerns bilayering will show that the (enhanced) stability of this first layer is decisive for the onset of multilayers.

## V Bilayer

We now consider the case where additionally PA chains are present (systems $`C`$ and $`D`$), so that we have a neutral polyelectrolyte complex (i.e., $`N_+N_m=N_{}N_m=200`$). Global electroneutrality is ensured by the counterions of the planar macroion as usual. For such parameters, the final equilibrium structure consists essentially of bilayers with sometimes the onset of a weakly stable third layer. Experimentally this would correspond to the process of the second polyelectrolyte layer formation (with system $`A`$ or $`B`$ as the initial state). We stress the fact that this process is fully reversible for the parameters investigated in our present study. In particular, we checked that the same final equilibrium configuration is obtained either by (i) starting from system $`A`$ or $`B`$ and then adding PAs or (ii) starting directly with the mixture of oppositely charged polyelectrolytes.

The profiles of the monomer density $`n_\pm (z)`$ at $`\chi _{vdw}=0`$ and $`\chi _{vdw}=5`$ are depicted in Fig. 5(a) and Fig. 5(b), respectively. The corresponding microstructures are sketched in Fig. 6. The density of PC monomers $`n_+(r)`$ near contact increases considerably with $`\chi _{vdw}`$ as expected. Interestingly, at $`\chi _{vdw}=0`$, a comparison with systems $`A`$ and $`B`$ (see Fig. 1) indicates that the adsorption of PC monomers is weaker when additional PAs are present. This effect was already observed with spherical substrates Messina\_Langmuir\_2003 , and the same mechanisms apply here to planar surfaces. More explicitly, the PC chain tends to build up a globular state (reminiscent of the classical bulk PE collapse Hayashi\_JCP\_2002 ) by getting complexed to the PA chain, as well illustrated in Fig. 6(a) and Fig. 6(b). Thereby, the mean monomer coordination number (or the mean number of monomer neighbors) gets higher which is both (i) entropically and (ii) energetically (at least from the PE complex viewpoint) favorable. This PC desorption is only appreciable at sufficiently low $`\chi _{vdw}`$ where the energy loss stemming from the PC desorption is well balanced (or even overcompensated depending generally on the parameters) by the energy gained in building a PC-PA globular structure. This “auto-globalization” is also enhanced by increasing $`N_m`$ as it should be \[compare Fig. 6(a) and (b)\]. Note also that there is a small second peak in $`n_+(z)`$ at $`z3.8a`$ \[see Fig. 5(a)\], which is rather the signature of a strong PC-PA globalization than a third PE layer. Besides, the peak in the PA density $`n_{}(z)`$ located at $`z=z^{}2.3a`$ \[see Fig. 5(a)\], which is relatively far from that of a compact bilayer where $`z^{}=1.5a`$, indicates the diffuse character of the bilayer at $`\chi _{vdw}=0`$.

At $`\chi _{vdw}=5`$, the $`n_+(z)`$-profiles are basically identical for $`N_m=10`$ and $`N_m=20`$. In contrast to $`\chi _{vdw}=0`$, $`n_+(z)`$ near contact is somewhat larger at $`\chi _{vdw}=5`$ and it is going to be explained later in the discussion of $`\overline{N}_+(z)`$. As far as the PA density $`n_{}(z)`$ is concerned, we see that the peak is roughly 2-3 times higher (depending on $`N_m`$) with $`\chi _{vdw}=5`$ than with $`\chi _{vdw}=0`$. Also, its position ($`z^{}1.5a`$) corresponds to that of a compact bilayer. A visual inspection of Fig. 6(c) and (d) confirms this feature. This again shows how important is the role of extra non-electrostatic attractive force for the stability of bilayers.

An intermediate conclusion can be drawn from the above findings and especially from the microstructures depicted in Fig. 6:

* True bilayering (i.e., flat and dense layers) can only occur at non-zero $`\chi _{vdw}`$, as already reported for spherical charged substrates Messina\_Langmuir\_2003  with large curvature.

An interesting common characteristic of the microstructures at $`\chi _{vdw}=0`$ and $`\chi _{vdw}=5`$ is the formation of small islands (along the substrate) made up of more or less flat (depending on $`\chi _{vdw}`$) PC-PA complexes, easily identifiable at $`N_m=20`$ \[see Fig. 6(b) and Fig. 6(d)\].

The fraction $`\overline{N}_+(z)/(N_+N_m)`$ of adsorbed monomers at $`\chi _{vdw}=0`$ and $`\chi _{vdw}=5`$ can be found in Fig. 7(a) and (b), respectively. A close look at Fig. 7(a) reveals a smaller PC monomer accumulation (at $`\chi _{vdw}=0`$) up to $`z3a`$ (independently of $`N_m`$) than in the case where PA chains were absent (compare with Fig. 2). This is fully consistent with the formation of PC-PA globules (relevant at $`\chi _{vdw}=0`$) leading to the effective PC desorption already discussed above. In parallel, this PC-PA globalization tends to cancel the effect of chain length $`N_m`$ on $`\overline{N}_+(z)`$. On the other hand, at $`\chi _{vdw}=5`$, the situation is qualitatively different where the presence of PAs now induces an increase of $`\overline{N}_+(z)`$ \[compare Fig. 7(b) and Fig. 2\]. This phenomenon can be explained by electrostatic correlation effects. Indeed, at $`\chi _{vdw}=5`$, the highly stable PC layer attracts more PA monomers than at $`\chi _{vdw}=0`$, and thereby, ”super” dipoles made of PC-PA monomer pairs build up, that are perceptible in Fig. 6(c) and Fig. 6(d). This leads to a strong attractive correlation interaction between the plate and those dipoles. In other terms the effect of finite $`\chi _{vdw}`$ is to (strongly) polarize the adsorbed charged chains. Note also that a $`\chi _{vdw}=5`$ the PC-PA globalization is much less favorable than at $`\chi _{vdw}=0`$ due to the higher cost of PC desorption energy in the former case. As a net effect there can be more adsorbed PC monomers compared to $`\chi _{vdw}=0`$. In that case of $`\chi _{vdw}=5`$, it is precisely this mechanism that tends to cancel the effect of $`N_m`$ on $`\overline{N}_+(z)`$. As far as the PA monomer fraction $`\overline{N}_{}(z)`$ is concerned, Fig. 7 shows that the adsorption of monomers is much weaker and more diffuse at $`\chi _{vdw}=0`$ than at $`\chi _{vdw}=5`$, as expected from Figs. 5 and 6.

The net fluid charge $`\sigma (z)`$ is reported in Fig. 8. In all cases, the planar macroion gets overcharged and undercharged as one gets away from its surface. That is we have to deal with charge oscillations. Our results clearly show that the amplitude of those oscillations is systematically larger at high $`N_m`$, as also observed without PAs (see Fig. 3). This is consistent with the idea that lateral electrostatic correlations are enhanced by increasing the valence of the polyions (here $`N_m`$). Nevertheless, as soon as oppositely charged polyions can interact, there is a subtle interplay between clustering and the lateral correlations of polyions that governs the degree of overcharging near the planar macroion. At $`\chi _{vdw}=5`$, we observe a significantly higher overcharging than without PAs (compare with Fig. 3). This is in agreement with the profiles of $`N_+(z)`$ discussed previously. However, the positions of the first peak ($`z^{}a`$ for $`\chi _{vdw}=5`$ and $`z^{}1.8a`$ for $`\chi _{vdw}=0`$) in $`\sigma (z)`$ remain nearly unchanged by the presence of PAs (compare Fig. 8 with Fig. 3).

## VI Multilayer

Presently, we consider the case where there are enough polyelectrolytes ($`N_+N_m=N_{}N_m=400`$) in the system to produce multilayers (systems $`E`$ and $`F`$). Hence, compared to systems $`C`$ and $`D`$, we have now doubled the polyelectrolyte concentration. Global electroneutrality is ensured by the counterions of the planar macroion as usual.

The density profiles of $`n_\pm (r)`$ for $`\chi _{vdw}=0`$ and $`\chi _{vdw}=5`$ are depicted in Fig. 9(a) and Fig. 9(b), respectively. The corresponding microstructures are sketched in Fig. 10. In general, the densities of PC and PA monomers are systematically larger than those found for systems $`C`$ and $`D`$ corresponding to a lower PE concentration. This effect is due to the fact that, at higher concentration of oppositely charged chains, the number of dipoles (i.e., PC-PA monomer pairs) are also larger, and from this it results larger plate-dipole correlations.

Even at $`\chi _{vdw}=0`$ with $`N_m=10`$, we can observe a non-negligible second peak in $`n_+(z)`$ (located at $`z3.8a`$) which is the signature of a third layer. This finding contrasts with what was observed at spherical substrates Messina\_Langmuir\_2003  (also with $`\chi _{vdw}=0`$, $`N_m=10`$, and with a similar macroion surface charge density), where not even a stable bilayer could build up. This radically different behavior can be accounted by geometrical arguments. Indeed, the potential of electrostatic interaction scales like $`1/r`$ in spherical geometry against $`z`$ in planar one. Hence, at sufficiently high curvature (as it was the case in Ref. Messina\_Langmuir\_2003  where $`N_ma/r_0>1`$ note\_curvature  with $`r_0`$ being the radius of the spherical macroion), qualitative differences from the planar case are then expected. However, the corresponding microstructure \[see Fig. 10(a)\] suggests a relatively large formation of PC-PA globules leading to a diffuse and porous multilayer. Always at $`\chi _{vdw}=0`$, but with longer chains ($`N_m=20`$), Fig. 9(a) shows that the degree of layering is higher as expected. This feature is well illustrated by Fig. 10(b), where the PA monomers are visibly more attracted to the planar macroion surface.

At $`\chi _{vdw}=5`$, the adsorption of monomers is drastically increased due to the enhanced stability of the first PC layer that, in turn, induces a larger adsorption of the subsequent PAs and PCs. Compared to $`\chi _{vdw}=0`$, all the peaks in $`n_\pm (z)`$ are shifted to smaller $`z`$, indicating a higher compaction. These higher ordering and compaction at $`\chi _{vdw}=5`$ can be visually checked in Fig. 10(c) and (d).

The net fluid charge $`\sigma (z)`$ is reported in Fig. 11. As expected charge oscillations are detected. However in all cases, the corresponding amplitude is decaying (as also found in Ref. Messina\_Langmuir\_2003  for spherical geometry). Compared to the bilayer situation (see Fig. 8), one remarks that the charge oscillations are now larger due to the enhanced “plate-dipole” correlations occurring at higher chain concentrations (as discussed above). On the other hand the positions of the extrema in the charge oscillations remain quasi unchanged.

## VII Concluding remarks

We first would like to briefly discuss our findings with some experimental examples. Our results concerning the first layer (i.e., single PC layer) show that an additional non-electrostatic force is needed to enhance its stability. Experimentally, this is achieved by choosing PCs with good ”anchoring” properties to a given substrate. In our model this was done by taking $`\chi _{vdw}=5`$. This being said, the case $`\chi _{vdw}=0`$ is from a fundamental point of view interesting, since it corresponds to a purely electrostatic regime.

Recently, Menchaca et al. found, by means of ”liquid-cell AFM”, that PE-complex grains appear at the first PE-layers Perez\_CollSurf\_2003 . This kind of structure (that we referred to as small islands - see Fig. 6) are confirmed by our simulations. Concomitantly, a significant roughness of the deposited bilayer was also detected in this experiment, which is directly linked to the presence of those grains. This microstructure seems also to be (indirectly) reported in other experiments using ellipsometry Harris\_Langmuir\_2000 , where it is found that the structure of the two first bilayers are more porous than that of later bilayers. This is also in qualitative agreement with our microstructures depicted in Fig. 6 and Fig. 10. However, more simulation data are needed to understand the PE structure beyond two bilayers.

The degree of charge inversion of the substrate can be indirectly obtained by measuring the $`\zeta `$-potential via electrophoresis Ladam\_Langmuir\_2000 . In their experiment, Ladam et al. Ladam\_Langmuir\_2000  observed that, after a few deposited PE layers note\_thickness , the $`\zeta `$-potential profile is symmetrically oscillating. This reveals a “stationary” regime where, successively, polycations and polyanions are adsorbed with the same strength. Unfortunately, it is not possible for us to investigate numerically this regime due to the highly prohibitive computation time required there. However, the charge oscillations observed in our $`\sigma (z)`$-profiles indicate that by increasing the amount of layers, one first increases the amplitude of these oscillations. This confirms at least the general experimental evidence of the non-stationary regime a the early stage of PE multilayering.

We also would like to mention the possible effect of image charges stemming from the dielectric discontinuity between the substrate (typically $`ϵ_r25`$) and the solvent (here $`ϵ_r80`$), as is the case under experimental conditions. It is expected that image forces become especially relevant for PE monolayering (i.e., when PCs solely are present) messina\_PE\_image . Indeed for multivalent ions, a strong self-image repulsion occurs and leads to a shifted density-profile $`n_+(z)`$ with a maximum located somewhat further than the contact region Messina\_image\_2002 ; Torrie\_JCP\_1982 . However for PE multilayering, due to the presence of oppositely charged PEs the effect of image forces is considerably reduced (especially sufficiently away from the wall) due to the (self-)screening of the image charges.

In summary, we have investigated by means of extensive MC simulations the equilibrium buildup of the few first layers adsorbed on a charged planar substrate. Two parameters were considered: (i) the chain length $`N_m`$ and (ii) the extra non-electrostatic short-range attraction (characterized here by $`\chi _{vdw}`$) between the planar macroion surface and the polycation chains.

For the bilayering, it was demonstrated that, within the electrostatic regime (i.e., $`\chi _{vdw}=0`$), significant PC-PA globules build up leading to a very “porous” and diffuse bilayer structure. The PC-PA globalization is enhanced with $`N_m`$. At sufficiently large $`\chi _{vdw}`$ (here $`\chi _{vdw}=5`$), the bilayer is much less diffuse and the oppositely charged chains are more polarized, leading to a high stability of the structure.

The same qualitatively applies to the case of the two-bilayer (i.e., four PE layers) adsorption. Within this regime of layering as investigated here (up to four layers), we also found a non-linear regime, where for instance the separation of the peaks in the monomer densities are not identical. This is in qualitative agreement with the finding of Ladam et al. Ladam\_Langmuir\_2000  where they reported a non-linear regime in the so-called “region I” corresponding to the PE multilayer-region close to the buffer note\_thickness . The effect of $`N_m`$ is to globally enhance the stability of the multilayer structure due to the higher electrostatic correlations and also due to entropic effects.

A future study should take into account the rigidity of the chain, which can drastically change the multilayer structure depending the stiffness. The formation of PE multilayers on cylindrical substrates seems also to be a promising research area, and to our knowledge it has never been investigated so far messina\_rod .

###### Acknowledgements.

The author thanks F. Caruso, H. Löwen, S. K. Mayya and E. Pérez for helpful and stimulating discussions.

TABLES

FIGURE CAPTIONS

1. Profiles of PC monomer-density $`n_+(z)`$ at different $`\chi _{vdw}`$ couplings (systems $`A`$ and $`B`$). The inset corresponds to $`\chi _{vdw}=5`$ where the two curves ($`N_m=10`$ and 20) are nearly indistinguishable.
2. Fraction $`\overline{N}_+(z)/(N_+N_m)`$ of adsorbed PC monomers at different $`\chi _{vdw}`$ couplings (systems $`A`$ and $`B`$).
3. Net fluid charge $`\sigma (z)`$ at different $`\chi _{vdw}`$ couplings (systems $`A`$ and $`B`$).
4. Typical equilibrium configurations for PC chains adsorbed onto an oppositely charged planar macroion (systems $`A`$ and $`B`$). The little counterions are omitted for clarity. (a) $`\chi _{vdw}=0`$, $`N_m=10`$ (b) $`\chi _{vdw}=0`$, $`N_m=20`$ (c) $`\chi _{vdw}=5`$, $`N_m=10`$ (d) $`\chi _{vdw}=5`$, $`N_m=20`$.
5. Profiles of monomer density $`n_\pm (z)`$ for oppositely charged polyelectrolytes (systems $`C`$ and $`D`$). (a) $`\chi _{vdw}=0`$. (b) $`\chi _{vdw}=5`$.
6. Typical equilibrium configurations for the adsorption of oppositely charged PE chains (systems $`C`$ and $`D`$) onto a planar macroion. The polycations are in white and the polyanions in red. The little ions are omitted for clarity. (a) $`\chi _{vdw}=0`$, $`N_m=10`$ (b) $`\chi _{vdw}=0`$, $`N_m=20`$ (c) $`\chi _{vdw}=5`$, $`N_m=10`$ (d) $`\chi _{vdw}=5`$, $`N_m=20`$.
7. Fraction of adsorbed monomers $`\overline{N}_\pm (z)`$ for oppositely charged polyelectrolytes (systems $`C`$ and $`D`$). (a) $`\chi _{vdw}=0`$. (b) $`\chi _{vdw}=5`$.
8. Net fluid charge $`\sigma (z)`$ at different $`\chi _{vdw}`$ couplings (systems $`C`$ and $`D`$).
9. Profiles of monomer density $`n_\pm (z)`$ for oppositely charged polyelectrolytes (systems $`E`$ and $`F`$). (a) $`\chi _{vdw}=0`$. (b) $`\chi _{vdw}=5`$.
10. Typical equilibrium configurations for the adsorption of oppositely charged PE chains (systems $`E`$ and $`F`$) onto a planar macroion. The polycations are in white and the polyanions in red. The little ions are omitted for clarity. (a) $`\chi _{vdw}=0`$, $`N_m=10`$ (b) $`\chi _{vdw}=0`$, $`N_m=20`$ (c) $`\chi _{vdw}=5`$, $`N_m=10`$ (d) $`\chi _{vdw}=5`$, $`N_m=20`$.
11. Net fluid charge $`\sigma (z)`$ at different $`\chi _{vdw}`$ couplings (systems $`E`$ and $`F`$).


text:
# Power-law expansion in 𝑘-essence cosmology

## 1 introduction

To address the outstanding theoretical challenges of modern cosmology, especially the so-called coincidence problem, which questions as why is it exactly now, the universe driven by some sort of dark energy is accelerating, several authors have introduced and studied the so-called $`k`$-essence models .

Originally, the $`k`$-essence, or $`k`$-inflation, was introduced in  in order to bridge phenomenologically the string theories with inflation (see Ref.  for a recent review). The main ingredient of the $`k`$-essence is a scalar field, with non standard higher order kinetic terms. Interestingly enough, and contrary to what one could have expected, these non-standard terms do not necessarily lead to acausal propagation of the $`k`$-field . Studying inflationary patterns with the $`k`$-fields the authors of  were able to show that $`k`$-field may drive an accelerated expansion of the universe starting from a generic initial conditions without an assistance of the usual potential terms.

In a different development , the $`k`$-essence was proposed as a dynamical solution to the coincidence problem. The basic idea of  is that $`k`$-essence could play a role of a dynamical attractor at the onset of matter domination period introducing cosmic acceleration at present time. Further study of $`k`$-essence was performed recently in . It was argued that in certain dynamical regimes the $`k`$-essence is equivalent to quintessence and it may prove difficult to distinguish between the two fields. In this paper we make a step further and show that the dynamically different $`k`$-theories can produce kinematically equivalent cosmological models.

The construction of cosmological models with tracker-like, or the attractor behaviour , where the $`k`$-essence either mimics the equation of state of the matter-radiation component, or drives towards acceleration, is relied heavily on the existence of $`k`$-essence solutions which, re-written in terms of energy density and pressure, represent, hydrodynamically, fluids with a constant baryotropic index (BI). These, in turn, give rise to a power-law behaviour of the scale factor when the underlying geometry is that of a spatially flat isotropic universe.

In the $`k`$-essence models studied earlier  one would usually consider solutions where, during the $`k`$-field driven expansion with the constant BI, two things happen: i) The scalar field evolves linearly with time and ii) The $`k`$-field potential is an inverse square in terms of the $`k`$-field. The property ii) follows directly from i). The assumption i), triggered probably by the simplicity of finding solutions in the case of linear $`k`$-field, although permits to consider different $`k`$-theories, is too restrictive with the form of the k-potential and the evolution of the field itself.

The main purpose of this paper is to show that in the case of the $`k`$-essence, one may find solutions with arbitrary potentials and non-linear scalar fields, but still have a constant BI. For the solutions we find, one can have a fixed evolution of the geometry, yet incredibly rich repertoire of scalar field behaviour and its $`k`$-potential. This kind of a degeneracy is quite problematic for the model building. Nevertheless, there seems to be a way to reduce this degeneration. For solutions with constant BI we find that the linear $`k`$-field model is isomorphic to a divergent sound speed model. In fact, the former presents a superior type of degeneracy as compared to the solutions obtained from the non-linear $`k`$-field model, therefore fully justifying our quest for a different type of $`k`$-field solutions. The isomorphism between the linear $`k`$-field model and the divergent sound speed model looks especially interesting in the light of the results of the recent publication  where the behaviour of the solutions near divergent sound speed was thoroughly investigated.

## 2 The general framework

We start with a general Lagrangian

$$=V(\varphi )F(x),x=g^{\mu \nu }\varphi _\mu \varphi _\nu ,$$
(1)
where $`\varphi `$ is the scalar field and $`\varphi _\mu =\varphi /x^\mu `$, and do not impose any conditions neither on $`V`$, nor $`F`$ at this stage. One may easily figure out the energy-momentum tensor for (1):

$$T_{\mu \nu }=V(\varphi )\left[2F_x\varphi _\mu \varphi _\nu g_{\mu \nu }F\right],F_x=\frac{dF}{dx}.$$
(2)
Identifying (2) with the energy-momentum tensor of a perfect fluid we have

$$\rho _\varphi =V(\varphi )[F2xF_x],p_\varphi ==V(\varphi )F.$$
(3)
As usual in this setting we assume a spatially flat homogeneous and isotropic spacetime with line element

$$ds^2=dt^2+a^2(t)\left[dx^2+dy^2+dz^2\right],$$
(4)
where $`a(t)`$ is the scale factor and the expansion rate is defined as $`H=\dot{a}/a`$. The Einstein field equations then reduce to

$$3H^2=\rho _\varphi ,2\dot{H}=\rho _\varphi +p_\varphi ,$$
(5)
and the conservation equation reads

$$\dot{\rho }_\varphi +3H(\rho _\varphi +p_\varphi )=0.$$
(6)
The field equation for the $`\varphi `$ field may be either obtained by substituting expressions (3) into the conservation equation (6), or varying directly the Lagrangian (1). Doing so, we get:

$$[F_x+2xF_{xx}]\ddot{\varphi }+3HF_x\dot{\varphi }+\frac{V^{}}{2V}[F2xF_x]=0,$$
(7)
where $`V^{}=dV/V\varphi `$. On the other hand, assuming a “formal” equation of state of the form $`p_\varphi =(\gamma 1)\rho _\varphi `$ for the k-essence and using Eqs. (3), (5) we obtain the BI $`\gamma `$

$$\gamma =\frac{2\dot{H}}{3H^2}=\frac{2xF_x}{F2xF_x}.$$
(8)
We now assume that the BI is a constant. This kinematically leads to a power-law scale factor $`a=a_0t^{2/3\gamma }`$.

The first question we ask is, how stable are the solutions with the constant BI $`\gamma =\gamma _0`$?. To answer this question, we allow $`\gamma `$ to vary with time. Differentiating the equation of state and using the conservation equation we find

$$\dot{\gamma }=3H\gamma (\gamma 1)+\frac{\dot{p}_\varphi }{\rho _\varphi },$$
(9)
which together with (3) and (8) lead to

$$\dot{\gamma }+\left[3H\gamma +\frac{\dot{V}}{V}+\frac{\dot{F}}{F}\right](1\gamma )=0.$$
(10)
We further check as to whether $`\gamma =\gamma _0`$ are solutions to this equation at all. Obviously, there are two different ways for this to happen: either $`\gamma _0=1`$, or generically, the following stationary condition holds

$$\frac{\dot{V}}{V}+\frac{\dot{F}}{F}=3H\gamma _0,$$
(11)
When the stationarity condition (11) holds, the potential $`V`$ and the function F are related by:

$$VF=\frac{4(1\gamma _0)a_0^{3\gamma _0}}{3\gamma _0^2a^{3\gamma _0}}.$$
(12)
Here, we have integrated (11) and inserted the solution into the Einstein equation (5) to fix the integration constant. For a positive potential $`V`$, the constrain (12) gives rise to two different theories depending on whether $`\gamma _0<1`$ or $`\gamma _0>1`$. In the case $`\gamma _0<1`$ we take the function $`F`$ to be positive, whereas in the case $`\gamma _0>1`$ we take it negative. We denote these as $`F^+`$ and $`F^{}`$ respectively.

We now assume that the stationarity condition (11) holds. So, the equation (10) reads

$$\dot{\gamma }+3H(\gamma \gamma _0)(1\gamma )=0.$$
(13)
Integrating, we find:

$$\gamma =\frac{\gamma _0+ca^{3(1\gamma _0)}}{1+ca^{3(1\gamma _0)}}.$$
(14)
Here $`c`$ is an integration constant. For the expanding universe and $`\gamma _0<1`$ we see that the solutions of (13) have the asymptotic limit $`\gamma _0`$. Therefore, the solutions with constant BI $`\gamma _0`$ are attractors in the case $`\gamma _0<1`$. This attractor behaviour holds even for superaccelerated universes  with $`\gamma _0<0`$.

The limit $`\gamma _01`$, should be considered apart, and the solution of the equation (13) is

$$\gamma =1\frac{1}{c+\mathrm{ln}a^3},$$
(15)
where $`c`$ is an integration constant. Hence, for an expanding universe the solution with $`\gamma _0=1`$ is stable as well . The $`\gamma _0=1`$ solutions separate stable from unstable regions in the phase space (for a positive expansion rate) as can be easily seen from the equation (13), and since $`\gamma _0=1`$ corresponds to dust, we conclude that the dust-like solutions define the border line between stable and unstable behavior. It is probably worthwhile to mention that the above stability analysis is simple and direct as compared to the study performed directly in the field variables using the solution $`\varphi t`$ as an input.

## 3 Power-law solutions

As from now we stick to the solutions with the constant BI $`\gamma `$. It follows then that the Einstein and the field equations (5), (7) have two different classes of solutions:

1) The solutions with constant $`x=x_0=\dot{\varphi }^2`$.

In this case for $`\gamma =const0`$ we have $`a=a_0t^{2/3\gamma }`$, the first term of the l.h.s. of the Eq. (7) vanishes, and the consistent solution of Eqs. (5), (7) becomes $`\varphi =\pm \sqrt{x_0}t`$ and

$$V=\frac{4x_0}{3\gamma ^2[F2xF_x]}\frac{1}{\varphi ^2},$$
(16)
with an arbitrary $`F`$ evaluated at $`x=x_0`$. We will not discuss these solutions further, since these were thoroughly investigated and exploited in model building in . The particular case with $`x_0=0`$ ($`\varphi =\varphi _0`$) must be solved apart and gives a de Sitter solution $`a=a_0\text{e}^{\sqrt{VF/3}t}`$ for arbitrary $`F`$ evaluated at $`x=0`$ and constant potential $`V`$.

2) Solutions with $`xconst`$.

In this case the conservation equation (6) can be readily integrated to find the first integral of the field equation (7)

$$VF_\gamma =\frac{\rho _0}{a^{3\gamma }}.$$
(17)
Comparing this expression with the constrain equation (5) we are lead to the relation (12) between the potential $`V`$ and the function $`F`$. Hence, the integration constants $`a_0`$ and $`\rho _0`$ are left fixed to $`\rho _0=4(1\gamma )a_0^{3\gamma }/3\gamma ^2`$.

We now look at (8) as a differential equation for $`F(x)`$. Its immediate general solution is

$$F_\gamma (x)=c(x)^{\frac{\gamma }{2(\gamma 1)}}.$$
(18)
Without any loss of generality, one may fix the integration constant $`c=\pm 1`$. The two corresponding families of solutions are then $`F_\gamma ^+`$ and $`F_\gamma ^{}`$ respectively. Inserting the last equation into (3) we get two possibilities

$$\rho _\varphi ^+=\frac{VF_\gamma ^+}{1\gamma },p_\varphi ^+=VF_\gamma ^+,\gamma <1,$$
(19)
and

$$\rho _\varphi ^{}=\frac{VF_\gamma ^{}}{1\gamma },p_\varphi ^{}=VF_\gamma ^{},\gamma >1,$$
(20)
where we have assumed that both the $`k`$-potential and the energy density are positive definite.

Inserting (18) into (12) one gets a relation of the form $`t^2V\dot{\varphi }^{\gamma /1\gamma }`$. Finally, the general relations connecting the field $`\varphi `$ and the potential $`V`$ follows:

$$t^{\frac{2\gamma }{\gamma }}=\frac{2\gamma }{\gamma }\left[\pm \frac{3\gamma ^2}{4(1\gamma )}\right]V^{\frac{\gamma 1}{\gamma }}𝑑\varphi ,\gamma 2,$$
(21)
$$\mathrm{ln}t=\sqrt{3}\sqrt{V}𝑑\varphi ,\gamma =2.$$
(22)
The $`(+)`$ branch in equation (21) corresponds to $`\gamma <1`$ while the $`()`$ branch to $`\gamma >1`$. For linear $`\varphi `$ the integral (21) is not defined and this situation corresponds to the first class (i) of the solutions. The relations (21) and (22) should be read as follows: given $`V(\varphi )`$, one may integrate and obtain $`t=t(\varphi )`$, invert and find $`\varphi =\varphi (t)`$. Then $`F(x)`$ is still given by (18). Note, that for a fixed $`\gamma `$ (fixed power of the scale factor) one has different potentials and different field evolutions, and consequently different $`k`$-theory. It looks as the $`k`$-essence theories have a considerable amount of freedom in choosing the theory, the potential and the scalar field behaviour, all describing the same kinematics of the universe. This sounds somewhat “fantastic” for these are not just simple field redefinitions, and all the theories with the different $`\varphi `$ and $`V`$ are dynamically different.

We now show how the power law solutions with the linear scalar field and the inverse square potential are related to the family of solutions with the divergent velocity of sound. We do so by constructing a one-to-one mapping between these solutions. We suggest that this might be the reason as to why the solutions with the linear scalar field run into trouble as discussed in a recent paper by Malquarty etal .

To do so we introduce what we call a “divergent” $`k`$\- essence Lagrangian with the kinetic energy proportional to the velocity. For such a theory one may take the function $`F^{\mathrm{}}`$ as

$$F^{\mathrm{}}=c+b\sqrt{x}.$$
(23)
It is easy to see that the above function leads immediately to a divergent sound velocity $`C_s^2=F_x^{\mathrm{}}/(F_x^{\mathrm{}}+2xF_{xx}^{\mathrm{}})`$ and to an inverse square potential by using the $`k`$-field equation (7). This does not constitute a major problem in itself, for one could have just avoided using this sort of a model. It follows, however, that the solutions of the models with linear scalar field and the inverse square potential discussed in  are isomorphic to those obtained in the divergent models.

To see the relation between the models we work with the power law solutions. Consider a typical model cosmology given by:

$$a=t^n,V=\frac{\beta }{\varphi ^2},\varphi =\varphi _0t,$$
(24)
obtained by evaluating $`F`$ and $`F_x`$ at $`x=x_0=\varphi _0^2`$. We further use $`f=F(\varphi _0^2)`$ and $`f^{}=F_x(\varphi _0^2)`$. Substituting these constants into the Friedmann and $`k`$-field equations (5), (7) we find that the index $`n`$ and the slope of the potential $`\beta `$ are given by

$$n=\frac{1}{3}\frac{f+2\varphi _0^2f^{}}{\varphi _0^2f^{}}=\frac{2}{3\gamma },\beta =\frac{n}{f^{}},$$
(25)
On the other hand, if in the divergent model we choose the constants $`c=3n^2\varphi _0^2/\beta `$ and $`b=2n\varphi _0/\beta `$ we obtain the same solution. Therefore all the power-law solutions obtained from the model with the linear scalar field and the inverse square potential map into the solutions of the divergent model with the same potential.

Moreover, the following reasoning underlines the highly degenerate character of the linear $`k`$-field solutions. Consider the series expansion of the function $`F(x)`$ around $`x=x_0`$. The background cosmology is completely determined by just the first two coefficients in the expansion $`(f,f^{})`$ and the value of $`\varphi _0`$ as seen from (25). Put in different words, the model is insensitive to keeping the first two coefficients in the expansion of the function $`F`$ and the same value of $`\varphi _0`$, but varying the rest of the higher order terms. Since, given the same value of $`\varphi _0`$, the first two terms in the $`F`$ expansion for the linear and divergent models coincide, they should be thought of as equivalent. Therefore, the models with the linear $`k`$-field and the inverse square potential possess a symmetry, or rather a degeneracy in the sense that all the solutions for which the first two terms of the expansion of the function $`F`$ around $`x=x_0`$ coincide are equivalent among themselves and also equivalent to the divergent model. This does not happen with the solutions with the non-linearly behaved scalar field, and suggests that physically the later are more acceptable, thus partially removing the degeneracy of the solutions.

A more subtle distinction between the models, to completely remove the degeneracy, would be probably seen by perturbing these solutions. This, however, is beyond the scope of the present paper.

## 4 Conclusions

In this paper we have studied particular solutions to the Einstein Equations coupled to $`k`$-essence. Imposing spatially flat isotropic geometry we have shown that different $`k`$-theory Lagrangians may lead to the same kinematical evolution of the universe.

We have seen, however, that the linear $`k`$-field power-law solutions possess an odd property of being isomorphic to a family of solutions with a divergent speed of sound generated by the function $`F^{\mathrm{}}`$. This relation of isomorphism induces, in fact, problems with the power-law solutions regardless of the model (function$`F`$) as long as the potential is inverse square and the field is linear, leading to consider different solutions. It has been recently argued  that the cosmological models based on linear $`k`$ field lead to serious problems. These problems are associated with the behaviour of the models in the divergent sound speed region. We believe that our findings relating the linear $`k`$-field models with the divergent models sheds new light on the reasons of the peculiar behavior of those models in the region of the divergent speed of sound .

This work was supported by the University of Buenos Aires under Project X223. A.F. acknowledges the support of the University of the Basque Country Grants 9/UPV00172.310-14456/2002, and The Spanish Science Ministry Grant 1/CI-CYT 00172. 310-0018-12205/2000.