Clémentine Courtès

Matthieu Boileau

Raphaël Côte

Paul-Antoine Antoine Hervieux

Giovanni Manfredi

Micromagnetic

Micromagnetic simulations of the size dependence of the Curie temperature in ferromagnetic nanowires and nanolayers

Keywords: Micromagnetism, Landau-Lifshitz-Gilbert equation, Curie temperature, nanowire, nanolayer, size effect

come

Introduction

Interest in nano-scale ferromagnetic objects (nanowires, nanolayers, etc.) has grown dramatically in recent years, as these objects are integrated into many current devices with the aim of storing, reading and writing digital information. Hence, there is a growing need for accurate numerical simulations capable of predicting the behaviour of such objects, and possibly predict new properties and features. Computational models of ferromagnetic materials can be roughly divided into two families. On the one hand, atomistic models such as VAMPIRE [START_REF] Evans | User Manual[END_REF] describe the magnetic interactions microscopically, at the natural atomic length scale of the material. Although in principle very accurate, they demand a considerable computational cost and are thus limited to relatively small systems.

In contrast, micromagnetism describes the structure and dynamics of a ferromagnetic object at an intermediate mesoscopic scale, averaging over a large number of atoms. Micromagnetic codes, based on the Landau-Lifshitz-Gilbert (LLG) equation, describe the dynamics of the average magnetic moment m(t, x), which is a continuous function of the spatial coordinate x. Notable projects in computational micromagnetism are the OOMMF project [START_REF] Donahue | OOMMF User's Guide, Version 1.0[END_REF] (Object Oriented MicroMagnetic Framework) for the development of a public micromagnetic program in C++, the mumax3 project [START_REF] Vansteenkiste | The design and verification of MuMax3[END_REF], a GPU-accelerated micromagnetic simulation program, or else the tetmag project, a 3D micromagnetic finite-element simulation software [START_REF] Hertel | tetmag[END_REF].

Here, we will adopt the micromagnetic approach (LLG equation) to study the influence of the temperature on some fundamental properties of both 1D (nanowires) and 2D (nanolayers) ferromagnetic nano-objects, for which thermal effects may become important. In order to model thermal fluctuations in the context of micromagnetics, in 1963 Brown [START_REF] Brown | Thermal fluctuations of a single-domain particle[END_REF] proposed to add a stochastic term to the LLG equation, in the form of a randomly fluctuating magnetic field with zero mean and variance that is proportional to the temperature. However, the Brown method suffers from a fundamental problem: while the LLG equation is valid at a mesoscopic scale and the corresponding magnetic moment m(t, x) represents an average over many atomic spins, thermal fluctuations occur at the atomic level. Hence, one is mixing two different levels of descriptions in the same LLG equation: the mesoscopic level for the deterministic terms, and the microscopic level for the stochastic terms. Then, if applied without further corrections, this procedure entails that temperature effects (for instance, the numerically-calculated Curie temperature) depend on the computational cell, which is obviously a spurious result.

Indeed, the computational cell size ∆x is noticeably larger than the physical lattice constant a eff . Thermal fluctuations occur at the length scale a eff , but are necessarily implemented at the scale ∆x in the micromagnetic codes. This induces an error in the computed properties, in particular near the Curie temperature T C , which can be overestimated by one order of magnitude or more [START_REF] Grinstein | Coarse Graining in Micromagnetics[END_REF].

In order to mitigate this spurious effect, several strategies have been proposed. In one approach [START_REF] Kirschner | Cell size corrections for nonzero-temperature micromagnetics[END_REF][START_REF] Kirschner | Relaxation times and cell size in nonzero-temperature micromagnetics[END_REF], the magnetization at saturation M s is scaled with the temperature and the computational cell size, using a Bloch-like law that is similar to the well-known temperature dependence of M s . By comparing the micromagnetic results to those obtained with an atomistic model, the authors of [START_REF] Kirschner | Cell size corrections for nonzero-temperature micromagnetics[END_REF][START_REF] Kirschner | Relaxation times and cell size in nonzero-temperature micromagnetics[END_REF] obtain a difference of less than 1% in the estimation of the equilibrium magnetization, for a temperature T = 0.38 T C and a computational cell size ∆x = 1.5 nm.

Another approach [START_REF] Grinstein | Coarse Graining in Micromagnetics[END_REF] consists in defining a rescaled temperature to take into account the fact that thermal averages (coarse graining) are performed on a computational cell that is larger than the lattice constant. In particular, the effective exchange constant varies with the size of the coarsegraining block, and this dependence should be taken into account. Grinstein and Koch [START_REF] Grinstein | Coarse Graining in Micromagnetics[END_REF] used renormalization-group techniques to unravel this dependence. In the same spirit, Hahn [START_REF] Hahn | Temperature in micromagnetism : cell size and scaling effects of the stochastic Landau-Lifshitz equation[END_REF] proposed a simple scaling law between the physical temperature and a "numerical" temperature to be used in the micromagnetic code, which depends on the ratio between the computational cell size and the lattice constant. This method was tested for nickel, cobalt and iron objects using the OOMMF code [START_REF] Donahue | OOMMF User's Guide, Version 1.0[END_REF], and yielded Curie temperature that were virtually independent on the computational cell and very close to the experimental values.

A possible limitation of the LLG approach is that the amplitude of the local magnetic moment |m(t, x)| remains constant in time, which is not necessarily true at high temperatures, notably near T C . Chubykalo-Fesenko et al. [START_REF] Chubykalo-Fesenko | Dynamic approach for micromagnetics close to the Curie temperature[END_REF] have investigated this issue using an atomistic time-dependent model and indeed they found that the modulus of the magnetization varies in time (see figure 1 in [START_REF] Chubykalo-Fesenko | Dynamic approach for micromagnetics close to the Curie temperature[END_REF]). However, this variation is limited to a dip during an initial transient, after which |m(t, x)| recovers approximately its initial value. As our results are obtained by taking time-averages at longer times, this variation should not be too significant. But indeed, when studying transient phenomena, it may be necessary to take this effect into account, for instance by using a Landau-Lifshitz-Bloch approach, as suggested in [START_REF] Chubykalo-Fesenko | Dynamic approach for micromagnetics close to the Curie temperature[END_REF][START_REF] Atxitia | Fundamentals and applications of the Landau-Lifshitz-Bloch equation[END_REF].

