text:
# A molecular dynamics study on the equilibrium magnetization properties and structure of ferrofluids

## I Introduction

Ferrofluids (dipolar magnetic fluids) are colloidal suspensions of ferromagnetic particles of about $`10nm`$ diameter dispersed in a carrier liquid rosensweig85a  that are usually stabilized against agglomeration by coating particles with long-chain molecules (sterically) or decorating them with charged groups (electrostatically). The small size of the particles favors magnetic mono-domains with a magnetic moment proportional to the volume of the magnetic grains. As a result, the particles interact with each other by the long-range anisotropic dipole-dipole potential as well as the short-range symmetric potentials, such as the steric repulsion, the electrostatic repulsion and the van der Waals attraction. The study of both the magnetization properties and structure of ferrofluids are of both fundamental and application interests.

A ferrofluid system of sufficiently low concentration behaves like an ideal paramagnetic gas. The interactions between the particles can be neglected, and the physical properties of the system are well described by the one-particle model shliomis74a . In this case the equilibrium magnetization is written in terms of the Langevin function $`(\alpha )=\mathrm{coth}(\alpha )1/\alpha `$

$$M_L=M_{sat}(\alpha ),M_{sat}=\frac{N}{V}\frac{m}{\mu _0}.$$
(1)
Here, $`M_{sat}`$ is the saturation magnetization of the fluid, $`\alpha =mH/kT`$ is the Langevin parameter, and $`m`$ is the magnetic moment, $`N/V`$ denotes the number density of the particles, $`H`$ is the magnetic field, $`T`$ the temperature, $`k`$ the Boltzmann constant, and $`\mu _0=4\pi \times 10^7H/m`$. Eq.(1) leads to Curie’s law for the initial susceptibility $`\chi _L`$ which depends linearly on the particle concentration

$$\chi _L=\frac{N}{V}\frac{m^2}{3\mu _0kT}.$$
(2)
However, experiments on concentrated ferrofluids have revealed an essential deviation from the Langevin formula pshenichnikov95a . The initial susceptibility increases with the particle concentration much faster than that predicted in Eq.(2). This deviation is clearly due to the particle interactions.

The influence of the particle interactions on the physical properties of ferrofluids is investigated in terms of two dimensionless parameters. One is the volume fraction of the particles

$$\varphi =\frac{N}{V}\frac{\pi \sigma ^3}{6},$$
(3)
with $`\sigma `$ the particle diameter. The other is the dipolar coupling constant

$$\lambda =\frac{m^2}{4\pi \mu _0kT\sigma ^3},$$
(4)
which relates the dipole-dipole interaction energy of two contacting particles to the thermal energy $`kT`$. The two parameters can be combined into the Langevin initial susceptibility $`\chi _L`$ in Eq.(2) as $`\chi _L=8\varphi \lambda `$.

A number of theoretical models allow the evaluation of the equilibrium magnetization of ferrofluids onsager36a ; grady83a ; wertheim71a ; morozov90a ; buyevich92a ; pshenichnikov96a ; ivanov01a ; huke00a . Among them, the mean-spherical model wertheim71a ; morozov90a , the thermodynamic perturbation model buyevich92a  and the modified variant of the effective field model pshenichnikov96a  have been shown to give good results for the magnetic properties of ferrofluids with low or moderate concentration of magnetic particles (up to $`1012\%`$). Taking $`\chi _L`$ as a universal parameter in weak magnetic fields, all the above approaches lead to the same expression for the initial susceptibility

$$\chi =\chi _L(1+\chi _L/3).$$
(5)
For ferrofluids of higher concentration, larger particle size, or at low temperature, Eq.(5) is insufficient ivanov01a . Higher order corrections in the dependence of $`\chi `$ or $`M`$ on $`\varphi `$, $`\lambda `$, or in general on $`\chi _L`$, turn out to be crucial.

An expression for $`\chi `$ up to cubic accuracy in $`\chi _L`$ can be obtained by the mean-spherical model wertheim71a , the Born-Mayer expansion method huke00a  and the statistical model based on the pair correlation function ivanov01a , which reads asivanov01a

$$\chi =\chi _L(1+\chi _L/3+\chi _L^2/144).$$
(6)
The cubic term makes a noticeable difference of 5% over Eq.(5) for $`\chi _L4.14`$. The significance of the third order contribution has been verified in comparison with experimental measurements on the temperature dependence of $`\chi `$ ivanov01a . However, in what we classify here as the strong coupling limit, $`\lambda >2`$, the enhanced particle aggregation is not sufficiently accounted for by the analytical studies, thus giving rise to considerable deviations from the predictions of Eqs.(5,6). Here the magnetic properties turn out to depend on $`\mathrm{\Phi }`$ and $`\lambda `$ separately rather than $`\chi _L`$ alone, a point which we will try to elucidate with our simulation study.

Previous simulation works on dipolar (ferro) fluids were mainly undertaken to explore the phase diagram. Large ranges of particle concentration and dipolar coupling strength have been examined. The results are found to depend substantially on the employed short range interactions, such as hard sphere weis93a ; weis93b ; levesque94a ; camp00a , soft sphere wei92a ; wei92b ; stevens94a ; stevens95a  and Lennard-Jones (LJ) potentials leeuwen93a ; stevens95b . As an example, the question whether a minimum amount of dispersive energy, i. e., the attractive van der Waals energy, is required to observe the liquid-vapor coexistence is still under discussion camp00a ; leeuwen93a ; stevens94a ; stevens95a ; stevens95b ; levin99a ; tlusty00a . Recent studies on the magnetic properties include the Monte Carlo simulations on the dispersions of interacting super-paramagnetic particles in a solid matrix chantrell00a  and on ferrofluids with a finite spherical boundary pshenichnikov00a . In order to compare with the macroscopic theories on ferrofluids, it is necessary to simulate the ferrofluid systems with truly periodic boundary conditions.

In this paper, the Langevin dynamics simulation method is used to study the equilibrium properties of ferrofluids. The long-range dipolar interactions are dealt with the Ewald summation. The initial susceptibility and the magnetization curves are calculated as a function of $`\varphi `$ and $`\lambda `$. The simulation results are directly compared with the theoretical models. The influence of particle aggregation is investigated systematically by performing a cluster analysis of the microstructure of the systems. The formation of aggregates is found to enhance the magnetization of ferrofluids at weak magnetic fields. This is confirmed by comparison with the calculated initial susceptibility for a ferro-solid, i. e., a ferrofluid with the same properties but with the particles’ translational degree of freedom frozen in at random positions.

The paper is organized as follows. We describe the simulation method in Sec.II. The simulation results on the initial susceptibility and the magnetization curves are given in Sec.III (A) and (B), respectively. The cluster analysis of the microstructure is presented in Sec.III (C). The initial susceptibility of the ferro-solid systems is calculated in Sec.III(D). We end with our conclusions in Sec.IV.

## II Simulation Method

The investigated ferrofluid systems consist of $`N`$ spherical particles of diameter $`\sigma `$ distributed in a cubic simulation box of side length $`L`$. Each particle has a permanent point dipole moment $`𝐦_𝐢`$ at its center. Using periodic boundary conditions in all spatial directions, the dipole-dipole interaction potential between particle $`i`$ and $`j`$ is given by

$$U_{ij}^{dip}=\frac{1}{4\pi \mu _0}\underset{𝐧^3}{}\left\{\frac{𝐦_i𝐦_j}{𝐫_{ij}+𝐧L^3}\frac{3[𝐦_i(𝐫_{ij}+𝐧L)][𝐦_j(𝐫_{ij}+𝐧L)]}{𝐫_{ij}+𝐧L^5}\right\},$$
(7)
where $`𝐫_{ij}=𝐫_i𝐫_j`$ is the displacement vector of the two particles. The sum extends over all simple cubic lattice points, $`𝐧=(n_x,n_y,n_z)`$ with $`n_x`$, $`n_y`$, $`n_z`$ integers.

In this work, we use the Ewald summation for dipolar systems to evaluate Eq.(7) effectively, which gives leeuw80a ; allen87a ; wang01a

$$U_{ij}^{dip}=U_{ij}^{(r)}+U_{ij}^{(k)}+U_{ij}^{(self)}+U_{ij}^{(surf)},$$
(8)
where the real-space $`U_{ij}^{(r)}`$, the k-space (reciprocal-space) $`U_{ij}^{(k)}`$, the self $`U_{ij}^{(self)}`$, and the surface $`U_{ij}^{(surf)}`$ contributions are respectively given by:

