Frank elasticity of composite colloidal nematics with anti-nematic order H. H. Wensink 1

1 Laboratoire de Physique des Solides, CNRS, Université Paris-Sud, Université Paris-Saclay, 91405 Orsay, France * (Dated: October 17, 2018)

Mixing colloid shapes with distinctly different anisotropy generates composite nematics in which the order of the individual components can be fundamentally different. In colloidal rod-disk mixtures or hybrid nematics composed of anisotropic colloids immersed in a thermotropic liquid crystal, one of the components may adopt so-called anti-nematic order while the other exhibits conventional nematic alignment. Focussing on simple models for hard rods and disks, we employ Onsager-Straley's second-virial theory to derive scaling expressions for the elastic moduli of rods and disks in both nematic and anti-nematic configurations and identify their explicit dependence on particle concentration and shape. We demonstrate that the splay, bend and twist elasticity of anti-nematically ordered particles scale logarithmically with the degree of anti-nematic order, with the bend-splay ratio for anti-nematic discotic nematics being far greater than for conventional nematic systems. The impact of surface anchoring on the elastic properties of hybrid nematics will also be discussed in detail. We further demonstrate that the elasticity of mixed uniaxial rod-disk nematics depends exquisitely on the shape of the components and we provide simple scaling expressions that could help engineer the elastic properties of composite nematic liquid crystals.

I. INTRODUCTION