In the present work, we adopt Hahn's method to model thermal effects [START_REF] Hahn | Temperature in micromagnetism : cell size and scaling effects of the stochastic Landau-Lifshitz equation[END_REF] and use it to study the dependence of the magnetization law (in particular the Curie temperature) with the size of the system under consideration. We will focus on two nano-objects, namely one-dimensional (1D) nanowires and 2D nanolayers. Theoretical considerations [START_REF] Fisher | Scaling Theory for Finite-Size Effects in the Critical Region[END_REF] indicate that the Curie temperature follows a power-law of the type:

T C (∞) -T C (d) T C (∞) = ξ 0 d λ ,
where d is the smallest size of the system, T C (∞) and T C (d) are the Curie temperatures of the bulk and of the finite system respectively, ξ 0 is the correlation length at zero temperature, and λ is the critical exponent. The main purpose of this work will be to validate the above law and obtain the exponent and the correlation length from micromagnetic simulations with thermal effects.

From an experimental point of view, several works considered this problem. Early studies on thin nickel films, cobalt films and Co 1 Ni 9 alloys yielded critical exponents λ = 1.25, 1.34 and 1.39, respectively [START_REF] Schneider | Curie temperature of ultrathin films of fcc-cobalt epitaxially grown on atomically flat cu(100) surfaces[END_REF][START_REF] Huang | Finite-size scaling behavior of ferromagnetic thin films[END_REF], while measurements on nickel films [START_REF] Li | Dimensional crossover in ultrathin Ni(111) films on W(110)[END_REF] revealed an exponent λ = 1.4. Later work on nickel nanowires [START_REF] Sun | Finite-size effects in nickel nanowire arrays[END_REF], with diameters ranging from 30nm to 500nm, yielded ξ 0 = 2.2 nm and λ = 0.94. More recent works [START_REF] Zhang | Thickness-Dependent Curie Temperatures of Ultrathin Magnetic Films: Effect of the Range of Spin-Spin Interactions[END_REF][START_REF] Wang | Finite-size scaling behavior and intrinsic critical exponents of nickel: Comparison with the three-dimensional Heisenberg model[END_REF] mention larger values for the exponent, up to λ = 2.8.

On the theory front, statistical estimations based on a first-principles-based Monte Carlo approach yielded a critical exponent λ = 1.47 for Pb(Zr 0.5 Ti 0.5 O 3 ) (PZT) thin films [START_REF] Almahmoud | Dependence of Curie temperature on the thickness of an ultrathin ferroelectric film[END_REF]. Similar values were obtained using other approaches [START_REF] Yang | Size and interface effects on critical temperatures of ferromagnetic, ferroelectric and superconductive nanocrystals[END_REF]. These exponents should be compared to the theoretical values predicted by the 3D Heisenberg and Ising models (respectively, λ = 1.4 and λ = 1.58). In contrast, an atomistic mean-field model [START_REF] Penny | Mean-field modelling of magnetic nanoparticles: The effect of particle size and shape on the Curie temperature[END_REF] yielded larger values, close to λ = 2 (range: λ ∈ [1.82 -2.17]), for magnetite nanoparticles of different shapes and sizes. As we shall see, our own work suggests a critical exponent close to λ ≈ 2, for both nanowires (1D) and nanolayers (2D). The observed correlation length is found to be ξ 0 ≈ 3 nm for nanowires and ξ 0 ≈ 1.6 nm for nanolayers.

We further note that our results are obtained using a time-dependent model, in contrast to statistical and Monte Carlo approaches used in other studies [START_REF] Almahmoud | Dependence of Curie temperature on the thickness of an ultrathin ferroelectric film[END_REF]. In other words, we solve the timedependent LLG equation with thermal effects and, once a fluctuating equilibrium is reached, we measure the relevant magnetic properties by performing ensemble averages over many realizations and/or time averages over a certain duration. This approach is less computationally expensive than fully atomistic simulations. In addition, it is amenable to investigating the dynamical properties of magnetism, such as the propagation of domain walls and other transient effects, which will make the object of future work. Here, we have used this time-dependence to show that statistical fluctuations explose near T C , confirming the presence of a phase transition at that temperature.

The present paper is organized as follows. Section 2 details the mathematical setting, namely the LLG equation at finite temperature. After recalling the various terms involved in the effective magnetic field in Sec. 2.1, the following subsections are devoted to the modelling of thermal effects through a stochastic magnetic field (Sec. 2.2) and to the implementation of the temperature scaling [START_REF] Hahn | Temperature in micromagnetism : cell size and scaling effects of the stochastic Landau-Lifshitz equation[END_REF] (Sec. 2.3). Section 3 contains the details of the numerical scheme. Section 4 is devoted to the validation of the numerical code with several test cases (Sec. 4.1), clearly proving that the Curie temperature does not depend on the computational cell size or the time step (Sec. 4.2). The validity of the Bloch law (at low temperatures, T ≪ T C ) and Curie law (at T ≲ T C ) are also tested (Sec. 4.3). Finally, Section 5 contains the main results of this work, namely the size dependence of the Curie temperature for two types of nano-objects: 1D nanowires and 2D nanolayers. Conclusions are drawn in Section 6.

Micromagnetic model at finite temperature

The present numerical study focuses on a ferromagnetic domain modeled either as an infinite nanowire along the e x axis, where (e x , e y , e z ) is the canonical orthonormal basis in R 3 , or an infinite nanolayer in the (e x , e y ) plane. For all times t ≥ 0 and positions x ∈ R 3 , let m(t, x) ∈ S 2 = {m ∈ R 3 , ∥m∥ = 1} be the magnetic moment vector field normalized to the saturation magnetization M s . Here, S 2 is the unit sphere. The precession dynamics of m(t, x) is described by the Landau-Lifshitz-Gilbert (LLG) equation:

∂m ∂t = -γ 0 m × H eff -γ 0 α m × (m × H eff ) , (1)
where H eff is the effective magnetic field. Here, γ 0 = γµ 0 > 0 is the scaled gyromagnetic ratio, with γ = e/2m (where e > 0 and m are the charge and mass of the electron, respectively), and α > 0 is the dimensionless damping constant [START_REF] Gilbert | A phenomenological theory of damping in ferromagnetic materials[END_REF]. Table 1 lists all the physical variables used in this work, their units, and their numerical values.

Effective magnetic field H eff

The effective magnetic field H eff results from the sum of the exchange field H exch and the anisotropy field H ani :

H eff = H exch + H ani .
The exchange field is due to the Heisenberg exchange interaction and is written as

H exch = 2A µ 0 M s ∆m, (2)
with A > 0 the exchange constant and µ 0 > 0 the vacuum permeability (see Table 1).

