Role of length-polydispersity on the phase behavior of freely-rotating hard-rectangle fluid

TL;DR

This study uses density functional theory to analyze how length-polydispersity influences phase transitions and stability in a fluid of freely-rotating hard rectangles, revealing effects on transition order and phase stability.

Contribution

It provides a detailed theoretical analysis of the impact of polydispersity on phase behavior and phase transitions in hard-rectangle fluids, including stability of the tetratic phase.

Findings

Polydispersity increases the first-order nature of the isotropic-nematic transition.

The tricritical point shifts to higher aspect ratios with increased polydispersity.

The tetratic phase is slightly destabilized as polydispersity increases.

Abstract

We used the Density Functional formalism, in particular the Scaled Particle Theory, applied to a length-polydisperse hard-rectangular fluid to study its phase behavior as a function of the mean particle aspect ratio ($\kappa_0$) and polydispersity ($\Delta_0$). The numerical solutions of the coexistence equations were calculated by transforming the original problem with infinite degrees of freedoms to a finite set of equations for the amplitudes of the Fourier expansion of the moments of the density profiles. We divided the study into two parts: The first one is devoted to the calculation of the phase diagrams in the packing fraction ($\eta_0$)- $\kappa_0$ plane for a fixed $\Delta_0$ and selecting parent distribution functions with exponential (the Schulz distribution) or Gaussian decays. In the second part we study the phase behavior in the $\eta_0$-$\Delta_0$ plane for fixed $\kappa_0$ while $\Delta_0$ is changed. We characterize in detail the orientational ordering of particles and the fractionation of different species between the coexisting phases. Also we study the character (second vs. first order) of the Isotropic-Nematic phase transition as a function of polydispersity. We particularly focused on the stability of the Tetratic phase as a function of $\kappa_0$ and $\Delta_0$. The Isotropic-Nematic transition becomes strongly of first order when polydispersity is increased: the coexisting gap widens and the location of the tricritical point moves to higher values of $\kappa_0$ while the Tetratic phase is slightly destabilized with respect to the Nematic one. The results obtained here can be tested in experiments on shaken monolayers of granular rods.

Equations

$$f_0(\kappa )=C_{\nu ,q}({\kappa }_0){\left(\frac{\kappa -1}{{\kappa }_0-1}\right)}^{\nu }\exp \left[-{\lambda }_{\nu ,q}{\left(\frac{\kappa -1}{{\kappa }_0-1}\right)}^q\right],$$

$$C_{\nu ,q}({\kappa }_0)=q{\left({\kappa }_0-1\right)}^{-1}\frac{{\Gamma }^{\nu +1}\left[(\nu +2)/q\right]}{{\Gamma }^{\nu +2}\left[(\nu +1)/q\right]},$$

$${\lambda }_{\nu ,q}^{1/q}=\frac{\Gamma \left[(\nu +2)/q\right]}{\Gamma \left[(\nu +1)/q\right]},$$

$$\Delta \equiv \sqrt{\frac{⟨{\kappa }^2⟩-⟨\kappa {⟩}^2}{⟨\kappa {⟩}^2}}=\left(1-\frac{1}{{\kappa }_0}\right){\Delta }_0,$$

$${\Delta }_0\equiv \sqrt{\frac{\Gamma \left[(\nu +1)/q\right]\Gamma \left[(\nu +3)/q\right]}{{\Gamma }^2\left[(\nu +2)/q\right]}-1},$$

$${\Phi }_{\mathrm{id}}[\rho ]\equiv \beta {\mathcal{F}}_{\mathrm{id}}[\rho ]/A=\int d\kappa \int d\varphi \rho (\kappa ,\varphi )\left[\log \rho (\kappa ,\varphi )-1\right]$$

$${\Phi }_{\mathrm{ex}}[\rho ]\equiv \beta {\mathcal{F}}_{\mathrm{id}}[\rho ]/A=-{\rho }_0^{(0)}\ln \left(1-{\rho }_0^{(1)}{\sigma }^2\right)$$

$$+\frac{S_0\left[\rho \right]}{1-{\rho }_0^{(1)}{\sigma }^2},$$

$$S_0\left[\rho \right]\equiv \int d{\kappa }_1\int d{\kappa }_2\int d{\varphi }_1\int d{\varphi }_2\rho ({\kappa }_1,{\varphi }_1)\rho ({\kappa }_2,{\varphi }_2)$$

$$\times A_0({\kappa }_1,{\kappa }_2,{\varphi }_1-{\varphi }_2),$$

$$A_0({\kappa }_1,{\kappa }_2,\varphi )=\frac{{\sigma }^2}{2}\{({\kappa }_1{\kappa }_2+1)|\sin \varphi |$$

$$+\left({\kappa }_1+{\kappa }_2\right)|\cos \varphi |\}.$$

$${\rho }_0^{(i)}\equiv {\int }_0^{\pi }d\varphi {\rho }^{(i)}(\varphi )={\int }_0^{\pi }d\varphi \left[{\int }_1^{\infty }d\kappa {\kappa }^i\rho (\kappa ,\varphi )\right],$$

$${\rho }_0{f}_0(\kappa )=\frac{1}{\pi }\sum _{\alpha }{\gamma }_{\alpha }{\int }_0^{\pi }d\varphi {\rho }^{(\alpha )}(\kappa ,\varphi ),$$

$${\rho }^{(\alpha )}(\kappa ,\varphi )=\frac{{\rho }_0{f}_0(\kappa )\exp \left[-\beta {\mu }_{\mathrm{ex}}^{(\alpha )}(\kappa ,\varphi )\right]}{{\pi }^{-1}{\sum }_{\tau }{\gamma }_{\tau }{\int }_0^{\pi }d{\varphi }^{\prime }\exp \left[-\beta {\mu }_{\mathrm{ex}}^{(\tau )}(\kappa ,{\varphi }^{\prime })\right]},$$

$${\rho }^{(i,\alpha )}(\varphi )=\frac{1}{\pi }\sum _{j=0}^N{\rho }_j^{(i,\alpha )}\cos (2j\varphi ),$$

$${\rho }_j^{(i,\alpha )}={\rho }_0{D}_{j0}{\int }_1^{\infty }d\kappa {\kappa }^i{f}_0(\kappa )\frac{T_j^{(\alpha )}(\kappa )}{{\sum }_{\tau }{\gamma }_{\tau }{T}_0^{(\tau )}(\kappa )},$$

$$T_j^{(\tau )}(\kappa )\equiv {\pi }^{-1}{\int }_0^{\pi }d\varphi \exp \left[-\beta {\mu }_{\mathrm{ex}}^{(\tau )}(\kappa ,\varphi )\right]\cos (2j\varphi ),$$