Understanding the structure-property relationship of nematic liquid crystals -fluids without long-ranged positional order composed of aligned non-isotropic molecules or nanoparticles transmitting long-ranged orientational order -relies strongly on knowledge of the nematoelastic properties and their dependence on thermodynamic properties such as temperature or density [START_REF] De Gennes | The Physics of Liquid Crystals[END_REF]. For a bulk nematic fluid in three spatial dimensions there are essentially three fundamental modes, namely splay, twist and bend, that describe the extent of thermally driven elastic deformations of the uniform director field [START_REF] Frank | [END_REF]. In liquid-crystal based devices, the electro-optic switching characteristics of the director field are largely controlled by the splay and bend elastic modes, while the twist mode plays a key role in determining the helical pitch of chiral liquid crystals [3,[START_REF] Hiltrop | Chirality in Liquid Crystals[END_REF].

Connecting the amplitude of the various elastic modes to the microscopic properties of the constituents remains a challenging problem in view of the orientationdependent interactions between the molecules as well as the various attractive and repulsive forces acting among them whose range and amplitude are often unknown. In common modelling practice, one usually resorts to simple coarse-grained potentials or emblematic particle shapes (e.g. a cylindrical rod, disk or ellipsoid) to arrive at a satisfactory microscopic description of the material, including its elastic properties [START_REF] Mederos | [END_REF][6][7][8]. For lyotropic systems consisting of rigid or semi-flexible nanoparticles in a structureless solvent, the role of the (effective) particle shape is routinely expressed in terms of the aspect ratio of the components which, along with the particle concentration, constitutes the main parameter controlling the liquid crystal structure and phase behavior [9][10][11][12].

Obtaining first-principle predictions of the elastic properties from theory or simulation for elongated colloidal particles is a challenging task in view of the mesoscopic long-wavelength character of director fluctuations which necessitates extensive simulation setups and considerable numerical effort [7,[13][14][15][16][17]. Moreover, strongly anisotropic cylindrical objects suffer from poor sampling statistics in Monte Carlo methods rendering equilibrium structure investigations numerically cumbersome [12,16].

Density functional theoretical approaches for common nematics of cylindrical mesogens have been applied quite fruitfully to quantifying the elastic moduli as a function of particle anisotropy, concentration or nematic order parameter [18][19][20][21]. However, extending these theoretical treatments to analyze composite nematics comprising several components of different size and shape has been mostly overlooked thus far. Naively, one could argue that the elasticity of these multicomponent nematics should simple follow from a linear superposition of the elastic contributions of the individual components, weighted by their respective mole fraction. While this approach could certainly be a plausible one for nematic compounds in which all components obey the same nematic symmetry, it is no longer a viable route to describing multi-component nemato-elasticity when different nematic symmetries are present. For instance, particles may, under certain circumstances, generate anti-nematic orientational order characterized by a uniaxial symmetry of particles aligned across a plane perpendicular to the nematic director [22]. This gives rise to a negative nematic order parameter, contrary to what is observed in most liquid crystals.

Examples where such order may occur are so-called hybrid colloidal-molecular liquid crystals (Fig. 1a) or nematics of mixed mineral colloids involving a component with a distinct anti-nematic orientational organization [23][24][START_REF] Smalyukh | Annu Rev Condens Matter Phys[END_REF]. In hybrid nematic systems, small concentrations of strongly anisotropic colloidal particles are immersed in a thermotropic liquid crystal and spontaneously adopt anti-nematic order driven by an interplay between surface-anchoring and elasticity-mediated interactions enforcing the colloids to align orthogonal to the molecular director [24,[START_REF] Smalyukh | Annu Rev Condens Matter Phys[END_REF]. Other examples are mixed nematics composed of purely repulsive rod-and disk-shaped colloids [START_REF] Van Der Kooij | [END_REF][27][28] where the anti-nematic order of one component is generated spontaneously through excluded-volume-minimizing steric interaction with the other component (see Fig. 1b). For single component systems, anti-nematic order has been observed in some soft-interaction models for clay particles [29], deformable dendrimers [30], clay platelets exposed to an external electric-field [31] and lacuna smectic phases composed of ring-shaped mesogens [32,33].

An important condition for stable anti-nematicity in the composite nematics mentioned above is that the mole fraction of the anti-nematic component should be low enough to ensure uniaxial order to be preferable over biaxial order in which both species are oriented along mutually orthogonal directors, thus imparting an orthorhombic symmetry onto the fluid [24,34].

In an effort to further our understanding of these composite nematic materials we wish to address the question as to how nemato-elasticity could be quantified in a multi-component nematic in which one or several components exhibit anti-nematic order. In this work, we use a tractable second-virial theory applied to a representative model system of hard rods and disks to systematically investigate the elastic moduli for both nematic and anti-nematic order. By building on Odijk's asymptotic analysis for the orientation distribution [18,35] we can work out the different elastic modes in explicit form and gauge their scaling with particle concentration and aspect ratio. While the results for the rod-and diskbased nematics are in qualitative agreement with earlier predictions, we find that the interrelation of the splay, twist and bend modes for the anti-nematic fluid is different from that of a common nematic fluid. Focussing on hybrid nematics, involving strongly anisotropic colloids immersed in a thermotropic nematic background, we argue that surface anchoring has a considerable impact on the elastic properties of the compound nematic material, mostly boosting the splay elasticity in case of needle-type inclusions while for disk-shaped particles the bend mode is significantly enhanced. We finish our analysis by addressing the nemato-elasticity of rod-disk mixtures where coupling terms describe the elastic response transmitted by steric interactions between the rods and the disks. We propose scaling results that enables one to predict in detail how doping a rod-based nematic with disks (or vice versa) will affect the elastic properties of the reference system.

The rest of this paper is organized as follows. We begin by describing the principal context of the model and introducing the starting expressions. The subsequent two Sections are devoted to a computation of the scaling expressions for the elastic moduli of nematic and antinematic fluids of disks which we compare to previously derived results for rods. In Section V we address the role of colloid surface anchoring on the elastic properties of hybrid molecular-colloidal nematic materials. In Section VI we set out to quantify the effect of interspecies interactions on the elastic properties of nematic phases composed of a rod-disk mixture. Finally, some concluding remarks are formulated in Section VII.

II. FREE ENERGY OF DIRECTOR DEFORMATION

Let us consider a bulk nematic fluid composed of anisotropic building blocks and define elementary director fluctuations of amplitude ε σ -1 (with σ the typical particle size) around the nematic director n = ẑ fixed along the z-direction of a Cartesian frame defined as R = {Xx, Y ŷ, Zẑ} [16]:

n1 (R) = n + εY ŷ n2 (R) = n + εX ŷ n3 (R) = n + εZŷ, (1)
corresponding to the splay, twist and bend modes, respectively. Following Onsager-Straley theory [16,18,36], we assume that the free energy due to elastic distortions of the director field n can be quantified by expanding Onsager's original free energy functional [9] for nematics with a weakly non-uniform director field:

δF el ∼ - ρ 2 4 dR ∆rΦ(∆r; û, û )( ∇û n)( ∇û n) ḟ = 1 2 V Kε 2 , (2)
Here, we adopt short-hand notation for the double orientational average (

•) ḟ = dû ḟ (û • n) dû ḟ (û • n)(•)
and director gradient ∇û n = (∆r • ∇)(n(R) • û). For notational compactness, we set the thermal energy k B T to unity without loss of generality. The principal variables in the above expressions are ḟ (x) = ∂f /∂x the derivative of the uniaxial orientational probability f (û • n) for the orientation vector û of each particle, with V is the system volume, ρ = N/V the particle density, and K the elastic modulus. The key input in the present second-virial approach is the Mayer function Φ = e -U -1 which contains all the information on the shape and interaction between the constituents through the pair potential U [9]. In this, study we assume the colloidal interactions to be hard so that for any particle pair at centre-of-mass distance ∆r and orientations û and û :

Φ(∆r; û, û ) = -1, if cores overlap 0, otherwise (3)
The principal model we consider is a hard cylinder with aspect ratio L r /D r → ∞ corresponding to infinitely long hard rods and the inverse limit L d /D d → 0 describing infinitely thin hard disks. These shapes are considered as emblematic models for a vast range of many lyotropic liquid crystals. Soft interactions arising from e.g, chargemediated interactions or any other type of coarse-grained dispersion forces, could, in principle, be included through some effective second-virial theory [37,38] but this seriously complicates the present analysis and goes beyond the scope of the current project. Plugging in the fluctuations ni (R) into Eq. ( 2) enables us to read off the microscopic expressions for the respective modes:

K 1 ∼ - ρ 2 2 d∆r(∆y) 2 Φ(∆r, û, û )u y u y ḟ K 2 ∼ - ρ 2 2 d∆r(∆x) 2 Φ(∆r, û, û )u y u y ḟ K 3 ∼ - ρ 2 2 d∆r(∆z) 2 Φ(∆r, û, û )u y u y ḟ , (4)
These expressions serve as our starting point for an explicit calculation of the moduli for different particle shapes (which determines the interaction kernel Φ) and orientational symmetry. The latter is governed by the orientation distribution f (û • n)) which in turn follows from an optimization of the orientational and excludedvolume entropies contributions to the free energy. This leads to the following self-consistency expression:

f (û • n) = N exp ρ d∆rΦ(∆r, û, û ) f = N exp [-ρ v ex (û, û ) f ] (5)
with N a normalization factor ensuring f f = 1. The Boltzmann probability is thus defined in terms of the ensemble-average of the excluded volume v ex (û, û ) between anisotropic particles which itself depends on f . To describe common nematic order, we shall employ a much simpler Gaussian form appropriate for situations where the nematic alignment of the particles is asymptotically strong [35]:

f G (θ) ∼ α 4π exp - 1 2 αθ 2 , ( α 1) (6)
in terms of a polar angle θ = cos -1 (û • n) and concentration-dependent variational parameter α. Since the orientational order we consider is strictly apolar, the distribution must be supplemented with its mirror form f G (π -θ) to reflect the equivalence between alignment along n and -n. Further on in this study, we shall define the equivalent distribution for the anti-nematic case.

Although the Gaussian distribution is not thermodynamically consistent [39] it enables us to render the orientational averages involved in the elastic moduli analytically tractable as we shall demonstrate in the remainder of this paper.

III. ELASTIC MODULI OF NEMATIC FLUIDS OF DISKS AND RODS

We begin by presenting an detailed calculation of the splay, twist and bend elastic modes for a simple uniaxial nematic fluid of infinitely thin rigid hard disks of diameter D. The spatial integral required for the elastic moduli takes the form of a generalized or weighted excluded volume:

M n (û, û ) = -d∆r(∆r • êi ) n Φ(∆r, û, û ) = d∆r dd (∆r dd • êi ) n , (7)
with n = 2 and êi (i = x, y, z) indicating the Cartesian base vectors. Most importantly, ∆r dd parameterizes the overlap zone between two infinitely flat disks (illustrated in Fig. 2) at fixed orientations û and û within a particlebased frame { ŵ, ŵ , v} with v = (û × û )/| sin γ| and ŵ( ) = û( ) × v:

∆r dd = - D d 2 t 1 ŵ- D d 2 t 2 ŵ + D d 2 t 3 [(1-t 2 1 ) 1 2 +(1-t 2 2 ) 1 2 ]v, (8)
with integration variables -1 ≤ t 1,2,3 ≤ 1. The Jacobian associated with the transformation from the Cartesian to the particle-based frame reads J =

D 3 d 8 | sin γ|[(1 -t 2 1 ) 1 2 + (1 -t 2 2 ) 1 2
]. The excluded volume between two hard disks simply follows from Eq. ( 7) taking n = 0 and reads:

M 0 (û, û ) = ∆r dd ∼ π 2 D 3 d | sin γ|. (9)
The kernels M 2 needed for the elastic moduli take on a more complicated orientational dependence, namely:

M (i) 2 (û, û ) = πD 5 d 192 | sin γ|[22(v•ê i ) 2 +7( ŵ•ê i ) 2 +7( ŵ •ê i ) 2 ].
(10) The remaining task is to perform a double weighted orientational averages in Eq. ( 4):

K 1 ∼ ρ 2 d 2 M (y) 2 (û, û )u y u y ḟ K 2 ∼ ρ 2 d 2 M (x) 2 (û, û )u y u y ḟ K 3 ∼ ρ 2 d 2 M (z) 2 (û, û )u y u y ḟ . ( 11
)
with ρ d the number concentration of disks. Proceeding toward an explicit calculation of the moduli we consider a simple Gaussian Ansatz Eq. ( 6) for the orientational fluctuations around the colloidal director n. The derivative needed for the computation of the elastic moduli simply reads ḟG = αf G . In the weak fluctuation limit θ 1 we write up to leading order û ∼ {θ sin φ, θ cos φ, 1}. If we further assume uniaxial nematics then f G does not depend on the azimuthal angle φ. Furthermore, we introduce a dimensionless disk concentration c d = ρ d D 3 d and express, from here on, all elastic modes in units k B T /D d with D d the cylinder diameter.

After performing tedious trigonometric manipulations we can express the elastic moduli in terms of the relevant angular variables. Retaining only the leading order terms for small polar angles we obtain a set of rather hefty expressions shown in the Appendix. Fortunately, the results can be greatly compactified by invoking a basic scaling relation for the Gaussian averages, namely

θ n /|γ| ḟG ∝ α -n 2 + 5 2 .
Combining this with the common quadratic dependency of the Gaussian variational parameter with concentration, α ∼ (4/π)(π 2 /16) 2 c 2 d , well-known from bulk nematics [10,35], we readily infer a number of simple scaling expressions for the discotic elastic moduli:

K 1 ∼ 0.303c 3 d K 2 = 2K 1 K 3 ∼ 0.234c d . (12)
The prefactors were computed numerically even though it may be possible to derive rational prefactors upon extensive mathematical analysis that we did not pursue.

The large splay-tot-bend ratio K 1 /K 3 1 (assuming c d > 1 in the nematic regime) is in agreement with earlier theoretical predictions [6,22,40,41] and experiments [42] on various thermotropic discotic nematics. The interrelation of the different moduli turns out to be quite different from those established for rod-based nematics. Based on an analogous analysis Odijk [18] derived the following results for thin hard rods:

K 1 ∼ 7 32 c r K 2 = 1 3 K 1 K 3 ∼ π 48 c 3 r . ( 13
)
where c r = ρ r L 2 r D r is a dimensionless rod concentration defined in terms of the principal rod dimensions. The results demonstrate that the principal elastic response of discotic nematics stems from the splay (K 1 ) and twist (K 2 ) modes, while bend elasticity (K 3 ) dominates rodbased nematic fluids. The discrepancy among the modes is expected to be considerable for crowded nematics since all dominant modes scale with the cube of the particle concentration while the minor contributions increase linearly with concentration. The basic interrelation of the elastic moduli for rod (N r+ ) and disk-based (N d+ ) nematics are summarized in Table I.

IV. ELASTIC MODULI FOR ANISOTROPIC COLLOIDS WITH SURFACE-ANCHORING STABILIZED ANTI-NEMATIC ORDER

We now proceed towards analyzing the case of antinematic order in which the particle orientation vectors point perpendicular to the principal director n = ẑ without exhibiting a preferred direction of alignment across the x-y plane.

Hybrid colloidal-molecular liquid crystals may be characterized by anti-nematic order that is spontaneously TABLE I. Overview of the elastic moduli ratio for nematic (+) and anti-nematic (-) fluids of hard rods (r) and disks (d) as predicted from second-virial theory in the limit of asymptotically strong (anti-)nematic order. The interrelation of the elastic moduli depends on particle concentration ρ or the effective anti-nematic field with strength W , specified in the main text. The field amplitude is tuned by the temperature of the molecular host phase in case of a hybrid colloidal-molecular liquid crystal or by the interspecies interaction in case of a mixed rod-disk nematic characterized by a ratio of particle dimensions q = Lr/D d (see Fig. 1). Also shown are the twist-splay and bend-splay ratios that stem from the surface-anchoring of the colloids (indicated by (s)). These contributions are only specified for the single component anti-nematic systems.

Phase

twist-splay:

K2/K1 K (s) 2 /K (s) 1 bend-splay: K3/K1 K (s) 3 /K (s) 1 N r+ rods 1/3 1 (∝ ρ 2 r ) N r-rods 1/3 1/3 1 (∝ 1/W ) 1 (∝ 1/W ) N d+ disks 2 1 (∝ 1/ρ 2 d ) N d-disks 3 3 44 7 4 N d+/r-mixed (q 1) ≈ 1 1 N d+/r-mixed (q 1) ≈ 1/3 < 1 N r+/d-mixed (q 1) 3 ≈ 4 N r+/d-mixed (q 1) 3 1(∝ q)
generated through the surface anchoring energy of the anisotropic colloids embedded in a thermotropic liquid crystalline solvent [23,24]. This is illustrated in Fig. 1.

Let us consider a simple Rapini-Papoular expression for the surface anchoring energy per colloid [43]:

U s = - 1 2 w 0 dû dS(n • n0 (S)) 2 f (û • n) (14)
with w 0 > 0 the surface anchoring coefficient, n the uniform director of the molecular host in which the colloid is embedded and n0 (S) the unit vector describing the direction of preferred surface alignment across the colloid surface S. Let us assume that the molecules of the thermotropic background prefer to align parallel to the face of the disk (planar anchoring) or perpendicular to the long axis of the rod (homeotropic anchoring), then we may write in both cases n0 (ξ) = cos ξv + sin ξ ŵ (0 < ξ < 2π) in terms of an orthonormal particle frame depicted in Fig. 2. It is then easily established that the surface anchoring energy of the rod or disk changes with the polar angle θ between the molecular director and the main particle orientation via:

U s ∼ -r/d sin 2 θ ( 15
)
with r = π 4 w 0 L r D r (rods) and d = π 4 w 0 D 2 d (disks). Clearly, the surface anchoring energy reaches a minimum at θ = π/2 when the particle orientation is perpendicular to the molecular director [23,24]. It is important to note that, throughout this work, we assume that the surface anchoring is weak and that nematic distortions generated by the colloidal inclusions are minimal [44,45], a condition that should be satisfied as long as the rod or plate thickness is very small, more specifically D r ≈ L d K/w 0 with K the typical elastic modulus of the molecular host. Taking order-of-magnitude estimates for K ∼ 10 -12 N and w 0 ∼ 10 -5 N/m one infers a surface anchoring extrapolation length of about K/w 0 ∼ 10 -1 µm which is much larger than the typical thickness of most strongly anisotropic mineral colloids [24,46].

Since the basic uniaxial symmetry is retained for antinematic order the basic director deformations are identical to those defined in Eq. ( 1) and we start from the microscopic definitions laid out in Sec. II to analyze the nemato-elastic properties.

A. Rod-based anti-nematics

Starting with hard rods, we define M 2 as the generalized excluded volume kernel for two hard cylinders with aspect ratio tending to infinity (L r /D r → ∞) [36,47]:

M (i) 2 (û, û ) = -d∆r(∆r • êi ) 2 Φ(∆r, û, û ) ∼ 1 6 L 4 r D r | sin γ|(u 2 i + (u i ) 2 ) (16)
It is easily checked that the zeroth order moment simply yields the excluded volume M 0 = 2L 2 r D r | sin γ| of the lozenge depicted in Fig. 2.

Defining a meridional angle ψ = π 2 -θ, parameterizing the unit vector û = {cos ψ sin φ, cos ψ cos φ, sin ψ} and taking the leading order contributions for small ψ 1 in the orientational kernels above we obtain explicitly:

K 1 ∼ c 2 r 32 | sin γ| cos ∆φ cos ψ cos ψ (cos 2 ψ + cos 2 ψ ) ḟU K 2 ∼ 1 3 K 1 K 3 ∼ c 2 r 24 | sin γ| cos ∆φ cos ψ cos ψ (sin 2 ψ + sin 2 ψ ) ḟU (17)
with sin 2 γ = 1-(cos ∆φ cos ψ cos ψ +sin ψ sin ψ ) 2 . The leading order contribution for small meridional angles can be obtained by expand the trigonometric functions up to second order in ψ, using circular coordinates ψ = R sin χ and ψ = R cos χ. The integration over the relative azimuthal angle can then be cast in terms of complete elliptic integrals E(1 -R 2 ) of which we retain only the dominating contribution for R 1. The mathematical details of this procedure are outlined in Ref. [48]. After some algebra, we obtain up to second order in R:

K 1 ∼ c 2 r 16π R 2 sin χ cos χ(4 ln 2 -2 -2 ln R) ḟU K 2 ∼ 1 3 K 1 K 3 ∼ c 2 r 24π R 4 sin χ cos χ(4 ln 2 -2 -2 ln R) ḟU (18)
The splay-twist ratio K 1 /K 2 = 3 is identical to the one predicted for conventional rod-based nematics, while the bend-splay ratio K 3 /K 1 turns out much smaller than unity, similar to a that of a discotic nematic.

Let us now attempt an explicit calculation of the averages by taking the anti-nematic orientation distribution function for the rods proposed by Mundoor et al. [24] to describe the anti-nematic order of rods in a hybrid molecular-colloidal liquid crystal:

f U (θ) = N exp [-W P 2 (cos θ)] (19)
in terms of a variational parameter W quantifying the degree of anti-nematic order:

-W = 5π 16 c r S r -r . (20)
and N a normalization constant. Furthermore, S r < 0 is the nematic order parameter being negative for the anti-nematic state. Crucially, r > 0 is the strength of a temperature-controlled effective external field imparted by the surface anchoring properties of the rods with the molecular host forcing the rods to align perpendicular to the molecular director (cf. Eq. ( 15)). In the experimental system considered in Ref. [24], the effective field originating from the surface-anchoring energy is quite strong ( r 1 and W 1) and so will be the degree of anti-nematic order. The rod vectors will therefore adopt a very small equatorial tilt angle ψ = π 2 -θ 1 which justifies the use of a simple Gaussian asymptotic of the form:

f U (ψ) ∼ 3W 2π 1/2 exp - 3 2 W ψ 2 , W 1 ḟU (ψ) ∼ ∂f U ∂ψ ∼ -3W ψf U (ψ) (21)
with the variational parameter W related to the nematic order parameter via S r ∼ (1 -W )/2W . The limit S r → -1 2 for W → ∞ is easily verified. The double orientational averages in Eq. ( 18) are given by:

(•) ḟU = 27W 3 2π 2π 0 dχ sin χ cos χ ∞ 0 dRR 3 (•)e -3 2 W R 2 (22)
and the remaining integrals can be solved in closed form with the help of the mathematical analysis reported in Ref. [48]. Rearranging terms we obtain for the elastic moduli of an anti-nematic phase of rods:

K 1 ∼ 1 16π c 2 r (ln W + C 1 ) K 2 ∼ 1 3 K 1 K 3 ∼ 1 12π c 2 r ln W + C 3 W (23)
with constants C 1 = γ E -7/2 + ln 24 ≈ 0.255269 and C 3 = γ E -23/6 + ln 24 ≈ -0.0780638 and γ E being Euler's constant. In good approximation we may assume that S r ≈ -1 2 so that W ∼ r + 5π 32 c r . This renders the elastic moduli fully explicit in terms of the rod concentration c r and field strength r 1.

B. Discotic anti-nematics

We may repeat the analysis for an anti-nematic arrangement of hard disks with their normals pointing perpendicular to the nematic director. Starting from the kernel Eq. ( 10), we insert the anti-nematic parameterization of the unit vectors of the particle frame, introduced above for the rod case, and expand up to leading order. The anti-nematic moduli for disks take on a simple form and turn out to differ only by a constant factor:

K 1 ∼ 7πc d 1536 cos ∆φ(sin ∆φ) 2 |γ| ḟU K 2 ∼ 3K 1 K 3 ∼ 44 7 K 1 ( 24
)
The angular average in the first expression can be worked out in analytical form using the coordinate transform introduced below Eq. ( 17) which enables us to cast the average over the azimuthal angle ∆φ in terms of elliptic integrals for which the leading contributions for small R are of logarithmic form. This leads to:

K 1 ∼ 7πc 2 d 768 R 2 sin χ cos χ(4 ln 2 -4 -2 ln R) ḟU (25)
Employing the double averaging Eq. ( 22) weighted by the orientational distribution Eq. ( 19), with W suitably redefined in terms of the order parameter and concentration for the disks -W = 5π 2 32 c d S dd where S d > 0. Working out the averages we find that the splay modulus exhibits a similar logarithmic scaling with field strength W as observed for the rods, namely

K 1 ∼ 7π 3 c 2 d 768 (ln W + C 1 ) (26)
with C 1 = γ E -11/2 + ln 24.

FIG. 3. A splay deformation of a thermotropic liquid crystal is suppressed by the presence of a thin rigid rod with homeotropic surface anchoring oriented perpendicular to the average nematic director, while a thin disk with planar surface anchoring counteracts a bend fluctuation of the host director.

The surface-anchoring effect will be proportional to the colloid size (rod length or disk diameter). The direction of preferred surface alignment is indicated by the vector n0. The director fluctuations n1,2,3 have been specified in Eq. ( 1).

The main characteristics of the predicted elastic modes gathered so far are summarized in Table I. The logarithmic form of the moduli implies a much weaker concentration (or temperature) sensitivity for anti-nematic systems than expected for conventional nematics. Also, while both splay and twist elasticity of antinematically ordered rods are (logarithmically) increasing functions of W , the bend contribution K 3 ∝ W -1 ln W exhibits a maximum at W = e ≈ 2.71 and drops at large W as the degree of anti-nematic order increases. The opposite trend is observed for the disks, where the bend elastic term increases monotonically with field strength W . Clearly, the bend-splay ratio of nematic materials can be boosted considerably by introducing anti-nematically ordered colloidal disks.

V. SURFACE-ANCHORING ELASTICITY FOR HYBRID COLLOIDAL-MOLECULAR NEMATICS

So far, our attention has been on describing antinematic elasticity transmitted through colloid-colloid interactions. However, in case of anisotropic colloidal particles immersed in a thermotropic liquid crystals the inclusion of even a single colloid will perturb the elastic properties of the embedding thermotropic liquid crystals. The mechanism is quite simple; a long-wavelength perturbation of the molecular director in the presence of a colloidal particle will locally alter the anchoring condition at the particle surface and induce a change in free energy. This is illustrated in Fig. 3. The associated free energy change will be proportional to the main size of the colloidal inclusion and can be identified with a surfaceanchoring mediated elasticity that we will quantify for both rods and platelets.

A. Rods with homeotropic surface anchoring

Starting with infinitely thin rods with homeotropic surface anchoring, the molecular director is preferentally aligned along the surface normal vector n0 = cos ξv + sin ξ ŵ the rod surface in terms of an orthonormal particle frame {û, v, ŵ} (see Fig. 2). The Rapini-Papoular expression Eq. ( 14) can be easily generalized to the case of a non-uniform molecular director field n(R) formulated in Eq. ( 1). Expanding the surface anchoring free energy up to O( 2 ) we can deduce the free energy change associated with infinitely weak long-wavelength molecular director deformations. The associated elastic contributions originating from the surface anchoring condition are then defined as

K (s) n = δF s / 1 2 2
V which read in explicit form:

K (s) n = - 1 4 w 0 ρ r D r dû 2π 0 dξ Lr/2 -Lr/2 dt (tû • ên ) 2 × (ŷ • n0 (ξ)) 2 -(ẑ • n0 (ξ)) 2 f (û • n) (27)
with ê1 = ŷ, ê2 = x and ê3 = ẑ as per the different modes [cf. Eq. ( 1)]. The expression above involves an integral over the rod contour t probing the protrusion of the rod in the direction ên along which the host director is distorted. The integrals are easily solved using the Gaussian orientational distribution for the anti-nematic configuration formulated above. The resulting expressions for the elastic modes take a particularly simple form:

K (s) 1 ∼ π 128 w 0 c r L r D r K (s) 2 = 1 3 K (s) 1 K (s) 3 ∼ π 48 w 0 c r L r D r 1 W (28)
Similar to rod-interaction driven elasticity the twist contribution is one third of the splay term, and that splay elasticity dominates the other two modes (since W > 0). An important distinction, however, is that the surface anchoring-mediated contributions, representing principally a single rod effect, increase linearly with rod concentration. Clearly, the surface-anchoring elasticity is most noticeable for the splay mode while practically negligible for the bend mode provided the rods are strongly antinematic (W 1).

B. Disks with planar surface anchoring

Turning now to the case of infinitely thin disks with planar surface anchoring we formulate the surfaceanchoring elasticity analogous to Eq. ( 27):

K (s) n = - 1 2 w 0 ρ d dû 2π 0 dξ D d /2 0 drr dr (rr • ên ) 2 × (ŷ • n0 (ξ)) 2 -(ẑ • n0 (ξ)) 2 f (û • n) (29)
with r = cos ξv + sin ξ ŵ parameterizing the disk surface (see Fig. 2). The resulting expressions (in units ) are again quite simple and read:

K (s) 1 ∼ π 2 1024 w 0 c d K (s) 2 = 3K (s) 1 K (s) 3 = 4K (s) 1 (30)
Since disks are much more isotropic objects than slender rods they tend protrude in all Cartesian directions. Consequently, their surface anchoring properties exert a considerable impact on all director deformations of the molecular host structure with the bend mode being most affected.

For either colloid shape, the importance of the correlation-driven versus surface-anchoring mediated elasticity depends quite sensitively on the anchoring strength w 0 , the main particle dimension and the concentration of inclusions which should be restricted to the regime c < 1 to keep the colloids from aligning inplane and preserve anti-nematicity. We anticipate that for hybrid nematics with strongly nematic colloidal order (S > 0 and c 1) the correlation-mediated elasticity will outweigh the one driven by surface anchoring effects.

VI. ELASTICITY GENERATED BY INTERSPECIES INTERACTIONS IN ROD-DISK NEMATICS

We conclude our analysis by investigating the elastic contributions arising from rod-disk interactions that characterize the mixed colloidal nematics illustrated in Fig. 1b. In these situations, the nemato-elasticity is primarily transmitted by the principal component which is assumed to be ordered nematically, with the anti-nematic dopant perturbing the elastic properties of the mixed phase through rod-disk cross-interactions that we will attempt to quantify in the following.

In keeping with our original model of infinitely thin hard uniaxial cylinders we identify the principal dimensions as the disk diameter D d and rod length L r , and parameterize the rod-disk (rd) excluded volume as follows:

∆r rd = D d 2 t 1 sin ξv + D d 2 t 1 cos ξ ŵ + L r 2 t 2 û ( 31
)
in terms of a particle-based orthogonal frame {û, v, ŵ} for the disk and {û , v , ŵ } for the rod, with integration variables 0 < t 1,2 < 1 and 0 < ξ < 2π. The Jacobian associated with the coordinate transformation reads J =

1 8 L 2 r D d t 1 |v • ( ŵ × û )| so that: M 0 (γ) = d∆r rd ∼ π 4 L r D 2 d | cos γ| (32)
yields the excluded volume between a rod and a disk with orientation vectors making an angle γ. As in the preceding calculations we have neglected all subleading correction terms arising from the particles' finite thicknesses which are notoriously difficult to quantify [49] but are deemed unimportant for sufficiently slender mineral colloids. The second-order kernel reads:

M (i) 2 (û, û ) ∼ π 192 D 5 d | cos γ| 4q 3 (u i ) 2 + 3q(v 2 i + w 2 i ) , (33)
(with i = x, y, z) in terms of the ratio of the principal particle dimensions q = L r /D d . As before, the moduli are given by [cf. Eq. ( 11)]:

K 1 ∼ ρ r ρ d 2 M (y) 2 (û, û )u y u y ḟU K 2 ∼ ρ r ρ d 2 M (x) 2 (û, û )u y u y ḟU K 3 ∼ ρ r ρ d 2 M (z)
2 (û, û )u y u y ḟU (34) in terms of the partial concentrations ρ r and ρ d . It is important to note that both components obey distinctly different orientation distributions where one species is ordered nematically with the other adopting an antinematic configuration, so that the double orientational averaging proceeds via unequal orientational distributions, i.e. ḟU (û) = ḟU (û ). In the following we shall demonstrate that these orientational averages do not need be computed in explicit form in order to gauge some basic scaling properties pertaining to the ratio of the elastic moduli.

The first system under consideration is labeled N d+/r- and consists of a nematic arrangement of disks mixed with an anti-nematic organization of rods (see Fig. 1). As in the preceding, we use k B T /D d as a force unit to render the scaling expressions for the elastic moduli dimensionless and consider normalized colloid concentrations c r and c d . The results are as follows:

K 1 ∼ c r c d 2 π 128 θ cos ∆φ| cos γ| q -1 sin 2 γ + q ḟU K 2 ∼ c r c d 2 π 128 θ cos ∆φ| cos γ| q -1 sin 2 γ + q 3 ḟU K 3 ∼ c r c d 2 π 96 qθ(ψ ) 2 cos ∆φ| cos γ| ḟU (35)
with θ 1 the polar angle of the disks and ψ 1 the meridional one for the rods (see Fig. 2). Furthermore, we note that | cos γ| 1 and sin 2 γ ∼ O(1). The second case, N r+/d-, is a rod-based nematic in which the disks are dispersed anti-nematically. In this situation the moduli take on the following form:

K 1 ∼ c r c d 2 π 512 θ cos ∆φ| cos γ| q -1 sin 2 γ ḟU K 2 = 3K 1 K 3 ∼ c r c d 2 π 128 θ cos ∆φ| cos γ| q -1 sin 2 γ + 2q ḟU (36)
Refraining from any further analysis, we easily infer the interrelation of the elastic moduli for the two mixed roddisk uniaxial nematics depicted in Fig. 1b. The results obtained so far are tabulated in Table I. It is apparent that the elasticity of the mixed nematics can be carefully tuned by varying the ratio q of the rod length versus disk diameter, as well as by changing the partial concentrations. A general conclusion one could draw from looking at the results for the mixed systems is that antinematically order disks have a much more significant impact on the elastic properties (in particular enhancing the bend mode) of rod-dominated nematics (N r+/d-) while rod-shaped inclusions turn out to only weakly influence the elasticity of a discotic nematic. Furthermore, the results above suggest an interesting non-monotonic trend with the rod-disk size ratio q. This happens for the splay and twist elasticity of the discotic nematic and the bend elasticity of the rod-based nematic, where the moduli exhibit a minimum around q ∼ 1, that is, for rods and disks of about equal particle dimensions. We wish to demonstrate that it is possible to specify the effective anti-nematic field strength of the roddisk nematics in terms of relevant component variables by considering the coupled rod-and disk orientation distributions for a mixed uniaxial nematic proposed in Ref. [50]:

f r (cos θ) ∼ exp 5 4 (ρ r v rr ex S r -2ρ d v rd ex S d )P 2 (cos θ) f d (cos θ) ∼ exp 5 4 (-2ρ r v rr ex S r + ρ d v rd ex S d )P 2 (cos θ) (37)
in terms of the isotropized excluded volume of two rods v rr ex = (π/4)L 2 r D r and a rod-disk pair v rd ex = (π/8)L r D 2 d . The order parameters are such that S r > 1 and S d < 0 for the rod-dominated nematic (N r+/d-) and S d > 1 and S r < 0 for the disk-dominated one (N d+/r-). Comparing with Eq. ( 19)) we immediately deduce that the effective anti-nematic field strength for the two mixed nematics in Fig. 1

read: r ∼ 5π 16 c d qS d , (N d+/r-) d ∼ 5π 16 c r q -2 r S r , (N r+/d-) (38)
with r = L r /D r 1 the length-to-width (aspect) ratio of the rods. We observe that, in line with expectation, the effective anti-nematic field imposed on one species imparted by the nematic alignment of the other is proportional to the degree of nematic order (note that the values of S r,d are positive and close to unity) and the concentrations c r,d 1 of the other component. Keeping in mind that validity of our asymptotic analysis requires that the anti-nematic field strength be much larger than the thermal energy k B T we infer that the condition r,d

1 imposes two important criteria on the size disparity of the components, namely q > 1 and q -2 r > 1. The latter criterion is easily satisfied for rods in the Onsager limit r → ∞ whereas the former imposes the rod length to be at least the disk diameter.

A much stricter criterion that we have not addressed here is the thermodynamic stability of the uniaxial mixed nematic with respect to the formation of a biaxial phase in which the rods and disks are each aligned along mutually orthogonal directors. Intuitively, for rod-disk mixtures, one would expect uniaxial order to be preferable over biaxial nematics for strongly asymmetric mixtures with a strong disparity in the excluded volume among the components [48,51,52], which is broadly satisfied through the aforementioned criteria. In most practical situations studied thus far, the biaxial nematic tends to be unstable with respect to demixing into two fractionated uniaxial nematics [27,34,53].

VII. CONCLUSION

Inspired by recent experimental advances in the fabrication of well-controlled composite nematics [24,54] involving a colloidal component with distinct anti-nematic order we have embarked on computing the elastic moduli associated with such anti-nematically ordered anisotropic particles. For simplicity we have kept our focus on simplified but emblematic models for lyotropics, namely hard cylindrical rods and hard disks with vanishing particle thickness and address the elastic properties through Onsager's second-virial theory, suitably extended by Straley [36] to capture the effect of weak deformations of the nematic director field. While the elastic moduli for common nematics are in full agreement with those reported in experimental, theoretical and simulation studies, the elastic moduli for the anti-nematic case had not yet been identified and unveil a remarkable logarithmic scaling with the degree of anti-nematic order. More importantly, the interrelation between splay, twist and bend moduli turns out to be quite different from that of conventional nematics composed of colloidal rods.

Most notably, we find that the bend elasticity of rods vanishes with increasing anti-nematic order while the bend modulus for disks increases as the disks adopt a more pronounced anti-nematic configuration. In addition to the elasticity imparted by colloidal interactions, we quantify the effects of colloid surface anchoring which alters the deformation energy of the director field of the molecular liquid crystal in which the colloids are immersed. These moduli are all linear in colloid concentration and increase with the colloid size and the surface anchoring amplitude. We find that the ratio of the twistsplay and bend-splay elastic modes induced by surface anchoring effects are qualitatively similar to those generated by colloid interactions.

Extending our treatment to mixed uniaxial rod-disk nematics we present a preliminary analysis of the elastic moduli for coupled rod-disk interactions and argue that bend fluctuations of rod-dominated nematics can be strongly suppressed by adding disk-shaped inclusions. A full analytical characterization of the elastic moduli for mixed rod-disk nematics in terms of the partial concentrations of the components, not attempted in this paper, seems a highly non-trivial task but should at least be feasible numerically using the preliminary scaling expressions for the moduli and the Gaussian orientational distributions proposed in this work.

It would be intriguing to compare our predictions with experimental measurements of the moduli for e.g. colloidal rod-disk mixtures or other composite nematic phases involving anti-nematic order of one or several components. For thermotropics, there are a number of papers reporting experimental observations of remarkable changes in elastic response upon adding colloidal dopants, often (but not always) leading to a strong reduction of the elastic resistance of the composite material [55][56][57][58]. To the best of our knowledge no such experimental data have been reported thus far for the lyotropic systems under scrutiny and we hope that our work will stimulate experimental and simulation efforts to characterize the intricate nemato-elasticity in mixed colloidal or hybrid molecular-colloidal nematics. We finally wish to point out that it should be feasible to extend the current second-virial analysis for biaxial colloidal composites for which there are twelve independent elastic moduli [59,60]. It would be intriguing to theoretically quantify those moduli along the lines of the present analysis and compare with experimental measurements in biaxial hybrid molecular-colloidal materials. Efforts in this direction are currently being undertaken.

APPENDIX

We show here the explicit angular dependence of the discotic elastic constants:

K 1 ∼ - c 2 d 2 - 7π 192
θ
with θ and θ the polar angles, ∆φ = φ -φ the relative azimuthal angle and γ the angle between the orientation vectors of the disks (see Fig. 2).

  FIG. 1. (a)Anti-nematic order in hybrid nematic materials generated through surface-anchoring of colloidal inclusions (thin rods or disks) immersed in a thermotropic nematic liquid crystal with molecular director nm. (b) Anti-nematic order in mixed colloidal nematics composed of rod-and disk-shaped particles. In the rod-dominated phase (left) the disk normals are aligned in an anti-nematic configuration with their normals pointing perpendicular to the overall director n. In the disk-dominated phase (right) the rod-shaped inclusions are aligned anti-nematically.

  FIG. 2. (a) Overview of the lab frame, nematic director n, and principal angles used in the present analysis. (b) Excluded volume between two thin hard rods at fixed mutual angle γ and aspect ratio Lr/Dr 1. The orthonormal particle frame {û, v, ŵ} is indicated by the red arrows. (c) Same for infinitely thin disks with diameter D d . For non-perpendicular orientations (γ = π/2) the rods generate a basic lozenge shape while for the disks a slanted sphero-cuboid is obtained.

ACKNOWLEDGEMENTI

  am grateful to Ivan Smalyukh of the University of Colorado (Boulder, USA) for fruitful discussions and for stimulating my interest in this topic.


Henri Bachau
email: henri.bachau@u-bordeaux.fr

Matabara Dieng

Nondipole effects in helium ionization by intense soft x-ray laser pulses

Keywords:

We consider the helium atom exposed to laser pulses in soft-x regime at intensities of the order of 10 17 -10 20 W/cm 2 . Two cases are investigated, first a pulse with a central frequency of 20 a.u. (∼ 544.2 eV) and a pulse with a central frequeny of 50 a.u. (∼ 1.36 keV). We solve the time-dependent Schrödinger equation (TDSE), nondipole (retardation) effects are included up to O(1/c) (c being the velocity of light).

We study the single photoionization of helium with a focus on nondipole effects. We calculate the total and partial photoelectron energy spectra, and the photoelectron angular distributions. The nondipole effects associated with the A • P and A 2 coupling terms in the Hamiltonian (A denoting the vector potential of the field and P the momentum operator) are thoroughly analyzed in the region of one-photon resonance.

We show that the nondipole contribution associated with A • P varies linearly with the intensity I while a three-photon transition involving A 2 increases like I 3 , in agreement with lowest order perturbation theory (LOPT). At high intensities these contributions compete and the nondipole transition becomes nonlinear.

The physical mechanisms associated with nondipole couplings are discussed, as well as their effects on photoelectron energy and angular distributions.

Introduction

The advent of free electron lasers delivering x-rays at high intensity [1,2] opens the way to the exploration of new regimes in physics, chemistry and bio-imaging. In atomic physics multiphoton processes have been identified in neon in soft x-ray regime at intensities of the order of 10 17 -10 18 W/cm 2 [3,4] (LCLS facility, USA) and in xenon above 5 keV at SACLA facility (Japan) [5]. High intensities are also obtained in the hard x-ray regime, to cite a recent application in atomic physics multiphoton multiple ionization of a high-Z atom has been reported at ultra-intense intensity (10 19 W/cm 2 ) and photon energies from 5.5 to 8.3 keV [6]. Many FEL sources are now in operation, upgrade or in construction throughout the world, their status as well as the scientific achievements are reported in recent reviews [7][8][9], the latter focusing more on technology and physics of FELs. Concerning now the target, light atoms, like hydrogen or helium, are ideal candidates to investigate the general mechanisms governing their interaction with light since they are accessible to an accurate theoretical description. But, due to the small photoionization cross section, the x-ray photoelectrons are hardly detectable.

Nevertheless, it must be noticed that bright-enough synchrotron sources allows now to study light atoms using x-ray photoelectron spectroscopy [10], and multiphoton processes have been detected in helium at photon energies of 42.8 and 38.4 eV at the FEL facility FLASH (Germany) [11]. FELs operating at high intensities open the way to pursue these investigations in x-ray regime.

In the field of photoelectron spectroscopy of atoms at short wavelengths, observation of nondipole effects in photoelectron angular distributions (PAD) were made in the late 1920s [12] in the hard x-ray photon energy range, while the first experimental studies below 5 keV date back to the late 1960s-early 1970 [13,14]. Thanks to the unique properties of the photons generated by synchrotron facilities, it has been possible to measure nondipole parameters in x-ray regime (like for example the nondipolar asymmetry parameters of PADs in Ar [15]), a tremendous amount of work has been devoted to this field these last fifty years.

Due to the relatively low photon fluence at these facilities the theoretical work has been developed in the framework of LOPT. The nondipole corrections are associated with the coupling terms A • P and A 2 (diamagnetic term) in the Hamiltonian. In dipole approximation (DA) A 2 does not contribute to light-matter dynamics while the coupling term A • P reduces to an electric-dipole (E1) transition. Beyond the DA the correction term A 2 , associated with a two-photon transition, can be neglected in the context of the early experiments above cited, the photoionization process involves one-photon absorption and it is dominated by the coupling term A • P. Using nonrelativistic, independent-particle wavefunctions and developing A • P in first-order Cooper [16,17] (see references therein and the Errata [18]) proposed a simple parametrization of the differential cross section, in the form

dσ nl dΩ = σ nl 4π [1 + βP 2 (cos θ) + (δ + γ cos 2 θ) sin θ cos φ] (1)
for linearly polarized light. In this expression, θ is the angle between the polarization vector and the electron mo- Experimental results obtained around FELs at high intensity stimulate the development of theoretical approaches able to describe nonlinear processes with intense and short

x-ray pulse. These studies raise the problem of the effects of the retardation terms, with in particular the role of A 2 . In this context two-photon ionization processes, like above-threshold ionization in x-ray regime, have been investigated within LOPT [19][20][21]. Non-perturbative methods have been also developed to investigate the relative contributions of A • P and A 2 in one-and two-photon ionization of hydrogen. For example in [22] photons energies ranging from 200 eV to 3 keV have been considered, at the atomic unit of intensity, i.e., I 0 = 3.51 × 10 16 W/cm 2 .

The conclusion of this study is that, if the contribution of A 2 is non-negligible, the nondipole correction mainly comes from the term A • P. Incidentally, we notice that the nondipole diamagnetic term A 2 plays a crucial role in nonlinear processes like stimulated Compton Scattering [23] and stimulated Raman Scattering [24]. Analytical approaches, theory and methodologies of the TDSE in the framework of ab initio calculations, with a particular focus on dynamical effects induced by the A 2 term, can be found in [25][26][27][28][29][30]. Here it is worth recalling that, even at intensities of the order of I 0 , LOPT applies in x-ray regime since the ponderomotive energy U p (U p = I/4ω 2 where I is the intensity and ω the photon energy) is much smaller than the photon energy. At the same time the value of the Keldysh parameter [31], given by I p /2U p (where I p is the atomic ionization potential) is much greater than 1, unfavorable for tunnel ionization. As a matter of fact, as the intensity increases beyond I 0 and/or the photon energy decreases from x-ray to xuv domain the validity of LOPT and DA must be carefully checked.

A recent theoretical study [30]

Theory

In this section we refer to the nonperturbative (TDSE) approach. We follow the numerical method previously developed in [22]. We recall only the main steps.

We start from the well-known expression of the nonrelativistic Hamiltonian operator for a one-active electron atom in an external electromagnetic field of vector potential A (in the Coulomb gauge):

H = 1 2 [P + A(r, t)] 2 + V (r). (2)
The atomic potential V (r) is a central one, it has the form

V (r) = - 1 r - 1 r (1 + λr)e -2λr . (3)
This potential is widely used in atomic and molecular physics (see [START_REF] Bachau | [END_REF] and references therein). It verifies the correct asymptotic conditions for He; V (r) ∼ r→0 -2/r and V (r) ∼ r→∞ -1/r. Here λ = 1.688, this parameter is chosen such that the absolute value of the lowest eigenenergy (here associated with the 1s state of the model atom) of the field free Hamiltonian reproduces the ionization potential of helium (∼ 0.9037 a.u.). Taking the z-axis (unit vector e z ) along the polarization vector and the wavevector κ along e x (unit vector of x-axis) as shown in Fig. 1, the vector potential reads;

A(r, t) = A(t -αx) e z (4)
with α the fine structure constant (α = 1/c). Considering that α ω ≤ 1 (or ω ≤ 3.73 keV) A(t -αx) is approximated by its truncated power series in x. Keeping only the firstorder correction in α one can show that the Hamiltonian can be written

H ≈ H a + H (1)
DA + H

RET + H

RET .

It is expressed as a sum of the atomic Hamiltonian H a ≡ P 2 /2 + V (r), of the interaction term in DA

H

(1)

DA = A(t) P z , (6)
and of two other terms, describing nondipole corrections;

H (1) RET = α F (t) x P z (7)
and

H (2) RET = α F (t) A(t) x (8)
where F (t) = -Ȧ(t). Within first-order perturbation theory, the terms A • P ≡ H

DA + H

RET and A 2 ≡ H and for nondipole coupling terms (∆l = 0, ±2; ∆m = ±1)

for H

RET and (∆l = ±1; ∆m = ±1) for H

RET .

The TDSE reads

i ∂ ∂t ψ(r, t) = Hψ(r, t). (9)
The above equation is resolved by using four different approximations for the Hamiltonian H. First we consider the DA with H ≡ H a + H

DA (see Eq. 5). Second we consider the Hamiltonian H ≡ H a + H

DA + H

RET which (only) includes the contributions associated with the coupling term

A • P. Then we use H ≡ H a + H (1)
DA + H

RET where the nondipole correction associated with A 2 is added to the

coupling term A • P in DA (H (1)
DA ). Finally the full expression ( 5) is considered;

H ≡ H a + H (1)
DA + H

RET + H

RET .

These four approaches allow to evaluate the relative contributions of DA and nondipole coupling terms. In the following they will be referred to as TDSE-DA, TDSE-AP, TDSE-DA-A 2 and TDSE-FULL, respectively.

A spectral method is used to resolve the TDSE. First, we define the amplitude A(t) of the field vector potential such that it is non-zero over the time interval (-T /2,T /2)

where T is the total pulse duration (i.e., the full extension of the pulse support), A(t) is chosen to have a cos 2 envelope. Under these conditions the laser bandwidth, defined as the FWHM of the square of the Fourier transform of the field amplitude, is given by 1.44×2π/T . As mentioned above the laser field is linearly polarized along the z-axis.

In order to calculate the time-dependent wave function ψ(r, t), solution of the TDSE, we use the expansion

ψ(r, t) = n, l, m e -iE nl t c (l,m) n (t) u nlm (r) (10)
in a discrete basis of H a eigenfunctions

u nlm (r) = X nl (r)/r Y lm (r/r). ( 11
)
The eigenvalues of H a , indexed as E nl , are determined together with the radial eigenfunctions X nl (r) by solving numerically the radial Schrödinger equation in a basis of B-spline functions [36]. The coupled differential equations satisfied by the coefficients c Once the wave function ψ(r, t) is known at the end of the pulse, the angle-integrated photoelectron energy spectrum (PES) and PADs can be extracted from it. Partial photoelectron energy spectrum (PPS), associated with channels of given l and m can be also extracted, in the following they will be referred to as PPS (l = l , m = m ). As we will see below, the PPSs associated with the different approximations of H allow to identify and to quantify the effects of nondipole coupling terms, as well as the contribution in DA. The procedure used to extract the angular distributions is presented in the appendix of [22], we do not repeat it here. The only difference with [22] is that the eigenfunctions of the field free Hamiltonian H a are not hydrogenic ones, since V (r) is not a pure Coulomb potential.

In the context of TDSE, this has not a major impact on calculations (for the extraction of angular distributions a phase shift must be calculated for each angular momentum l in addition to the Coulomb one). Nevertheless it is worth noticing that the presence of a screened Coulomb potential has a major impact in the context of perturbation theory (PT) since analytic approaches, valid in the case of hydrogenic systems (see [20,21] and references therein), cannot be used.

Before presenting the results we briefly discuss the limitations of the one-active electron model. It only represents singly excited (and ionized) states of He and not the doubly-excited states and ionization-excitation. The latter channels are energetically open here, through Raman-type processes or direct photon absorption. First we note that the photoionization of He is dominated by the channel leaving the He + ion in 1s state at the wavelengths considered here [32]. Given the laser bandwidths (about 0.18 a.u. for the case of frequency ω = 20 a.u. and 0.33 a.u.

for ω = 50 a.u.) the ionization-excitation channels should lead to structures in the electron spectrum well separated in energy from the main photoionization peak. Regarding now the excitation of bound states of He it relies on stimulated Raman scattering, we have demonstrated in previous work [24] that the dominant contribution is associated with the diamagnetic term H

RET in first order of PT, it populates the 1snp series with m = ±1. The other contribution is associated with the DA term H

(1) DA in second-order of PT, leading to the population of 1sns and 1snd series with m = 0 (1sns dominating over 1snd).

In this context the contribution of H

RET is negligible. We have checked that these trends hold in this work. Similar processes could populate doubly excited series as well the associated ionization-excitation channels in the low energy continuum. This is an interesting subject in itself but due to the laser bandwidths considered here, which are much smaller than the He ionization threshold, these transitions are not expected to play a significant role. A two-color scheme [23] should be more appropriate to study electronic double excitation. Incidentally we mention that the selection rules given above for DA and nondipole correction terms are strictly valid for the one-active electron model (see [START_REF] Bethe | Quantum Mechanics of One-and Two-Electron Atoms (Pluenum/Rosetta Edition[END_REF] and [START_REF] Bransden | Physics of Atom and Molecules[END_REF] for the rules in DA and retardation couplings in the general case). 3

). This broadening is due to the slow numerical convergence of the TDSE-FULL calculations with l in expansion (10). shown in the figure). Therefore, even if the channels with l = 1, 2 dominate at the end of the pulse (at t = T /2) in the region of peaks 1, the integration of the TDSE-FULL requires a much larger number of angular momenta l. As expected, higher is the intensity, larger is the value of l max required to reach the numerical convergence of the calculations [29]. Considering that -l ≤ m ≤ l this rapidly leads to an important increase of the basis set and a dramatic growth of the computational demand.