$`U_{ij}^{(r)}`$ $`=`$ $`{\displaystyle \frac{1}{4\pi \mu _0}}{\displaystyle \underset{𝐧^3}{}}\left\{(𝐦_i𝐦_j)B(𝐫_{ij}+𝐧)[𝐦_i(𝐫_{ij}+𝐧)][𝐦_j(𝐫_{ij}+𝐧)]C(𝐫_{ij}+𝐧)\right\},`$ (9)
$`U_{ij}^{(k)}`$ $`=`$ $`{\displaystyle \frac{1}{4\pi \mu _0L^3}}{\displaystyle \underset{𝐤^3,𝐤0}{}}{\displaystyle \frac{4\pi }{k^2}}\mathrm{exp}[(\pi k/\kappa L)^2](𝐦_i𝐤)(𝐦_j𝐤)\mathrm{exp}(2\pi i𝐤𝐫_{ij}/L),`$ (10)
$`U_{ij}^{(self)}`$ $`=`$ $`{\displaystyle \frac{1}{4\pi \mu _0}}{\displaystyle \frac{2\kappa ^3}{3\sqrt{\pi }}}(m_i^2+m_j^2),`$ (11)
$`U_{ij}^{(surf)}`$ $`=`$ $`{\displaystyle \frac{1}{4\pi \mu _0}}{\displaystyle \frac{4\pi }{(2\mu _{BC}+1)L^3}}𝐦_i𝐦_j.`$ (12)
with

$`B(r)`$ $`=`$ $`[\text{erfc}(\kappa r)+(2\kappa r/\sqrt{\pi })\mathrm{exp}(\kappa ^2r^2)]/r^3,`$ (13)
$`C(r)`$ $`=`$ $`[3\text{erfc}(\kappa r)+(2\kappa r/\sqrt{\pi })(3+2\kappa ^2r^2)\mathrm{exp}(\kappa ^2r^2)]/r^5.`$ (14)
Here $`\text{erfc}(x):=2\pi ^{1/2}_x^{\mathrm{}}\mathrm{exp}(t^2)dt`$ denotes the complementary error function. The inverse length $`\kappa `$ is the splitting parameter of the Ewald summation. Eq.(12) assumes that the set of the periodic replications of the simulation box tends in a spherical way towards an infinite cluster and that the medium outside this sphere is a uniform medium with magnetic permeability $`\mu _{BC}`$. In this work, we use the metallic boundary condition $`\mu _{BC}=\mathrm{}`$. The surface term vanishes in this case and demagnetization effects do not occur. Thus the applied external magnetic field $`𝐇`$ coincides exactly with the internal field. The influence of different boundary conditions on the simulation results will be the topic of a subsequent publication wang02b . The related formulas for the dipolar forces and torques can be found in a previous paper wang01a . The theoretical estimates of the cutoff errors in the Ewald summation derived there are used to determine the optimal values for the Ewald parameters which enable us to minimize the overall computational time at a predefined accuracy.

The short range interaction potential between the particles is chosen to mimic the steric mechanism for stabilizing the colloidal suspension. As mentioned in Sec.I, it is still not settled which kind of potential is most suitable to model real ferrofluids. It might actually depend on whether the solution is electrically or sterically stabilized. In this work we adopt the Lennard-Jones potential

$$U_{ij}^{LJ}=4\epsilon \left[\left(\frac{\sigma }{r_{ij}}\right)^{12}\left(\frac{\sigma }{r_{ij}}\right)^6C(R_c)\right],$$
(15)
where $`C(R_c)=(\sigma /R_c)^{12}(\sigma /R_c)^6`$ with a cutoff radius of $`R_c=2^{1/6}\sigma `$. In this way the particles have a purely repulsive interaction force which smoothly decays to zero at $`R_c`$. This cutoff range is smaller than the one usually used in the soft sphere potential which is commonly adopted as $`U_{SS}(r_{ij})=4\epsilon (\sigma /r_{ij})^{12}`$ wei92a ; wei92b ; stevens94a ; stevens95a , and thus closer to the hard sphere potential as employed in most other theoretical calculations ivanov01a ; huke00a . In this context the analytical calculations of Ivanov ivanov01a  reveal that at least for the spatially homogeneous systems the influence of different short range repulsive potentials seems to be marginal.

The translational and rotational Langevin equations of motion of particle $`i`$ are given by allen87a ; murashov00a

$`_i\dot{𝐯_i}`$ $`=`$ $`𝐅_i\mathrm{\Gamma }_T𝐯_i+𝝃_i^T,`$ (16)
$`𝐈_i\dot{𝝎_i}`$ $`=`$ $`𝝉_i\mathrm{\Gamma }_R𝝎_i+𝝃_i^R,`$ (17)
where $`_i`$ and $`𝐈_i`$ are the mass and inertia tensor of the particle, $`\mathrm{\Gamma }_T`$ and $`\mathrm{\Gamma }_R`$ are the translational and rotational friction constants, respectively. The Gaussian random force and torque have the properties that the first moments vanish

$`<\xi _{i\alpha }^T(t)>`$ $`=`$ $`0,`$ (18)
$`<\xi _{i\alpha }^R(t)>`$ $`=`$ $`0,`$ (19)
while the second moments satisfy

$`<𝝃_{i\alpha }^T(t)𝝃_{j\beta }^T(t^{})>=6kT\mathrm{\Gamma }_T\delta _{ij}\delta _{\alpha \beta }\delta (tt^{}),`$ (20)
$`<𝝃_{i\alpha }^R(t)𝝃_{j\beta }^R(t^{})>=6kT\mathrm{\Gamma }_R\delta _{ij}\delta _{\alpha \beta }\delta (tt^{}),`$ (21)
where $`\alpha `$ and $`\beta `$ denote the $`x,y,z`$ components in the Cartesian coordinates, respectively. Introducing the dipolar and short range terms into Eqs.(16) and (17), the dimensionless equations of motion can now be written as:

$`\dot{𝐯_i^{}}={\displaystyle \underset{ji}{}}(𝐅_{ij}^{dip}+𝐅_{ij}^{LJ})\mathrm{\Gamma }_T^{}𝐯_i^{}+𝝃_i^T,`$ (22)
$`𝐈_i^{}\dot{𝝎_i^{}}={\displaystyle \underset{ji}{}}𝝉_{ij}^{dip}+𝐦_i^{}\times 𝐇^{}\mathrm{\Gamma }_R^{}𝝎_i^{}+𝝃_i^R,`$ (23)
where the variables are given in dimensionless form reduced by the following units: length $`r^{}=r/\sigma `$, dipole moment $`m^2=m^2/4\pi \mu _0ϵ\sigma ^3`$, moment of inertia $`I^{}=I/(\sigma ^2)`$, time $`t^{}=t(ϵ/\sigma ^2)^{1/2}`$, the friction constants $`\mathrm{\Gamma }_T^{}=\mathrm{\Gamma }_T(\sigma ^2/ϵ)^{1/2}`$ and $`\mathrm{\Gamma }_R^{}=\mathrm{\Gamma }_R/(\sigma ^2ϵ)^{1/2}`$, magnetic field $`H^{}=H(4\pi \mu _0\sigma ^3/ϵ)^{1/2}`$ as well as temperature $`T^{}=kT/ϵ`$. The values of the dimensionless friction constants do not affect the equilibrium properties. Here we adopt $`\mathrm{\Gamma }_T^{}=10.0`$ and $`\mathrm{\Gamma }_R^{}=3.0`$, because in our simulations these values had been found to speed up the equilibration time. These values are also close to that used in Ref.murashov00a . A value of $`I^{}=0.4`$ is used for the dimensionless moment of inertia, corresponding to that of a rigid sphere with diameter $`\sigma `$.

The simulations were performed at constant temperature $`T^{}=1`$. The orientational coordinates of the particles were expressed in terms of quaternion parameters and the equations of motions were integrated with a leap-frog algorithm allen87a . A reduced time step $`\mathrm{\Delta }t^{}=0.002`$ was employed in all simulations. The runs were started from initial configurations with random particle positions and dipole moment orientations. For each case, the system was at first equilibrated for a period of 50 000 time steps. The magnetization and structural properties were then calculated from the data for another period of at least 200 000 time steps. Error bars for the simulation results were determined by dividing the simulation runs into blocks of 10 000 time steps and calculating an estimate for the standard deviation of the mean kusalik90a .

## III Results and Discussion

Results were obtained from simulations on systems with $`N=1000`$ particles. The number density of the particles was mostly taken in the range of $`\rho ^{}=N/V^{}0.45`$, corresponding to volume fractions $`\varphi =\rho ^{}\pi /60.236`$. This covers the typical particle concentration ranges which are obtainable in real ferrofluids. A few more runs were also carried out at higher concentrations in case of smaller $`\lambda `$ in order to reach a sufficient large $`\chi _L`$. The values of the dipolar coupling constant $`\lambda `$ were mainly chosen to be $`1,2,3`$ and $`4`$. To get an estimate of realistic $`\lambda `$-values, we consider a ferrofluid sample consisting of monodispersed magnetite particles. The bulk magnetization of this material is $`M_d=4.46\times 10^5`$ A/m. With a magnetic core diameter of the particles of $`\sigma =10`$ nm, the value of $`\lambda `$ amounts to about $`1.3`$ at room temperature $`T=300K`$. For larger particles with $`\sigma =13`$ nm, $`\lambda `$ rises up to $`2.9`$ under the same conditions. Note that $`\lambda `$ may further be enlarged by decreasing the temperature or using cobalt as working material ($`M_d=14.01\times 10^5`$ A/m), for which the previous $`\lambda `$-values increase by a factor of around $`3.14`$. We therefore also carried out some simulations for $`\lambda =8`$.