$$\beta {\mu }_{\mathrm{ex}}^{(\alpha )}(\kappa ,\varphi )=-\ln \left[1-{\rho }_0^{(1,\alpha )}{\sigma }^2\right]+\frac{S_1^{(\alpha )}(\kappa ,\varphi )}{1-{\rho }_0^{(1,\alpha )}{\sigma }^2}$$

$$+\beta p^{(\alpha )}\kappa {\sigma }^2,$$

$$S_1^{(\alpha )}(\kappa ,\varphi )\equiv \frac{2{\sigma }^2}{\pi }\{\left({\rho }_0^{(1,\alpha )}+{\rho }_0^{(0,\alpha )}\right)\left(\kappa +1\right)$$

$$-\sum _{j\ge 1}\frac{\left({\rho }_j^{(1,\alpha )}+(-1{)}^j{\rho }_j^{(0,\alpha )}\right)\left(\kappa +(-1{)}^j\right)}{4j^2-1}\cos (2j\varphi )\}.$$

$$\beta p^{(\alpha )}=\frac{{\rho }_0^{(0,\alpha )}}{1-{\rho }_0^{(1,\alpha )}{\sigma }^2}+\frac{S_0^{(\alpha )}[\rho ]}{{\left[1-{\rho }_0^{(1,\alpha )}{\sigma }^2\right]}^2},$$

$$S_0^{(\alpha )}[\rho ]=\frac{{\sigma }^2}{\pi }[{\left({\rho }_0^{(1,\alpha )}+{\rho }_0^{(0,\alpha )}\right)}^2$$

$$-\frac{1}{2}\sum _{j\ge 1}\frac{{\left({\rho }_j^{(1,\alpha )}+(-1{)}^j{\rho }_j^{(0,\alpha )}\right)}^2}{4j^2-1}]$$

$$p^{(\alpha )}\left(\{{\rho }_j^{(i,\alpha )}\},{\rho }_0\right)=p^{(\beta )}\left(\{{\rho }_j^{(i,\beta )}\},{\rho }_0\right),$$

$${\eta }_0^{(\mathrm{I}-\mathrm{N})}={\left\{1+\frac{2}{3\pi {\kappa }_0}({\kappa }_0-1{)}^2\left(1+{\Delta }_0^2\right)\right\}}^{-1},$$

$${\eta }_0^{(\mathrm{I}-\mathrm{T})}={\left\{1+\frac{2}{15\pi {\kappa }_0}\left[({\kappa }_0+1{)}^2+({\kappa }_0-1{)}^2{\Delta }_0^2\right]\right\}}^{-1}.$$

$${\kappa }_0^{\ast }=\frac{3+2{\Delta }_0^2+\sqrt{5+4{\Delta }_0^2}}{2(1+{\Delta }_0^2)}.$$

$${\eta }_0^{(\mathrm{N}-\mathrm{T})}=\{1+\frac{2}{3\pi {\kappa }_0}[({\kappa }_0-1{)}^2\left(1+{\Delta }_0^2\right)$$

$${+{\chi }_2\left({\eta }_0^{(N-T)}\right)]\}}^{-1},$$

$${\chi }_2({\eta }_0)=\int d\kappa f_0(\kappa )(\kappa -1{)}^2\frac{T_2^{(\mathrm{T})}(\kappa )}{T_0^{(\mathrm{T})}(\kappa )},$$

$${\eta }_j^{(0,\alpha )}={\rho }_j^{(0,\alpha )}{\kappa }_0{\sigma }^2,{\eta }_j^{(1,\alpha )}={\rho }_j^{(1,\alpha )}{\sigma }^2.$$

$$h^{(\alpha )}(\kappa ,\varphi )\equiv \frac{{\rho }^{(\alpha )}(\kappa ,\varphi )}{{\int }_0^{\pi }d{\varphi }^{\prime }{\rho }^{(\alpha )}(\kappa ,{\varphi }^{\prime })}.$$

$$h^{(\alpha )}(\varphi )={\left(\pi {\rho }_0^{(0,\alpha )}\right)}^{-1}\sum _j{\rho }_j^{(0,\alpha )}\cos (2j\varphi ).$$

$$Q_j^{(\alpha )}={\int }_0^{\pi }d\varphi h^{(\alpha )}(\varphi )\cos (2j\varphi ).$$

$$Q_j^{(\alpha )}=\frac{{\rho }_j^{(0,\alpha )}}{2{\rho }_0^{(0,\alpha )}}.$$

$$f^{(\mathrm{N})}(\kappa )=f_0(\kappa )e^{\beta {\mu }_{\mathrm{ex}}^{(\mathrm{I})}(\kappa )}{T}_0^{(\mathrm{N})}(\kappa ),$$

$$f^{(\mathrm{I})}(\kappa )=f_0(\kappa )\frac{e^{-\beta {\mu }_{\mathrm{ex}}^{(\mathrm{I})}(\kappa )}}{T_0^{(\mathrm{N})}(\kappa )},$$

$$f^{(\mathrm{N})}(\kappa )=f_0(\kappa )\frac{T_0^{(\mathrm{N})}(\kappa )}{T_0^{(\mathrm{T})}(\kappa )},$$

Full text

hard rectangles, polydispersity, phase transitions

I Introduction

The size polydispersity is an important feature present in experiments conducted on colloidal suspensions of anisotropic, rod or plate-like, particles Kooij1 ; Kooij2 ; Petukhov ; Woolston ; Vroege ; Pol ; Verhoeff ; Sun ; Byelov ; Leferik . The process of synthesis of mineral particles can not avoid the presence of polydispersity in sizes and only a sequence of fractionation steps can reduce it considerably. In experimental situations the number of species with different sizes is so large that we can consider a continuous size polydispersity as a reasonable approximation. However, the inclusion of a continuous size-polydispersity in theoretical models considerably complicates its numerical implementation. This constitutes the reason why the moments method was developed as a powerful theoretical tool to approximately solve the equations resulting from the two-phase coexistence calculations Clarke ; Sollich ; Speranza1 ; Speranza2 . The polydispersity is present not only in colloids but also in biological Fraden and granular Nguyen systems. The polymerization of F-actine filaments confined in quasi two-dimensional geometries constitutes a particular realization of a two-dimensional liquid-crystal system made of polydisperse flexible rods Viamontes .