The anisotropy field is due to the existence of preferred directions in the magneto-crystalline structure of the material. Throughout the following, two cases of anisotropy field will be considered: uniaxial anisotropy (for cobalt systems) and cubic anisotropy (for nickel and iron systems), whose expressions are given below:

H ani, uniaxial = 2K µ 0 M s (e x • m)e x , (3a)
H ani, cubic = - 2K µ 0 M s (i,j,k)∈I (e j • m) 2 + (e k • m) 2 + (e j • m) 2 (e k • m) 2 (e i • m)e i , (3b)
where I = {(x, y, z), (y, z, x), (z, x, y)} and K > 0 is the anisotropy constant (assumed identical in all three directions of space for the cubic case), see Table 1.

In the following, we will assume that the ferromagnetic domain is not subjected to any external magnetic field, so that no Zeeman energy is present. The demagnetizing field (due to the magnetic field generated by the nanowire or nanolayer itself) and the dipolar interactions are also not taken into account.

Universal constants

γ gyromagnetic ratio 1.76 × 10 11 rad s -1 T -1 µ 0 vacuum permeability 4π × 10 -7 NA -2 γ 0 rescaled gyromagnetic ratio γµ 0 mA -1 s -1 k B Boltzmann constant 1.38 × 10 -23 JK -1 Magnetic bulk parameters Cobalt Iron Nickel A exchange constant 3 × 10 -11 Jm -1 2.1 × 10 -11 Jm -1 9 × 10 -12 Jm -1 K anisotropy constant 5.2 × 10 5 Jm -3 4.8 × 10 4 Jm -3 -5.7 × 10 3 J.m -3 M s
saturation magnetization 1.4 × 10 6 Am -1 1.7 × 10 6 Am -1 4.9 × 10

Thermal fluctuations

The deterministic LLG equation considered above is valid in the zero-temperature regime. However, thermal effects obviously influence the magnetic properties, first and foremost by cancelling out the spontaneous magnetization of a ferromagnetic material above a certain critical temperature (Curie temperature T C ). The material then goes from a ferromagnetic to a paramagnetic state when T C is crossed.

In order to model thermal fluctuations, an additional random field is added to the effective magnetic field, following the idea of W.F. Brown [START_REF] Brown | Thermal fluctuations of a single-domain particle[END_REF], so that we have:

H eff = H exch + H ani + H stocha .
The random thermal field H stocha is an isotropic Gaussian white-noise vector process of variance ν 2 ∈ R. More precisely, H stocha may be written as:

H stocha (t)dt = ν dW (t), where W (t) = W 1 (t), W 2 (t), W 3 (t)
T is a classical time-continuous Wiener process.

Wiener process. The main properties of the stochastic field W (t) (and so H stocha (t)) are listed below, denoting ⟨•⟩ = E(•) the statistical average [START_REF] Klughertz | Autoresonant control of the magnetization switching in single-domain nanoparticules[END_REF][START_REF] Aquino | Midpoint numerical technique for stochastic Landau-Lifshitz-Gilbert dynamics[END_REF][START_REF] Ragusa | Full micromagnetic numerical simulations of thermal fluctuations[END_REF]:

• Space homogeneity: W only depends on t and not on x,

• Continuous time random process: W (t), ∀t ≥ 0,

• Null mean: ⟨W (t)⟩ = 0, ∀t ≥ 0,

• Decorrelated spatial components and vanishingly small autocorrelation time:

⟨W i (t)W j (t ′ )⟩ = δ ij δ(t -t ′ ), ∀i, j ∈ {1, 2, 3}
(the indices of spatial components) and t, t ′ ∈ R + . Here, δ(•) is the Dirac distribution and δ ij is the Kronecker symbol.

Variance. The standard deviation ν (and thus the variance ν 2 ) is directly related to the temperature T thanks to the following relation, obtained by the fluctuation-dissipation theorem [START_REF] Brown | Thermal fluctuations of a single-domain particle[END_REF][START_REF] Ragusa | Full micromagnetic numerical simulations of thermal fluctuations[END_REF] adapted here to the expression of Eq. ( 1):

ν 2 = 2αk B T γ 0 µ 0 M s V , (4)
where k B is the Boltzmann constant (see Table 1) and V stands for a characteristic volume that depends on the internal crystalline structure of the material (which can be a Face-Centered Cubic lattice (FCC) as in the case of nickel, a Body Centered Cubic lattice (BCC) as in iron, or a Hexagonal Close-Packed lattice (HCP) for cobalt). Following the notation of [START_REF] Hahn | Temperature in micromagnetism : cell size and scaling effects of the stochastic Landau-Lifshitz equation[END_REF], the shorter lattice distance is called the lattice constant a eff and the corresponding characteristic volume is V = a 3 eff .

Stochastic Landau-Lifshitz-Gilbert equation. Consequently, the Landau-Lifshitz-Gilbert equation ( 1) is modified to take the form of a stochastic partial differential equation (SPDE)

dm = -γ 0 m × (H eff dt + νdW ) -γ 0 αm × [m × (H eff dt + νdW )] . (5)
All physical constants used in the forthcoming simulations are summarized in Table 1.

In order to preserve the constraint that the magnetic moment m be on the unit sphere, i.e. ∥m(t, x)∥ = 1, ∀x ∈ R 3 and t ≥ 0, we interpret the above stochastic LLG equation in the Stratonovich sense; see [START_REF] Labbé | Stochastic modelling of thermal effects on a ferromagnetic nano particle[END_REF] for the issues posed by using Itô's approach.

Temperature scaling with the computational cell size

The variance of the stochastic magnetic field used to model thermal effects depends not only on the temperature T , but also on the characteristic volume V , see Eq. ( 4). In principle, this volume is related to the lattice constant, i.e. V = a 3 eff , but in practice a eff is much smaller than the spatial step ∆x used in the simulations. However, if one takes instead V = ∆x 3 , the simulation results will depend on the computational cell size, which is not an acceptable situation. At a fundamental level, this is due to the fact that the effective exchange constant varies with the size of the coarse-graining block [START_REF] Grinstein | Coarse Graining in Micromagnetics[END_REF].

In order to suppress this unwanted numerical effect, we follow the procedure recently suggested by Hahn [START_REF] Hahn | Temperature in micromagnetism : cell size and scaling effects of the stochastic Landau-Lifshitz equation[END_REF], which relies on a scaling of the temperature with the spatial step ∆x. The argument goes as follows: near the Curie temperature, the ferromagnetic behaviour disappears because the energy density of the thermal fluctuations k B T /V , which favor random orientations of the spin, becomes similar to the energy density of the exchange interactions A|∇m| 2 , which favor magnetic order. Hence, we write:

k B T a 3 ∼ A|m| 2 a 2 , ( 6
)
where a is a characteristic length that can be either the lattice constant a eff or the computational cell size ∆x. From Eq. ( 6) it is clear that, in order for the magnetic moment to be independent on the averaging volume a 3 , the temperature must scale as T ∼ a.

Therefore, we define a "numerical" temperature as

T num = ∆x a eff T.
According to the above considerations, taking a volume V = ∆x 3 together with the temperature T num should give results that are independent on ∆x and identical to those obtained using the "physical" volume V = a 3 eff and the real temperature T . This approach was recently tested by Hahn [START_REF] Hahn | Temperature in micromagnetism : cell size and scaling effects of the stochastic Landau-Lifshitz equation[END_REF], who indeed observed near independence of T C on the cell size up to ∆x = 4 nm for ferromagnetic thin films. Following this procedure, the numerical variance of the fluctuating field is defined as

ν 2 num = 2αk B ∆x a eff T γ 0 µ 0 M s ∆x 3 , (7)
which replaces Eq. ( 4) in the simulations. This expression may also be interpreted as meaning that the numerical characteristic volume to be taken into account is V = a eff ∆x 2 .

Computational method

The stochastic LLG equation ( 5) is solved numerically using a Python code. 1 The simulations are performed on the time interval [0, t f ], with t f the final time, and on a 3D finite domain [0,

L x ] × [0, L y ] × [0, L z ].
To mimic a 1D nanowire, this numerical domain is taken with a small square crosssection in the (e y , e z ) plane: [0,

L x ] × [0, d] × [0, d], with d ≪ L x .
To mimic a 2D nanolayer, the 3D domain is taken with a small thickness in the e z direction: [0,

L x ] × [0, L y ] × [0, d], with d ≪ L x , L y .
Figure 1 illustrates those two geometries. All numerical parameters are listed in Table 2. Discretization. Numerically, let ∆t > 0 and ∆x = ∆y = ∆z > 0 be the time step and the space steps in each space direction, respectively. We define time-discrete and space-discrete points with

N = ⌊ t f ∆t ⌋, J x = ⌊ Lx ∆x ⌋, J y = ⌊ Ly ∆x ⌋, J z = ⌊ Lz ∆x ⌋ t n = n∆t, 0 ≤ n ≤ N and (x i , y j , z k ) = (i∆x, j∆x, k∆x), 0 ≤ i ≤ J x , 0 ≤ j ≤ J y , 0 ≤ k ≤ J z .
The numerical solution m n i,j,k approximates the exact one m(t n , x i , y j , z k ) on each discrete point.

As the LLG equation ( 1) is valid at mesoscopic -and not atomistic -length scales, the spatial steps ∆x, ∆y, ∆z are each much larger than the lattice constant a eff , and the continuous magnetic moment vector filed m(t, x) actually represents an average of the atomic spins over a volume ∆x∆y∆z.

= L z = d nm) ject nanolayer (L x = L y = 600 nm, L z = d nm) J x , J y , J z number of mesh points nanowire (J x = 6000, J y = J z = ⌊d⌋) nanolayer (J x = J y = 600, J z = ⌊d⌋) i, j, k space indices 0 ≤ i ≤ J x , 0 ≤ j ≤ J y , 0 ≤ k ≤ J z
Numerical variables m 0 initialization of the magnetic moments

  1 0 0   for all points (x i , y j , z k ) m 1 (t) spatial average of one realization 1 L x L y L z [0,Lx]×[0,Ly]×[0,Lz] m 1 (t,
1 t f -τ t f τ 1 L x L y L z [0,Lx]×[0,Ly]×[0,Lz] m 1 (t, x)dxdt t conv convergence time for m 1 (t) to reach its plateau state M tot inf t∈[0,t f ] |m 1 (t) -M tot | < 0.1
Table 2: Numerical parameters used in the Python code Time evolution: Heun method. For consistency with the continuous problem, the stochastic LLG equation ( 5) must be discretized using a numerical method whose solution converges to the Stratonovich continuous solution. For this purpose, the modified Heun method [START_REF] Rüemelin | Numerical Treatment of Stochastic Differential Equations[END_REF] is chosen for the time integration and a second-order finite difference method is used for the discretization of the Laplacian operator.

Following [START_REF] Klughertz | Ultrafast magnetization dynamics in magnetic nanoparticles[END_REF], for further simplicity, we rewrite Eq. ( 5) as

dm = F (m, t)dt + j∈{1,2,3} G j (m) ν dW j (t), (8)
with

F (m, t) = -γ 0 m × H eff -γ 0 α m × (m × H eff )
being the deterministic part and

G j =   G 1,j G 2,j G 3,j 
 the factor term of the stochastic process, with

G i,j = γ 0 m k ϵ i,j,k -γ 0 α(m i m j -δ ij ),
where ϵ i,j,k is the Levi-Civita symbol.

After initializing the magnetization at the initial time t = 0:

m 0 i,j,k =   m 1 m 2 m 3   0 i,j,k = m(0, x i , y j , z k ),
the (stochastic) Heun method consists in the following steps to go from the time step n to the time step n + 1:

• Generate a random vector R n ∈ R 3 according to a reduced centered normal distribution, using a pseudo-random number generator. Define ∆W n = √ ∆tR n ;

• Compute ν num , the numerical version of the standard deviation ν (more details are given below);

• Define κ 1 = F (m n , t n ) and s 1,j = G j (m n ); • Define κ 2 = F (m n + ∆tκ 1 + j∈{1,2,3} s 1,j ν num ∆W n j , t n + ∆t) and s 2,j = G j (m n + ∆tκ 1 + j∈{1,2,3} s 1,j ν num ∆W n j ); • Update m n+1 = m n + 1 2 κ 1 + 1 2 κ 2 ∆t + j∈{1,2,3} 1
2 s 1,j + 1 2 s 2,j ν num ∆W n j ;
• Renormalize the magnetic moment: m n+1 = m n+1 ∥m n+1 ∥ , so that it remains on the unit sphere S 2 .

A mid-point numerical technique would be also possible alternative to the Heun method, see [START_REF] Aquino | Midpoint numerical technique for stochastic Landau-Lifshitz-Gilbert dynamics[END_REF][START_REF] Ragusa | Full micromagnetic numerical simulations of thermal fluctuations[END_REF].

The choice of a time-explicit discretization of the Laplacian operator induces a restrictive condition on the time step ∆t and the space step ∆x to ensure the stability of the scheme: ∆t ≲ ∆x 2 /2 (Courant-Friedrichs-Lewy condition).

At the domain boundaries, we take the usual Neumann condition: ∂m/∂n = 0, where n is the outgoing normal vector.

Numerical code validation

In the forthcoming simulations, three types of ferromagnetic materials will be considered: (i) nickel with a FCC lattice and cubic anisotropy, (ii) iron with a BCC lattice and a cubic anisotropy, and (iii) cobalt with a HCP lattice and uniaxial anisotropy. For each case, we shall study both 1D nanowires and 2D nanolayers. All physical parameters are listed in Table 1.

Here, we will perform several tests to validate our numerical code. In the following sections, we will focus on the scaling of the Curie temperature with the size of the system, for each type of nanoobject. We recall that our code relies on the time-dependent LLG equation. We solve numerically this equation with a fully magnetized initial condition, where m(t = 0) is uniform and directed along the e x axis, and a given temperature. Then, we wait that magnetic moment relaxes to a lower value under the effect of the temperature and determine its value by averaging over a sufficiently long period of time.

Test-case details

Except for the following Sec. 4.2 -where the independence of the results on the space and time discretizations are tested on a 3D cube -all numerical simulations are preformed on three ferromagnetic materials (cobalt, iron and nickel) and two geometries (see Figure 1):

• 1D nanowire with small square cross-section in the (e y , e z ) plane: [0, The numerical parameters are always fixed to (see Table 2): ∆t = 10 fs, t f = 50 ps, N = 5000, ∆x = ∆y = ∆z = 1 nm.

L x ] × [0, d] × [0, d] with L x =
The initialization of the magnetic moments m 0 is chosen uniform and directed along the e x axis in all test cases:

m 0 i,j,k =   1 0 0   for all i, j, k.
With the choice of initialization, the magnetization is initially equal to 1 (when the magnetic moments are all aligned and have a norm equal to 1), then falls to zero at the Curie temperature (when the magnetic moments are randomly aligned due to the thermal noise). Since Eq. ( 5) is stochastic, the average of the magnetic moments m over the entire ferromagnetic domain may differ from one realization to another, so this average should be calculated over several realizations. Since the initialization is oriented along the e x axis and in the absence of any external magnetic field, the average of m along this direction, i.e. m 1 , is enough to characterize the magnetic state of the system. Hence, we define the total (in space) magnetization

M tot = 1 L x L y L z [0,Lx]×[0,Ly]×[0,Lz] m 1 (t f , x)dx , (9)
with ⟨•⟩ denoting the statistical average over many realizations. In order to simplify this calculation, we assume that the stochastic process of Eq. ( 5) is ergodic, so that the statistical average may be replaced by the time average for a single, sufficiently long realization. In practice, we plot several realizations, look at the time τ after which the transient regime gives way to the stationary regime, and finally take the time average from this transient time τ up to the final time t f . Thus, the total magnetization is now defined as

M tot = 1 t f -τ t f τ 1 L x L y L z [0,Lx]×[0,Ly]×[0,Lz] m 1 (t, x)dxdt = 1 t f -τ t f τ m 1 (t)dt. ( 10
)
Figure 2 illustrates the spatial average of the x component of the magnetic moment m 1 (t) = 10)).