In all simulations, we fixed the root mean square (rms) absolute errors in the dipolar forces to $`\mathrm{\Delta }F^{dip}10^4m^2/4\pi \mu _0\sigma ^4`$, which corresponds to $`10^4/6`$ of the attractive force between two contacting particles with dipole moments in parallel alignment. The optimal values of the Ewald parameters had been determined separately for each system wang01a . As an example, for a sample with $`\rho ^{}=0.1`$ we get $`r_c=0.425L`$, $`\kappa =6.78/L`$ and $`K_c=6\times 2\pi /L`$. While $`\rho ^{}=0.3`$ yields $`r_c=0.46L`$, $`\kappa =6.84/L`$ and $`K_c=7\times 2\pi /L`$ to be the optimum values, respectively.

### III.1 Initial susceptibility

The influence of the particle interactions on the magnetic properties is most evident at weak fields, so we start with the investigation of the initial susceptibility $`\chi `$. The value of $`\chi `$ is determined from the linear relationship of $`𝐌=\chi 𝐇`$ at $`H0`$. Simulations are thus performed to get the $`M(H)`$ curves at weak fields. The equilibrium magnetization $`𝐌`$ is evaluated by the prescription

$$𝐌=\frac{1}{\tau }_{t_0}^{t_0+\tau }(\frac{1}{\mu _0V}\underset{i=1}{\overset{N}{}}𝐦_i)𝑑t,$$
(24)
where $`t_0`$ is the time needed for the system to establish the equilibrium state.

To study the weak field magnetization, we increase the Langevin parameter $`\alpha `$ from $`0`$ to $`0.5`$. As expected, the calculation of $`M`$ at $`\alpha =0`$ does not show the existence of a spontaneous magnetization, because $`\varphi `$ is still not high enough. Other works stevens95a  have shown that the transition to a ferromagnetic state occurs at $`\varphi 0.44`$ for $`\lambda =4`$ . For higher values of $`\lambda `$, the critical concentration is around $`\varphi 0.35`$ for $`\lambda =6.25`$ weis93b ; stevens95a  and $`\varphi 0.31`$ for $`\lambda =9`$ wei92b . These conditions are clearly not encountered in our simulations.

The dimensionless magnetization $`M^{}(=M/\sqrt{4\pi ϵ/\mu _0\sigma ^3})`$ of the system with $`\lambda =3`$ is shown in Fig.1 as a function of $`\alpha `$ and $`\varphi `$. The results for other values of $`\lambda `$ display a similar behavior. It can be seen that $`M^{}`$ at first increases linearly with $`\alpha `$, then at higher field values sub-linearly. We introduce $`\alpha _{lin}`$ as a definition for the maximum value of $`\alpha `$ where $`M^{}`$ still can be considered to obey the linear relationship. Evidently, at a fixed value of $`\lambda `$ this $`\alpha _{lin}`$-value decreases with growing $`\varphi `$. The data in Fig.1 show that when $`\varphi `$ is increased from $`0.052`$ to $`0.209`$, $`\alpha _{lin}`$ decreases from $`0.45`$ to about $`0.25`$. On the other hand, at a fixed particle concentration, rising $`\lambda `$ also implies a reduction of $`\alpha _{lin}`$. This demonstrates that the range over which $`M`$ is proportional to $`H`$ (or equivalently $`\alpha `$) shrinks with the increase of $`\chi _L`$. For a quantitative evaluation of the initial susceptibility we need to determine the linear region for each case separately. This is accomplished with the help of a linear regression fit with error weighting press92a . The basic idea of that method is to minimize the so-called $`\xi `$-square merit function

$$\xi ^2=\underset{i=1}{\overset{N_d}{}}(\frac{M_iM(0)\chi H_i}{\mathrm{\Delta }M_i})^2$$
(25)
so as to get $`\chi `$ and its error bar. Here $`\mathrm{\Delta }M_i`$ is the error bar related to the simulation data of $`M_i`$. In Fig.1 the data were taken with a step size of $`\mathrm{\Delta }\alpha =0.05`$ giving $`11`$ data points on each curve. Since no spontaneous magnetization occurs, the initial data point is always at $`M(0)=0.0`$. For each curve, a series of linear fittings was generated by increasing $`N_d`$ in Eq.(25) from $`5`$ to $`11`$. The fit which gives the minimum value of $`\xi ^2`$ per degree of freedom is taken to be the optimum linear fitting to the data. The associated $`N_d`$-value is referred as $`N_{d,opt}`$. The slope $`\chi `$ and its error bar $`\mathrm{\Delta }\chi `$ determine the initial susceptibility. The value of $`\alpha _{lin}`$ is roughly given by the product of $`\mathrm{\Delta }\alpha `$ with $`N_{d,opt}1`$. The results of this data fit including the error estimates are compiled in Table.I. The degree of reliability of the results is reflected by the high goodness-probability $`Q`$ press92a . The resulting straight lines $`M=\chi H`$, or in dimensionless form $`M^{}=\chi \alpha /4\pi \sqrt{\lambda /T^{}}`$, are also displayed in Fig.1 to show how they fit to the simulation data.

Fig.2(a) shows the calculated initial susceptibility $`\chi `$ as a function of $`\lambda `$ and $`\varphi `$, along with the analytical results from Eqs.(5) and (6). Since the Langevin susceptibility $`\chi _L`$ has been taken as the universal parameter in these theoretical predictions buyevich92a ; pshenichnikov96a ; ivanov01a , we also map the data of $`\chi `$ on $`\chi _L`$ in Fig.2(b) to show the scaling behavior more directly. In addition to the results in Fig.2(a), some more data was produced at higher concentrations (up to $`\varphi =0.419`$) for the case of $`\lambda =1`$ in order to get a larger range of $`\chi _L`$. We did not combine such high concentrations with larger values of $`\lambda `$, since the close proximity to the magnetic phase separation would make the discussion more complex, rendering the results more sensitive to the specific short-range interaction.

For the weak interaction range $`\lambda 2`$ we find perfect agreement with Eqs.(5) and (6) for all investigated volume fractions. Actually, because for these coupling parameters we can reach only up to $`\chi _L<4`$, both analytical predictions yield equally good results for our data and are indistinguishable within the error bars. For $`\lambda =3`$ we see small but significant deviations at low values of $`\varphi `$, meaning that already at low values of $`\chi _L`$ neither Eqs.(5) nor (6) are accurate. However, at high volume fractions ($`\chi _L>3`$), the third order equation (6) is compatible within error bars. For the numerical values of the data see Table 1, where also the values of the analytical results for different orders are given for comparison. Even more drastic deviations from the theoretical predictions are found for $`\lambda =4`$ at low concentrations, in the range $`\chi _L6`$. We observe again the trend, that at higher volume fractions the third order equation (6) becomes compatible with our data. We can conclude that the parameter $`\lambda `$, rather than $`\chi _L`$, becomes the dominating control parameter for $`\lambda >2`$. This effect can be attributed directly to the microstructure formation in the system. When $`\lambda 2`$, what we consider to be the weak interaction limit, the dipolar potential is comparable to the thermal fluctuation. Then the interparticle correlations are weak, and the present theoretical models work well. However, going to stronger interactions with $`\lambda >2`$, leads to a considerable amount of particle aggregation. The resulting clusters amplify the sensitivity of the system to weak external magnetic fields and thus lead to an enhanced magnetic susceptibility. We will come back to this point in Sec. III.3.

The influence of the volume fraction seems to be weak in the investigated parameter regime. Noteworthy is that the third order equation seems to describe the available data well at larger volume fraction and high coupling parameter $`\lambda `$, even better than it does at intermediate values. On the basis of our available data we cannot decide if this finding is a mere coincidence. To our knowledge there is no theoretical reason to believe that the third order equation should become reliable again at these high values of $`\chi _L`$, where higher order corrections of the form $`\varphi ^m\lambda ^n`$ with $`m`$, $`n`$ arbitrary, are expected to become important. Corrections linear in $`\varphi `$, but up to ninth order in $`\lambda `$ have been assembled for the Born-Mayer expansion method in Ref. huke00a, . These additional terms are all positive and therefore necessarily increase the susceptibility. We have found that they improve the theoretical prediction of $`\chi `$ only at very small concentrations ($`\varphi <0.05`$) for $`\lambda >2`$, but they dramatically overestimate the simulation results at larger values of $`\varphi `$.