Several experimental Kooij1 ; Kooij2 ; Petukhov ; Woolston ; Vroege ; Pol ; Verhoeff ; Sun ; Byelov ; Leferik and theoretical Sluckin ; Bates1 ; Bates2 ; M-R1 ; M-R2 ; Wensink ; Speranza3 ; Richter ; Green2 ; Pol2 ; M-R3 ; Green1 ; M-R4 ; M-R5 ; Escobedo2 works have studied the effect of polydispersity in particle sizes on the phase behavior and orientational ordering properties of colloidal suspensions made of anisotropic particles. Some general trends found in these studies concerning to the isotropic (I)-nematic (N) phase transition can be summarized as follows: (i) the polydispersity dramatically widens the region of I-N two-phase coexistence, i.e. the density gap between the cloud and shadow points is considerably enlarged making this transition strongly of first order. (ii) There exists strong fractionation in sizes between coexisting phases from which small and large particles usually populate the I and N phases respectively. The above mentioned experimental and theoretical works were conducted on three-dimensional systems, however the effect of polydispersity on two-dimensional liquid-crystals has been scarcely explored Ioffe ; Tavares ; Lopez ; Almarza ; Constantin . An experimental realization of a quasi 2D liquid-crystal system, apart from that already pointed out before, can be obtained by vertically shaking a monolayer of granular rods. Recent experiments on these monolayers showed the presence of liquid-crystal textures as stationary states Narayan1 ; Galanis1 ; Galanis2 ; Heras .

Theoretical studies on 2D monodisperse hard rods have shown that the transition from the high density N phase (without long range orientational order) to the low-density I phase is continuous, via a Kosterlitz-Thouless disclination unbinding type mechanism Kosterlitz , rather than first order Bates3 ; Lagomarsino . However, recent studies proved that when particle interactions are of a certain type the I-N transition becomes of first order in 2D Aernout ; Vink . A recent controversy resulted on the universality class (Ising vs. that corresponding to $q=1$ Pott s-type models) of the I-N transition of hard Zwanzig self-assembled rods on a lattice Lopez ; Almarza : A mean-field theory shows that these polydisperse self-assembled rods exhibit a first order I-N transition while MC results discard this fact Linares . Finally, a recent experiment on two-dimensional polydisperse hard rods (where magnetic nanorods were strongly confined between layers of a lamellar phase) have shown a first order I-N transition Constantin .

The particular shape of two-dimensional rods is also determinant to stabilize the N or tetratic (T) phases at low aspect ratios $\kappa$. In the latter the angular distribution function is invariant under $\pi/2$ rotations. Hard ellipses (HE) Cuesta1 ; Cuesta2 ; Xu ; Bautista and hard discorectangles (HDR) Bates3 only show the usual I-N transition followed by transitions to the plastic or orientationally ordered crystals. Hard rectangles (HR), on the other hand, also show fluid or crystalline phases with tetratic symmetry at low $\kappa_{0}$ Schlacken ; Frenkel ; M-R6 ; Donev ; M-R7 ; Triplett . The T phase of HR was theoretically predicted by mean-field DF studies Schlacken ; M-R6 . MC simulations on hard squares (HS) showed the presence of strong tetratic correlations of quasi-long-range order Frenkel and simulations on HR of $\kappa_{0}=2$ showed a liquid with the same T correlations and no N order Donev . The solid phase in the latter system was observed with a nonperiodic T ordering and having the structure of a random tiling of the square lattice with dimers of HR randomly oriented Donev . Experiments on a monolayer of disks standing on edges observed the conventional Kosterlitz-Thouless transition from I to N with almost smectic behavior at high density. But on the isotropic side of I-N transition an unusual regime of short-range T correlations dominates over N ones Chaikin . Recent experiments on monolayers of hard microscale square platelets showed a phase transition between an hexagonal rotator crystal and a rhombic crystal as the packing fraction is increased Zhao . The absence of T ordering in this system was further explained resorting to the roundness of the square corners Escobedo . If the roundness is sufficiently small the particles behave like perfect squares and the T phase is recovered Escobedo . Shaken monolayers of granular cubes Narayan2 and cylinders Narayan1 ; Heras ; G-P exhibit stationary textures with strong T ordering even for aspect ratios of cylinders as large as 7.

All these studies show the profound effect that particle shapes and pair interactions have on the symmetry of the orientationally ordered phases. Even the character (continuous vs. first order) can be influenced by them. The mean-field SPT for HR predicts a first-order I-N (or T-N) transition located between two tricritical points. This is a peculiar feature of HR because the I-N transition for HDR and HE following the same theory is always of second order M-R6 . This difference can be explained by the peculiar form of the excluded volume between two HR as it will be discussed in Sec. II.2.

The main purpose of the present work is to study how the polydispersity affects the phase behavior of HR. To this purpose we will extend the mean-field SPT formalism from its multicomponent version (already described in Ref. M-R8 ; Heras2 ) to the length-polydisperse case. The minor length of rectangles is considered constant while the major one is polydisperse according to a Schulz-type (with exponential decay) probability density distribution or also according to a distribution with a Gaussian decay. We calculated the phase diagrams for a fixed $\Delta_{0}$ while $\kappa_{0}$ varies, and also for certain fixed values of $\kappa_{0}$ while $\Delta_{0}$ changes. We measured the degrees of orientational ordering of particles and their fractionation between the coexisting phases. The addition of polydispersity dramatically increases the interval of aspect ratios where the I-N is of first order and the I-N tricritical point moves to higher values of $\kappa_{0}$. Also we prove that the stability of the T phase is not severely affected by polydispersity: the area of the phase diagram where the T phase is stable has not a severe decrement when the polydispersity is increased from zero to its maximum allowed value.

II Model and Theory

In this section we present the theoretical tools used to study the phase behavior of a length-polydisperse HR fluid. In Sec. II.1 we formulate the HR model and define the family of length-polydisperse probability distributions we have used. In Sec. II.2 we generalize the SPT formalism from the multicomponent HR mixture M-R8 to a continuous polydisperse fluid. Further we explicitly obtain the set of equations used in Sec. II.3 to calculate the two-phase coexistence giving also details on their numerical implementation. A concise summary of the I-N, I-T and T-N bifurcation analysis will be given in Sec. II.4. Finally, in Sec. II.5, we define the main variables and functions used to measure the orientational ordering and fractionation between coexisting phases.

II.1 Polydisperse HR model

Our model consists on a collection of freely-rotating HR. They move and rotate in 2D and cannot overlap. We take the minor edge-length of rectangles, $\sigma$, to be constant while its mayor length, $L$ (with $\sigma<L<\infty$), defining the main particle axis, is considered continuously polydisperse according to a fixed, so called parent, probability density distribution function:

[TABLE]

defined in terms of the particle aspect ratio $\kappa\equiv L/\sigma$. This function is normalized to unity, $\int_{1}^{\infty}d\kappa f_{0}(\kappa)=1$, and its first moment, $\langle\kappa\rangle\equiv\int_{1}^{\infty}d\kappa\kappa f_{0}(\kappa)=\kappa_{0}$, is selected to be equal to $\kappa_{0}$ (the mean aspect ratio). The polydisperse coefficient, or mean square deviation, can be calculated as

[TABLE]

