X Tellier
email: xavier.tellier@enpc.fr

C Douthe

O Baverel

L Hauswirth

Designing funicular grids with planar quads using i-liwein surfaces

Keywords: Architectural geometry, funicular structures, gridshells, cable-nets, structural design, Monge-Ampère equation, discrete differential geometry

The design of curved structural building envelopes is challenging as it requires to account for a wide variety of constraints. In particular, the shape must be mechanically efficient, fabricable, and fit the site geometry. This article shows how a family of surfaces, called isotropic Linear Weingarten surfaces (nicknamed i-liwien), may fulfil all these constraints together, and be used as an intuitive design tool. We start by showing that these shapes are funicular for a uniform vertical load, and that principal stress lines form a conjugate net. This allows in particular for the design of gridshells with planar faces and low bending moments or for the design of self-stressed cable nets cladded by planar glass panels. We then propose a discrete model based on recent advances in discrete differential geometry. We use this model to propose a generation method from boundary curves with two additional control parameters. We demonstrate the shape potential on several examples. The application to shells, membranes, cable-nets and gridshells are discussed.

Curved structural envelopes combine mechanical efficiency, architectural potential and a wide design space to fulfil functional and environment goals. However, curvature complicates manufacturing, design and assembly, such that cost and construction time are usually significantly higher than for parallelepipedic structures. Furthermore, curvature does not guarantee structural efficiency, but may actually result in the opposite if mechanical considerations are forgotten.

Manufacturing complexity may be reduced by making sure that the structure fulfils geometrical properties specific to the fabrication method. Similarly, mechanical performance may be obtained by structural optimization or form-finding of funicular shapes, which are shapes that can resist a given load with pure axial forces. However, considering the other constraints of a project such as functional requirements and aesthetic, accounting for both fabricability and mechanics leads to highly constrained design problems. In this problem, not only the shape -the surface described by the structure -is important, but also the structural pattern, a term by which we refer to both the orientation and the connectivity of the structural components.

The high complexity of this problem can be a barrier to the use of curved structures on projects, a barrier that happens at two levels for a designer: firstly, what algorithms work and how to use them? Secondly, how to apprehend intuitively the possible shapes?

It is from this perspective that we focus in this article on isotropic Linear Weingarten, that we will nickname i-liwein for brevity, as they happen to have the following properties:

• They offer a rich design space, that contains some popular architectural shapes such as parabolic hyperboloids. The shapes are, we think, smooth and harmonious (see for example Figure 1 and Figure 2).

• They can be covered by a quadrangular pattern with planar faces, such that the edges constitute a funicular beam structure under uniform vertical load -thus combining mechanical performance with a geometrical property which is highly convenient for fabrication of cladding panels. This pattern appears from the so-called isotropic principal curvature lines, and tend to form smooth networks.

• Their behaviour and generation is intuitive and well suited to architectural context: They can be controlled from closed boundary curves, letting then to the designer two degrees of freedom -the first one corresponding to the "inflation" of the shape, the second to the "bumpiness" of its crest line (Figure 1).

Overview of paper

After a literature review in subsection 1.1, section 2 presents i-liwein surfaces and the basics of isotropic geometry -the geometrical point of view from which they emerge. In section 3, we show that i-liwein surfaces have an interesting mechanical property: they are funicular under a uniform vertical load, and the principal membrane stress directions then form a conjugate net. This makes them well-suited to design mechanically efficient support structures covered with planar quadrangles.

In section 4, we show how they can be generated on target boundaries. The generation is equivalent to solving a Dirichlet problem of the Monge-Ampere equation. This method allows to control the shape by its boundary curve and by two degrees of freedom. Then, in section 5, we highlight some design applications to membranes (Figure 2), cable-nets (Figure 3), gridshells (Figure 1) and funicular shells. Section 6 presents two case studies that investigate the mechanical performance of i-liwein shells for various values of the parameter b, considering in particular non-uniform loads. We finish in section 7 by a discussion of design applications and performances. Before getting started on i-liwein surfaces, we will be reviewing here the design methods for double-curvature structures that account for both mechanical efficiency and fabricability. We subdivide the literature based on how the topology and shape are dealt with, as this is the main difference of our approach with the core of the literature on the topic.

From topology to shape

Many design methods take topology as an input, and explore the possible shapes. For example, hanging models, which have been used abundantly for shell and membrane design by masterbuilders such as Heinz Isler, Frei Otto or Antonio Gaudi, give funicular shapes that correspond to a specific net configuration. Using a hanging cable net with constant edge length gives a shape that is both funicular and realizable as an elastic gridshell -a strategy used to find the form of the Mannheim gridshell (Happold and Liddel 1975). The numerical form-finding methods, such as force density (Schek 1974), particle-spring method, or the thrust network analysis TNA [START_REF] Block | Thrust Network Analysis. Exploring Three-Dimensional Equilibrium[END_REF], also take pattern as an input. The orientation and connectivity of the starting pattern may change drastically the shape resulting from these methods, as highlighted in [START_REF] Ramm | Shell Structures -A Sensitive Interrelation between Physics and Numerics[END_REF]. Similarly, orientation and connectivity are key parameters when optimizing a geometry to combine mechanical efficiency and fabrication properties such as face planarity, as is for example done in [START_REF] Adriaenssens | Finding the Form of an Irregular Meshed Steel and Glass Shell Based on Construction Constraints[END_REF], [START_REF] Feng | Shape Optimization Method of Free-Form Cable-Braced Grid Shells Based on the Translational Surfaces Technique[END_REF], [START_REF] Mesnil | Fabrication-Aware Shape Parametrisation for the Structural Optimisation of Shell Structures[END_REF] or [START_REF] Tang | Form-Finding with Polyhedral Meshes Made Simple[END_REF].

From shape to topology

Another common design approach is to take a shape as an input, and to find an optimal pattern on it. For structural performance, optimal pattern for a quad gridshell or a waffle shell is commonly searched by aligning the pattern with principal stress nets. Similarly, considering fabrication alone is done by aligning the pattern with specific surface curves or parametrizations, such as conjugate nets [START_REF] Vouga | Design of Self-Supporting Surfaces[END_REF][START_REF] Zadravec | Designing Quad-Dominant Meshes with Planar Faces[END_REF], geodesics, asymptotic curves [START_REF] Schling | Designing Grid Structures Using Asymptotic Curve Networks[END_REF], conformal patterns or principal curvature lines [START_REF] Bo | Circular Arc Structures[END_REF]Tellier et al. 2019b;[START_REF] Wang | A Note on Planar Hexagonal Meshes[END_REF].

Structural performance and fabricability can also be searched together with this approach. For example, [START_REF] Schiftner | Statics-Sensitive Layout of Planar Quadrilateral Meshes[END_REF] propose to pattern a surface following maximum stress lines and their conjugate directions. [START_REF] Oval | Topology Finding of Patterns for Structural Design[END_REF][START_REF] Oval | Topology Finding of Structural Patterns[END_REF] propose a method to explore various topologies and compare their performances.

From boundaries to shape and topology

The two above approaches may be limited when trying to obtain too many properties from the geometry at the same time. For example, [START_REF] Sun | Discrete Curvature Theories and Applications[END_REF]) studied how to design a surface fitting a target boundaries in which principal stress lines coincide with principal curvature lines in a "topology to shape" approach. The limited success of this approach led the team to relax entirely the topology: in (Pellis and Pottmann 2018), the alignment of principal stresses and curvature is obtained by a new constraint on triangular meshes. The topology of the structural pattern is not imposed beforehand; it instead appears in a subsequent remeshing step. This article will use the same approach: given our mechanical and fabrication constraints (funicularity and face planarity), we derive both the optimal topology and the optimal shape from the boundaries. This is permitted by a formulation of these properties which is independent of surface parametrization.

This article extend the work presented at the AAG2020 conference in the continuation of the work presented in [START_REF] Tellier | Discrete CMC Surfaces for Doubly-Curved Building Envelopes[END_REF](Tellier et al. , 2019a[START_REF] Tellier | Linear-Weingarten Membranes with Funicular Boundaries[END_REF][START_REF] Tellier | Funicularity of Conics[END_REF], in which we looked at the generation of constant mean curvature and Linear Weingarten surfaces (not isotropic): on these surfaces, principal stresses and curvature are aligned under a uniform normal pressure load. The present work differs in the sense that we will now consider vertical uniform loads, which model better self-weight and snow loads.

i-liwein surfaces

Isotropic geometry

Isotropic Linear Weingarten surfaces (i-liwein) are objects of isotropic geometry. Isotropic geometry is a way to describe the geometry of shapes with a specific set of tools, in which the vertical direction plays a particular role. A thorough treatment can be found in [START_REF] Sachs | Isotrope Geometrie Des Raumes[END_REF]. [START_REF] Pottmann | Discrete Surfaces in Isotropic Geometry[END_REF] gives a concise introduction.

In isotropic geometry, the notions of distance, angles and curvature differ from the traditional Euclidian geometry. They are all computed relatively to a vertical projection on the horizontal plane. The notion of surface curvature is also described differently. This difference is best understood by looking at a planar curve 𝑧 = 𝑓(𝑥) (Figure 4). There are two ways to describe the inclination of the tangent of such a curve:

• The angle it makes with a reference axis (e.g. a horizontal axis);

• The slope 𝑑𝑧/𝑑𝑥.

With the first one (the angle), the curvature is naturally understood as the variation of the angle of the tangent as we move along the curve (𝑑𝛼/𝑑𝑠, with 𝑑𝑠 = √𝑑𝑥 2 + 𝑑𝑧 2 ), which is the inverse of the radius of the tangent circle: this is the Euclidian curvature. With the second one (the slope), the curvature is best defined as the variation of slope with respect to 𝑥 (𝑑²𝑧/𝑑𝑥²): this is the isotropic curvature. We remark that the slope of a vertical line is infinite: isotropic geometry is not the right tool to describe curves with vertical tangents or surfaces with vertical tangent planes. It is then natural to describe surfaces as height fields:

(𝑥, 𝑦) ↦ (𝑥, 𝑦, 𝑓(𝑥, 𝑦))

(1)

We will mostly use this way of parametrizing surfaces throughout this paper. The role of the curvature tensor is played by the hessian of 𝑓:

∇ 2 𝑓 = [ 𝜕 𝑥𝑥 𝑓 𝜕 𝑥𝑦 𝑓 𝜕 𝑥𝑦 𝑓 𝜕 𝑦𝑦 𝑓 ]
∇ 2 𝑓 is symmetric, so it admits two orthogonal eigenvectors in the 𝑥𝑦 plane, which are called the i-principal directions. The maximum and minimum eigenvalues, 𝑓 1 and 𝑓 2 , are called i-principal curvatures. By integrating the i-principal directions, one obtains the i-principal curvature lines.

They differ from the Euclidian principal curvature lines. However, they still form a conjugate net, which is a smooth equivalent of planar quad meshes [START_REF] Sauer | Differenzengeometrie[END_REF]. Furthermore, their vertical projection on the 𝑥𝑦 plane is an orthogonal net.

i-liwein surfaces

Similarly, to Euclidian geometry, the Gauss and mean curvature can be defined as the determinant and half-trace of the curvature tensor (an i superscript is used to differentiate them from the traditional Euclidian curvatures):

𝐾 𝑖 = 𝜕 𝑥𝑥 𝑓𝜕 𝑦𝑦 𝑓 -𝜕 𝑥𝑦 𝑓 2 = 𝑓 1 𝑓 2 𝐻 𝑖 = 1 2 (𝜕 𝑥𝑥 𝑓 + 𝜕 𝑦𝑦 𝑓) = 1 2 Δ𝑓 = 1 2 (𝑓 1 + 𝑓 2 )
i-liwein can then be defined as the surfaces verifying at each point:

𝑎𝐻 𝑖 + 𝑏𝐾 𝑖 = 𝑐 (2)
for given constant real coefficients a, b and c.

Examples

It turns out that some well-known families of surfaces satisfy equation ( 2), including popular shapes in architecture.

Paraboloids

Paraboloids are a very common family of shapes for shells (see Figure 5). With a proper choice of axes and origin, their equation reads:

𝑧 = 𝛼𝑥 2 + 𝛽𝑦²
If 𝛼𝛽 < 0, the surface is called a hyperbolic paraboloid (or « hypar »), and contains two families of straight lines. If 𝛼𝛽 > 0, the surface is called an elliptic paraboloid. The Hessian of the height field reads:

∇ 2 𝑓 = [ 2𝛼 0 0 2𝛽 ]
Paraboloids are isotropic LW surfaces, as their mean and Gaussian curvature are both constant: There is actually one degree of freedom to define the (𝑎, 𝑏, 𝑐) coefficients of the LW relationship (more precisely, one more degree of freedom than the trivial scaling of a, b and c). If we fix the trivial scaling of (𝑎, 𝑏, 𝑐) by setting either 𝑐 = 0 or 𝑐 = 1, equation ( 1) is verified for any value of 𝑎 as long as:

𝐻 𝑖 = 𝛼 + 𝛽 𝐾 𝑖 = 4𝛼𝛽

𝑏 = c -𝑎(𝛼 + 𝛽) 4𝛼𝛽

(3)

Pelikan surfaces

Membranes for which the membrane stress tensor, when projected in the horizontal plane, is constant, are called Pelikan surfaces. They were introduced in [START_REF] Pelikan | Membrane Structures[END_REF]). These surfaces have been generated by dynamic relaxation [START_REF] Hincz | The Effect of the Apporximations Used during Generation of Membrane Cutting Pattern[END_REF], and have been used for example to design of membranes and concrete shells, such as the Keramion ceramics museum in Frechen (Figure 6), for which a tensile membrane was used as a formwork for a concrete shell (www.keramion.de).