Note that the simulations presented in Ref.pshenichnikov00a,  give rise to a susceptibility $`\chi `$ that decreases with growing $`\lambda `$ if $`\chi _L`$ is fixed, which is opposite to our findings in Fig.2(b). We suspect that this trend may be induced by finite size effects, since a finite spherical simulation volume with only a few hundred particles had been used in their simulations (without a periodic replication to properly account for the long range dipolar interaction as is done in the present work). We shall clarify this point in a subsequent publicationwang02b .

### III.2 Magnetization curves

In this section, we extend the calculation of $`M`$ to a larger range of the applied external field strength $`\alpha `$ in order to get the full magnetization curves. For moderately concentrated ferrofluids, the ‘modified mean-field model’ pshenichnikov96a  describes the magnetization curve in terms of the Langevin magnetization $`M_L`$ \[Eq.(1)\] as follows

$$M(H)=M_L(H_e),$$
(26)
where the effective field acting on an individual ferro-particle is given by

$$H_e=H+M_L(H)/3,$$
(27)
rather than $`H_e=H+M/3`$ as done in the Weiss model. The resulting expression for the initial susceptibility is Eq.(5). A higher order description is obtained by replacing $`H_e`$ in Eq.(27) by

$$H_e=H+\frac{M_L(H)}{3}+\frac{1}{144}M_L(H)\frac{dM_L(H)}{dH},$$
(28)
which in the weak field limit reduces to Eq.(6) for $`\chi `$ ivanov01a . The effective Langevin parameters corresponding to Eqs.(27) and (28) can be written as $`\alpha _e=\alpha +\chi _L(\alpha )`$ and $`\alpha _e=\alpha +\chi _L(\alpha )+\frac{\chi _L^2}{16}(\alpha )\frac{d(\alpha )}{d\alpha }`$, respectively. This reflects the $`\chi _L`$-dependence of these two theoretical predictions. The Born-Mayer expansion described $`\chi `$ in a different way, where its dependence on $`\varphi `$ and $`\lambda `$ are given separately huke00a .

We investigated the magnetization curves for two different cases. In the first case $`\chi _L`$ is fixed to a small value of $`1.256`$ where the predictions using Eqs.(27) and (28) practically coincide. Then $`\lambda `$ is taken to be $`1,2,4`$ and $`8`$, respectively. To keep $`\chi _L`$ constant, the volume fraction $`\varphi `$ had to be adapted from $`0.157`$ down to $`0.0196`$. In the second case, we take a larger value of $`\chi _L=5.026`$ which makes the difference between the two theoretical predictions distinguishable by a few percentage. Only $`\lambda =3`$ and $`4`$ are discussed. The corresponding volume fractions are $`\varphi =0.209`$ and $`0.157`$, respectively. In the latter case, we did not include the results for $`\lambda 2`$ due to the limitation of the particle concentrations and for $`\lambda =8`$ due to the extremely slow dynamics which prohibited us from reaching equilibrium distributions within our computational limits.

The calculated magnetization curves are shown in Fig.3 together with the theoretical predictions using Eqs.(26-28). For a more explicit demonstration of the particle interaction the inlays display the difference between the simulation data and the Langevin model $`(MM_L)/M_{sat}`$. The characteristic maximum of this difference is found to be at intermediate field strength of $`\alpha 1`$, where it may reach more than $`35\%`$. The results in Fig.3(a) reveal again the essential role of the dipolar coupling constant $`\lambda `$. When $`\lambda `$ equals $`1`$ and $`2`$, the simulation and theoretical curves agree with each other very well. But for larger $`\lambda `$ the calculated magnetization $`M`$ rises with $`\alpha `$ much faster than the theoretical predictions, especially in the weak field regime. This reflects the same result already observed in Fig.2, that $`\chi `$ increases with $`\lambda `$ at low and moderate concentrations if $`\chi _L`$ is held fixed. Due to this effect, the peak positions of the increment $`(MM_L)/M_{sat}`$ shift to smaller values of $`\alpha `$ with the increase of $`\lambda `$. In Fig.3(b), the volume fraction of the particles are relatively high. The maximum absolute difference between the two analytical predictions occurs at $`\alpha =0.62`$ where the result given by using Eq.(28) is about $`4.3\%`$ larger than that by Eq.(27). Our simulation results are found to be consistent with both of them within the range of error bar. A quantitative analysis of the data indicates that Eq.(28) gives a better description in the weak field regime. This consequently leads to a better agreement between the simulation results on $`\chi `$ with the prediction of Eq.(6). However, both data sets for $`\lambda =3`$ and $`\lambda =4`$ show deviations from the theoretical curves at larger values of $`\alpha `$ which cannot be explained by statistical errors. Also here we infer that higher order corrections are needed to more accurately describe the magnetization curves.

### III.3 Microstructure analysis

The microstructure and its relation to the magnetic properties of ferrofluids are investigated by performing a cluster analysis. In conventional cluster analysis, agglomerates are usually defined on the base of the spatial proximity between the particles or by means of an energy criterion. For ferrofluids the latter one is more favorable, because the anisotropic dipolar interaction implies that two neighboring particles can form a stable bonding only if their dipole moments are roughly aligned in head-to-tail orientation.

The energetically based cluster analysis of ferrofluids has been done in several different ways weis93a ; levesque94a ; stevens95a ; pshenichnikov00a ; coverdale98a . In this work, we adopt the definition of Refs. weis93a, ; levesque94a,  and stevens95a,  where two particles are considered to be bound if their dipolar potential energy is less than a predetermined value $`U_{bond}`$. The influence of using different threshold values for $`U_{bond}`$ has been investigated in a recent work tavares99a . The values of $`1.4\lambda kT`$ weis93a ; levesque94a  and $`1.5\lambda kT`$ stevens95a  were shown both to give good results. In our calculations we have varied $`U_{bond}`$ from $`1.4\lambda kT`$ to $`1.6\lambda kT`$. Only slight quantitative changes of the results have been found without a qualitative difference. The results shown here have been evaluated for $`U_{bond}=1.5\lambda kT`$, i. e., $`75\%`$ of the contact energy of two perfectly co-aligned dipolar particles. The above cluster definition implies that aggregates are mainly linear chains and rings. Thicker particle agglomerates or branched chains, however, occur very rarely in our simulations.

The size of a cluster is defined as the number of particles belonging to it. We first set up the connectivity matrix by means of our cluster definition, and then use it to evaluate the sizes and numbers of clusters in a given configuration. The average cluster size is defined by

$$S_{avg}=<\underset{s}{}sn_s/\underset{s}{}n_s>.$$
(29)
where $`n_s`$ is the number of clusters having size $`s`$, and the triangular brackets denote the time average, or equivalently, the average over the configuration space.

Fig.4 plots $`S_{avg}`$ as a function of $`\lambda `$ and $`\varphi `$ for the case of zero field ($`\alpha =0`$). Clearly, increasing the dipolar coupling $`\lambda `$ leads to a larger cluster size $`S_{avg}`$. However, in the weak coupling limit $`\lambda =1`$ or $`2`$, $`S_{avg}`$ exceeds $`1`$ only slightly, meaning that there are very few dimers present. Inter-particle correlations seem to be weak and we already noted that the simulation results agree well with the theoretical models in these cases.

As $`\lambda `$ is increased, the larger values of $`S_{avg}`$ indicate the formation of more clusters. To elucidate this, we show snapshots of the three-dimensional configurations in Fig.5 for the case of $`\varphi =0.052`$ ($`\rho ^{}=0.1`$) at different $`\lambda `$ values. The increase of the number and size of clusters with $`\lambda `$ can be seen clearly in these figures. For $`\lambda =3`$ we find $`S_{avg}=1.123`$ at $`\varphi =0.052`$ . About $`10\%`$ of the particles are organized in dimers and less than $`1\%`$ in trimers. At the same density, $`S_{avg}`$ increases to $`1.481`$ for $`\lambda =4`$. Correspondingly, $`22\%`$ of the particles are in dimers, $`7\%`$ in trimers and $`2\%`$ in quadrumers. The density dependence of these numbers is shown in Fig.6 where the average percentage of particles in $`n`$-mers are plotted as a function of $`\varphi `$ for $`\lambda =3`$ and $`4`$. It can be seen that more than $`60\%`$ of the particles are still in monomers and the clusters with $`s5`$ are rarely found in these cases. However, when $`\lambda `$ equals to $`8`$, all the particles get involved in the formation of long chain-like structures. The average chain length is larger than $`13\sigma `$ at $`\varphi =0.052`$. Due to thermal fluctuations, these chains are quite flexible. They are continuously breaking and recombining with others weis93a ; levesque94a . Even closed ring structures can be observed occasionally. In this strong coupling case, the chains have a close analogy to living polymers levesque94a ; tavares99a .