where $\Gamma(x)$ is the standard Gamma function. The parameters $\nu$ and $q$, used to define the family of functions (1), control the width of $f_{0}(\kappa)$, with $q$ being the one that dictates the decay of the tail: for $q=1$ (corresponding to the Schulz distribution function) we obtain an exponential decay while for $q=2$ a Gaussian tail is obtained. In the following we will use $\Delta_{0}$ as a measure of the degree of length-polydispersity.

II.2 SPT for length-polydisperse HR

We use the SPT formalism for a multicomponent HR fluid M-R8 extended to the continuously length-polydisperse limit: $\rho_{\nu}(\phi)\to\rho(\kappa,\phi)$, where $\rho_{\nu}(\phi)$ corresponds to the angular density profile of species $\nu$ of the multicomponent mixture while $\rho(\kappa,\phi)$ denotes the continuously distributed angular density profile of species with aspect ratio $\kappa$. All sums over species in equations used to define the free-energy density in Ref. M-R8 are now substituted by integrals over the aspect ratio: $\sum_{\nu}\to\int d\kappa$. The resulting free-energy density in reduced thermal units, $\Phi[\rho]\equiv\Phi_{\rm id}[\rho]+\Phi_{\rm ex}[\rho]$, splitted in its ideal and excess contributions, $\Phi_{\rm id}[\rho]$ and $\Phi_{\rm ex}[\rho]$ respectively, is a functional of $\rho(\kappa,\phi)$:

[TABLE]

In (6) and (7), $\beta\equiv 1/k_{B}T$, is the Boltzmann factor, ${\cal F}_{\rm id,ex}[\rho]$ correspond to the ideal and excess parts of the free-energy density functional while $A$ is the total area of the system. We have also defined the $i$th moment of the integrated, over $\kappa$ and $\phi$, density profile:

[TABLE]

where $i=1,2$ and the integration with respect to $\phi$ was taken between 0 and $\pi$ (instead of 2$\pi$) due to the head-tail symmetry of rectangles. Note that these moments can be obtained from the $i$th moment angular profile $\rho^{(i)}(\phi)$ as defined in Eq. (10). The magnitude $A_{0}(\kappa_{1},\kappa_{2},\phi)$ is directly related to the excluded area, $A_{\rm excl}(\kappa_{1},\kappa_{2},\phi)$, between two HR of aspect ratios $\kappa_{1}$ and $\kappa_{2}$ with a relative angle between their main axes equal to $\phi$. The relation between both magnitudes is $2A_{0}(\kappa_{1},\kappa_{2},\phi)=A_{\rm excl}(\kappa_{1},\kappa_{2},\phi)-\left(\kappa_{1}+\kappa_{2}\right)\sigma^{2}$. In Fig. 1 (a) we schematically show the excluded volume between two HR of $\kappa=2$ while in (b) we show the function $A_{\rm excl}(2,2,\phi)$. The secondary minimum at $\phi=\pi/2$ is responsible for the stability of the T phase at low aspect ratios which is a peculiar property almost unique to the rectangular shape. The most common functional form of the excluded area is exemplified for hard discorectangles (HDR) in (b) for comparison (note the presence of a maximum at $\pi/2$ instead of a local minimum).

II.3 Coexistence calculations and numerical schemes

We now proceed to present the set of equations we used to calculate the phase coexistence between two phases with different orientational symmetries (I, N or T phases). The total free-energy density of a phase-separated system with a fraction $\gamma_{\alpha}$ of the total volume occupied by phase $\alpha$ ($\alpha=\rm I,N,T$) can be expressed as $\Phi^{(\rm t)}[\rho]=\sum_{\alpha}\gamma_{\alpha}\Phi^{(\alpha)}[\rho]$ (obviously $\sum_{\alpha}\gamma_{\alpha}=1$).

From now on the coexistence between (I,T) and N phases occupying fractions of the total volume $\gamma_{I,T}=\gamma_{0}$ and $\gamma_{\rm N}=1-\gamma_{0}$ respectively will be denoted as (I,T)${}_{\gamma_{0}}$-N${}_{1-\gamma_{0}}$. For example I0-N1 coexistence implies that the I phase occupies a vanishing small volume (the so called shadow phase) while the N phase takes up the total volume (the so called cloud phase).

The lever rule guarantees the conservation of the total number of species of aspect ratio $\kappa$:

[TABLE]

where $\rho_{0}=N/A$ (with $N$ the total number of particles) is the total number density of the system. Using (11) as a constraint in the minimization of the total free-energy density $\Phi^{(\rm t)}[\rho]$ with respect to $\rho^{(\alpha)}(\kappa,\phi)$ gives us the following set of equations:

[TABLE]

where $\alpha\in\{\rm I,T,N\}$ and the excess chemical potential of the coexisting $\alpha$-phase is defined as $\displaystyle{\beta\mu_{\rm ex}^{(\alpha)}(\kappa,\phi)=\frac{\delta\Phi^{(\alpha)}_{\rm ex}[\rho]}{\delta\rho^{(\alpha)}(\kappa,\phi)}}$. From now on we will use Greek letters $\alpha,\tau,\dots$ as superscripts to label coexisting phases while the Latin letters $i,j,\dots$ will be used to label the moments.

The numerical procedure used to solve Eqs. (12) is based on the reduction of this infinite set of equations (note that $\kappa$ and $\phi$ are continuous variables) in a finite set. To this purpose we firstly introduce a truncated Fourier expansion of the moment profiles:

[TABLE]

with $\{\rho_{j}^{(i,\alpha)}\}$ the Fourier amplitudes. Secondly we substitute this expansion into $\beta\mu_{\rm ex}^{(\alpha)}(\kappa,\phi)$, multiply (12) by $\kappa^{i}\cos(2j\phi)$ and integrate over $\kappa$ and $\phi$ to obtain a closed set of equations:

[TABLE]

where $j=1,\dots,N,\ i=0,1$ and $D_{j0}=2/(1+\delta_{j0})$ with $\delta_{j0}$ the Kronecker delta. The success of the present strategy is based on the dependence of the excess part of the chemical potential on the angular density profile $\rho(\kappa,\phi)$ only through its moments $\{\rho_{j}^{(i,\alpha)}\}$:

[TABLE]

where $\beta p^{(\alpha)}$, defined in Eq. (18), is the pressure of phase $\alpha$. Thus we finally obtain a closed set of equations (14) to be solved for the Fourier amplitudes $\{\rho_{j}^{(i,\alpha)}\}$. This set in turn guarantees the equality between chemical potentials of each species in both coexisting phases. An equation corresponding to the equality between pressures at each phase should be added to get the mechanical equilibrium, and consequently the pressure in Eq. (16) can be dropped during the numerical implementation of (14). A fixed point iteration method was used to solve (14). To perform the numerical integrals over $\phi$ and $\kappa$ in Eqs. (14) and (15) we used, after some convenient changes of variables, a Gauss-Legendre and Gauss-Laguerre quadratures respectively. Once the numerical convergence of Eqs. (14) was reached for a fixed $\rho_{0}$, we used the equality