Harmonic functions as height fields

Scalar harmonic functions of the plane are functions for which the Laplacian is null at any point: Δ𝑓 = 0. They are therefore field heights of i-minimal surfaces, i.e. surfaces verifying 𝐻 𝑖 = 0.

There is a rich mathematical theory behind these surfaces. One method to generate a rich variety of analytical harmonic functions is via complex holomorphic functions [START_REF] Cartan | Elementary Theory of Analytic Functions of One or Several Complex Variables[END_REF]. A holomorphic function of a region 𝐷 ⊂ ℂ is a function 𝑓: 𝐷 → ℂ which is conform, i.e. it transforms an infinitesimal square into another infinitesimal square. Most of the complex functions that one could think about that do not contain the complex conjugate 𝑧̅ are holomorphic. This is the case in particular of the following functions:

• Polynomials. Ex : 𝑓(𝑧) = 𝑧 7 -42𝑖𝑧 , 𝑧 ∈ ℂ

• Rational functions (fractions of polynomials). Ex :

𝑓(𝑧) = 𝑧+1 𝑧-1 • Exponential functions. Ex : 𝑓(𝑧) = 𝑒 𝑖𝑧 +𝑒 -𝑖𝑧 2 = cos 𝑧
• Linear combination of harmonic functions. Ex : 𝑓(𝑧) = 𝑒 𝑧 + 1/𝑧

The link with i-minimal surfaces appears with an important theorem about holomorphic functions: Their real part and imaginary part are both harmonic. Figure 7 shows several i-minimal surfaces generated using this property:

• Top left: 𝑓(𝑧) = 𝑅𝑒(log 𝑧) where 𝑅𝑒(𝑧) is the real part of a complex number 𝑧.

• Top middle: Monkey saddle of order 3: 𝑓(𝑧) = 𝑅𝑒(𝑧 3 ) = 𝑥 3 -3𝑥𝑦 2

• Top right: 𝑓(𝑧) = 𝑅𝑒(1/‖𝑧‖ 2 )

• Bottom left: 𝑓(𝑥, 𝑦) = cos 𝑥 exp(𝑦 -1) + cos(𝑥 + 𝜋) exp(-𝑦 -1)

• Bottom right: 𝑓(𝑧) = 𝑅𝑒(tan -1 (𝑧))

𝐾 = 𝜕 𝑥𝑥 𝑓𝜕 𝑦𝑦 𝑓 -𝜕 𝑥𝑦 𝑓 2 (1 + 𝜕 𝑥 𝑓 2 + 𝜕 𝑦 𝑓 2 ) 2 = 𝐾 𝑖 (1 + ‖∇𝑓‖ 2 ) 2
Therefore, any developable surface (𝐾 = 0) that can be described as a height field is i-liwein (𝐾 𝑖 = 0).

Note:

The mean curvature can also be expressed as a function of the isotropic mean curvature:

𝐻 = 𝐻 𝑖 + 1 2 ∇𝑓. ∇ 2 𝑓 ̂. ∇𝑓 (1 + ‖∇𝑓‖ 2 ) 3/2
From that expression, it appears clearly that there is no equivalence between minimal and iminimal.

Surfaces of revolution and helices

(Ogrenmis 2016) solved the ordinary differential equation fulfilled by the meridian of an i-liwein surface of revolution, and gives the analytical solution. There are two degrees of freedom to design the meridian for a given set of coefficients (𝑎, 𝑏, 𝑐) -except for some degenerate coefficients for which there is only one. (Yoon and Lee 2016) performed the same work for i-liwein surfaces with helicoidal symmetry.

3 Mechanics of i-liwein surfaces

Funicularity

This section shows that i-liwein surfaces enjoy a mechanical property very similar to their Euclidian counterparts. This property is best expressed using the tensor 𝐍 𝐢 of projected membrane stresses in the 𝑥𝑦 plane.

Proposition 1

Let us consider a shell with membrane-compatible support conditions and with an i-liwein reference surface fulfilling equation ( 2) (a𝐻 𝑖 + 𝑏𝐾 𝑖 = 𝑐).

If 𝑐 ≠ 0, the shell is funicular for uniform vertical loads 𝑤𝑒 𝑧 (𝑤 is in kN/m 2 ). The projected membrane stresses are given by:

𝐍 𝐢 = 𝑤 2𝑐 (𝑎 𝐈 + 𝑏 ∇ 2 𝑓 ̂)
Where 𝐈 the identity tensor of the 𝑥𝑦 plane, and ∇ 2 𝑓 ̂ is the tensor obtained by rotating the hessian ∇ 2 𝑓 by +90° around the 𝑧 axis:

∇ 2 𝑓 ̂= [ 𝜕 𝑦𝑦 𝑓 -𝜕 𝑥𝑦 𝑓 -𝜕 𝑥𝑦 𝑓 𝜕 𝑥𝑥 𝑓 ]
If 𝑐 = 0, the shell admits a one parameter family of self-stress fields given by:

𝐍 𝐢 = 𝜆(𝑎 𝐈 + 𝑏∇ 2 𝑓 ̂) , 𝜆 ∈ ℝ

Corollary

Under a uniform vertical load, principal membrane stresses on an i-liwein shell form a conjugate net and coincide with i-curvature directions. Discretization by stress lines hence give quadrangles which are near-planar.

Proof:

Membrane equations can be written in Cartesian coordinates, as detailed for example in [START_REF] Ventsel | Thin Plates and Shells[END_REF]. The projected membrane stress field 𝐍 𝐢 is at equilibrium under a vertical load 𝐺𝑒 𝑧 if it verifies:

𝑑𝑖𝑣 𝐍 𝐢 = 𝟎 𝑑𝑖𝑣(𝐍 𝐢 . ∇𝑓) = 𝐺
These equations correspond respectively to the horizontal and vertical equilibrium. The fact that 𝑑𝑖𝑣 𝐍 𝐢 = 0 can be used to simplify the vertical equilibrium equation in the form:

𝑡𝑟(𝐍 𝐢 . ∇ 2 𝑓) = 𝐺
The divergence of the stress fields from the proposition are:

𝑑𝑖𝑣 𝐍 𝐢 = 𝑤 2𝑐 (𝑎 𝑑𝑖𝑣 𝐈 + 𝑏 𝑑𝑖𝑣(∇ 2 𝑓 ̂)
Obviously, 𝑑𝑖𝑣 𝐈 = 0. Furthermore, by expressing ∇ 2 𝑓 ̂ in 𝑥, 𝑦 coordinates:

𝑑𝑖𝑣(∇ 2 𝑓 ̂) = 𝑑𝑖𝑣 [ 𝜕 𝑦𝑦 𝑓 -𝜕 𝑥𝑦 𝑓 -𝜕 𝑥𝑦 𝑓 𝜕 𝑥𝑥 𝑓 ] = 0
The proposed stress fields are therefore at horizontal equilibrium.

The vertical equilibrium is obtained by expressing ∇ 2 𝑓 and ∇ 2 𝑓 ̂ in the i-principal curvature directions:

If 𝑐 ≠ 0: 𝑡𝑟(𝐍 𝐢 . ∇ 2 𝑓) = 𝑤 2𝑐 𝑡𝑟 ((𝑎 [ 1 0 0 1 ] + 𝑏 [ 𝑓 2 0 0 𝑓 1 ]) [ 𝑓 1 0 0 𝑓 2 ]) = 𝑤 2𝑐 (𝑎(𝑓 1 + 𝑓 2 ) + 2𝑏𝑓 1 𝑓 2 ) = 𝑤 2𝑐 (2𝑎𝐻 𝑖 + 2𝑏𝐾 𝑖 ) = 𝑤
The same reasoning applies to the case 𝑐 = 0. The corollary stems from the fact that eigenvectors of ∇ 2 𝑓 are also eigenvectors of ∇ 2 𝑓 ̂, which in turn are eigenvectors of 𝑎 𝐈 + 𝑏∇ 2 𝑓 ̂.

The corollary is proved by remarking that the eigenvectors of 𝐍 𝐢 = 𝜆(𝑎 𝐈 + 𝑏∇ 2 𝑓 ̂) are the eigenvectors of ∇ 2 𝑓 ̂ :

𝑤𝑖𝑡ℎ 𝑘 ∈ ℝ , 𝜆 (𝑎 𝐈 + 𝑏∇ 2 𝑓 ̂)𝒖 = 𝑘𝒖 ⟺ ∇ 2 𝑓 ̂𝒖 = ( 𝑘/𝜆 -𝑎 𝑏 ) 𝒖
Since ∇ 2 𝑓 is symmetric, its eigenvectors are orthogonal. As ∇ 2 𝑓 ̂ is obtained by rotating ∇ 2 𝑓 inplane by 90°, eigenvectors of ∇ 2 𝑓 ̂ and ∇ 2 𝑓 coincide. We can conclude that the eigendirections of 𝐍 𝐢 coincide with i-principal curvature directions of the surface 𝑓. ∎

Link with Airy stress function

Planar stress tensors at equilibrium without external loads can be expressed with the hessian of a scalar function Φ, called the Airy stress function. When a shell is subject to purely vertical forces, such a function can be used to describe the projected stress tensor 𝑵 𝒊 :

𝑵 𝒊 = ∇ 2 Φ ̂
In our case, the Airy stress function Φ is :

Φ(x, y) = 𝑎 2 (𝑥 2 + 𝑦 2 ) + 𝑏𝑓(𝑥, 𝑦)
The Airy stress function is hence a linear combination of the heights of the isotropic sphere and of the surface height itself, which is a remarkable property.

Application to paraboloids

The above result can be applied to paraboloids with an interesting consequence: they are funicular for two different load cases. Given a paraboloid 𝑧 = 𝛼𝑥 2 + 𝛽𝑦² we have two possible choices for parameter 𝑐, 𝑐 = 0 or 𝑐 = 1:

With 𝑐 = 1, we can interpret the shape as funicular under uniform vertical. Using equation (3), there is a one-parameter family of admissible stress fields :

𝐍 𝐢 = 𝑤 2 (𝑎 𝐼 + 𝑏 ∇ 2 𝑓 ̂) = 𝑤 2 (𝑎 [ 1 0 0 1 ] + ( 1 -𝑎(𝛼 + 𝛽) 4𝛼𝛽 ) [ 2𝛽 0 0 2𝛼 ]) , 𝑎 ∈ ℝ
With 𝑐 = 0, we find another family of admissible stress fields, which corresponds to self-stressed states:

𝑵 𝒊 = 𝑎 𝑰 + 𝑏𝛻 2 𝑓 ̂= 𝑎 ( [ 1 0 0 1 ] - (𝛼 + 𝛽) 2𝛼𝛽 [ 𝛽 0 0 𝛼 ]) , 𝑎 ∈ ℝ
The fact that paraboloids are at the same time self-stress shapes and funicular shaped under uniform load is remarkable. For example, if a hypar membrane or cable-net (case 𝛼𝛽 < 0) is used as a formwork for casting a concrete shell, the membrane deflection under the dead load of the fresh concrete is expected to be very low compared to other arbitrary membrane shapes:

deflections will be mostly due to material stretching rather than to a deflection to get into a funicular shape. We remark that if 𝛼 + 𝛽 = 𝐻 𝑖 ≠ 0, warp and weft of the membranes need to be aligned with the horizontal 𝑥 and 𝑦 axes -otherwise the stress field is not admissible.

Patching paraboloids

Rows of paraboloids can be assembled to form a C 1 surface that is still i-liwein. Indeed, let us consider the paraboloids of Figure 8 which form a row on the 𝑥 axis. Blue and white patches can be described respectively by equations 𝑧 = 𝛼 𝑏𝑙𝑢𝑒 𝑥 2 + 𝛽𝑦² and 𝑧 = 𝛼 𝑤ℎ𝑖𝑡𝑒 𝑥 2 + 𝛽𝑦² , 𝛽 being identical to insure continuity. If we set:

𝑎 = 1 ; 𝑏 = -1 4𝛽 ; 𝑐 = 𝛽
Then, for both blue patches:

𝑎𝐻 𝑖 + 𝑏𝐾 𝑖 = (𝛼 𝑏𝑙𝑢𝑒 + 𝛽) + -1 4𝛽 4𝛼 𝑏𝑙𝑢𝑒 𝛽 = 𝛽 = 𝑐
The same results hold for the white patches. The surface is therefore i-liwein. The fact that it is not C 2 does not alter the mechanical properties: we can find a continuous membrane stress field across the whole surface under uniform vertical load. It does not either alter the panelling by planar quads, because the i-principal curvatures are aligned with the 𝑥 and 𝑦 axes for all patches. where 𝑔 is an arbitrary real continuous function and 𝛽 a real number.

Invariance by x, y or z scaling

One can easily prove that the vertical scaling conserve the i-liwein property. However, a horizontal scaling (for example a 1D-scaling of 𝑥 direction) transforms equation ( 2) into the form:

𝑎 𝑥 𝜕 𝑥𝑥 𝑓 + 𝑎 𝑦 𝜕 𝑦𝑦 𝑓 + 𝑏 0 𝐾 𝑖 = 𝑐 0
The same reasoning as in 3.1 shows that these surfaces are also funicular under a uniform vertical load. In particular, scaling an i-CMC surface (𝐻 𝑖 = 𝑐𝑡𝑒) yields surfaces with constant projected stress tensors, i.e. a general Pelikan surfaces.

Extension to surfaces funicular to self-weight

Proposition 1 considers a uniform vertical load. This is often a good approximation of dead load if the slope of an envelope is small, and it is always a good model of symmetric snow loads.

However, in shell structure of constant thickness with high slopes, this approximation might not be precise. The result of proposition 1 can actually be adapted for this setting to search for surfaces with constant projected stresses:

Proposition 2

A shell of constant thickness with reference surface 𝑧 = 𝑓(𝑥, 𝑦) such that:

𝑡𝑟(∇ 2 𝑓) √1 + ‖∇f‖² = 𝑐 , 𝑐 ∈ ℝ
is self-supported: It can resist its own weight with pure membrane stresses, such that the projected membrane stress tensor in the horizontal plane is isotropic and constant.

Proof:

Self-weight is a vertical load of amplitude:

𝐺(𝑥, 𝑦) = 𝜇 √1 + ‖∇f‖² , 𝜇 ∈ ℝ
Where G is expressed for example in kN per m² of projected area in the xy plane. Let us consider the projected stress tensor 𝐍 𝐢 = 𝜇 𝑐 𝐈. It verifies 𝑑𝑖𝑣 𝐍 𝐢 = 𝟎. Furthermore:

𝑡𝑟(𝐍 𝐢 . ∇ 2 𝑓) = 𝑡𝑟 ( 𝜇 𝑐 ∇ 2 𝑓) = 𝜇 √1 + ‖∇f‖² = 𝐺
𝐍 𝐢 is therefore at equilibirum under the load 𝐺. ∎

This equation opens the way to the search of analytical self-supporting grids. One could for example search for a family of funicular surfaces of revolution: Their meridians could be found by solving an ordinary differential equation, with one variable. One could also look for funicular surfaces of translations. These surfaces will not be addressed further in this article.

Recovering membrane stresses

The membrane stress tensor 𝐍 can be recovered from the projected tensor 𝐍 𝐢 following equations given for example in [START_REF] Frey | Analyse Des Structures et Milieux Continus -Coques[END_REF]. Let us call:

• 𝒆 𝒙 and 𝒆 𝒚 the base vectors of the plane 𝑥𝑦;

• 𝒆 𝒙 𝑺 and 𝒆 𝒚 𝑺 the projections of the base vectors 𝒆 𝒙 and 𝒆 𝒚 in the tangent plane of the surface. These vectors are in general not orthogonal; Note: Despite the fact that 𝑁 𝑥𝑦 = 𝑁 𝑥𝑦 𝑖 , principal directions of 𝐍 𝐢 (for which shear is null) are not the projections of the principal directions of 𝐍. This is due to the fact that the 𝑥 and 𝑦 axes do not project in general to orthogonal axes on the membrane tangent plane.

Generation of i-LW surfaces from boundary curves

In this section, a generation method for elliptical i-liwein surfaces (verifying 𝑎 2 + 4𝑏𝑐 > 0) from a boundary curves is presented. The Monge Ampere equation has many physical applications, ranging from optimal transport problems [START_REF] Benamou | Numerical Solution of the Optimal Transportation Problem Using the Monge -Ampère Equation[END_REF][START_REF] Villani | Topics in Optimal Transportation[END_REF] to reflector design [START_REF] Wu | Freeform Illumination Design : A Nonlinear Boundary Problem for the Elliptic Monge -Ampére Equation[END_REF]. It also has many applications in differential geometry, allowing for example the proof of the existence of particular types of surface [START_REF] Trudinger | The Monge-Ampère Equation and Its Geometric Applications[END_REF]. Numerical resolution has been addressed quite recently by research. A review of is given in [START_REF] Neilan | The Monge-Ampère Equation[END_REF].

This problem is particularly complex to solve for the following reasons:

• It is a fully nonlinear problem: the term with highest derivative is not linear;

• Solutions have low regularity and are not unique;

• Weak solutions are not based on variational principles. A significant portion of the literature is devoted to viscosity solutions, which may have weaker regularity;

• Existence and uniqueness of solutions is subject to a constraint on the convexity of the solution and domain -which is a strong constraint for application to architecture.

Many resolution methods are based on finite differences and a quad grid discretization of the plane, [START_REF] Benamou | Two Numerical Methods for the Elliptic Monge-Ampère Equation[END_REF]Dean andGlowinski 2003, 2004;[START_REF] Oberman | Wide Stencil Finite Difference Schemes for the Elliptic Monge-Ampère Equation and Functions of the Eigenvalues of the Hessian[END_REF]).