When the external field is switched on, the tendency of the dipole moments to co-align with the field direction enhances the aggregation probability. Accordingly the average cluster size is expected to increase with the field strength. Fig.7 shows the field-dependence of $`S_{avg}`$ for the systems studied in Fig.3(a). The increment of $`S_{avg}`$ with $`\alpha `$ is only very weak for $`\lambda 2`$, but it is more evident at $`\lambda =4`$. As exhibited in Fig.8, the structure character of the system varies from randomly distributed short clusters to string-like alignments along the magnetic flux with the increase of $`\alpha `$. Consequently the cluster analysis in Fig.9 reflects that more and more particles merge into larger aggregates. Cluster as long as $`s=15`$ could be observed at high $`\alpha `$. But these string-like structures in Fig.8(d) are much more unstable than the chains formed at $`\lambda =8`$ \[see Fig.5(d)\]. They posses a smaller contact energy so that thermal fluctuations can easily break long clusters into smaller segments. This effect is field-independent. In fact, most of the strings observed in Fig.8(d) consist of a few continuously breaking and recombining segments. As a result, the average cluster size $`S_{avg}`$ does not present a dramatic change as expected from the snapshot of the configurations. When $`\lambda =8`$, however, the long-chains had already been formed even at zero field. The increase of the field strength just arrange these chains to have a better alignment along the field direction. But the average chain length does not change too much.

The results in Fig.7 reveal that the magnetic field dependence of $`S_{avg}`$ remains weak at small $`\alpha `$. It is therefore legal to invoke the zero-field cluster analysis in order to interpret the outcomes related to the initial susceptibility (Fig.2). Let us focus to the cases of $`\lambda =3`$ and $`4`$ at low concentrations where the deviation between theory and simulation are most pronounced. In these cases the small clusters are relatively stable and spatially well separated from each other. For a rough approximation, we treat the clusters as rigid particles thus considering our system as a poly-dispersed sample. We define the effective dipole moment of each cluster according to

$$𝐦_s=\underset{i=1}{\overset{s}{}}𝐦_i.$$
(30)
The effective Langevin susceptibility $`\chi _L^{eff}`$ of such a poly-dispersed system can be computed as

$$\chi _L^{eff}=\frac{1}{V}\frac{<_{allclusters}m_s^2>}{3\mu _0kT}.$$
(31)
Since most of the clusters are short linear $`n`$-mers, the value of $`\chi _L^{eff}`$ is larger than $`\chi _L`$ \- the appropriate value for the mono-dispersed system. Replacing $`\chi _L`$ by $`\chi _L^{eff}`$ in Eq.(6) we get new modified theoretical curves. In Fig.10 it can be seen that this simple operating prescription considerably improves the agreement between simulation and theory, at least at low concentrations. This proves that the formation of small clusters is responsible for the enhancement of the initial susceptibility.

At particle concentration above $`\varphi 0.1`$, the modified theoretical curves start to overestimate the simulated $`\chi `$-values. In these cases, the increased angular correlations between the particles are not just limited to the two next nearest neighbors and thus make the above definition of clusters ambiguous levesque94a . Although we keep on using the energy criterion to obtain the cluster informations for Fig.4 and 6, these values are more likely reflecting the local orientational order of the dipole moments. The clusters we find are not stable entities involving well-defined particle numbers as in the low density regime. This can be partially seen from the $`S_{avg}(\varphi )`$ behavior in Fig.4. When $`\lambda 2`$, $`S_{avg}`$ slightly increases with $`\varphi `$. But it saturates at $`\lambda =4`$ and even decreases at $`\lambda =8`$ for $`\varphi 0.1`$. Next we examine the radial distribution function $`g(r)`$ of the system (Fig.11), which gives the probability of finding a pair of particles a distance $`r`$ apart, relative to the probability expected for a completely random distribution at the same densityallen87a . We find that the height of the first and second peaks, respectively, associated with dimer and trimers, decrease with the increase of $`\varphi `$. This implies that the number of clusters does not increase as rapidly as the density stevens95a . In other words, upon rising $`\varphi `$ the spatial distribution of the particle positions re-approaches to that of a homogeneous system. Therefore we expect our suggested modification with Eq.(31) only to be valid at low volume fractions below $`\varphi =0.1`$.

The formation of long chain-like structures in the strong coupling regime makes the systematic investigation of the magnetization properties at $`\lambda =8`$ difficult. At $`\varphi =0.0196`$, there are very few overlaps between the elongated chains. They have enough free volume to re-arrange in response to the external field. The magnetization of the system increases very fast in the weak field regime, as shown in Fig.3(a). The occasional formation of a few closed ring structures is not found to affect noticeably the global magnetization behavior. On the other hand, when $`\varphi `$ is enlarged, the chains become strongly entangled, like in a semi-flexible polymer solution. This slows down the dynamics of the system dramatically. An extremely long simulation time is required to reach the equilibrium state in the weak field regime. For example, at a density of $`\varphi =0.0785`$, the equilibrium value of $`𝐌`$ at $`\alpha =0.5`$ is obtained only after 1.5 million simulation time steps. Simultaneously, since the number of particles is fixed in our simulation, the increase of $`\varphi `$ makes the chain length become comparable to the size of the simulation box. This unavoidably induces finite size effects. For these reasons, we were unable to provide systematic investigations at $`\lambda =8`$ at higher concentrations (or equivalently larger $`\chi _L`$), but only for $`\chi _L=1.256`$ shown in Fig.3(a). These numerical difficulties could well be related to experimental observations which find a dramatic increase of viscosities at similar dipolar couplings odenbachpc .

### III.4 Comparison with ferro-solids

To further investigate the influence of the cluster formation process we now turn to ferro-solid systems. By this we denote a suspension of dipolar particles which have their spatial positions frozen in, but are permitted to rotate freely. In this way the aggregation related effects in ferrofluids can easily be separated.

We prepared the initial configurations of the ferro-solid system by placing the particles randomly in the simulation box with random dipole moment orientations. For each given $`\lambda `$ and $`\varphi `$, several different initial configurations were generated and simulated. They were found to give very similar results. Computer time limitations prohibited us to average over a large sequence of quenched configurations, which would be desirable to obtain a good ensemble average for a ferro-glass system. However, we have investigated this case only to study the qualitative influence of freezing the positional degrees of freedom. Accordingly, only the rotational equations of motion \[Eq.(17)\] of the particles need to be integrated in the simulations. Fig.12 displays the simulation results for the initial susceptibility at $`\lambda =3`$ and $`4`$, along with the theoretical predictions for the fluid system. In ferro-solids $`\chi `$ also exhibits a nonlinear increase with the volume fraction $`\varphi `$. Comparing our data with the results of the ferrofluid system in Fig.2, we see that the suppression of clusters leads to a slower increase of the magnetization at weak fields, and hence to a lower $`\chi `$. It is interesting to note that the discrepancy between $`\chi _{solid}`$ and $`\chi _{fluid}`$ goes up when $`\lambda `$ or $`\varphi `$ are increased. This shows that the positional correlations are also very important at high density, and couple strongly to the dipolar rotational degrees of freedom, even if they do not evidently form stable clusters. We consider the fact that the measured values of $`\chi _{solid}`$ coincide with the theoretical curves for the fluid system, Eqs.(5) and (6), almost exactly at low concentration for an accidental correspondence. At higher volume fractions the ferro-solid displays a susceptibility even below the theoretical prediction. This effect is more evident for larger $`\lambda `$ at higher $`\varphi `$. Note, that in Ref. pshenichnikov00a,  a similar solidified ferro-colloid was investigated. There $`\chi _{solid}`$ was found to increase with $`\chi _L`$, but the values are larger than that for the associated ferrofluid. Again we suspect that this due to finite size effects of their sample.

We point out that the particles in our ferro-solid systems are supposed to be super-paramagnetic. This means that they have a negligible crystallographic magnetic anisotropy energy. There is no energy barrier which prevents the magnetic dipoles from a free rotation. That makes the present systems different from that one studied in Ref.chantrell00a,  where the initial susceptibility is found to decrease with the increase of the particle concentrations.

## IV Conclusion

We have investigated in detail the initial susceptibility, magnetization curve and microstructure of ferrofluids in various particle concentration and dipolar coupling strength ranges by means of molecular dynamics simulations. The ferrofluid was modeled as a soft sphere system by means of a purely repulsive Lennard-Jones potential, and the magnetic cores interacting with a dipolar coupling strength $`\lambda `$ from $`1`$ to $`8`$. The investigated volume fractions were in the range of $`\varphi 0.24`$ mostly. Our simulation method used the Ewald summation with metallic boundary conditions to deal with the long-range dipolar interactions, and took explicitly into account the translational and rotational degrees of freedom of the dipolar particles. The temperature was kept constant by means of a Langevin thermostat for all degrees of freedom. Our simulation results for the susceptibility show good agreement with the theoretical predictions in the weak coupling limit $`\lambda 2`$. For higher dipolar couplings, however, we find systematic deviations from the theoretical predictions at low and intermediate concentrations. The detailed cluster analysis demonstrates that this is due to the particle aggregation at higher $`\lambda `$. The formation of clusters tend to increase the magnetization at weak field regime and consequently induce a larger initial susceptibility. At high densities, the spatial distribution of the particles starts to homogenize again, and the significance of clusters goes down. In fact we speculate that at higher concentrations there are almost no stable clusters any more, because they can disintegrate and reform with other particles rapidly. In this regime we find that the third order equation in $`\chi _L`$ is compatible with our data.