In the following we focus on the region of peak 1,

where the PPS (l = 2, m = ±1) dominates over the others nondipole channels. In Fig. 5 we present the PPS (l = 2, m = 1) for TDSE-FULL, TDSE-DA-A 2 and TDSE-AP calculations, as defined below the Eq. ( 9). First the right figure shows results for I = 10 18 W/cm 2 , we clearly see that the main contributions to the nondipole effect is related to the coupling term

A • P (included in TDSE-AP) while the contribution of A 2 (included in TDSE-DA-A 2 )
is several orders of magnitude smaller. We recall that the selection rule for the term H

RET is (∆l = 0, ±2; ∆m = 1), thus the PPS (l = 2, m = 1) is mainly due the term H W/cm 2 match. We reconsider now the angular distribution expressed in Eq. ( 1). We have noticed in the introduction that it takes into account only the A • P coupling term in first-order. For a one-electron s-p transition β = 2 and δ = 0 [17], therefore the PAD differential probability reads

dP 1s dΩ ∝ (3 cos 2 θ + γ cos 2 θ sin θ cos φ). (12)
We recall that γ is the nondipole angular parameter. As expected the DA distribution is flat (γ 0 in DA) in Fig. 6, the TDSE-FULL PADs clearly show the influence of the nondipole effects. There is a very good fit to the analytic formula (12) at both intensities, it is obtained with γ ≈ 0.51, in agreement with the value of γ calculated by Amusia et al [37] beyond the framework of one-electron approximation (see Fig. 1 in [37]). It is worth noticing that in [37] electron correlations have been found to be unimportant in the calculation of γ in the case of He 1s shell, this explains the agreement between our one-active electron model atom and the above-cited multielectronic approach. We notice the agreement obtained for I = 10 19 W/cm 2 , although the contribution of A 2 is non-negligible at this intensity. In order to check the influence of the latter term we have calculated the PAD with TDSE-DA-A 2 (always integrating the differential distribution over W/cm 2 . The others laser parameters are the same as in Fig. 2.

For I = 10 18 W/cm 2 the angular distribution has been scaled such that the angular distributions in DA match for both intensities. The figure also shows the angular distribution calculated with the formula (12), with γ = 0.51. The PAD is shown in Fig. 8, as in the precedent case the polar angle is fixed at θ = π/4. It is obtained by integrating the triply differential distribution in the region of one-photon resonance, i.e., over the energy range 48-50 a.u. The PAD fits the analytic formula (12) with γ ≈ 0.99, a value in agreement with the above-mentioned theoretical work [37]. As in the precedent case, the diamagnetic term A 2 has a little impact on the PAD.

Conclusions

We have presented a study of nondipole effects in the ion-

Authors contributions

All the authors were involved in the preparation of the manuscript.

They have read and approved the final manuscript.

Fig. 1 .

 1 photon direction

  RETinduce one and two-photon transitions, respectively. The selection rule is (∆l = ±1; ∆m = 0) in DA (term H (1) DA )

  from t = -T /2 to t = T /2 using an explicit time-adaptive Runge-Kutta method, with initial conditions (at t = -T /2) corresponding to He in its ground state, i.e., here the 1s state. In a typical calculation we use angular momenta l = 0 -10, with |m| ≤ l, and a basis of 2500 B-spline functions of order k = 7, distributed linearly inside a box of length b = 600 a.u. for the case ω = 20 a.u. and b = 300 a.u. for ω = 50 a.u. The box radius b is chosen such that the probability to find the photoelectron outside the box at the end of the pulse is negligible. The convergence of the calculations is checked by varying the box size and/or the number of basis functions and angular momenta, we will return to this point later.

3 Results 3 . 1 Fig. 2 .

 3312 Fig. 2. PES in DA, calculated at various intensities (see the inset). The central frequency ω = 20 a.u. and the pulse duration T ≈ 50.26 a.u.

Fig. 3 . 19 W/cm 2 .

 3192 Fig. 3. PESs calculated in DA (TDSE-DA) and with full Hamiltonian (TDSE-FULL). PPS (l = 2, m = 1; m = -1) (sum of the PPS (l = 2, m = 1) and PPS = 2, m = -1)) calculated with TDSE-FULL. The central frequency ω = 20 a.u., the pulse duration T ≈ 50.26 a.u. and the intensity I = 10 19 W/cm 2 .