These are less suited for irregular contours, as are likely to occur in an architectural context. Also, these methods only search for a convex solution to the problem, which we will refer to as type II, in section 4.3.1.

(Feng and Jensen 2017) introduced a finite-element formulation based on a triangular discretization. It was used in [START_REF] Jensen | Numerical Solution Of The Simple Monge-Ampère Equation With Nonconvex Dirichlet Data On Nonconvex Domains[END_REF] to explore some solutions of the Monge-Ampere equation on non-convex domains. As this method is relatively complex, we will use the fact that, in our context, 𝑔 = 𝑐𝑡𝑒, to propose a simpler generation method based on a quadratic discrete operator on triangular meshes from discrete differential geometry.

Discrete i-liwein surfaces

Discrete curvature models

We discretize i-liwein surfaces by using the vertex-wise discrete model for 𝐻 𝑖 and 𝐾 𝑖 for unstructured triangular meshes proposed in [START_REF] Pottmann | Discrete Surfaces in Isotropic Geometry[END_REF]. The calculation steps to compute the curvature of a vertex 𝑣 in a mesh M are the following:

i.

Compute the metric dual of the planes of the faces adjacent to vertex 𝑣. The dual of the plane of a face 𝑖 is a point, which, once projected vertically onto the 𝑥𝑦 plane, yields a point 𝑚 𝑖 . The points 𝑚 𝑖 form a closed polygon 𝐿 𝑚 . The area of this polygon, 𝐴(𝐿 𝑚 ) can be calculated as:

𝐴(𝐿 𝑚 ) = 1 2 ∑ det(𝑚 𝑖 , 𝑚 𝑖+1 )
ii.

Project the 1-star of 𝑣 vertically on the Maxwell paraboloid Σ. Compute the metric dual of the projected faces: they form a closed polygon. Project it onto the 𝑥𝑦 plane, thus forming the polygon 𝑠 1 𝑠 2 … called 𝐿 𝑠 , with oriented area A(𝐿 𝑠 ).

iii.

Compute the mixed area of polygons 𝐿 𝑀 and 𝐿 𝑆 , defined as:

𝐴(𝐿 𝑀 , 𝐿 𝑆 ) = 1 4 ∑(det(𝑚 𝑖 , 𝑠 𝑖+1 ) + det (𝑠 𝑖 , 𝑚 𝑖+1 ))
iv.

The discrete Gaussian and mean curvature are then computed as:

𝐾 𝑖 (𝑣) = 𝐴(𝐿 𝑚 ) 𝐴(𝐿 𝑠 ) , 𝐻 𝑖 (𝑣) = 𝐴(𝐿 𝑚 , 𝐿 𝑠 ) 𝐴(𝐿 𝑠 )
Notes:

• With this model, vertices should not have too high 𝑥 or 𝑦 coordinates: The projection to the isphere gets points to an altitude 𝑥² + 𝑦², working tolerance can be reached easily.

• 𝐻 𝑖 (𝑣) could be calculated in a more simple way by taking a discrete Laplacian of the 𝑧 coordinates of the vertices. The chosen discretization is more complex, but more consistent with the discrete model of 𝐾 𝑖 .

Assessment of discrete model

Figure 9 shows the values given by these two discrete models for two triangular meshes inscribed in paraboloids -for which the Gaussian and mean curvature are known analytically. By inscribed, we mean that the mesh vertices are exactly on the surfaces. We observe that the discrete curvatures are very close to the ones of the smooth paraboloids for the synclastic paraboloid.

However, high discrepancies are observed for 𝐾 𝑖 for the hypar, especially for nodes of valence different from 6, and in areas of high slope. This problem of approximation of Gaussian curvature for inscribed meshes was identified in [START_REF] Borrelli | On the Angular Defect of Triangulations and the Pointwise Approximation of Curvatures[END_REF] in Euclidian geometry.

We observe the same issue in isotropic geometry. However, we shall see in section 4.5 that these discrepancies are not an issue in our method because they are due to the fact that we use inscribed meshes. From these discrete models, i-liwein triangular meshes can naturally be defined as meshes for which a linear combination of these two curvatures is constant for each vertex:

∀ 𝑣𝑒𝑟𝑡𝑒𝑥 𝑣, 𝑎𝐻 𝑖 (𝑣) + 𝑏𝐾 𝑖 (𝑣) = 𝑐, 𝑎, 𝑏, 𝑐 ∈ ℝ

Preliminaries on smooth i-liwein surfaces

Type I and II solutions

Important geometrical insight on i-liwein surfaces can be obtained from the curvature diagram, i.e. a diagram plotting the value of 𝐻 𝑖 and 𝐾 𝑖 at each surface point (Figure 10). For i-liwein surfaces, this diagram is included in a straight line. We observe that the isotropic Gaussian and mean curvature verify the inequality:

𝐻 𝑖 2 -𝐾 𝑖 = ( 𝑓 1 -𝑓 2 2 ) 2 ≥ 0
There is hence a domain of the diagram which is impossible geometrically (in red). If 𝑎 2 + 4𝑏𝑐 ≥ 0, the line defined by equation (2) crosses this forbidden domain. The curvature diagram is hence included in the union of two rays (semi-infinite lines). Now, in a smooth surface, curvatures vary continuously: the curvature diagram is a one-connected domain. As a result, the curvature diagram of an i-liwein surface with 𝑎 2 + 4𝑏𝑐 ≥ 0 cannot be on both rays. We will call type I (resp.

type II) the i-liwein surfaces having their curvature diagram entirely in the left (resp. right) ray.

From a design perspective, type I surfaces are more interesting formally than type II, because the latter tend to look quite spherical. (𝑥 2 + 𝑦 2 ):

4𝑏 2 𝐾 𝑖 (𝑔) = 𝑎 2 + 4𝑏𝑐
Solving equation ( 2) with 𝑎 2 + 4𝑏𝑐 > 0 is therefore equivalent to computing the geometry of a surface with positive constant i-Gaussian curvature. This problem is elliptic, and is hence well suited to a generation from a closed boundary.

Generation method

We introduce an iterative method to solve the equation, with the following simple pseudo-code:

i. Triangulate domain ii. Create initial solution (e.g. z=0 or minimal)

iii. Iterate the following steps until displacements are below tolerance:

a. Solve the equation 𝑎𝐻 𝑖 + 𝑏𝐾 𝑖 = 𝑐 for each vertex individually, assuming that neighbours are fixed.

b. Move each vertex to the calculated position.

Local resolution

The key step of the algorithm is the resolution with fixed neighbours (step iii.a). The equation is nonlinear, but we remark that, as the centre vertex of a 1-star is moved along a vertical axis:

• 𝐴(𝐿 𝑠 ) does not change;

• 𝐴(𝐿 𝑚 ) varies quadratically;

• 𝐴(𝐿 𝑚 , 𝐿 𝑠 ) varies linearly.

Hence, 𝐻 𝑖 and 𝐾 𝑖 vary respectively linearly and quadratically with a vertex height. Therefore, 𝑎𝐻 𝑖 + 𝐵𝐾 𝑖 -𝑐 is quadratic in 𝑧. Thus, there exists real coefficients 𝛼, 𝛽 and 𝛾 such that:

𝑎𝐻 𝑖 + 𝑏𝐾 𝑖 -𝑐 = 𝑃(𝑧) = 𝛼𝑧 2 + 𝛽𝑧 + 𝛾
The resolution in 𝑧 is then straightforward, we use the following method:

Case 𝑏 = 0

If 𝑏 = 0, the equation is linear. This case corresponds to i-CMC surfaces. 𝛼, 𝛽 and 𝛾 are given by: 𝛼 = 0 ; 𝛾 = 𝑃(0) ; 𝛽 = 𝑃(ℎ) -𝛾 ℎ Where ℎ is an arbitrary length, taken for example as the average edge length of the 1-star boundary. Then:

𝑧 = -𝛾/𝛽 Case 𝑏 ≠ 0
In that case, the equation is quadratic. We calculate coefficients 𝛼, 𝛽 and 𝛾 as follows:

𝛾 = 𝑃(0) ; 𝛽 = 𝑃(ℎ) -𝑃(-ℎ) 2ℎ
; 𝛼 = 𝑃(ℎ) -𝛽ℎ -𝛾 ℎ² We then compute the discriminant ∆= 𝛽 2 -4𝛼𝛾. We choose the height z as follows:

• If ∆> 0 : For a type I solution, we choose the root solution 𝑧 = -𝛽-√∆

2𝑎

; for a type II solution, we choose 𝑧 = -𝛽+√∆ 2𝑎

• If ∆= 0, we obviously return the only solution;

• If ∆< 0, we output the critical point, which minimizes |𝑎𝐻 + 𝑏𝐾 -𝑐|: 𝑧 = -𝛽 2𝛼

Validation

The generation method is validated using again paraboloids. By generating a surface with the analytical boundaries of a paraboloid, one can obtain an i-liwein mesh very close to the paraboloid. Examples are given in Figure 11. Curvatures, Gaussian and mean, both converge well to their expected value. Results are given for 1000 iterations To finish our discussion about vertex-wise approximation of curvatures, related to [START_REF] Borrelli | On the Angular Defect of Triangulations and the Pointwise Approximation of Curvatures[END_REF], we remark that that the discrete scheme of 𝐾 𝑖 and 𝐻 𝑖 , which were giving poor results for a mesh inscribed in an i-liwein surface, can actually converge very well to the smooth curvatures of the underlying surface (whatever the vertex valence) if the mesh approximates the surface without being inscribed in it, i.e. if its vertices are not constrained to be on the surface.

Cases of non-existence of solutions

We do not have theoretical results on the convergence of our method. For convex boundaries (i.e. boundaries for which the projection in the 𝑥𝑦 plane defines a convex domain of the plane), existence of solutions to our form of the Monge-Ampere is guaranteed, and our method proves to be fairly robust. However, there is at the moment no existence result regarding non-convex boundaries. Practically, instability may be encountered on non-convex boundaries: if parameter 𝑏 is pushed to a certain value, and if |𝑐| is low, then the vertex height may diverge. This may happen even at low values of 𝑏 if the domain has holes with a highly curved concave boundary.

The search of existence of solutions for non-convex boundaries for given parameters (𝑎, 𝑏, 𝑐) is beyond the scope of this paper.

Designing with i-liwein surfaces

The generation algorithm presented in section 4 and the analytical primitives shown in section 2 may be used to design various typologies of structures which simultaneously have planar panels and are funicular or self-stressed.

Membranes and cable nets

With 𝑐 = 0, surfaces 𝑎𝐻 𝑖 + 𝑏𝐾 𝑖 = 0 correspond to the equilibrium shapes of self-stressed membranes or cable nets. The shape of i-liwein membranes can be understood by comparing iminimal with minimal surfaces. Minimal surfaces correspond to shapes of soap films, they have been a great source of inspiration for the design of textile membranes in architecture [START_REF] Bach | IL18 Forming Bubbles[END_REF]. Compared with minimal surfaces, i-minimal surfaces:

• Cannot be self-stressed with constant isotropic membrane stresses. They still correspond to self-stressed membranes (proposition 1), but the constant isotropic quantity is their projected stresses. This implies that membrane stresses are higher in areas of high slope;

• Are restricted to shapes that are field heights;

• Can join boundaries with significant height difference, which cannot be interpolated by a single minimal surface.

The differences are illustrated in Figure 12. i-minimal and minimal surfaces are created on the same three boundaries, with increasing height differences between the large oval on the ground and the two loops at the centre. In the top row, the height difference is quite low, such that the difference between the two surfaces is hardly noticeable. In the middle row, the height difference is increased: the difference between the two surfaces is clear. Necking is observed in the minimal surface around the center loops, but not for the i-minimal surface, which is actually just a scaled copy (in the z-direction) of the one above. In the bottom row, the height difference is further increased. There does not exist anymore a stable minimal surface fitting the boundary curves. i-liwein membranes do not require that the isotropic direction is vertical: whatever the direction, they correspond to a self-stressed membrane geometry. Therefore, on a given boundary, the isotropic direction can be tilted to explore shape variations. An example is shown in Figure 13, where three i-minimal surfaces are generated on the same boundaries (two coaxial circles), but with different isotropic directions (indicated by red arrows). In each case, the stresses projected in the plane perpendicular to the red arrows are constant. Figures 12 to 14 illustrate how i-liwein surfaces may be used to design membranes. They were generated from the workflow of section 4. Figure 14 shows a cable net designed from a monkey saddle (such as the one shown in Figure 7). The pattern and the surface were derived from the analytical equation (given in 2.3.3).

Gridshells

For a gridshell, i-principal curvature directions of i-liwein surfaces combine both fabrication and mechanical properties. An application is illustrated in Figure 1 and Figure 15. An i-liwein surface is generated on a boundary composed of closed boundary curves -the boundary of Figure 15 mimics the one of the Chadstone gridshell. An i-principal curvature net is computed, and is used to define the position and orientation of the beams of a gridshell. Since the i-principal net forms a conjugate curve network, mesh faces are almost planar (an optimization could make them exactly planar with small vertex displacements). Furthermore, beams are aligned with principal stress directions under uniform vertical pressure (snow load), so they have a relevant orientation.

Finally, we can notice that pattern drawn by i-curvature lines is smooth, and gives faces with approximately uniform size.

If one were to search by optimization of vertex positions the geometry of a gridshell combining planar quads and beams aligned with principal stresses, a main question would be: what is the right mesh combinatorics to start with? In i-liwein surfaces, this combinatorics appears naturally from the curvature field.

Funicular shells

i-liwein surfaces can also be used to design funicular shells. The method allows a quick exploration of a 2-parameter family of funicular surfaces fitting a target boundary. For a concrete shell, the fact that i-principal curvature lines form a conjugate net could be used to simplify the fabrication of formwork. For example, the surface can be covered with developable strips so that thin flexible panels can be used as formwork. The effect of parameters b and c on the surface geometry are the following: Parameter c allows to more or less inflate the surface (with a pressure acting vertically), while parameter b allows the control of the height difference between saddles and hill tops. This is illustrated in Figure 16, in which the reference surface of Figure 15 (which is i-CMC) is deformed. The left picture shows the case b<0: the surface looks like domes connected by low corridors. The right picture shows the case b>0: the surface now looks level, a more unified space is defined underneath.

6 Mechanical behavior: two case studies

In this section, we perform two case studies to investigate first the effect of parameters a, b and c on the mechanical behavior of a concrete shell, then the efficiency of i-curvature lines for defining the layout of a gridshell structure.

Application to concrete shell rooves

We compare three i-liwein surfaces (Figure 17) with the same boundary (that resembles an eight or a peanut), and same area (430m²), but different values of b. The constant area is ensured by adjusting the value of c while b varies (see table 1). With b<0, the height difference between the two domes and the valley is higher, while with b>0, the crest line is more uniform (see figure 15 and table 1). All three surfaces present two portions with positive Gaussian curvature on the side (the domes) and one portion with negative Gaussian curvature in the middle (the valley).

Practically, the reference mesh is the same triangulated mesh with 1852 faces for the three surfaces. Supposing now that these three surfaces are compared in terms of mechanical performances in order to build a concrete shell roof. The rise over span ratio of the surfaces remaining between 1/7 and 1/5, the shells cannot be considered as low shells and their behaviour should be similar and dominated by membrane actions with limited bending stresses. The thickness of the shell is uniform and equal to 10cm, so approximately 1/100 of the average span (10.5m). As the area is constant over the three surfaces, the quantity of material used for the three shells is the same.

Differences in the structural response will thus be uniquely related to the shape of the shell.

Two support conditions are investigated: i) membrane supports with zero displacements in the tangent plane and ii) fixed supports with no translations and no rotations along the edge. The first one corresponds to the theoretical membrane model where the shell develop forces only within the surface. The second one corresponds to more practical concrete shell roof boundary conditions and will develop bending moments around the edge.

Two loading conditions are then applied to the shell: i) uniform vertical load on the whole structure, ii) uniform vertical load on one fourth of the structure. The shells are all funicular for the first load and should resist without bending, while they should exhibit some under the second load.

Two kinds of analysis are conducted: i) standard static linear analysis and ii) linear buckling analysis. The first one is used to retrieve principal stress pattern, as well as displacement measures; the second one critical buckling loads. Key information for comparison is not absolute but relative, the intensity of load and Young modulus being somewhat arbitrary (here 25kN/m² and 31GPa respectively). All calculations are conducted with finite element software Karamba3D.

Alignment of i-principal stresses with i-principal curvatures

For the three surfaces, one first compares the layout of i-principal stresses under uniform vertical load with i-principal curvature lines. The formers are calculated by simple integration of the membrane stress field in the shell by Karamba3d, the later from the principal curvature lines obtained by EvoluteTools Tmap. Both networks are calculated on a scaled down geometry by a factor of 10 in the z-direction so that principal network and i-principal network can be assimilated. As one can see on figure 16 for the iCMC surface on membrane supports, i-principal stresses and i-principal curvature coincide very well: the surface is hence funicular in the sense of isotropic geometry.

Buckling capacity

Concerning mechanical performances, several criteria were tested: maximum displacement, strain energy and buckling capacity which all lead to the same conclusions regarding the influence of the shell shape. For conciseness, it was thus chosen to focus on the buckling capacity which is often the sizing criterion for shell structures. Critical factors for the three surfaces under both loadings with both boundary conditions have been reported in Table 2. One remarks that, for all cases, the capacity of the three shells compares well. Indeed, relative differences (maxmin)/average range from 10% to 27% with membrane supports and fall down to a range of 4% to 6% with fixed edges. One can thus conclude that the shape has less influence than the boundary conditions. Once appropriated boundary conditions are given to the shell (here fixed supports), then the shape of the i-liwein shell has few influence on the capacity. This is true for the funicular load as well as for the non-symmetric loading on one quarter of the structure. This small example demonstrates that the framework of i-liwein surfaces provides to the designer a whole family of shapes with constant area (and by there, constant volume of material) and similar mechanical performances. This family of shapes can be explored easily with help of two parameters, b and c (which are not independent if the area has to remain constant) in order to satisfy other criteria than mechanical performance, such as aesthetic or functionality.

Application to steel gridshell structures

The interest of i-liwein surfaces for structural engineers is further investigated through a gridshell variant of the iCMC surface studied in previous section. To obtain the layout, the following procedure was applied:

Form the desired surface by adjusting parameters b and c for a given boundary;

Scale down the surface and find the network of principal curvatures;

Project this network vertically on the surface.

Build the gridshell on this network orienting cross-section according the bisecting plane of adjacent faces.

The generated grids have by construction quasi-flat quadrangular panels because they follow a conjugate network on the reference smooth surfaces. In the present case (the iCMC surface studied in previous section), the maximum out-of plane default (calculated as the distance between quad diagonals relatively to the length of the largest diagonal) is hence below 1%. It can be easily reduced by a short optimisation for planarity of the node positions. Indeed, authorising a displacement of maximum 2cm around the initial position on the mesh allow to fit standard construction tolerances for glass panels (i.e. out of plane default below 0.5%). Theoretically, non-trivial offsets with torsion free nodes are guaranteed on gridshells built on curvature lines, not on i-curvature lines. However, the geodesic torsion of i-curvature lines is still small and they can build some reasonable inputs for optimisation. In the present case after planarization of faces, optimisation is computed but torsion at nodes is evaluated considering that the local beam planes are set according to the bisecting plane of adjacent panels (i.e. normal to the surface). One can hence build four axes at nodes by intersecting the beam planes and compute the angles between node axes which gives an estimation of torsion at nodes. Here the maximum is of 10 mrad, which corresponds to 1 mm for a 10 cm beam and is almost acceptable as it is.