For the magnetization curves we find that both, the second and third order theoretical curves, describe our data well for couplings $`\lambda 2`$. For higher values of $`\lambda `$ we find again systematic deviations. The system shows at coupling strength of $`\lambda =8`$ an almost gel-like behavior due to entangled chain structures. We expect this effect to become experimentally relevant in Cobalt based ferrofluids, where $`\lambda `$ is large.

The influence of the aggregation phenomena was isolated by studying ferro-solids where the particles are randomly placed inside the simulation volume, and only rotations are allowed to occur. The observed initial susceptibilities are lower compared to the ferrofluid systems due to the suppression of cluster formation. The deviations to the fluid system grow with increasing dipolar coupling strength and increasing volume fraction.

## Acknowledgments

We thank B. Huke and S. Odenbach for helpful discussions. Financial support from the DFG under grant No. HO 1108/8-2 is greatly appreciated.


format: text

text:
# Examples of Anosov diffeomorphisms

## 1. Introduction

Anosov diffeomorphisms play an important and beautiful role in dynamics as the notion represents the most perfect kind of global hyperbolic behavior, giving examples of structurally stable dynamical systems (see $`\mathrm{\S }`$2 for a precise definition). In , S. Smale raised the problem of classifying the compact manifolds (up to homeomorphism) which admit an Anosov diffeomorphism. At this moment, the only known examples are of algebraic nature, namely Anosov automorphisms of nilmanifolds and infranilmanifolds. It is conjectured that any Anosov diffeomorphism is topologically conjugate to an Anosov automorphism of an infranilmanifold. J. Franks  and A. Manning  proved the conjecture for Anosov diffeomorphisms on infranilmanifolds themselves.

All this certainly highlights the problem of classifying all nilmanifolds which admit Anosov automorphisms, which are easily seen in correspondence with hyperbolic automorphisms of nilpotent Lie algebras over $``$. Nevertheless, not too much is known on the question since it is not so easy for an automorphism of a (real) nilpotent Lie algebra being hyperbolic and unimodular at the same time. We propose the following

###### Definition 1.1.

A rational Lie algebra $`𝔫`$ (i.e. with structure constants in $``$) of dimension $`d`$ is said to be Anosov if it admits a hyperbolic automorphism $`A`$ (i.e. all their eigenvalues have absolute value different from $`1`$) such that $`[A]_\beta \mathrm{G}L(d,)`$ for some basis $`\beta `$ of $`𝔫`$, where $`[A]_\beta `$ denotes the matrix of $`A`$ with respect to $`\beta `$.

It is well known that any Anosov Lie algebra is necessarily nilpotent, and it is easy to see that the classification of nilmanifolds which admit an Anosov automorphism is essentially equivalent to that of Anosov Lie algebras (see $`\mathrm{\S }`$2). If $`𝔫`$ is a rational Lie algebra, we call the real Lie algebra $`𝔫`$ the real completion of $`𝔫`$.

Curiously enough, the only explicit examples of Anosov Lie algebras in the literature so far are an infinite family of $`6`$-dimensional rational Lie algebras all having real completion isomorphic to $`𝔥_3𝔥_3`$ (see ) –where $`𝔥_3`$ is the 3-dimensional Heisenberg Lie algebra–, the free $`k`$-step nilpotent Lie algebras on $`n`$ generators with $`k<n`$ (see , and also  for a different approach) and certain $`k`$-step nilpotent Lie algebras of dimension $`d+\left(\genfrac{}{}{0pt}{}{d}{2}\right)+\mathrm{}+\left(\genfrac{}{}{0pt}{}{d}{k}\right)`$ with $`dk^2`$ (see ). For the known examples of infranilmanifolds which are not nilmanifolds and admit Anosov automorphisms we refer to .

This motivates the following natural questions:

* The minimal possible dimension for an Anosov non-abelian Lie algebra is $`6`$ (see ). Is there an Anosov Lie algebra (without a nonzero abelian factor) in every dimension $`d6`$?.
* It is proved in \[12, Theorem 3.2\] that any Anosov $`k`$-step nilpotent Lie algebra $`𝔫`$ has $`dim𝔫2k+2`$. Is there for any $`k2`$ an Anosov $`𝔫`$ with $`dim𝔫=2k+2`$?.
* Does there exist a dimension $`d`$ for which there are infinitely many Anosov Lie algebras such that their real completions are pairwise non-isomorphic?.
* Is it feasible a more or less explicit classification of all Anosov Lie algebras (or at lest of their real completions) up to isomorphism?.
* Does every Anosov automorphism on a nilmanifold come from the Anosov action of some lattice in a semisimple Lie group? (see ).

The purpose of this paper is to give a procedure to construct explicit examples of Anosov Lie algebras, which is inspired on the famous $`6`$-dimensional example given by S. Smale in  (due to A. Borel), but from the point of view of L. Auslander and J. Scheuneman in \[1, Section 4\]. We prove that $`𝔫𝔫`$ is Anosov for any graded rational nilpotent Lie algebra $`𝔫`$. The construction is very direct and simple, but it produces examples enough to answer positively questions (i) (for $`d`$ even), (ii) and (iii), and negatively questions (iv) and (v). Also, two examples of Anosov Lie algebras which are not of the form $`𝔫𝔫`$ are given. One of them has dimension $`9`$, and so by considering direct sums we get examples for question (i) for any $`d7,11,13`$.

We note that any Anosov Lie algebra is necessarily graded as they always admit semisimple hyperbolic automorphisms (see ).

Acknowledgements. I wish to express my gratitude to F. Grunewald for pointing me out Remark 3.5. I am also grateful to D. Fisher and R. Miatello for helpful observations.

## 2. Anosov diffeomorphisms on nilmanifolds

A diffeomorphism $`f`$ of a compact differentiable manifold $`M`$ is called Anosov if the tangent bundle $`\mathrm{T}M`$ admits a continuous invariant splitting $`\mathrm{T}M=E^+E^{}`$ such that $`\mathrm{d}f`$ expands $`E^+`$ and contracts $`E^{}`$ exponentially, that is, there exist constants $`0<c`$ and $`0<\lambda <1`$ such that

$$\mathrm{d}f^n(X)c\lambda ^nX,XE^{},\mathrm{d}f^n(Y)c\lambda ^nY,YE^+,$$

for all $`n`$. The condition is independent of the Riemannian metric.

Let $`N`$ be a real simply connected nilpotent Lie group with Lie algebra $`𝔫`$. Let $`\phi `$ be a hyperbolic automorphism of $`N`$, that is, all the eigenvalues of its derivative $`A=(\mathrm{d}\phi )_e:𝔫𝔫`$ have absolute value different from $`1`$. If $`\phi (\mathrm{\Gamma })=\mathrm{\Gamma }`$ for some lattice $`\mathrm{\Gamma }`$ of $`N`$ then $`\phi `$ defines an Anosov diffeomorphism on the nilmanifold $`N/\mathrm{\Gamma }`$, which shall be called an Anosov automorphism. The subspaces $`E^+`$ and $`E^{}`$ are obtained by left translation of the eigenspaces of eigenvalues of $`A`$ of absolute value greater than $`1`$ and less than $`1`$, respectively.

We now review the well-known relationship with Definition 1.1 (see  for more detailed expositions). The lattice $`\mathrm{\Gamma }`$ contains, as a subgroup of finite index, some full lattice $`\mathrm{\Gamma }_1`$, that is, $`\mathrm{log}(\mathrm{\Gamma }_1)`$ is a $``$-Lie subalgebra of $`𝔫`$ (i.e. with integer structure constants). Thus $`\phi ^m(\mathrm{\Gamma }_1)=\mathrm{\Gamma }_1`$ for some $`m`$ and hence $`A^m`$ is an hyperbolic automorphism of the $``$-Lie algebra $`\mathrm{log}(\mathrm{\Gamma }_1)`$. This implies that the rational completion $`\mathrm{log}(\mathrm{\Gamma }_1)`$ is an Anosov Lie algebra.

Conversely, let $`A`$ be an automorphism of a rational Lie algebra $`𝔫`$ satisfying the conditions in Definition 1.1. Then there exists a $``$-Lie subalgebra $`𝔫_1`$ of $`𝔫`$ such that $`\mathrm{\Gamma }=\mathrm{exp}(𝔫_1)`$ is a lattice in $`N`$, where $`N`$ is the simply connected nilpotent Lie group with Lie algebra $`𝔫`$, the real completion of $`𝔫`$. Now for some $`m`$, $`A^m(𝔫_1)=𝔫_1`$ and hence $`A^m`$ determines an Anosov automorphism of the nilmanifold $`N/\mathrm{\Gamma }`$. The second condition on $`A`$ in Definition 1.1 is equivalent to the fact that the characteristic polynomial of $`A`$ has integer coefficients and constant term equal to $`\pm 1`$ (see ).