[TABLE]

numerically solved using the Brent’s zero-finding method, to find the value of $\rho_{0}$ at equilibrium.

II.4 Bifurcation analysis

To find the onset of N and T ordering from the I phase we performed a bifurcation analysis of Eqs. (14) particularized for $\gamma_{\rm I}=1$ with respect to the small Fourier amplitudes $\{\rho^{(i,\alpha)}_{j}\}$ ($j=1$ for $\alpha=$N while $j=2$ for $\alpha=$T). The resulting packing fractions, $\eta_{0}\equiv\rho_{0}\kappa_{0}\sigma^{2}$, at which the I phase is destabilized with respect to N and T phases (the so called spinodal instability) result in the following analytical expressions:

[TABLE]

For small (large) mean aspect ratios $\kappa_{0}$ it is expected that the I phase destabilizes firstly with respect to the T (N) phase. Consequently there should exist a crossing point $\kappa_{0}^{*}$ between both, I-T and I-N, spinodals. This point can be calculated from the equality $\eta_{0}^{(\rm I-N)}=\eta_{0}^{(\rm I-T)}$ which results in

[TABLE]

From this equation we obtain, for the one-component fluid ($\Delta_{0}=0$), the value $\kappa_{0}^{*}=1+\varphi$ with $\varphi=(1+\sqrt{5})/2\approx 1.618$ the golden ratio. For the highest polydispersity, $\Delta_{0}=1$ (reached for the Schulz distribution), we obtain $\kappa_{0}^{*}=2$. Thus we can extract as a preliminary conclusion that the stability of T phase slightly decreases with polydispersity.

For small aspect ratios, $\kappa_{0}\sim 1$, it is also expected an stable T phase up to the density at which a T-N transition takes place. To calculate the T-N spinodal we have performed a bifurcation analysis of Eqs. (14), particularized for $\gamma_{\rm T}=1$, with respect to the small Fourier components $\{\rho_{1}^{(i,\rm N)}\}$. The obtained result for the packing fraction at the N-T bifurcation can be calculated from

[TABLE]

where we have defined

[TABLE]

with the functions $T_{j}^{(T)}(\kappa)$, already defined in (15), calculated once the Fourier amplitudes $\{\rho_{2j}^{(i,\rm T)}\}$ at the equilibrium T phase were found for a given $\eta_{0}$. We have solved iteratively the nonlinear integral equation (24) with respect to $\eta_{0}^{(\rm N-T)}$ to find the N-T spinodal instability.

We plot the three found spinodals $\eta_{0}^{(\rm I-N)}(\kappa_{0})$, $\eta_{0}^{(\rm I-T)}(\kappa_{0})$ and $\eta_{0}^{(\rm N-T)}(\kappa_{0})$ in Fig. 2 for $\Delta_{0}=0$ and 1. We can see that the I-N spinodal is severely affected by polydispersity (is significantly below its one-component counterpart). Another conclusion we can extract from the figure is that the stability of T phase only slightly decreases with polydispersity as we have already point out before by comparing the values of the crossing points $\kappa_{0}^{*}$. This destabilization can be explained by the presence of long rods (those which destroy the T symmetry) with aspect ratios much larger than $\kappa_{0}$. The long rods have a profound effect on the orientational ordering properties of the fluid, usually favoring the N ordering. However the decrease in the T phase stability is not so strong as expected because the mean-square deviation, $\Delta$, of the function $f_{0}(\kappa)$ decreases with the mean aspect ratio $\kappa_{0}$ [see Eq. (4) ]: when $\kappa_{0}\to 1$ we obtain $\Delta\to 0$ for any $\Delta_{0}$.

II.5 Measuring the orientational ordering and fractionation effects

To present the results in Sec. III it is convenient to define, apart from the already introduced packing fraction $\eta_{0}$, the following dimension-less variables:

[TABLE]

To measure the orientational ordering of rectangles of aspect ratio $\kappa$ at the coexisting $\alpha$-phase we use the angular probability distribution function:

[TABLE]

The integration of $h^{(\alpha)}(\kappa,\phi)$ over $\kappa$ gives us the averaged angular distribution function which, resorting to the Fourier expansion (13), can be easily calculated using the expression:

[TABLE]

The N ($j=1$) and and T ($j=2$) order parameters are defined through

[TABLE]

Inserting (28) into (29) we obtain

[TABLE]

To measure the fractionation between the cloud and shadow coexisting phases we should bear in mind that the cloud phase always has the parent distribution function $f_{0}(\kappa)$. Thus we will concentrate only on the shadow-phase distribution functions. For the I1-N0 and I0-N1 coexistences the distributions corresponding to the N and I shadow phases are

[TABLE]

while those corresponding to N and T shadow phases of the T1-N0 and T0-N1 coexistences are

[TABLE]

All these distributions will be normalized as $\tilde{f}^{(\alpha)}(\kappa)\equiv f^{(\alpha)}(\kappa)/\int d\kappa f^{(\alpha)}(\kappa)$ for plotting. We will use the averaged over $\tilde{f}^{(\alpha)}(\kappa)$ aspect ratio,

[TABLE]

to quantify the fractionation. Note that it is equal to $\kappa_{0}$ when averaged with respect to $f_{0}(\kappa)$.

III Results

In this section we show the results as obtained from the numerical implementation of SPT-DF for length-polydisperse HR. We divide this section in different subsections each one devoted to discuss the thermodynamic and orientational properties of polydisperse HR. In Sec. III.1 we firstly describe the resulting phase diagrams using two different parent distribution functions (those with exponential and Gaussian-like decays). Further we proceed, in Sec. III.2, to present a quantitative analysis of the orientational ordering of particles at coexistence. To this purpose we use the angular distribution functions and order parameters. The fractionation of species with different aspect ratios between coexisting phases as a function of the mean aspect ratio and polydispersity is studied in Sec. III.3 and III.4 respectively. Finally, in Sec. III.5 we study the effect of polydispersity on the location of the I-N tricritical point.

III.1 Phase diagrams for fixed polydispersity $\Delta_{0}$