Once the grid is found, it is tested for funicularity and submitted to a uniform projected load (UPL): each node is loaded according to its tributary area in the top view with 1kN/m². It is also loaded with the same load applied on one quarter of the surface. All connexions between members are fixed, as well as all supports on the ground. A rectangular hollow section of 100x50x8 is attributed to all members orienting the cross-section according to the bisecting plane of adjacent faces. Extremal values of internal forces are shown in table 3. Under uniform projected load, the structure is hence under pure compression: no tension arises and the average bending stresses represent less than 1% of the total stresses. One can thus conclude that the gridshell is funicular of UPL, verifying hence the discrete version of the theorem demonstrated earlier for the smooth geometry. Under non symmetric loading (UPL on one quarter of the grid), bending forces are multiplied by 40 in the weak direction and torsion compared to UPL. These moments correspond to shear forces in the equivalent continuous membrane which were almost equal to zero under UPL because the beams are aligned with iprincipal stress directions. In the strong axis, the increase is large but smaller because under UPL some bending moments were already induced by edge effects near the supports. However, stresses due to normal forces still represent 87% of the total stress in average (more than 80% for 80% of beams) which demonstrates that this kind of funicular gridshell is also very efficient under non-symmetrical loads.

7 Discussion

Application to design and comparison with other methods

Our proposed framework has the following limitations. Firstly, although i-liwein are efficient for uniform vertical load, a final design would require checking asymmetrical loads and buckling. They might not be optimal then. Secondly, the proposed method does not give access to the entire design space of surfaces of funicular planar quads.

However, a strength of the proposed method is the speed at which a designer can explore a rich design space in a project. Compared to the popular "topology to shape" methods discussed in section 1.1, designers need not input trial topologies, a step which may by time consuming for complex boundaries. Even though this preliminary shape might not be the final one, the designer knows a shape that works, and can still deform it afterwards to account for other constraints.

Considering both strengths and limitations, we believe our method is particularly suited for conceptual design phases, especially in situations in which one wants to explore curved options such as gridshells or membranes with very limited available time.

Implementation details

The algorithm described and used in sections 4 and 5 was implemented by scripting the software Rhino/Grasshopper. All the numerical images of this section were generated from that script. The linear algebra package Math.net numerics is used to diagonalize curvature and stress matrices.

The i-Principal net in Figure 1 or 14 are computed in three steps:

1. The surface is scaled down in z direction such that its height is significantly lower than its width. In this state, principal and i-principal directions nearly coincide.

2. A principal net is computed, using commercial software Evolute Tmap.

3. The net and the surface are scaled back in the z direction to fit the original geometry.

We use the fact that i-principal directions are invariant under z scaling (the eigenvectors of ∇ 2 𝑓 are the eigenvectors of λ∇ 2 𝑓 (with λ ∈ ℝ)).

In the current setup, a designer can interactively deform a surface by controlling parameters. It takes about 10 seconds to deform a mesh with 500 vertices to a near LW surface on a 3.5GHz

CPU with 16Go of RAM, and get a good idea of the final shape. Full convergence (displacements below 1/1000 of surface dimensions) takes about a minute, but could be made much faster by not displaying intermediate calculation steps and using a second order solver, like L-BFGS. In particular, we remark that if 𝑏 = 0, our algorithm does not benefit from the linearity of the problem. Faster resolution methods are a topic of future work.

Elastic deformations

i-liwein surfaces an admissible membrane stress field under uniform vertical load. They hence correspond to equilibrium shapes of deformed configurations. However, there is no guaranty that an i-liwein will stay i-liwein under elastic deformation, or even that i-liwein surfaces can be obtained exactly by elastic deformation of an initially stress-free material configuration.

The study of this problem could be addressed in future work. However, these are not issues in the considered applications:

• In the case of a funicular shell, displacements are very small. If the initial state is i-liwein, the deformed state will be nearly i-liwein. The analytical stress fields can then be used for limit analysis, for example for checking that a funicular surface can be found within the middle third of the cross section of a masonry shell [START_REF] Heyman | The Stone Skeleton: Structural Engineering of Masonry Architecture[END_REF]). In the case of a concrete shell, the surface also does not need to be exactly i-liwein, as it can accommodate some bending.

• In the case of a tensile membrane, our approach is similar to popular methods such as force densities or dynamic relaxation method: Both methods also yield a deformed configuration at equilibrium. Potential incompatibilities are usually not a problem, because patterning and pretensioning is an approximate task anyway: tension is applied in warp and weft directions to membrane patches, so the actual deformed configuration is always a bit different from the target surface. In worst case, this leads to small waves on the surface.

Buckling

As stresses are known analytically, there might be ways to assess buckling loads under uniform load analytically. This could be the topic of future work.

Conclusion

This article showed the relevance of isotropic Linear Weingarten ("i-liwein") for architectural design, with applications to many structural typologies: gridshells, membranes, shells and cable nets. It was proven that these surfaces remarkably combine mechanical and fabrication 32 properties: they may be covered by a pattern of planar quadrangles which is also funicular under a uniform load or self-stressed. A generation method based on discrete differential geometry was proposed. i-liwein surfaces appear to offer an interesting design space, with a control from boundary curves which is well suited for many architectural projects. Parameters 𝑎, 𝑏 and 𝑐 offer a global shape control, which is particularly interesting for editing complex shapes. The low amount of pre-processing required for the designer makes this method particularly suited in the context of conceptual design explorations.

On a final note, it is remarkable to observe how constant mean-curvature surfaces and Linear

Weingarten surfaces are funicular for a normal pressure or a vertical pressure depending on the type of geometry used to define the curvature: Euclidian or isotropic.

Figure 1 :

 1 Figure 1: Design workflow for i-liwein surfaces. A designer can control boundary curves and 2 parameters

Figure 2 :

 2 Figure 2: isotropic Linear Weingarten (i-liwein) tensile membranes generated on the same boundary curves (in black)

Figure 3 :

 3 Figure 3: A self-stressed cable-net with planar quad glass panels supported at cable crossings 1.1 Previous work

Figure 4 :

 4 Figure 4: Description of curvature in Euclidian and isotropic geometry

Figure 5 :

 5 Figure 5: These three structure are composed of paraboloids, which are a special types of isotropic LW surfaces. Top left: El Oceanográfico (F. Candela), Valencia, Spain (©Felipe Gabaldón); Top right: Teepott restaurant, Warnemünde, Germany (©An-D) ; Bottom: Courtyard roof Industriepalast, Leipzig (reproduced with permission from (Glymph et al. 2004))

Figure 6 :

 6 Figure 6: The concrete shell of the Keramion museum in Frechen, Germany, is a Pelikan surface -a special type of i-LW surface also referred to as i-minimal. ©Klaas Vermaas

Figure 7 :

 7 Figure 7: Analytical i-minimal surfaces -each surface can be realized as a tension membrane 2.3.4 Developable surfaces Developable surfaces are characterized by 𝐾 = 0. For a height field surface defined as in equation (1), 𝐾 and 𝐾 𝑖 are related by (do Carmo 1976):

Figure 8 :

 8 Figure 8: C 1 i-liwein surface composed of paraboloids This reasoning can be actually extended to an arbitrary number of patches 𝑧 = 𝛼 𝑖 𝑥 2 + 𝛽𝑦 2 , 𝑖 = 1,2 … , 𝑛, or even to a translation surface generated by translating a parabola of the 𝑦𝑧 plane along a curve of the 𝑥𝑧 plane, i.e. a surface with equation: 𝑧 = 𝑔(𝑥) + 𝛽𝑦 2

4. 1 ;

 1 The Monge-Ampere equation As i-liwein surfaces are necessarily height fields, they can be generated by solving the partial differential equation (2) in a function 𝑓 of 𝑧. This is a special case of the Monge-Ampère equation as it can be written in the form (Courant and Hilbert 1962): 𝐸(𝜕 𝑥𝑥 𝑓𝜕 𝑦𝑦 𝑓 -𝜕 𝑥𝑦 𝑓 2 ) + 𝐴 𝜕 𝑥𝑥 𝑓 + 2𝐵 𝜕 𝑥𝑦 𝑓 + 𝐶 𝜕 𝑦𝑦 𝑓 + 𝐷 = 𝑔(𝑥, 𝑦) With 𝐸 = 𝑏 ; 𝐴 = 𝐶 = 𝑎 2 𝐵 = 0 ; 𝐷 = -𝑐 ; 𝑔(𝑥, 𝑦) = 0

Figure 9 :

 9 Figure 9: Evaluation of discrete models for isotropic Gaussian and mean curvatures

Figure

  Figure 10: Curvature diagram

Figure 11 :

 11 Figure 11: Validation of generation method -Construction of paraboloids from their boundaries as i-liwein surfaces Note:

Figure 12 :

 12 Figure 12: Comparison of the shape of i-minimal and minimal membranes on three sets of boundary curves, with increasing height difference.

Figure 13 :

 13 Figure 13: i-minimal membranes constructed on the same boundary, but with different isotropic directionsThe effect of parameter 𝑏 on the geometry can be visualised in Figure2. It makes the surface tend towards the upper or lower portion of the convex envelope of the boundary curves.

Figure 14 shows

 14 Figure 14 shows Pelikan surfaces obtained by scaling the boundary in the x or y direction of an iminimal surface (middle), creating an i-minimal surface on it, and scaling everything back to the original boundary. This result in with increased transverse (left figure) and longitudinal (right) tension.

Figure 14 :

 14 Figure 14: Middle: i-Minimal surface. Left and right: Pelikan surfaces with respectively increased transverse and longitudinal tension

Figure 15 :

 15 Figure 15: Funicular gridshell with planar quads, constructed by discretizing an i-CMC surface along i-principal curvature directions on the boundary of the Chadstone gridshell (right picture ©Seele)

Figure 16 :

 16 Figure 16: i-LW surfaces generated on the same boundaries as the gridshell of Figure 15. Left: b>0. Right: b<0, close to the limit a²+4bc=0 (this is a close-up of Figure 1)

Figure 17 :

 17 Figure 17: Influence of parameter b on the shape of i-liwein surfaces with same boundary curve and same area. Top: Longitudinal cut on the three shapes. Bottom: Perspective view (left: b=1.4, middle: b=0 iCMC and right: b=-1.7).

Figure 16 :

 16 Figure 16: Comparison of i-curvature lines (left) and i-stress line (right) for an iCMC surface under UPL. The alignment was tested for fixed boundary conditions and, as expected, the perturbations induced by the boundaries affect mainly the neighbourhood of the boundary. It remains very good in the central part of shell.

Figure 16 :

 16 Figure 16: Planarity of quad panels (coloured quad from green 0% to red 0.5% out-of plane default) and torsion at node (coloured spheres from green 0 mrad to red 15 mrad angular torsion at nodes) of a gridshell built on icurvature lines of an i-CMC surface.

•

  𝑁 𝑥𝑥 , 𝑁 𝑥𝑦 , 𝑁 𝑦𝑦 , the components of the membrane tensor 𝐍 in the basis (𝒆 𝒙 𝑺 , 𝒆 𝒚 𝑺 ). 𝑁 𝑥𝑥 𝒆 𝒙 𝑺 + 𝑁 𝑥𝑦 𝒆 𝒚

		𝑺
	(respectively 𝑁 𝑥𝑦 𝒆 𝒙 𝑺 + 𝑁 𝑦𝑦 𝒆 𝒚 𝑺 ) is then the membrane stress vectors acting on a facet with
	normal 𝒆 𝒙 𝑺 (resp. 𝒆 𝒚 𝑺 );
	The coefficients of the membrane stress tensor are then:
	𝑁 𝑥𝑥 = 𝑁 𝑥𝑥 𝑖 cos 𝛽 cos 𝛼	; 𝑁 𝑥𝑦 = 𝑁 𝑥𝑦 𝑖 ; 𝑁 𝑦𝑦 = 𝑁 𝑦𝑦 𝑖 cos 𝛼 cos 𝛽

• 𝛼 and 𝛽 are the slope angles in the 𝑥 and 𝑦 directions: tan 𝛼 = 𝜕 𝑥 𝑓 ; tan 𝛽 = 𝜕 𝑦 𝑓

Table 1 : Geometric characteristics of the three i-liwein surfaces

 1

	a=0.5	C	Min	Max	Max	Min	% flat area
			Height of	Height of	Gaussian	Gaussian	(|K| < 0.0016)
			saddle	dome	curvature	curvature
	b=0	0.083	1.51	2.54	-15.2	6.6	11.5%
	(iCMC)
	b=1.7	0.076	1.13	2.61	-16.4	5.7	9.8%
	b=-1.4	0.086	1.79	2.46	-15.0	7.8	14.7%
	•

Table 2 : Geometric characteristics of the three surfaces

 2

	supports	Loaded area	b=0	b=1.7	b=-1.4	Max-
						min/aver.
	Membrane	Whole	317	288	300	10%
	supports	quarter	324	252	335	27%
	Fixed supports Whole	477	448	464	6%
		quarter	568	552	547	4%