In practice, the Curie temperature is determined numerically by plotting M tot as a function of the temperature T , and defining T C as the temperature for which M tot is smaller than a certain threshold, fixed to 0.1:

T C := argmin T M tot (T ) < 0.1.

Dependence on the numerical parameters

Here, we show that our results do not depend on either the time step ∆t or the computational cell size ∆x.

Figure 3 shows the total magnetization M tot as a function of the temperature for different cell sizes, going from 1 nm to 5 nm, for a cubic object of dimensions 50 nm × 50 nm × 50 nm, for cobalt, iron and nickel. All other parameters are identical (∆t = 10 fs, t f = 50 ps, N = 5000 and τ = 40.0 ps). The results are indeed independent on ∆x, and the computed Curie temperatures are very close to the experimental values for the bulk materials (see Table 1). A slight discrepancy starts occurring for iron at ∆x = 5 nm.

Figure 4 shows the dependence of the numerical results with respect to the time step ∆t, again for cubic nano-objects of side 50 nm, with computational cell size ∆x = 1 nm. The final time is always t f = 50 ps, so that the number of time steps is N = 5 × 10 3 for ∆t = 10 fs (blue curve), N = 10 4 for ∆t = 5 fs (red curve), and N = 2 × 10 4 for ∆t = 2.5 fs (green curve). According to the value of ∆t, the time-averaged M tot in Eq. ( 10) includes the last 1000, 2000 or 4000 time iterations. The computed Curie temperature varies only slightly with the time step, and appears to have converged for ∆t = 5 fs. 10), as a function of the temperature, for cubic cobalt (left panel), iron (right panel) and nickel (middle panel) nano-objects with dimensions 50 nm × 50 nm × 50 nm. The different symbols and colors stand for different computational cell sizes ∆x, going from 1 nm to 5 nm. The black vertical dash-dotted lines represent the bulk Curie temperatures as given in Table 1.

Magnetization curve: Bloch's law and Curie's law

In this section, we check that the numerically calculated magnetization M tot (T ) satisfies the standard Bloch's and Curie's laws, respectively at low temperatures T ≪ T C and near the Curie temperature, T ≲ T C . As we have already verified that the spatial and temporal steps do not influence the result, ∆t and ∆x will be fixed as specified in Table 2, i.e. ∆t = 10 fs and ∆x = 1 nm.

Figure 5 illustrates the behavior of the magnetization curve M tot (T ) for the two geometric configurations considered here (nanowires and nanolayers), for iron (red triangles), cobalt (blue circles) and nickel (green crosses). Results are in broad agreement with the expected magnetization curves, and the computed Curie temperatures are close to the experimental values found in the literature for bulk materials [START_REF] Kittel | Introduction to solid state physics[END_REF], see Table 1. 10), as a function of temperature, for a cobalt (left), iron (right) or nickel (middle) nanowire with dimensions 50 nm×50 nm×50 nm. The different symbols stand for different time steps ∆t = 2.5 fs (green crosses), 5 fs (red circles), and 10 fs (blue crosses). The computational grid size is ∆x = 1 nm. The black vertical dash-dotted lines represent the bulk Curie temperatures as given in Table 1.

Bloch's law states that total magnetization M tot (T ), for low temperatures, behaves as follows:

M tot (T ) ∼ T →0 1 - T T C 3/2
, which can be rewritten as: log(1 -M tot ) ∼ 3 2 log (T /T C ). Figure 6 checks this behaviour on a log-log scale, for a nanowire (left) and a nanolayer (right), with sizes corresponding to the two extreme cases in Table 3, i.e. d = 11nm and 41 nm. The theoretical 3/2 slope is represented as a solid black line and matches the numerical results quite well.

Next, we check the Curie law, valid near T C :

M tot (T ) ∼ T → < T C 1 - T T C 1/2
, which may be rewritten as: log(M tot ) ∼ 1 2 log 1 -T T C . Figure 7 shows the behavior of the magne-Figure 5: Total magnetization M tot [START_REF] Gilbert | A phenomenological theory of damping in ferromagnetic materials[END_REF] with respect to temperature with numerical parameters detailed in Table 2. The Curie temperature corresponds to the first temperature at which magnetization falls to zero. Simulation results are represented by dots, the solid curves are an interpolation based on cubic splines. tization M tot as a function of 1 -T /T C in logarithmic scale, for the same cases as those of Figure 6. Again, the numerical results match rather well the theoretical 1/2 slope. 5 Finite-size effects on the Curie temperature

This section is devoted to the study of the influence of size effect on the magnetization curve and Curie temperature for nanometric objects, both in 1D (nanowires) and 2D (nanolayers). The objective is to vary the cross-section of the nanowire or the thickness of the nanolayer and study the variations induced in the Curie temperature. Throughout this section, numerical parameters are chosen as listed in Table 2.

Size effects on the magnetization curve

For all tested geometries, the computed Curie temperatures are close enough to the experimental values reported in the literature. This is achieved thanks to the scaling of the fluctuating thermal field as detailed in Sec. 2.3. Nevertheless, we observe small variations with the size parameter d, which corresponds to the side of the square cross-section of a nanowire or the thickness of a nanolayer. Figure 8 (for nanowires) and Figure 9 (for nanolayers) show the magnetization curves M tot (T ), and illustrate how the Curie temperature increases with increasing size d. The computed values of T C are summarized in Table 3. We also observe greater variability in M tot (T ) for nanowires than for nanolayers. Size effects then appear to be stronger for lower-dimensional structures.

Power-law scaling of the Curie temperature

Theoretical considerations [START_REF] Fisher | Scaling Theory for Finite-Size Effects in the Critical Region[END_REF] indicate that the Curie temperature should vary with the size d following a power-law of the type: where λ is the critical exponent, ξ 0 is the correlation length at zero temperature, T C (∞) is the Curie temperature of the bulk, and T C (d) the Curie temperature of a finite nano-object of size d. This type of power-law has been observed in many experiments [START_REF] Schneider | Curie temperature of ultrathin films of fcc-cobalt epitaxially grown on atomically flat cu(100) surfaces[END_REF][START_REF] Huang | Finite-size scaling behavior of ferromagnetic thin films[END_REF][START_REF] Li | Dimensional crossover in ultrathin Ni(111) films on W(110)[END_REF][START_REF] Sun | Finite-size effects in nickel nanowire arrays[END_REF][START_REF] Zhang | Thickness-Dependent Curie Temperatures of Ultrathin Magnetic Films: Effect of the Range of Spin-Spin Interactions[END_REF][START_REF] Wang | Finite-size scaling behavior and intrinsic critical exponents of nickel: Comparison with the three-dimensional Heisenberg model[END_REF] and numerical simulations [START_REF] Almahmoud | Dependence of Curie temperature on the thickness of an ultrathin ferroelectric film[END_REF][START_REF] Penny | Mean-field modelling of magnetic nanoparticles: The effect of particle size and shape on the Curie temperature[END_REF]. Experimental works yielded correlation length ξ 0 of the order of a few In the analysis of our simulation results, the Curie temperature of the bulk T C (∞) is in fact replaced by the Curie temperature of the largest structure that we consider T C (d max ), that is d max = 41 nm, see Table 3. Figure 10 (for nanowires) and Figure 11 (for nanolayers) illustrate this powerlaw behaviour for the three materials considered here, both on a linear scale (left panels) and on a logarithmic scale (right panels). Blue circles correspond to the numerical results T C (d) extracted from Table 3. The red solid lines correspond to the theoretical power-law [START_REF] Grinstein | Coarse Graining in Micromagnetics[END_REF]. The exponent λ is deduced by fitting a log-log straight line through the numerical points using a least-square method, and then ξ 0 is obtained by finding the intercept λ log ξ 0 of this line with the vertical axis. The last two or three points further to the right deviate from the power-law, and were therefore discarded in the fitting procedure.

T C (∞) -T C (d) T C (∞) = ξ 0 d λ , (11)
The computed values of the correlation length ξ 0 and the critical exponent λ are reported on each figure. For nanowires, ξ 0 ranges from 2.77 nm to 3.02 nm, while for nanolayers it varies from 1.55 nm to 1.82 nm. In both cases, the smallest value is found to be for iron, whereas the largest values is for nickel. These values are broadly in good agreement with those observed in the experiments.

As to the critical exponent, we obtain λ = 2.12 -2.14 (a range that is probably within the numerical uncertainty) for all three materials in the case of nanowires. For nanolayers, we obtain λ = 1.90, 1.92 and 2.00, for cobalt, iron and nickel, respectively. It is remarkable that all exponents are close to λ = 2, which is also the value observed in atomistic mean-field simulations [START_REF] Penny | Mean-field modelling of magnetic nanoparticles: The effect of particle size and shape on the Curie temperature[END_REF].

Nonequilibrium properties

Although we have used our code to study steady-state phenomena, such as the dependence of the magnetization curve and Curie temperature with size, our approach is fundamentally time- dependent. Indeed, to obtain our results we solved numerically the time-dependent LLG equation (with thermal effects) and then deduced steady-state quantities through time and/or ensemble averages. However, the code could be used to investigate nonequilibrium properties, such as transient behaviors, the motion of domain walls, instabilities, etc. Over more commonly used statistical methods, this dynamical approach has the further advantage of being well-suited to study the fluctuations that appear in the vicinity of phase transitions. These fluctuations play a key role in finite systems, particularly when their dimensionality is low. This vast realm is left to future investigations. Here, we conclude our work by searching for a signature of the phase transition occurring at T C in the dynamical behaviour of the magnetization. We do this by looking at the convergence time towards the equilibrium state.