We numerically solved the set of equations (14) together with the mechanical equilibrium condition (20) selecting a parent distribution function (1) with $q=1$ and a polydisperse coefficient $\Delta_{0}=0.408$, corresponding to the value $\nu=5$ in Eq. (1). For small mean aspect ratios, $\kappa_{0}\sim 1$, we calculated the T-N coexistence from Eqs. (14) and (20) fixing the odd Fourier coefficients $\rho_{2j+1}^{(i,\rm T)}$ to zero. This simplification is justified by the symmetry of the T phase: $\rho^{(\rm T)}(\kappa,\phi)=\rho^{(\rm T)}(\kappa,\phi+\pi/2)$. For medium and large aspect ratios we solved the same set of equations taking into account all the Fourier coefficients $\rho_{j}^{(i,\rm N)}$ to look for the I-N coexistence. We have found that the I-T transition is always of second order. The total packing fractions $\eta_{0}$ of the I1 and N1 (coexisting with their shadow counterparts) phases are shown in Fig. 3. Together with these coexisting binodals we also plot the I-N, I-T and N-T spinodals as calculated in Sec. II.4. For comparison we also plot in (a) the phase diagram corresponding to the one-component HR fluid ($\Delta_{0}=0$) calculated in Ref. M-R6 . The most prominent effects of polydispersity consist of: (i) the dramatic broadening of the I-N coexisting region, (ii) the displacement of the I-N tricritical point to higher values of $\kappa_{0}$, and (iii) the lowering of the second-order I-N transition density (the spinodal curve). The widening of the two-phase coexistence region is a well known fact in many studies conducted on three-dimensional polydisperse hard-rods Speranza1 ; Speranza2 ; Speranza3 ; Wensink . As it is well known the I-N transition of hard rods in 3D is of first order for any aspect ratio. In 2D the mean-field DF predicts a second order I-N transition for most particle shapes, HR being an important exception: The I-N and T-N transitions are of first order for aspect ratios between two (N-T and I-N) tricritical points. The broadening of the I-N coexisting region found here is related, in analogy to 3D, to a demixing-like mechanism: the N and I phases are preferentially populated by long and small rods respectively (the so called fractionation effect). We can see from Fig. 3 that the T phase slightly destabilizes with polydispersity as we have already discussed in Sec. II.4.

In Fig. 3 (b) we show a zoom of the phase diagram around the left T-N tricritical point (solid circle). Note that the T1-N0 binodal (the lower dot-dashed curve) ends up at an aspect ratio value $\kappa_{0}$ (lower square) slightly below that corresponding to the end of the T0-N1 binodal (upper square). This result suggests that the system could exhibit a triple I-T-N coexistence for aspect ratios located between the two vertical dotted lines.

The phase diagram corresponding to the selected parent distribution function with $q=2$ (that with a Gaussian decay) is shown in Fig. 4 for the same $\Delta_{0}=0.408$. For comparison we show in the same figure the case $q=1$. We observe that the cloud binodals corresponding to the I1-N0 coexistences for $q=1$ and 2 coincide. On the other hand the cloud binodal of the of the I0-N1 coexistence for $q=2$ is located, in particular for large aspect ratios, below that of $q=1$, reducing in this way the coexisting gap and consequently moving the I-N tricritical point from $\kappa_{0}^{(\rm t)}\approx 7.9$ to $\kappa_{0}^{(\rm t)}\approx 7.05$. We can conclude from this result that the position of the tricritical point depends not only of the first two moments of the parent distribution function (both are the same for $q=1$ and 2), but also on higher moments which are certainly different for the Schulz and Gaussian-like distributions. Through a bifurcation analysis of the free-energy with respect to the Fourier components of an incipient N phase performed around the tricritical point we can obtain an analytic equation relating the position of this point with these moments. This analysis was done for binary mixtures of HR in Ref. M-R8 and we plan to extend it to the polydisperse case in a future. The reduction of the coexisting gap with $q$ reflects the lesser relevance of long rods on the phase behavior for distributions with strong decay.

III.2 Orientational ordering

In this section we discuss the orientational ordering of rectangles at equilibrium. To this purpose we firstly show the coexisting angular distribution functions $h^{(\alpha)}(\phi)$ as calculated from Eq. (28) corresponding to an stable I-N phase separation for $\kappa_{0}=3.2$, $q=1$ and $\Delta_{0}=0.408$. The cloud and shadow N distributions are plotted in Fig. 5 (a) with solid and dashed lines respectively. The sharp peaks at 0 and $\pi$ show the strong N ordering present in the fluid, while a rather small oscillation around $\pi/2$ indicates the closeness of the system to the region of T phase stability ($1\leq\kappa_{0}\lesssim 2$). Also, the cloud and shadow N distributions are very similar, with the latter showing an slightly higher ordering of particles. We show these distributions for the same values of $q$ and $\Delta_{0}$ but this time selecting $\kappa_{0}=1.9$ where the T-N coexistence is stable. The results are shown in Fig. 5 (b) and (c) for the N and T phases respectively. Now the N phase develops a well defined central peak at $\pi/2$ reflecting a high proportion of small rectangles with T-like ordering. Again the shadow N phase exhibits slightly higher degree of N ordering and consequently a lesser proportion of small rectangles with T-like configurations. The coexisting T distributions [see panel (c)] has a periodicity of $\pi/2$ (the central peak is of the same height as those located at 0 and $\pi$) reflecting the T symmetry: the system is invariant under rotations of $\pi/2$. The cloud T phase has a higher degree of ordering than its shadow counterpart.

A better understanding of the degree of particle ordering at I-N coexistence can be reached through the calculation of the angular distribution function of species of a certain aspect ratio $\kappa$, $h^{(\alpha)}(\kappa,\phi)$, as defined by Eq. (27). It is expected that the N (T) ordering increases (decreases) with $\kappa$, an assertion well supported by the N distributions shown in Fig. 6 and calculated for the parameters $q=1$, $\Delta_{0}=0.408$, $\kappa_{0}=3$, and for three different values of $\kappa$: 1.5 (dotted), 3 (dashed) and 5 (solid). The small species, being the less ordered, present a high proportion of T-like configurations while the large ones exhibit a high degree of uniaxial N ordering. Again the cloud distributions [panel (a)] reflect a higher ordering than their shadow counterparts [panel (b)], a fact already pointed out before when we discussed the behavior of $h^{(\alpha)}(\phi)$ for $\kappa_{0}=3.2$.