Table 3 : Internal forces in the gridshell built on the i-curvature lines of the iCMC surface

 3

		UPL		UPL on one quarter	ratio
	N [kN]	-11.5	0.02	-40.7	22	4
	Mt [kNm] -0.012	0.012	-0.51	0.46	43	38
	My [kNm] -0.17	0.13	-4.7	4.1	27	30
	Mz [kNm] -0.03	0.04	-1.5	1.5	44	37

Acknowledgement

We would like to warmly thank Pr. István Sajtos for bringing Pelikan surfaces to our attention, as it initiated this project on isotropic linear Weingarten surfaces. We would also like to warmly thank Laurent Monasse for fruitful discussions.


B Rousseau
email: benoit.rousseau@univ-nantes.fr

Y Favennec

J Vicente

F Enguehard

IMPACT OF THE MESO-TEXTURE ON THE EXCTINCTION BEHAVIOUR OF ARCHITECTED POROUS CERAMICS: FROM SPATIAL CORRELATION CHARACTERIZATIONS UP TO (NON) BEERIAN RADIATIVE BEHAVIOUR EVALUATION

Keywords: Architected refractory ceramics, correlated media, Radiative Distribution Function Identification, non beerian behaviour

So far the firm link between the spatial arrangement of consolidated porous ceramics and their ability to attenuate the thermal radiation along its direction of propagation remains a topic of debate. This is particularly true for architected ceramics with regular topologies, obtained by additive manufacturing. To cope with this issue, numerical 3D geometries with different arrangement of cells (regular and irregular) will be generated with genMat an homemade software. Then, the two-point correlation functions of each architecture will be computed to precise the type of spatial correlations prevailing in the medium. The analysis of the extinction cumulative distribution function curves obtained with the Radiative Distribution Function Identification method allows us to check whether or not if the attenuated thermal radiation follows an exponential decay as it the case for Beerian media. Effect of the cell size, cell numbers and cell geometries will be discussed. Finally, preliminar practical considerations will be given for making the distinction between Beerian and non Beerian 3D architectures.

INTRODUCTION

In the field of carbon-free heat generation and waste heat recovery at high temperatures (T>1000°C), there is a growing interest in the design of compact and long-life high temperature energy systems (HTES) such as gas-to-gas heat exchangers [START_REF] Bayaniahangar | Toward extreme high-temperature supercritical CO2 power cycles: Leakage characterization of ceramic 3D-printed heat exchangers[END_REF], volumetric solar receivers [START_REF] Heisel | Digital design and 3D printing of innovative SiC architectures for high temperature volumetric solar receivers[END_REF] (VSR) and radiant tube inserts [START_REF] Pelanconi | Thermal design, optimization and additive manufacturing of ceramic regular structures to maximize the radiative heat transfer[END_REF], among others. These HTES can be described as radiative convective heat exchangers [START_REF] Luque | Exploiting volumetric effects in novel additively manufactured open solar receivers[END_REF]. Advanced HTES, intensively studied today at lab scale, are mostly based on architected refractory ceramics (p ∼75-95%, cell size ∼0.5-10 mm) which can be depicted by a continuous ligament network (the solid phase) delimiting open or closed cells, periodically or pseudo-randomly arranged in space [START_REF] Pelanconi | Cellular ceramic architectures produced by hybrid additive manufacturing: A review on the evolution of their design[END_REF]. The whole 3D architecture exhibits therefore outstanding properties such as high strength-to-weight ratio, high specific surface area, good flow-mixing capacity. According to the intrinsic chemical and physical properties of the material used in the elaboration of the solid phase, the 3D architectures can display additionally high thermal shock resistance, good effective thermal conductivity, high thermal emissivity and high resistance to chemical corrosion. These latter properties can be reached in particular when silicon carbide (SiC)-based compounds are considered [START_REF] Heisel | Digital design and 3D printing of innovative SiC architectures for high temperature volumetric solar receivers[END_REF][START_REF] Pelanconi | Application of ceramic lattice structures to design compact, high temperature heat exchangers: material and architecture selection[END_REF]. The rapid development of additive manufacturing (AM) opens today the gates to elaborate a large range of 3D geometries going from regular lattice structures with different types of unit cells (cube, honeycomb, Kelvin, Weaire-Phelan, octet, star, among others), triply periodic minimal surface-based structures (D-surface, gyroid, neovius, among others) up to more classical irregular strut-based structures (open-cell foams, voronoi lattices). Multi-lattice structures can be also designed by combining several stacks constituted of a given type of cell [START_REF] Avila-Marin | Experimental study of innovative periodic cellular structures as air volumetric absorbers[END_REF]. Such multi-lattice structures can also exhibit a gradual cell size configuration [START_REF] Pelanconi | Cellular ceramic architectures produced by hybrid additive manufacturing: A review on the evolution of their design[END_REF][START_REF] Avila-Marin | Experimental study of innovative periodic cellular structures as air volumetric absorbers[END_REF] either axially, along the main axis of heat transfer, or radially [START_REF] Pelanconi | Thermal design, optimization and additive manufacturing of ceramic regular structures to maximize the radiative heat transfer[END_REF]. In order to propose a topology optimization process aiming at the definition of the best 3D architectures for a prescribed thermal objective, the fine modelling of their thermal behaviour is here essential [START_REF] Pelanconi | Thermal design, optimization and additive manufacturing of ceramic regular structures to maximize the radiative heat transfer[END_REF][START_REF] Ouchtout | Finite element framework for modeling conducto-radiative transfers within heterogeneous media at both discrete and continuous scales[END_REF]. In particular, one of the main challenges, is to be able to exactly describe the role of thermal and/or solar radiation in the heat balance of the system for both unsteady or steady regimes. To date, two main routes can be used to determine the temperature and/or heat flux fields within 3D structures [START_REF] Ouchtout | Finite element framework for modeling conducto-radiative transfers within heterogeneous media at both discrete and continuous scales[END_REF][START_REF] Avila-Marin | Modelling strategies for porous structures as solar receivers in central receiver systems: A review[END_REF][START_REF] Petrasch | Discrete vs. continuum-scale simulation of radiative transfer in semitransparent two-phase media[END_REF][START_REF] Xia | Integrated simulation of continuous-scale and discrete-scale radiative transfer in metal foams[END_REF]. The first class of methods, at the discrete mesoscopic scale, allows a treatment of the radiative heat exchange directly at the scale of the individual constituents (here the struts and walls constituting the basic cells). The discrete scale method gives results that are more accurate but at the price of high computational costs especially when the volume of the numerical domain imposed by the 3D architecture is higher or close to the representative elementary volume. Several parallelized approaches have been developed for that purpose : stabilized vectorial Finite Element schemes [START_REF] Ouchtout | Finite element framework for modeling conducto-radiative transfers within heterogeneous media at both discrete and continuous scales[END_REF], cell centered Finite Volume models associated to ray tracing [START_REF] Kumar | A numerical method based on domain decomposition to solve coupled conduction-radiation physics using parallel computing within large porous media[END_REF][START_REF] Perraudin | Numerical quantification of coupling effects for radiation-conduction heat transfer in participating macroporous media: Investigation of a model geometry[END_REF], hybrid scheme combining motion of Brownian walkers in the solid opaque phase and ray tracing in the fluid phase [START_REF] Vignoles | A hybrid random walk method for the simulation of coupled conduction and linearized radiation transfer at local scale in porous media with opaque solid phases[END_REF], and Monte Carlo strategy based on integral formulation [START_REF] Caliot | Combined conductive-radiative heat transfer analysis in complex geometry using the Monte Carlo method[END_REF]. The second class of method, at the continuous macroscopic scale, requires to rigorously solve the radiative transfer equation as long as the assumptions (randomness, homogeneity and continuity) that govern its statement are valid. Let us precise that an architected ceramic can be viewed as an equivalent domain absorbing, emitting and scattering the thermal radiation, these phenomena being induced by light-matter interaction at the scale of the ligaments. If the equivalent medium is optically thick, the radiative transfer equation can be reduced to a more simple diffusion equation by applying for example the spherical harmonics method at the order 1 (P1 approximation) and, even more crudely, the Rosseland approximation [START_REF] Avila-Marin | Modelling strategies for porous structures as solar receivers in central receiver systems: A review[END_REF]. Whatever the precision sought, these approaches can be only applied if the volume radiative properties (extinction, absorption and scattering coefficients, scattering phase function) are identified from an equivalent semitransparent continuous Beerian medium. The term Beerian stands, here, for an extinction law characterised by an exponential function of the optical thickness. From a radiative viewpoint, the medium must be, in principle, constituted of individual objects (particles, grains, fibers) spatially uncorrelated. This statement changes when spatial correlations occur meaning that the equivalent semitransparent medium becomes not easily homogeneisable. In this case, the extinction coefficient conventionally obtained by adding the individual contribution of the scatterers must be successively corrected by a scaling factor computed from the radial distribution function which depicts, here, the spatial correlations [START_REF] Drolen | Independent and dependent scattering in packed-sphere systems[END_REF]. This leads to the appearance of either super-exponential or sub-exponential extinction coefficients respectively for negative and positive correlated medium [START_REF] Shaw | Super-exponential extinction of radiation in a negatively correlated random medium[END_REF]. When the solid matter is connected, as it the case for architected ceramics, few studies have established a firm link between the spatial organisation of the scatterers and their radiative properties. In general, the samples are considered as Beerian media [START_REF] Pelanconi | Thermal design, optimization and additive manufacturing of ceramic regular structures to maximize the radiative heat transfer[END_REF][START_REF] Avila-Marin | Experimental study of innovative periodic cellular structures as air volumetric absorbers[END_REF][START_REF] Tseng | Effect of foam properties on radiative properties of opencell silicon carbide foams[END_REF] without a detailed investigation of their ability to transport radiation having been carried out beforehand. By using the radiative distribution function identification method, Taine et al. have shown that a dispersion of overlapping transparent spheres in an opaque continuous hosting medium deviates from the Beerian behavior when its volume fraction, 𝑓𝑓 𝑣𝑣 , becomes important (𝑓𝑓 𝑣𝑣 > 0.4) [START_REF] Taine | Generalized radiative transfer equation for porous medium upscaling: Application to the radiative Fourier law[END_REF]. Such a synthetic medium can be used to represent a statistically isotropic foam. Guevelou et al. computed the extinction cumulative distribution functions of virtual SiC-based irregular samples (cubic open-cell foams with domain size, 𝐿𝐿 = 4.8-mm) for variable nominal pore diameters (𝑑𝑑 𝑛𝑛𝑛𝑛𝑛𝑛 ∈ [0.25-1.1 mm]) and fixed porosity (𝑝𝑝=0.7) and concluded that they complied with a Beerian behavior [START_REF] Guevelou | Evolution of the homogenized volumetric radiative properties of a family of α-SiC foams with growing nominal pore diameter[END_REF]. More recently, Ouchtout et al. have computed the temperatures fields within a SiC-based regular cubic lattice (𝑝𝑝 = 0.9, 𝑑𝑑 𝑛𝑛𝑛𝑛𝑛𝑛 = 0.15 mm, 𝐿𝐿 = 0.3 and 0.65 mm) by solving a conductive-radiative numerical scheme with stabilized vectorial finite element methodologies applied respectively at both the continuous scale and the discrete scale. Both the temperature profiles are closely identical [START_REF] Ouchtout | Finite element framework for modeling conducto-radiative transfers within heterogeneous media at both discrete and continuous scales[END_REF]. Since the continuous approach assumes for the radiative part that the Beer law is valid, it can be concluded that this cubic lattice exhibits a Beerian behavior. Such a survey indicates that a study investigating the radiative behavior of synthetic lattices as a function of variable textural features (spatial repartition, shape and size of cells, size of the computing domain) is still relevant. To go one step further with these issues, a set of virtual regular and irregular SiC-based structures (cubic cell structures, Kelvin cell structures, gyroids, and open-cell foams) with prescribed textural features ( 𝑑𝑑 𝑛𝑛𝑛𝑛𝑛𝑛 ∈ [0.5-10 mm], 𝐿𝐿 ∈ [10 -100 mm], 𝑝𝑝 ∼ 0.92, regular and irregular) will be designed with the genMat software (C++, Qt) [START_REF] Guévelou | Representative elementary volumes required to characterize the normal spectral emittance of silicon carbide foams used as volumetric solar absorbers[END_REF]. Their two-point auto-correlation functions for the void phase will be systematically computed.

IHTC-17 | ID: 532

Then, the radiative distribution function identification method will be applied to compute the extinction cumulative distribution functions of the whole set of numerical samples. This will allow us to provide preliminar practical criteria for defining the radiative behaviors (non Beerian or Beerian) of the samples.

NUMERICAL METHODOLOGIES