In Figure 2, we observed that the convergence time to the plateau state depends on the temperature. More precisely, we define the convergence time as the first time at which the x component of the spatially-averaged magnetic moment m 1 (t) becomes sufficiently close to the total magnetization M tot (which corresponds to its plateau):

t conv = inf t∈[0,t f ] |m 1 (t) -M tot | < ε, (12)
with the tolerance parameter ε fixed to 0.1. Figure 12 illustrates schematically this convergence time t conv , and also shows the transient time τ , which is the time used to compute the average magnetization, see Eq. [START_REF] Gilbert | A phenomenological theory of damping in ferromagnetic materials[END_REF].

Figure 13 shows the convergence time, as a function of the temperature, for nanowires of size 6000 nm×d nm ×d nm and nanolayers of size 600 nm×600 nm ×d nm, for two values of d. First, we note that the convergence time t conv is always smaller than the transient time τ taken to compute the average M tot (t conv has a maximum value around 20 ps, which is always smaller than τ = 40 ps). This confirms that convergence to the plateau state takes place before τ and hence that the calculation of M tot is correct. This figure also shows that the convergence time peaks near the Curie temperature [START_REF] Maurat | Thermal properties of open-shell metal clusters[END_REF]. This is a dynamical signature that large fluctuations near T C translate into longer times taken by the system to relax to its equilibrium state. The rapid increase of the relaxation time close to T C is known as critical slowing down [START_REF] Peczak | Dynamical critical behavior of the three-dimensional Heisenberg model[END_REF][START_REF] Chubykalo-Fesenko | Dynamic approach for micromagnetics close to the Curie temperature[END_REF], an effect which is characteristic of second-order phase transitions. Moreover, Figure 13 illustrates the effect of size on this phase transition. The structures with the smallest size (d = 11 nm, blue curves) have longer convergence times t conv and larger widths than the larger structures (d = 41 nm, orange curves). Also, the convergence times and the widths are larger for 1D nanowires than for 2D nanolayers. One can deduce that fluctuations near the transition temperature T C are stronger for smaller and lower-dimensional structures.

Conclusion

In this work, we developed a computational code (written in Python) that solves the LLG equation of micromagnetism, including finite-temperature effects. By adopting an appropriate temperature scaling with the computational cell size ∆x [START_REF] Hahn | Temperature in micromagnetism : cell size and scaling effects of the stochastic Landau-Lifshitz equation[END_REF], it was possible to recover magnetization curves and Curie temperatures that match closely the experimental ones for cobalt, nickel and iron nanoobjects. Compared to fully atomistic simulations, our micromagnetic approach has the advantage of being less computationally expensive, as it relies on a mesoscopic description at a scale ℓ much larger than the lattice constant: ℓ ≫ ∆x ≫ a eff . It is also more easily amenable to dynamical simulations to study, for instance, the motion of domain walls or other time-dependent phenomena.

Using this accurate computational tool, we investigated the size-dependence of the Curie temperature for 1D nanowires and 2D nanolayers, by varying the smallest size d of the system. We confirmed that the difference between the computed finite-size T C and the bulk T C follows a powerlaw of the type: (ξ 0 /d) λ , where ξ 0 is the correlation length at zero temperature and λ is a critical exponent. We obtained values of ξ 0 in the nanometer range, also in accordance with other simulations and experiments. The computed critical exponent was found to be close to λ = 2 for all considered materials and geometries, which is the expected result for a mean-field approach, but slightly larger than the values observed experimentally. Finally, the time-dependent model developed here represents an effective tool for studying thermal fluctuations near the ferromagnetic phase transition.

All in all, the behaviour of 1D and 2D ferromagnetic nano-objects, as a function of both temperature and size, was recovered with good precision and a relatively low computational cost in comparison to fully atomistic simulations. This computational tool may therefore be applied in the future to more complex configurations, involving for instance 3D structures and dynamical effects.

Figure 1 :

 1 Figure 1: Illustration of the two generic geometries coresponding to a 1D nanowire (left) and a 2D nanolayer (right)

  6000 nm and 11 nm ≤ d ≤ 41 nm; • 2D nanolayer with small thickness in the e z axis: [0, L x ] × [0, L y ] × [0, d] with L x = L y = 600 nm and 11 nm ≤ d ≤ 41 nm.

Figure 2 :

 2 Figure 2: Spatially averaged x component of the magnetic moment m 1 (t) for one statistical realization, as a function of time t for a ferromagnetic nanowire. Colors correspond to different temperatures, as indicated on the figure. Top panels: cobalt; middle panels: iron; bottom panels: nickel. The left column corresponds to nanowires with dimensions (in nm): 1680 × 11 × 11, the right column to nanowires with dimensions (in nm): 120 × 41 × 41.

Figure 3 :

 3 Figure 3: Total magnetization M tot , from Eq. (10), as a function of the temperature, for cubic cobalt (left panel), iron (right panel) and nickel (middle panel) nano-objects with dimensions 50 nm × 50 nm × 50 nm. The different symbols and colors stand for different computational cell sizes ∆x, going from 1 nm to 5 nm. The black vertical dash-dotted lines represent the bulk Curie temperatures as given in Table1.

Figure 4 :

 4 Figure 4: Total magnetization M tot , from Eq. (10), as a function of temperature, for a cobalt (left), iron (right) or nickel (middle) nanowire with dimensions 50 nm×50 nm×50 nm. The different symbols stand for different time steps ∆t = 2.5 fs (green crosses), 5 fs (red circles), and 10 fs (blue crosses). The computational grid size is ∆x = 1 nm. The black vertical dash-dotted lines represent the bulk Curie temperatures as given in Table1.

Figure 6 :

 6 Figure 6: Bloch's law. Behavior of 1 -M tot (T ) as a function of T /T C in logarithmic scale. The materials are represented by different colors: iron (red), cobalt (blue) and nickel (green). The left panel corresponds to nanowires and the right panel to nanolayers, with sizes d = 11 nm (crosses) and d = 41 nm (dots). Each curve was multiplied by a multiplicative factor for easier reading. The black solid lines represent the theoretical 3/2 slope.

Figure 7 :

 7 Figure 7: Curie's law. Behavior of M tot (T ) as a function of 1 -T /T C in logarithmic scale. The materials are represented by different colors: iron (red), cobalt (blue) and nickel (green). The left panel corresponds to nanowires and the right panel to nanolayers, with sizes d = 11 nm (crosses) and d = 41 nm (dots). Each curve was multiplied by a multiplicative factor for easier reading. The black solid lines represent the theoretical 1/2 slope.