Two lattices $`\mathrm{\Gamma }_1,\mathrm{\Gamma }_2`$ of $`N`$ are called commensurable if $`\mathrm{\Gamma }_1\mathrm{\Gamma }_2`$ has finite index in both $`\mathrm{\Gamma }_1`$ and $`\mathrm{\Gamma }_2`$, or equivalently, the Lie algebras of their rational Mal’cev completions are isomorphic as rational Lie algebras. It is proved in  and  that if $`\mathrm{\Gamma }_1`$ and $`\mathrm{\Gamma }_2`$ are commensurable then $`N/\mathrm{\Gamma }_1`$ admits an Anosov automorphism if and only if $`N/\mathrm{\Gamma }_2`$ does. This fact is the main reason for defining the notion of Anosov in the context of rational Lie algebras.

Recall that if $`𝔫_1`$ and $`𝔫_2`$ are two non-isomorphic Anosov Lie algebras, then the corresponding Anosov diffeomorphisms are not topologically conjugate since the nilmanifolds $`N_1/\mathrm{\Gamma }_1`$ and $`N_2/\mathrm{\Gamma }_2`$ can never be homeomorphic. Indeed, $`\mathrm{\Pi }_1(N_i/\mathrm{\Gamma }_i)=\mathrm{\Gamma }_i`$ and $`𝔫_i`$ is the Lie algebra of the rational Mal’cev completion of $`\mathrm{\Gamma }_i`$.

Let $`𝔫`$ be a real nilpotent Lie algebra of dimension $`d`$. We note that $`𝔫`$ has a rational form which is Anosov (or equivalently, $`N`$ is the simply connected cover of a nilmanifold admitting an Anosov automorphism) if and only if there exists a $``$-basis $`\beta `$ of $`𝔫`$ and a hyperbolic $`A\mathrm{Aut}(𝔫)`$ such that $`[A]_\beta \mathrm{G}L(d,)`$, where $`[A]_\beta `$ is the matrix of $`A`$ in terms of $`\beta `$.

## 3. Construction of examples

A rational Lie algebra $`𝔫`$ is said to be graded if there exist (rational) subspaces $`𝔫_i`$ of $`𝔫`$ such that

$$𝔫=𝔫_1𝔫_2\mathrm{}𝔫_k\text{and}[𝔫_i,𝔫_j]𝔫_{i+j}.$$

Equivalently, $`𝔫`$ is graded when there are nonzero (rational) subspaces $`𝔫_{d_1},\mathrm{},𝔫_{d_r}`$, $`d_1<\mathrm{}<d_r`$, such that $`𝔫=𝔫_{d_1}\mathrm{}𝔫_{d_r}`$ and if $`0[𝔫_{d_i},𝔫_{d_j}]`$ then $`d_i+d_j=d_k`$ for some $`k`$ and $`[𝔫_{d_i},𝔫_{d_j}]n_{d_k}`$. Recall that any graded Lie algebra is necessarily nilpotent.

###### Theorem 3.1.

Let $`𝔫`$ be a graded rational nilpotent Lie algebra. Then the direct sum $`\stackrel{~}{𝔫}=𝔫𝔫`$ is Anosov.

###### Proof.

Let $`\{X_1,\mathrm{},X_d\}`$ be a $``$-basis of $`𝔫`$ compatible with the gradation $`𝔫=𝔫_{d_1}\mathrm{}𝔫_{d_r}`$, that is, each $`X_i𝔫_{d_j}`$ for some $`j`$. We denote by $`\{Y_i\}`$ a copy of the basis $`\{X_i\}`$, so that $`\{X_1,\mathrm{},X_d,Y_1,\mathrm{},Y_d\}`$ is a basis of $`\stackrel{~}{𝔫}`$ and

(1)
$$[X_i,X_j]=\underset{k=1}{\overset{d}{}}m_{ij}^kX_k,[Y_i,Y_j]=\underset{k=1}{\overset{d}{}}m_{ij}^kY_k,m_{ij}^k.$$

Every nonzero $`\lambda `$ defines an automorphism $`A_\lambda `$ of $`𝔫`$ by

$$A_\lambda |_{𝔫_{d_i}}=\lambda ^{d_i}I,$$

and also an automorphism $`\stackrel{~}{A}_\lambda `$ of $`\stackrel{~}{𝔫}`$ by

(2)
$$\stackrel{~}{A}_\lambda =\left[\begin{array}{cc}A_\lambda & \\ & A_\lambda ^1\end{array}\right].$$

If $`a`$, $`a2`$, then the roots $`\lambda ,\lambda ^1`$ of $`x^22ax+1[x]`$ are

(3)
$$\lambda =a+(a^21)^{\frac{1}{2}},\lambda ^1=a(a^21)^{\frac{1}{2}},$$

and hence $`0<\lambda ^1<1<\lambda `$. Consider the new basis of $`\stackrel{~}{𝔫}`$ defined by

$$\beta =\{X_1+Y_1,(a^21)^{\frac{1}{2}}(X_1Y_1),\mathrm{},X_d+Y_d,(a^21)^{\frac{1}{2}}(X_dY_d)\}.$$

It follows from (1) that

$$\begin{array}{c}[X_i+Y_i,X_j+Y_j]=\underset{k=1}{\overset{d}{}}m_{ij}^k(X_k+Y_k),\hfill \\ [X_i+Y_i,(a^21)^{\frac{1}{2}}(X_jY_j)]=\underset{k=1}{\overset{d}{}}m_{ij}^k(a^21)^{\frac{1}{2}}(X_kY_k),\hfill \\ [(a^21)^{\frac{1}{2}}(X_iY_i),(a^21)^{\frac{1}{2}}(X_jY_j)]=\underset{k=1}{\overset{d}{}}(a^21)m_{ij}^k(X_k+Y_k).\hfill \end{array}$$

This implies that $`\beta `$ is also a $``$-basis. On the other hand, it is easy to see by using (3) that if $`T`$ is the transformation whose matrix in terms of the basis $`\{X_i,Y_i\}`$ equals

$$\left[\begin{array}{cc}\lambda & 0\\ 0& \lambda ^1\end{array}\right],$$

then the matrix of $`T`$ with respect to the basis $`\{X_i+Y_i,(a^21)^{\frac{1}{2}}(X_iY_i)\}`$ is given by

$$B=\left[\begin{array}{cc}a& a^21\\ 1& a\end{array}\right]\mathrm{S}L(2,).$$

Therefore the matrix of $`\stackrel{~}{A}_\lambda \mathrm{Aut}(\stackrel{~}{𝔫})`$ (see (2)) in terms of the basis $`\beta `$ is given by

$$\left[\stackrel{~}{A}_\lambda \right]_\beta =\left[\begin{array}{ccccccc}B^{d_1}& & & & & & \\ & \mathrm{}& & & & & \\ & & B^{d_1}& & & & \\ & & & \mathrm{}& & & \\ & & & & B^{d_r}& & \\ & & & & & \mathrm{}& \\ & & & & & & B^{d_r}\end{array}\right]\mathrm{S}L(2d,).$$

Thus $`\stackrel{~}{A}_\lambda `$ determines a hyperbolic automorphism of $`\stackrel{~}{𝔫}`$ whose matrix in terms of $`\beta `$ is in $`\mathrm{S}L(2d,)`$, and hence $`\stackrel{~}{𝔫}`$ is Anosov. ∎

###### Remark 3.2.

In spite of the simply connected cover of the examples given by Theorem 3.1 are all of the form $`N\times N`$, the nilmanifolds $`(N\times N)/\mathrm{\Gamma }`$ which admit Anosov automorphisms are not a direct product of nilmanifolds when $`N`$ is irreducible.

We will now answer questions (i)-(v) in the Introduction. Theorem 3.1 will be applied in all cases without any further reference.

Question (iv). We first note that any two-step nilpotent Lie algebra is graded. Therefore, the classification of all Anosov nilpotent Lie algebras contains the classification of all rational two-step nilpotent Lie algebras. This is equivalent to the classification of all alternating bilinear maps

$$\mu :^n\times ^n^m$$

under the equivalence relation $`\mu \lambda `$ if and only if $`\mu (gX,gY)=h\lambda (X,Y)`$ for some $`g\mathrm{G}L(n,)`$, $`h\mathrm{G}L(m,)`$, which is considered for $`m>2`$ a wild problem (see ). Moreover, the classification up to isomorphism of the real two-step nilpotent Lie algebras admitting a rational form is also completely open.

Question (iii). We can answer this question positively for $`d=14`$ in the following way. For each $`0<t<1`$ consider the $`7`$-dimensional nilpotent Lie algebra $`𝔫_t`$ with Lie bracket $`\mu _t`$ defined by