3D structures generation

In this work, a homemade 3D structures generator, called genMat (C++, Qt) [START_REF] Guévelou | Representative elementary volumes required to characterize the normal spectral emittance of silicon carbide foams used as volumetric solar absorbers[END_REF], has been used to provide the whole set of numerical samples. The generator allows the numerical elaboration of regular and irregular 3D structures with prescribed textural features (porosity, nominal pore diameter). It was initially settled for designing open-cell foams with realistic textural features [START_REF] Guevelou | Evolution of the homogenized volumetric radiative properties of a family of α-SiC foams with growing nominal pore diameter[END_REF]. For irregular open cell foams, four main steps must be followed. Firstly, the cell centers are positioned at the apexes of a regular tetrahedral network. A random perturbation, whom the value is lower than the side tetrahedron, is performed for slightly displacing the centers. Secondly, a fast marching-cube process is implemented with the aim to reproduce a realistic bubbling process. A distance map, where each voxel knows the value of the space separating it from the nearest neighboring center, is created. This set-up permits the simultaneous growth of the entire cells and avoids eventual overlaying. A watershed algorithm is afterwards applied on the distance map for the segmentation of the void phase and the solid phase by a Voronoï tessellation process. As such, a thin continuous skeleton of the 3D structures can be obtained. Thirdly, each strut is shaped with a spherical structuring element of variable side. Finally, in order to accurately tune the final porosity, the fast-marching algorithm is again used to compute a distance map within the void phase from the solid-void interface. The resulting distance map is then converted into grey level for which the adequate binarization threshold provide the targeted porosity. For the regular structures, a simplified four-steps process can be followed where centers of all cells are firstly created. Secondly, the nodes, edges, faces of all cells are defined, and a structuration loop adds solid voxels between each apex, forming a thin skeleton. Thirdly, as for irregular foams, the ligaments can be shaped with a spherical structuring element. Lastly, the desired porosity is reached by creating, as for the irregular foams, a distance map converted in grey levels. Once the 3D binarized structures are generated, a robust marching cube process allows to mesh the fluid/solid interface. GenMat proposes an export module for recording surfacic meshes. Let us precise that for triply periodic minimal surface-based structures, a Matlab script is used to produce the 3D meshed geometries, which are afterwards imported in genMat.

Textural properties characterization

genMat contains also modules which provide the volume fraction of the void phase i.e. the porosity, 𝑝𝑝, and the volumetric surface, 𝐴𝐴 of each numerical sample. For this work, a supplementary routine has been added for computing their two-point correlation functions by following the process described in the reference book of Torquato [START_REF] Torquato | Random Heterogeneous Materials Microstructure and Macroscopic Properties[END_REF]. Briefly, by using a given binarized 3D image, its phase or indicator function, 𝑍𝑍 𝑓𝑓 (𝒙𝒙) of the void or fluid space can be determined. Here 𝒙𝒙 is the position vector. 𝑍𝑍 𝑓𝑓 (𝒙𝒙) is equal to 1 if 𝒙𝒙 belongs to the void phase and, otherwise, to 0. The two-point correlation function, 𝑆𝑆 2 𝑓𝑓 (𝒙𝒙 𝟏𝟏 , 𝒙𝒙 𝟐𝟐 ), is defined by:

𝑆𝑆 2 𝑓𝑓 (𝒙𝒙 𝟏𝟏 , 𝒙𝒙 𝟐𝟐 ) = 〈𝑍𝑍 𝑓𝑓 (𝒙𝒙 𝟏𝟏 )𝑍𝑍 𝑓𝑓 (𝒙𝒙 𝟐𝟐 )〉 (1)
The angular brackets represent the statiscal average. This second-order descriptor can be interpreted as the probability that two randomly points 𝒙𝒙 𝟏𝟏 and 𝒙𝒙 𝟐𝟐 both belong to the fluid phase. If the medium is statically homogenous and isotropic, 𝑆𝑆 2 𝑓𝑓 depend only on the segment of length 𝑟𝑟 = |𝒓𝒓|, such as :

IHTC-17 | ID: 532 𝑆𝑆 2 𝑓𝑓 (𝒙𝒙 𝟏𝟏 , 𝒙𝒙 𝟐𝟐 ) = 𝑆𝑆 2 𝑓𝑓 ( |𝒓𝒓|) (2)
where the displacement vector 𝒓𝒓 = 𝒙𝒙 𝟐𝟐 -𝒙𝒙 𝟏𝟏 . Generally, for 𝒓𝒓 = 0, 𝑆𝑆 2 𝑓𝑓 (0) = 𝑝𝑝 and in the absence of any long-range order,

lim |𝒓𝒓|→∞ 𝑆𝑆 2 𝑓𝑓 ( |𝒓𝒓|) → 𝑝𝑝 2 (3)
The analysis of the shape of 𝑆𝑆 2 𝑓𝑓 ( |𝒓𝒓|) provides useful information on how a natural or synthetic porous medium is structured [START_REF] Torquato | Random Heterogeneous Materials Microstructure and Macroscopic Properties[END_REF]. In the case of artificial media, many studies taking into account 𝑆𝑆 2 𝑓𝑓 ( |𝒓𝒓|) have focused on packed beds. Thus, for a non-penetrable monodisperse sphere model of packed bed, 𝑆𝑆 2 𝑓𝑓 ( |𝒓𝒓|)

decreases quasi lineraly at short |𝒓𝒓| and then exhibit an oscillary behavior up to the asymptotic value 𝑝𝑝 2 . The distance after which 𝑆𝑆 2 𝑓𝑓 ( |𝒓𝒓|) becomes asymptotic is the correlation length. On the contrary, when spheres can overlap, after the linearly decrease previously evoked, the oscillations disapear. Let us underline that for a perfect regular arrangement of spheres, the amplitudes and the periods of the oscillations remain identical for growing values of |𝒓𝒓|. 𝑆𝑆 2 𝑓𝑓 ( |𝒓𝒓|) can be mathematicaly described from an exponentially decaying function and more sophiticated expresions are available to reproduce damped oscillation. This brief overwiew can be used to analysis the textural behavior of any porous media.

Radiative distribution function identification method

In this work, the radiative distribution function identification (RDFI) method, originally proposed by Tancrez and Taine [START_REF] Tancrez | Direct identification of absorption and scattering coefficients and phase function of a porous medium by a Monte Carlo technique[END_REF], will be applied to determine the volumetric radiative properties of the 3D architectures numerically generated. This method had been implemented in the genMat framework [START_REF] Guevelou | Evolution of the homogenized volumetric radiative properties of a family of α-SiC foams with growing nominal pore diameter[END_REF]. More precisely, the RDFI is a collision-based Monte Carlo Ray Tracing method assuming that for each local interaction between a ray and the solid matter, the geometrical optic approximation is valid. It means that the mean characteristic size of the ligaments and of the cells constituting the generated 3D structures must be much larger than the incident wavelength i.e. the Mie size parameter is greater than 1. The RDFI method consists in launching a huge number of rays from 𝑀𝑀 points that are uniformly distributed in the fluid phase surrounding the ligaments. Their directions are randomly chosen within the entire 4𝜋𝜋 steradian sphere. Once a ray hits a solid phase element at the point 𝐼𝐼, it is considered as totally extinguished and the path length to collision, 𝑠𝑠 𝑒𝑒 = |𝑀𝑀𝐼𝐼| is computed. For a large number of rays, 𝑁𝑁 𝑅𝑅 , a set of path length 𝑠𝑠 𝑒𝑒,𝑗𝑗 values are obtained from which the corresponding normalized extinction distribution function, 𝐹𝐹 𝑒𝑒 (𝑠𝑠) is computed :

𝐹𝐹 𝑒𝑒 (𝑠𝑠) = 1 𝑁𝑁 𝑟𝑟 ∑ 𝛿𝛿�𝑠𝑠 -𝑠𝑠 𝑒𝑒,𝑗𝑗 � 𝑁𝑁 𝑅𝑅 𝑗𝑗=1 (4)
where 𝑠𝑠 is the path length variable and 𝛿𝛿 the Dirac delta distribution. The extinction cumulative distribution function 𝐺𝐺 𝑒𝑒 (𝑠𝑠) is then obtained by integrating 𝐹𝐹 𝑒𝑒 (𝑠𝑠) :

𝐺𝐺 𝑒𝑒 (𝑠𝑠) = ∫ 𝐹𝐹 𝑒𝑒 (𝑠𝑠 ′ )𝑠𝑠 ′ 𝑠𝑠 0 (5)
To decide whether the extinguished ray is absorbed or scattered, a random number, 𝜂𝜂, uniformly distributed over [0,1] is selected and compared to the local reflectivity, 𝜌𝜌, that can be given by the Fresnel law (if the interacting fluid/solid interface is considered as optically smooth). If 𝜂𝜂 > 𝜌𝜌 the ray is considered as absorbed and otherwise scattered. This allows to determine the absorption and the scattering cumulative probability with a similar numerical process than for 𝐺𝐺 𝑒𝑒 (𝑠𝑠). If the studied IHTC-17 | ID: 532 equivalent medium is Beerian, the quantity �1 -𝐺𝐺 𝑒𝑒 (𝑠𝑠)� shows a typical Beer-Lambert exponential decay behaviour expressed by:

1 -𝐺𝐺 𝑒𝑒 (𝑠𝑠) = 𝑒𝑒 -𝛽𝛽𝑠𝑠 = 𝐼𝐼(𝑠𝑠) 𝐼𝐼 0 (6)
where 𝐼𝐼(𝑠𝑠) is the radiant intensity at the distance 𝑠𝑠 and 𝐼𝐼 0 is the impinging intensity at 𝑠𝑠 = 0 from a collimated beam. In that situation, the extinction coefficient 𝛽𝛽 of the equivalent medium can be evaluated by a fit of the numerically evaluated �1 -𝐺𝐺 𝑒𝑒 (𝑠𝑠)� function.

RESULTS AND DISCUSSION

In the following, the room-temperature complex refractive index of silicon carbide at the wavelength of 4 µm, is affected to the solid skeleton of each numerical 3D structure [START_REF] Guévelou | Representative elementary volumes required to characterize the normal spectral emittance of silicon carbide foams used as volumetric solar absorbers[END_REF]. It ensures that their ligaments can be considered as optically thick, their lowest diameter being fixed at 400 µm. The computational domains have a cubic geometry with a side, 𝐿𝐿, of 50 mm. For a sake a clarity, the acronyms CCS, KCS, GS and OCFS represent respectively cubic cell structures, Kelvin cell structures, gyroid structures and open-cell foam structures. To complete their denomination, one adds an information for the 𝑑𝑑 𝑛𝑛𝑛𝑛𝑛𝑛 , the number of cell, 𝑁𝑁 𝑐𝑐 , and the size of the domain. CCS_10_125_50 depicts a cubic cell structures with 𝑑𝑑 𝑛𝑛𝑛𝑛𝑛𝑛 = 10 mm, 𝑁𝑁 𝑐𝑐 = 125 and 𝐿𝐿 = 50 mm.

Influence of the spatial correlations for regular and irregular structures