Fig. 4 . 3 .

 43 Fig. 4. PES calculated with various basis set expansion in l (see the inset). The laser parameters are the same as in Fig. 3.

( 1 ) 1 ) 2 ) 2 )

 1122 scheme explained above. In parallel with the resolution of the TDSE-AP we have performed time-dependent calculations where H (1) RET is treated in first-order PT, they are in good agreement with the TDSE-AP PPS (l = 2, m = 1) for both intensities, this confirms that it is a one-photon transition in the perturbative regime. If we consider now the TDSE-DA-A 2 contribution, we have identified two pathways involving one-and two-photon resonances. First a one-photon dipole coupling from 1s to the channel (l = 1, m = 0), followed by a Compton-like two-photon transition from the channel (l = 1, m = 0) to the PPS (l = 2, m = ±1) around the one-photon resonance (at energy ∼ 19.1 a.u.). The latter transition is induced by the nondipole coupling term H (2) RET (associated with A 2 ) according to the selection rule (∆l = ±1; ∆m ± 1), it is a two-photon absorption-emission process. The second process relies on a direct two-photon coupling (through H (2) RET ) from 1s to the channel (l = 1, m = ±1) at energy ∼ 39.1 a.u., followed by a dipole transition to the channel (l = 2, m = ±1) at energy ∼ 19.1 a.u. The TDSE does not allow to discriminate between these two pathways, but the second one involves a continuum-continuum dipole coupling with high energy continua (in the region of peak 1 and peak 2), this type of coupling is usually very small and we can therefore surmise that the first pathway dominates. Both