$$\begin{array}{ccc}\mu _t(X_1,X_2)=(1t)^{\frac{1}{2}}X_3,\hfill & & \mu _t(X_2,X_3)=X_5,\hfill \\ \mu _t(X_1,X_3)=X_4,\hfill & & \mu _t(X_2,X_4)=X_6,\hfill \\ \mu _t(X_1,X_4)=t^{\frac{1}{2}}X_5,\hfill & & \mu _t(X_2,X_5)=t^{\frac{1}{2}}X_7,\hfill \\ \mu _t(X_1,X_5)=X_6,\hfill & & \mu _t(X_3,X_4)=(1t)^{\frac{1}{2}}X_7.\hfill \\ \mu _t(X_1,X_6)=X_7,\hfill & & \end{array}$$

Each Lie algebra $`\mu _t`$ is isomorphic to the Lie algebra $`𝔤_t`$, denoted by $`1,2,3,4,5,7_I:t`$ in \[16, pp.494\] (see also the curve $`\stackrel{~}{g}(0,t,1,0,1,0,0,0)`$ in \[11, 5.2.3\]). The isomorphism $`\phi _t:𝔤_t𝔫_t`$ is given by

$$\begin{array}{c}\phi _tX_1=X_1,\phi _tX_2=t^{\frac{1}{2}}X_2,\phi _tX_i=t^{\frac{1}{2}}(1t)^{\frac{1}{2}}X_i,i=3,4,\hfill \\ \\ \phi _tX_j=t(1t)^{\frac{1}{2}}X_j,j=5,6,7.\hfill \end{array}$$

This proves that $`𝔫_t`$ ($`0<t<1`$) is a curve of pairwise non-isomorphic real Lie algebras. It is easy to check that $`\{X_1,\mathrm{},X_7\}`$ is a $``$-basis of $`𝔫_{t_k}`$ for every

$$t_k=\frac{4k^2}{(k^2+1)^2},k,$$

and that the subspaces $`𝔫_i=X_i`$ determines a gradation $`𝔫_{t_k}=𝔫_1\mathrm{}𝔫_7`$ for any $`k`$. Thus $`\stackrel{~}{𝔫}_k=𝔫_{t_k}𝔫_{t_k}`$ is Anosov for all $`k`$, and $`\stackrel{~}{𝔫}_k\stackrel{~}{𝔫}_k^{}`$ if and only if $`k=k^{}`$, since the function $`f(x)=\frac{4x^2}{(x^2+1)^2}`$ is strictly decreasing in $`[1,\mathrm{})`$.

We can also get examples for this question from any pairwise non-isomorphic infinite sequence $`\{𝔫_k\}`$ of real two-step nilpotent Lie algebras all admitting a rational form. There is an explicit example of such a sequence of dimension $`10`$ in , and so $`d=20`$.

Question (ii). One can easily answer affirmatively this question by considering for each $`k2`$ the following $`k`$-step nilpotent Lie algebra $`𝔫`$ of dimension $`k+1`$:

$$[X_1,X_2]=X_3,[X_1,X_3]=X_4,\mathrm{},[X_1,X_k]=X_{k+1},$$

for which $`\{X_1,\mathrm{},X_{k+1}\}`$ is clearly a $``$-basis, and $`𝔫`$ is graded by $`𝔫=X_1\mathrm{}X_{k+1}`$. Thus $`\stackrel{~}{𝔫}=𝔫𝔫`$ is Anosov, $`(2k+2)`$-dimensional and $`k`$-step nilpotent. Actually, any $`d`$-dimensional graded rational Lie algebra which is filiform (i.e. $`(d1)`$-step nilpotent) provides an example as required in this question.

Question (v). It is easy to prove that $`\mathrm{Aut}(𝔫𝔫)`$ is solvable for any $`𝔫`$ defined in question (ii) with $`k3`$ (and for any filiform example with $`d4`$ as well). Thus $`\mathrm{Aut}(𝔫𝔫)`$ does not contain any semisimple Lie subgroup, which implies that the Anosov automorphism on $`(N\times N)/\mathrm{\Gamma }`$ determined by $`\stackrel{~}{A}_\lambda `$ (see (2)) does not come from the Anosov action of any lattice in a semisimple Lie group (compare with ).

Question (i). We have required examples without a nonzero abelian factor in order to avoid Lie algebras of the form $`𝔫^r`$ with $`𝔫`$ Anosov. For $`r2`$, the corresponding nilmanifolds $`N/\mathrm{\Gamma }\times T^r`$, where $`T^r`$ is the $`r`$-dimensional torus, clearly admit Anosov automorphisms. The examples considered in question (ii) also give a positive answer to this question for any dimension $`d6`$ which is even. Thus, by using the $`9`$-dimensional Anosov Lie algebra given in Example 3.4, we can consider direct sums and get examples of Anosov Lie algebras without an abelian factor for any dimension $`d15`$. Thus the question only remains open for $`d=7,11,13`$. It should be noted that the examples of odd dimensions $`d15`$ are all direct products of nilmanifolds (compare with Remark 3.2).

We finally want to give two examples of Anosov two-step nilpotent Lie algebras which are not of the form $`𝔫𝔫`$. We note that the construction of the first one is definitely similar in spirit to that in the proof of Theorem 3.1, from where it can be deduced that such a construction might be just a particular case of a much more general procedure.

###### Example 3.3.

Let $`𝔫`$ be the two-step nilpotent Lie algebra with basis

$$\{X_1,X_2,X_3,Y_1,Y_2,Y_3,Z_1,Z_2\}$$

and Lie bracket defined by

$$[X_1,X_2]=Z_1,[X_1,X_3]=Z_2,[Y_1,Y_2]=Z_1,[Y_1,Y_2]=Z_2.$$

Consider $`A\mathrm{Aut}(𝔫)`$ given by

$$\begin{array}{c}AX_1=\lambda X_1,AX_2=\lambda X_2,AX_3=\lambda ^3X_3,\hfill \\ \\ AY_1=\lambda ^1Y_1,AY_2=\lambda ^3Y_2,AY_3=\lambda ^1Y_3,\hfill \\ \\ AZ_1=\lambda ^2Z_1,AZ_2=\lambda ^2Z_2,\hfill \end{array}$$

where $`\lambda `$ is as in (3). It is easy to see that

$$\begin{array}{cc}\hfill \beta =& \{X_1+Y_1,(a^21)^{\frac{1}{2}}(X_1Y_1),X_2+Y_3,(a^21)^{\frac{1}{2}}(X_2Y_3),X_3+Y_2,\hfill \\ & (a^21)^{\frac{1}{2}}(X_3Y_2),Z_1+Z_2,(a^21)^{\frac{1}{2}}(Z_1Z_2)\}\hfill \end{array}$$

is a $``$-basis of $`𝔫`$ and that if $`B=\left[\begin{array}{cc}a& a^21\\ 1& a\end{array}\right]`$ then

$$\left[A\right]_\beta =\left[\begin{array}{cccc}B& & & \\ & B& & \\ & & B^3& \\ & & & B^2\end{array}\right]\mathrm{S}L(8,),$$

showing that $`𝔫`$ is Anosov.

The following example is just a modification of the free $`2`$-step nilpotent Lie algebra on $`3`$ generators (see ).

###### Example 3.4.

Consider the $`(3r+3)`$-dimensional $``$-Lie algebra

$$𝔫=^3\mathrm{}^3\mathrm{\Lambda }^2^3$$

with Lie bracket defined by

$$[v_1+\mathrm{}+v_r,w_1+\mathrm{}+w_r]=v_1w_1+\mathrm{}+v_rw_r.$$

If $`A\mathrm{G}L(3,)`$ has eigenvalues $`\lambda _1,\lambda _2,\lambda _3`$ and their products $`\lambda _i\lambda _j`$ all of absolute value different from $`1`$ then it easy to see that

$$\left[\begin{array}{cccc}A& & & \\ & \mathrm{}& & \\ & & A& \\ & & & \mathrm{\Lambda }^2A\end{array}\right]\mathrm{G}L(3r+3,),$$

is an hyperbolic automorphism of $`𝔫`$. For instance, one can take

$$A=\left[\begin{array}{ccc}1& 1& 1\\ 1& 2& 2\\ 1& 2& 3\end{array}\right]$$

(see ). Thus $`𝔫`$ is Anosov of dimension $`3r+3`$.

###### Remark 3.5.

It was pointed me out by F. Grunewald that the construction given in Theorem 3.1 can be generalized by using the ring of integers in any totally real number field $`F`$ over $``$ (we have used $`F=(a^21)^{\frac{1}{2}}`$). The proof is essentially the same, it only uses Dirichlet Unit Theorem (compare with \[9, Proposition 3.2\] and ). If the degree of $`F`$ is $`m`$ then we obtain Anosov Lie algebras of the form $`𝔫\mathrm{}𝔫`$ ($`m`$ times) for any graded rational nilpotent Lie algebra $`𝔫`$.