Four 3D structures (𝐿𝐿 = 50 mm, 𝑑𝑑 𝑛𝑛𝑛𝑛𝑛𝑛 ∼ 10 mm, voxel number = 500 3 , voxel size = 100 µm) are generated with genMat. 1 gives their respective porosity and volumetric surface after having performed the final porosity tuning process as explained in section 2.1. Fig. 1 depicts the synthetic 3D images where open channels can be obviously observed for CCS_10_125_50 and KCS_10_125_50. The two-point correlations, computed for all the directions, are plotted in Fig. 2a). For CCS_10_125_50, KCS_10_125_50 and GS_10_125_50, the curves show similar features namely a rapid quite linear decay at short 𝑟𝑟 followed by several oscillations, at larger 𝑟𝑟, with the periodicity imposed by the spatial arrangement of their respective unit cells. This is reflective of spatial correlations between the constituents of the structures [START_REF] Torquato | Random Heterogeneous Materials Microstructure and Macroscopic Properties[END_REF]. Note that for KCS_10_125_50 and GS_10_125_50, the amplitude of the oscillations are less intense than for CCS_10_125_50, meaning that the effect of the spatial periodicity for the range of 𝑟𝑟 investigated is less pronounced. It can be noticed also that the slopes of the linear decay, as expected, is all the more sharp than 𝐴𝐴 is important (see Table 1). For OCF_10_125_50, the shape of 𝑆𝑆 2 𝑓𝑓 ( |𝒓𝒓|) is rather different with a first sharp decay at the shorts 𝑟𝑟 followed by a slower decreasing at higher 𝑟𝑟. Slight oscillations can be seen. This indicates that the ligaments of OFC_10_125_50 appears as being less correlated than the others structures. Fig. 2b) shows extinction cumulative distribution functions 𝑙𝑙𝑙𝑙�1 -𝐺𝐺 𝑒𝑒 (𝑠𝑠)� as a function of 𝑠𝑠, obtained for 𝑁𝑁 𝑅𝑅 = 10 6 . Such a value of 𝑁𝑁 𝑅𝑅 guarantees to results statistical flutuations less than 0.1. Special care has been taken to define a suitable shooting zone for computing 𝑙𝑙𝑙𝑙�1 -𝐺𝐺 𝑒𝑒 (𝑠𝑠)� i.e., not too small to be representative, and not too large to limit outgoing rays. Only a central zone of 400 3 voxels is considered for each computation. Two behaviours can be observed in Fig. 2b). CCS_10_125_50 and KCS_10_125_50 present curves which deviate clearly from what is expected for a Beerian behaviour, whereas GS and OCFS exhibit a more linear feature. In the insight part of Fig. 2b) the same curves are given in function of the optical thickness 𝛽𝛽 𝑛𝑛𝑛𝑛 (𝑠𝑠 ′ -𝑠𝑠), 𝛽𝛽 𝑛𝑛𝑛𝑛 being the value of the extinction of the optically thin limit (𝛽𝛽 𝑛𝑛𝑛𝑛 = 𝐴𝐴 4𝑝𝑝 ⁄ ). They display a similar radiative behaviour than those described by Shaw et al. for negatively correlated particulate media [START_REF] Shaw | Super-exponential extinction of radiation in a negatively correlated random medium[END_REF]. It is worth of notify that the extinction cumulative distribution function curve of OCFS_10_125_50 shows the more linear feature knowing that this 3D structure is surely the less spatially correlated. In the opposite, CCS_10_125_50 have the stronger non-Beerian feature whereas it is the more spatially correlated. Higher values of 𝛽𝛽 𝑛𝑛𝑛𝑛 that can be numerically settled, for similar values of 𝑝𝑝 and 𝐿𝐿 by an higher value of 𝐴𝐴, is likely to confer a Beerian behaviour to the 3D structures. In other words, it turns to decrease 𝑑𝑑 𝑛𝑛𝑛𝑛𝑛𝑛 once 𝑝𝑝 and 𝐿𝐿 are imposed. To go one step further, each extinction cumulative distribution function curve has been nicely fitted thanks to a Levengberg-Marquardt scheme by the following expression:

𝑦𝑦(𝑠𝑠) = -𝑎𝑎�𝑒𝑒 𝑏𝑏𝑠𝑠 -1� (7)
In Gomez-Garcia et al, the authors proposed a quite similar expression for analysing the radiation propagation within a stack of SiC-based square grids [START_REF] Gomez-Garcia | Numerical analysis of radiation propagation in a multi-layer volumetric solar absorber composed of a stack of square grids[END_REF]. of 𝑎𝑎 and the lowest value of 𝑏𝑏 are required to confer a Beerian behaviour for the numerical 3D architectures. Let us now introduce a criterion accuracy 𝜀𝜀 with 𝑏𝑏 = 𝜀𝜀 𝐿𝐿 ⁄ and 𝐿𝐿 = 50 mm. As an example, for 𝜀𝜀 = 0.1, one obtains 𝑏𝑏 = 2 m -1 that yields 𝑦𝑦(𝑠𝑠) ∼ 2𝑎𝑎𝑠𝑠 for 0 < 𝑠𝑠 <𝐿𝐿 when the exponential term is approximated by its Taylor expansion. In this case, 𝛽𝛽∼2𝑎𝑎 knowing that 𝑎𝑎 is obtained through the latter fitting curve process. It suggests that obtaining values 𝑎𝑎 close to threshold value of 2 m -1 must be concomitant with a Beerian behaviour. To highlight this quantitative criterion, two numerical samples (KFCS, OFCS) with 𝑑𝑑 𝑛𝑛𝑛𝑛𝑛𝑛 = 5 mm and 𝐿𝐿 = 50 mm are generated and the fits, by using Eq. 7, of their linear extinction cumulative distribution functions give respectively 𝑏𝑏 = 7.3 m -1 and 6.8 m -1 (in comparison for 𝑑𝑑 𝑛𝑛𝑛𝑛𝑛𝑛 = 10 mm, 𝑏𝑏 = 28.7 m -1 and 20 m -1 at is reported in Table 1). By taking a cruder criterion (𝜀𝜀 = 0.4, giving 𝑏𝑏 = 8 m -1 ), one can consider, in the following, that all extinction cumulative distribution functions for which the fitting curve process with Eq. 7 conduct to 𝑏𝑏 < 8 m -1 should exhibit therefore a Beerian behaviour.

Influence of the cell number for a regular structure

From the previous section, it can be argued that the conjunction of weak spatial correlations and low cell size allows to confer to a 3D structure a Beerian behaviour. Let us now consider the effect of a growing number of cells with 𝑑𝑑 𝑛𝑛𝑛𝑛𝑛𝑛 equals here to 5 mm. Table 2 lists the main textural parameters used to design the numerical samples by genMat with 2 3 , 4 3, 6 3 , 8 3 and 10 3 cells. The numerical samples are described in Fig. 3 Table 2 Their two-point autocorrelation functions are given as a function of 𝑟𝑟 𝐿𝐿 ⁄ in Fig. 4a). The amplitudes of the oscillations tend to decrease as 𝐿𝐿 becomes large. It means the influence of the lattice periodicity weakens when the size of the computational domain is important. Moreover, the slopes of the curves increase when 𝑁𝑁 𝑐𝑐 grows up that is coherent with the computation of 𝐴𝐴 with genMat. Corollary, as it is reported in Fig4b), the structures present a non Beerian feature less and less pronounced when 𝑁𝑁 𝑐𝑐 comes to be large. A similar trend can be found for numerical samples of the KFS family (except for 𝑁𝑁 𝑐𝑐 = 8) with same characteristic size. The non Beerian feature turns to be less important when the spatial correlations imposed by the 3D architecture vanish at large 𝑟𝑟. By performing a fit process, based on Eq. ( 7), the curves of Fig4b) are perfectly reproduced. For each curves the residuals are less than 0.05. As given in Table 2, it can be notified that for an increasing value of 𝑁𝑁 𝑐𝑐 , 𝑎𝑎 grow up and 𝑏𝑏 decline, the IHTC-17 | ID: 532 product 𝑎𝑎𝑏𝑏 becoming more important. But according to the criterion proposed in the last section, it would be necessary to have larger 𝐿𝐿 (> 50 mm) to tend towards the Beerian case.

Influence of the domain size for regular structures with same number of cells

Let us now investigate the effect of the size of the computational domain for Kelvin-cell structures with a same number of cells (𝑁𝑁 𝑐𝑐 = 1000) but with a decreasing domain size (𝐿𝐿 is going from 50 up to 5 mm). 𝑁𝑁 𝑐𝑐 = 1000 ensures that the medium presents the best Beerian feature. This means that 𝑑𝑑 𝑛𝑛𝑛𝑛𝑛𝑛 falls down from 5 to 0.5 mm. In parallel, 𝐴𝐴 grows up softly when 𝐿𝐿 decreases from 50 to 30 mm to reach then more abruptly very high values when 𝐿𝐿 drops at 5 mm (see Table 3). In Fig. 5a), the two-point autocorrelation functions for the fluid phase of the whole set of samples can be nicely superimposed when they are given as a function of 𝑟𝑟/𝐿𝐿. The homothetic process tends to preserve the spatial correlations prevailing in KFS_5_1000_50. The extinction cumulative distribution function curves depicted in Fig. 5b) appears to be quite similar that is consistent with the preservation of the spatial arrangement of the network during the homothetic process. The curve fitting process reveals that the product 𝑎𝑎𝑏𝑏 drops with the same rate than 𝛽𝛽 𝑛𝑛𝑛𝑛 . However, whereas b decreases from 𝐿𝐿 = 5 mm up to 𝐿𝐿 = 50 mm, a shows an inhomogeneous evolution. A similar homothetic scheme applied to a set of 125 Kelvin cells forming a cubic domain with decreasing 𝐿𝐿 (from 50 mm down to 5 mm) confirms the latter evolution namely the extinction behaviour is preserved. In this case, the extinction cumulative distribution function curves can be fitted with Eq. 7 given b = 34.2 m -1 that is typical of a non-Beerian feature. From Table 3, it appears that KFS_5_1000_50 is the sample that can show the best Beerian. In doing so, its extinction coefficient will be 8𝑎𝑎 namely 137.6 m -1 . Let us recall that only for such a case, the application of a continuous scale method, based on the Beerian framework, for solving the radiative transfer equation is required for any heat transfers issues.

CONCLUSION

This work combines the numerical generation of highly porous 3D architectures (CCS, KFS, GS and OFCS) with a continuous solid networks and the determination of their respective two-point correlation functions and extinction cumulative distribution functions. All the numerical samples are spatially correlated with a degree of correlations which reduces when the structures become irregular. In turn, weak spatial correlations favour a Beerian behaviour which are all so more important that the volumetric surface is important. For non Beerian 3D structures, which are mainly periodic architected structures, a practical exponential growth law is given to describe their extinction behaviour as a function of the variable path length. A preliminary criterion is provided to make a discrimination between Beerian and non Beerian architectures from the fitting curve process of their extinction cumulative distribution functions. For this class of architected materials, more works are now required to establish a firm link between their extinction behaviours and their first-and second-order textural features.

ACKNOWLEGDEMENT

The authors thank the French Research Network called GDR TAMARYS and funded by CNRS for providing motivating and high-level framework for scientific exchanges.

NOMENCLATURE

IHTC- 17 | ID: 532 Fig. 1 Fig. 2 a

 1753212 Fig.1 3D images of 4 porous structures with 𝐿𝐿 = 50 mm, 𝑝𝑝 ∼0.93 and 𝑑𝑑 𝑛𝑛𝑛𝑛𝑛𝑛 ∼10 mm (1) CCS_10_125_50 (2) KCS_10_125_50, (3) GS_10_125_50 (4) OCFS_10_125_50

Fig. 3

 3 Fig.3 3D images of (1) CCS_5_8_15 (2) CCS_5_64_25, (3) CCS_5_216_35, (4) CCS_5_512_45, (5) CCS_5_1000_55

Fig4

  Fig4 a) Two-point autocorrelation function: (1) CCS_5_8_10, (2) CCS_5_64_20, (3) CCS_5_216_30, (4) CCS_5_512_40, (5) CCS_5_1000_50 b) Extinction cumulative distribution function: (1) CCS_5_8_10, (2) CCS_5_64_20, (3) CCS_5_216_30, (4) CCS_5_512_40, (5) CCS_5_1000_50

Fig 5

 5 Fig 5 a) Two-point autocorrelation function: (1) KFS_5_1000_50, (2) KFS_4_1000_40, (3) KFS_3_1000_30, (4) KFS_2_1000_20, (5) KFS_1_1000_10, (6) KFS_0.5_1000_5 b) Extinction cumulative distribution function: (1) KFS_5_1000_50, (2) KFS_4_1000_40, (3) KFS_3_1000_30, (4) KFS_2_1000_20, (5) KFS_1_1000_10, (6) KFS_0.5_1000_5

Table 1

 1

	Textural and radiative parameters
	Samples	𝑝𝑝	𝐴𝐴 (m -1 )	𝛽𝛽 𝑛𝑛𝑛𝑛 (m -1 )	a	b (m -1 )
	CCS_10_125_50	0.93	138.2	36.9	1.20	43.1
	KFS_10_125_50	0.93	256.7	68.7	2.79	28.7
	GS_10_125_50	0.95	330.0	86.7	3.76	28.4
	OCFS_10_125_50	0.93	297.2	79.9	6.55	20.0
	Table

Table 1

 1

lists values of 𝑎𝑎 and 𝑏𝑏 for CCS_10_125_50, KCS_10_125_50, GS_10_125_50 and OFCS_10_125_50. Going from CCS_10_125_50 up to OFCS_10_125_50, it can be evidenced from the Table 1 that the highest value a) b) IHTC-17 | ID: 532

Table 3

 3 Textural and radiative parameters

	Samples	𝑝𝑝	𝐿𝐿 (mm)	𝐴𝐴 (m -1 )	𝛽𝛽 𝑛𝑛𝑛𝑛 (m -1 )	a	b (m -1 )
	KFS_0.5_1000_5	0.91	5	5049	1387.10	8.39	176.9
	KFS_1_1000_10	0.91	10	2888	793.10	13.6	65.5
	KFS_2 _1000_20	0.91	20	1309	359.62	10.1	39.2
	KFS_3 _1000_30	0.91	30	1007	273.64	22.7	13.6
	KFS_4 _1000_40	0.91	40	690	185.48	18.6	11.7
	KFS_5_1000_50	0.91	50	546	148.37	17.2	9.9

Proceedings of the 17th International Heat Transfer Conference, IHTC-