2 Fig. 5 .Fig. 2 .

 252 Fig. 5. PPS (l = 2, m = 1) calculated with full Hamiltonian (TDSE-FULL) and including the nondipole correction associated with A • P (TDSE-AP) and A 2 (TDSE-DA-A 2 ). The intensity I = 10 19 W/cm 2 (left figure) and I = 10 18 W/cm 2 (right figure). The others laser parameters are the same as in Fig. 2.

Fig. 6 .

 6 Fig. 6. PAD vs the azimuthal angle φ, the polar angle θ = 45 degrees. The intensities are I = 10 19 W/cm 2 and I = 10 18

  a.u. (150 opt. periods) ω=50 a.u.

Fig. 7 .Fig. 8 .

 78 Fig. 7. PES and PPS (l = 2, m = 1) calculated with the full Hamiltonian (TDSE-FULL). PPSs (l = 2, m = 1) are also calculated considering nondipole correction associated with A • P (TDSE-AP) and A 2 (TDSE-DA-A 2 ). The central frequency ω = 50 a.u., the pulse duration T ≈ 18.85 and the intensity I = 10 20 W/cm 2 .

1 ) 1 )

 11 photoionization is dominated by the dipole transition, it populates the channel (l = 1, m = 0) at energies close to 19.1 a.u. Below I = 10 19 W/cm 2 the dominant nondipole contribution is associated with A • P (H (1) RET in Eq. 5), it is about three orders of magnitude smaller than the dipole contribution and it populates the channel (l = 2, m = ±1) through one-photon absorption from 1s. At the intensity I = 10 18 W/cm 2 , both the dipole and the dominant nondipole transitions increase linearly with the intensity, in agreement with PT. At I = 10 19 W/cm 2 , the nondipole term A 2 comes into play. Combined with a dipole coupling, it populates the channel (l = 2, m = ±1) in the region of one-photon resonance (at energy ∼ 19.1 a.u.) and it is of the same order of magnitude as the nondipole contribution associated with A • P. This process is nonlinear, three photons are involved in different combinations where two photons are absorbed and one is emitted. The transition rate varies like I 3 , in agreement with LOPT. Regarding now the PAD it fits well the parametrization of Cooper [16, 17], based on a first-order treatment of the coupling term A • P in PT. The diamagnetic term A 2 does not affect the PAD associated with the energy integrated one-photon resonant peak at the intensities considered. The general trends of the distributions calculated at ω = 20 a.u. are also observed at ω = 50 a.u., but at higher intensities. At I = 10 20 W/cm 2 and T = 18.85 a.u. (150 optical cycles) the nondipole contributions from A • P and A 2 are of the same order of magnitude and they vary with the intensity like I and I 3 , respectively, in