The global ordering along the coexisting binodal and spinodal curves of the phase diagram, shown in Fig. 3, is calculated through the order parameters $Q_{1,2}^{(\alpha)}$ [see Eq. (30)]. These are plotted in Fig. 7 as functions of $\kappa_{0}$. For large values of $\kappa_{0}$, not close enough to the I-N tricritical point, $\kappa_{0}^{t_{2}}$, the N order parameter $Q_{1}^{(\rm N)}$ along the I-N coexistence is relatively high and always above the T order parameter $Q_{2}^{(\rm N)}$. Obviously at the I-N tricritical point they both become zero. As $\kappa_{0}$ decreases from $\kappa_{0}^{t_{2}}$ the T and N order parameters increase, both reaching their maxima at slightly different values close to $\kappa_{0}=3$. Just below $\kappa_{0}=3$ a cross-over between $Q_{1}^{(\rm N)}$ and $Q_{2}^{(\rm N)}$ takes place. A slight decrease of $\kappa_{0}$ from this cross-over gives $Q_{2}^{(\rm N)}>Q_{1}^{(\rm N)}$, a clear sign of the presence of a secondary peak at $\pi/2$ in the angular distribution function. Then the I-N transition is substituted by the T-N one and $Q_{1}^{(\rm N)}$ decreases to zero at the T-N tricritical point $\kappa_{0}^{t_{1}}$. The order parameter $Q_{2}^{(\rm N)}$ also decreases in a lesser extent up to the value shown with a grey circle (the T-N tricritical point). From this point $Q_{2}^{(\rm T)}$ increases along the T-N spinodal up to its maximum allowed value 1 as $\kappa_{0}\to 1$ while it decreases along the T-N coexistence binodal as $\kappa_{0}$ is increased up to $\kappa_{0}^{c}$, the value corresponding to the end-critical point, where the I-T spinodal and the cloud binodal of the T1-N0 coexistence meet. The difference, $\Delta\kappa\equiv\kappa_{0}^{c}-\kappa_{0}^{t_{1}}$, between the aspect ratios corresponding to those of the end-critical and tricritical points [see the square symbols in (a) and (b)] depends on which of the phases is considered to be the cloud one. When $\gamma_{\rm N}=1$ [panel (a)], i.e. for the (I,T)0-N1 coexistence, this difference is clearly larger than that corresponding to $\gamma_{\rm N}=0$, the (I,T)1-N0 coexistence [panel (b)]. This behaviour is related to the already discussed fact on the possible existence of a triple I-N-T coexistence [see Fig. 3 (b)] for values of mean aspect ratios located between these points.

III.3 Fractionation as a function of $\kappa_{0}$

This section is devoted to explore in detail the fractionation of different species between the coexisting phases. To this purpose we have fixed the values of the mean aspect ratio $\kappa_{0}$ to one of the values of the following set: $\{3,5,7\}$. Also we have selected the polydisperse coefficient of the Schulz-type ($q=1$) distribution function to be $\Delta_{0}=0.408$ as before. Further, we numerically solved the Eqs. (14) and (20) for chemical and mechanical equilibrium between cloud-shadow coexisting phases. As a result we found the normalized length-distributions, $\tilde{f}^{(\alpha)}(\kappa)$, obtained from (32)-(34) for the shadow I, T or N phases. We should bear in mind that the cloud phases always follow the distribution $f_{0}(\kappa)$. All these distributions are plotted in Fig. 8. We can appreciate a clear fractionation of small/long species between cloud-shadow phases. The small ones preferentially populate the I or T phases while the N phase is rich in long species: see how the maxima of the I/T distributions are clearly located left of those of the N distributions. Also the latter has a slower decay, indicating a higher fraction of long species. As we can see from the figure the fractionation is much more dramatic for $\kappa_{0}=3$ which it is clearly correlated with the value at which the N phase has the highest order parameter $Q_{1}^{(\rm N)}$ (see Fig. 7). As we move away from this value in the directions of both tricritical points $\kappa_{0}^{t_{i}}$ the system exhibits fractionation but in a lesser extent due to the weaker character of the (I,T)-N phase transition.

Another important consequence of the fractionation is related to the values of the unit-less moments of the distribution functions, $\eta_{0}^{(i,\alpha)}$, at coexistence. In Fig. 9 we plot these moments as a function of $\gamma_{\rm I}$ for $q=1$, $\kappa_{0}=3$ and $\Delta_{0}=0.408$. The inequalities $\eta_{0}^{(0,\rm I)}>\eta_{0}^{(0,\rm N)}$ and $\eta_{0}^{(1,\rm N)}>\eta_{0}^{(1,\rm I)}$ are always fulfilled which constitutes a direct consequence of how different are the shapes of $\tilde{f}^{(\rm I)}(\kappa)$ and $\tilde{f}^{(\rm N)}(\kappa)$ as discussed previously. Integrating the sharpest, strongly decayed function $\kappa^{i}\tilde{f}^{(\rm I)}(\kappa)$, over $\kappa$ for $i=0$ ($i=1$) gives us a value of $\eta_{0}^{(0,\rm I)}$ ($\eta_{0}^{(1,\rm I)}$) greater (less) than that obtained from the integration of the slower decaying function $\kappa^{i}\tilde{f}^{(\rm N)}(\kappa)$.

To finish this section we present in Fig. 10 how $\langle\kappa\rangle_{\tilde{f}^{(\alpha)}(\kappa)}$ (the averaged aspect ratio with respect to the coexisting shadow distribution $\tilde{f}^{(\alpha)}(\kappa)$) evolves as a function of $\kappa_{0}$ for $\Delta_{0}=0.408$. As expected, this average, which directly measures the fractionation, becomes equal to $\kappa_{0}$ at both tricritical points (solid circles) while it reaches its maximum (minimum) value at $\kappa_{0}\sim 3$ when averaged with respect to N0 (I0) distributions. This qualitative behavior is similar for both, the Schulz ($q=1$) and Gaussian-tailed ($q=2$) distributions. However the latter provokes less fractionation.

III.4 Fractionation as a function of $\Delta_{0}$

To finish the discussion on fractionation we present in Fig. 11 the results regarding the behaviour of the dimension-less zeroth moments, $\eta_{0}^{(0,\alpha)}$, of coexisting cloud and shadow phases as a function of polydispersity for fixed $\kappa_{0}=3$ [(a)] and 7 [(b)] and selecting $f_{0}(\kappa)$ to be of Schulz type ($q=1$). We observe a dramatic widening of the coexistence region, i.e. the gap between the moments corresponding to I1 and N1 phases (solid lines) greatly increases with $\Delta_{0}$ (for $\kappa_{0}=7$ the same occurs but beyond the tricritical point). Also, for $\kappa_{0}=3$ and zero polydispersity the moments of the shadow I and N phases (dashed lines) obviously coincide with those of the cloud phases and they compare as usually: $\eta_{0}^{(0,\rm I)}<\eta_{0}^{(0,\rm N)}$. However as $\Delta_{0}$ increases they exhibits a cross-over at $\Delta_{0}\approx 0.1$ and the former relation inverts: $\eta_{0}^{(0,\rm I)}>\eta_{0}^{(0,\rm N)}$. This in turns means that the I shadow distribution function $f^{(\rm I)}(\kappa)$ (that which it is not normalized) exhibits, as a result of strong fractionation, a sharper peak located at relatively small values of $\kappa$. On the other hand the N shadow distribution, $f^{(\rm N)}(\kappa)$, shows a more smeared peak located at larger values of $\kappa$ (see Fig. 8). These differences in the shapes of distributions result in $\int d\kappa f^{(\rm I)}(\kappa)>\int d\kappa f^{(\rm N)}(\kappa)$. Interestingly this cross-over does not exist for $\kappa_{0}=7$: At the tricritical point the moments are equal while beyond it the fractionation mechanism always gives $\eta_{0}^{(0,\rm I)}>\eta_{0}^{(0,\rm N)}$. A more clear measure of fractionation, given by $\langle\kappa\rangle_{\tilde{f}^{\alpha}}$, is plotted as a function of $\Delta_{0}$ along the same shadow curves in Fig. 12. For $\kappa_{0}=3$ (7) the aspect ratio, averaged with respect to the shadow N (I) distribution, can reach values as large (low) as $2\kappa_{0}$ ($\kappa_{0}/2$) for high enough $\Delta_{0}$, a clear sign of the presence of strong fractionation.