Figure 8 :

 8 Figure 8: Total magnetization M tot (T ) as a function of the temperature for nanowires of different materials, with cross-section sizes going from d = 11 nm to d = 41 nm. Top panels: cobalt; middle panels: iron; bottom panels: nickel. The right column is a zoom near the Curie temperature.

Figure 9 :

 9 Figure 9: Total magnetization M tot (T ) as a function of the temperature for nanolayers of different materials, with thiknesses going from d = 11 nm to d = 41 nm. Top panels: cobalt; middle panels: iron; bottom panels: nickel. The right column is a zoom near the Curie temperature.

Figure 10 :

 10 Figure 10: Behavior of the Curie temperature as a function of the cross-section size d, for nanowires of different materials: cobalt (top panels), iron (middle panels) and nickel (bottom panels). Numerical results are shown as blue dots, while the solid blue and red lines represent the theoretical power-law of Eq. (11). The left column shows results on a linear scale, while the right column on a log-log scale. The correlation length ξ 0 and critical exponent λ can be read on each figure of the right column as: (ξ 0 /d) λ .

Figure 11 :

 11 Figure 11: Behavior of the Curie temperature as a function of the thickness d, for nanolayers of different materials: cobalt (top panels), iron (middle panels) and nickel (bottom panels). Numerical results are shown as blue dots, while the solid blue and red lines represent the theoretical power-law of Eq. (11). The left column shows results on a linear scale, while the right column on a log-log scale. The correlation length ξ 0 and critical exponent λ can be read on each figure of the right column as: (ξ 0 /d) λ .

Figure 12 :

 12 Figure 12: Evolution of the spatially-averaged x component of the magnetic moment m 1 (t) for a nickel nanowire of cross-section d = 11 nm, for three different temperatures, T = 510 K, T = 587 K and T = 650 K. The corresponding convergence times (t conv = 1.9 ps, 6.2 ps, and 11.2 ps) are shown on the figure. We also show the transient time τ = 40 ps, which is the same for all simulations (time used to average the total magnetization M tot ).

Figure 13 :

 13 Figure 13: Convergence time t conv , from Eq. (12), as a function of the temperature for cobalt (top panels), iron (middle panels) and nickel (bottom panels). Left column: nanowire geometry with 6000 nm×d nm×d nm; right column: nanolayer geometry with 600 nm×600 nm×d nm. Two values of d are considered: d = 11 nm (blue solid curve) and d = 41 nm (orange dash curve). The vertical lines correspond to the numerical Curie temperatures reported in Table 3, for d = 11 nm (blue vertical line) and d = 41 nm (orange vertical line).

Table 1 :

 1

	5 Am -1

Top: Values of the universal constants used in this work. Bottom: Magnetic parameters for bulk cobalt, iron, and nickel. Sources:

[START_REF] Kittel | Introduction to solid state physics[END_REF][START_REF] Davey | Precision Measurements of the Lattice Constants of Twelve Common Metals[END_REF][START_REF] Ono | Determination of lattice parameters in Hcp cobalt by using X-Ray bond's method[END_REF]

.

Table 3 :

 3 ).Nanowire L x × d × d with L x = 6000 nm Nanolayer L x × L y × d with L x = L y = 600 nm Curie temperatures (in Kelvin) obtained from the numerical simulations, for nanowires (left) and nanolayers (right) of different sizes and different materials.

	d [nm] Cobalt	Iron	Nickel	d [nm] Cobalt	Iron	Nickel
	41	1427.26	1090.05	579.62	41	1433.91	1095.33	582.39
	36	1423.74	1087.57	578.13	36	1432.20	1093.41	581.67
	31	1421.56	1085.81	577.23	31	1430.68	1092.13	581.01
	26	1414.62	1082.05	573.87	26	1428.68	1091.49	579.81
	21	1406.83	1076.85	571.05	21	1422.79	1087.65	578.07
	20	1401.23	1072.13	569.31	20	1421.75	1087.25	577.71
	19	1397.04	1073.01	568.05	19	1421.27	1087.49	577.41
	18	1395.52	1069.65	566.91	18	1419.85	1085.01	576.33
	17	1391.72	1069.17	564.45	17	1416.14	1083.41	575.61
	16	1385.55	1062.61	561.69	16	1415.10	1082.05	574.71
	15	1382.13	1062.45	560.67	15	1413.20	1080.53	573.63
	14	1375.19	1053.65	558.15	14	1410.35	1080.13	572.73
	13	1360.18	1045.89	552.69	13	1406.55	1076.37	570.87
	12	1352.11	1043.65	549.33	12	1400.37	1072.21	568.95
	11	1342.04	1036.53	544.83	11	1398.09	1071.33	567.33

Table 3 ,

 3 for d = 11 nm (blue vertical line) and d = 41 nm (orange vertical line). France matthieu.boileau@math.unistra.fr raphaël côte Université de Strasbourg & University of Strasbourg Institute for Advanced Study CNRS, IRMA UMR 7501, 7 rue René Descartes, 67084 Strasbourg, France raphael.cote@unistra.fr paul-antoine hervieux Université de Strasbourg CNRS, Institut de Physique et Chimie des Matériaux de Strasbourg, UMR 7504 F-67000 Strasbourg, France paul-antoine.hervieux@ipcms.unistra.fr giovanni manfredi Université de Strasbourg CNRS, Institut de Physique et Chimie des Matériaux de Strasbourg, UMR 7504 F-67000 Strasbourg, France giovanni.manfredi@ipcms.unistra.fr

	clémentine courtès
	Université de Strasbourg & Inria Nancy Grand Est
	CNRS, IRMA UMR 7501,
	7 rue René Descartes,
	67084 Strasbourg, France
	clementine.courtes@unistra.fr
	matthieu boileau
	CNRS & Inria Nancy Grand Est
	CNRS, IRMA UMR 7501,
	67084 Strasbourg,

Available at: https://gitlab.math.unistra.fr/llg3d/llg3d

Acknowledgements

The authors gratefully acknowledge partial support by the ANR project MOSICOF ANR-21-CE40-0004. R. Côte acknowledges support from the University of Strasbourg Institute for Advanced Study (USIAS) for a Fellowship within the French national programme "Investment for the future" (IdEx-Unistra). The autors would like to thank Jérôme Lelong and Bertrand Dupé for a suggestion that helped improve the computer code, and Guillaume Ferriere, Ludovic Godard-Cadillac and Yannick Privat for the enriching discussions throughout the project. We also thank Riccardo Hertel for his thorough reading of the manuscript and many insightful comments.