III.5 Tricritical points

The last study we have carried out concerns the location of the I-N tricritical point as a function of polydispersity. To this purpose we fixed the value of $\kappa_{0}$ and calculated the I0-N1 coexistence for high enough values of $\Delta_{0}$ [selected in such a way to guarantee a first order I-N transition, and thus to find a nontrivial numerical solution of Eqs. (14)]. Then we gradually decrease the value of $\Delta_{0}$ and used as new initial guesses those found in the previously converged iterations. This process is repeated up to that value of $\Delta_{0}$ for which the order parameter $Q_{1}^{(\rm N)}$ is negligibly small. Finally an accurate extrapolation (using a cubic-spline fitting) of $Q_{1}^{(\rm N)}$ to zero gives us the value $\Delta_{0}^{*}$ at which the first order transition becomes of second order, i.e. the position of the tricritical point. This process was carried out for a set of aspect ratios and the results are shown in Fig. 13. The open circles show our selected values of $\kappa_{0}$ and the curve is an Akima spline fitting to guide the eyes. This curve has certain credibility only as an interpolation fitting. However we have decided to use the same fitting to conjecture the value of $\Delta_{0}^{*}$ in the limit $\kappa_{0}\to\infty$. This result might indicate that there exists a terminal polydispersity, around the value of $0.7$, beyond which the I-N transition becomes of first order for any $\kappa_{0}$. We will certainly solve this conjecture by performing a bifurcation analysis around the tricritical point in such a way to find an analytic expression for $\Delta_{0}^{*}$ as a function not only of $\kappa_{0}$ but also of higher moments of the distribution $f_{0}(\kappa_{0})$. We will carry out this study in the future and we are planning to generalize it for any particle shape (polydisperse hard discorectangles, hard ellipses, etc.). In such a way we will be able to shed some light on the character (second vs. first) of the I-N transition that mean-field theories predict for different polydisperse particle shapes.

IV Conclusions

We have extended the SPT of HR fluid from its multicomponent version to consider a continuous length-polydispersity. Using this formalism we derived a simplified coexistence equations by using the Fourier expansions of the angular moment profiles. The numerical solution of the obtained set of equations for the Fourier components and for different values of the mean aspect ratio and polydispersity allowed us to obtain the phase diagram of the system. We have not taken any functional ansatz to parameterize the density profile. Thus our calculations are exact up to the errors associated with the truncation of the Fourier series. The main results can be summarized as follows: (i) The I-N transition becomes stronger with polydispersity, i.e. the coexistence gap becomes wider and the tricritical point moves to higher values of $\kappa_{0}$. (ii) The stability of the T phase slightly decreases with $\Delta_{0}$. (iii) There exists strong fractionation between the coexisting phases which becomes stronger at $\kappa_{0}\sim 3$ and decreases as we approximate the T-N or I-N tricritical points. The fractionation for a Schulz-type parent distribution function is stronger than that obtained for a Gaussian-tailed distribution. As usually the N phase is rich in long particles while the I or T phases are highly populated by the smaller particles. (iv) The orientational ordering of coexisting species of different $\kappa$ is in general different: the small species, those with $\kappa\gtrsim 1$, have a low orientational ordering and tend to align in T-like configurations while long species contribute to stabilize the uniaxial N. This effect is more pronounced for $\kappa_{0}\sim 3$. (v) The locations of tricritical points are severely affected by polydispersity and our results might suggest the existence of a terminal polydispersity beyond which the I-N transition becomes of first order. We require further calculations to definitively settle out this point.

A recent experiment on magnetic polydisperse nanorods which are strongly confined between lamellar layers (an experimental realization of a quasi-2D hard-rod fluid) shows a first order I-N transition with the N director being parallel to the plane of the layers Constantin . As we show in the present work when the polydispersity is high enough the character of the I-N transition could change from second to first order. Thus the polydispersity could explain the experimental results obtained in Ref. Constantin .

We have not considered the presence of non-uniform phases (like smectic or crystal) in the present study. For low polydispersity they become stable at high densities. However when polydispersity increases they rapidly destabilize so the results shown here, specially those obtained for high enough $\Delta_{0}$, could be qualitatively confirmed by experimental realizations of a 2D hard-rectangular fluid. One of them could be mechanically vibrated monolayers of granular cylinders. Recent experiments on such systems have confirmed the presence of strong T, N Narayan1 ; Heras and smectic G-P ordering of rods, which strongly depends on the aspect ratio and packing fraction. Our group will conduct experiments on shaken monolayers of granular rods to study the effect of length-polydispersity on the stability of those liquid-crystal textures already found in the one-component case G-P .

Acknowledgements.

Financial support from MINECO (Spain) under grant FIS2015-66523-P is acknowledged.

Appendix A The I1-N0 coexistence equations

For $q=1$ and $\gamma_{\rm I}=1$ the integral with respect to $\kappa$ in the set of equations (14) corresponding to the I-N coexistence calculations can be performed analytically resulting in the following set of equations for the variables $\eta_{j}^{(i,\rm N)}$:

[TABLE]

where we have defined the functions

[TABLE]

with $s_{j}^{(\rm N)}\equiv\eta_{j}^{(1,\rm N)}+(-1)^{j}\kappa_{0}^{-1}\eta_{j}^{(0,\rm N)}$, and $g_{j}\equiv\left(4j^{2}-1\right)^{-1}$. For the cloud I coexisting phase we have $\eta_{0}^{(i,\rm I)}=\eta_{0}$.

For $q=2$ and $\gamma_{\rm I}=1$ we obtain

[TABLE]

where the function $R_{0}(\phi)$ is the same as (37), while

[TABLE]

$D_{\mu}(x)$ is the parabolic cylinder function. We used a Fortran subroutine provided in Ref. parabolic to numerically evaluate $D_{\mu}(x)$.