Torben Prill
email: prill@itwm.fhg.de

Dominique Jeulin

François Willot

Juan Balach

Flavio Soldera

Prediction of Effective Properties of Porous Carbon Electrodes from a Parametric 3D Random Morphological Model

Keywords: Porous electrodes, Double-layer Capacitor, FIB-SEM nanotomography, Stochastic modeling

Pore structures have a major impact on the transport and electrical properties of electrochemical devices, such as batteries and electric double-layer capacitors (EDLCs). In this work we are concerned with the prediction of the electrical conductivity, ion diffusivity and volumetric capacitance of EDLC electrodes, manufactured from hierarchically porous carbons. To investigate the dependence of the effective properties on the pore structures, we use a structurally resolved parametric model of a random medium. Our approach starts from 3D FIB-SEM imaging, combined with automatic segmentation. Then, a random set model is fitted to the segmented structures and the effective transport properties are predicted using full field simulations by iterations of FFT on 3D pore space images and calculations based on the geometric properties of the structure model. A parameter study of the model is used to investigate the sensitivity of the effective conductivity and diffusivity to changes in the model parameters. Finally, we investigate the volumetric capacitance of the EDLC electrodes with a geometric model, make a comparison with experimental measurements and do a parameter study to suggest improved microstructures.

Introduction

Porous carbon materials are widely used as electrode materials in energy storage devices, such as electrical double-layer capacitors (EDLCs) [START_REF] Conway | Electrochemical Supercapacitors: Scientific Fundamentals and Technological Applications[END_REF]. These devices are used as alternatives or in combination with batteries. Yet, in general, they have a lower specific capacitance per volume and per weight and a higher specific power. In EDLCs, the electric energy is stored in a thin layer on the surface of the porous electrode which makes micro-or nanoporous materials with their large specific surface area especially suitable, since this results in a very high capacitance. Yet, the performance, measure in specific power and capacitance depends on the specific pore shape and size of the electrode materials. Hence, to predict the performance of an electrode material, knowledge and modeling of the morphology of the pore space is necessary, aside from the electrochemical modeling.

To this end, two nanoporous carbon-based materials, used in electrodes of EDLCs, are investigated. The two samples are imaged with FIB-SEM and the images are segmented using a new segmentation algorithm using mathematical morphology as in [START_REF] Prill | Morphological segmentation of FIB-SEM data of highly porous media[END_REF]. Using the segmented microstructures, a stochastic model is defined for a two-phase heterogeneous material. The first part of the modeling consists of defining a random set model depending on a set of free parameters. Then the best fit parameters are identified by matching the morphological characteristics of the observed microstructures to model realizations. The second part of the paper consists of a model for the effective transport properties of the modeled microstructures. Using the model, the specific conductivity and diffusion resistance of both materials are predicted by simulations using FFT-based methods. Then, we performed a parameter study to estimate the sensitivity of the effective properties with respect to the model parameters. Finally, we use geometric properties of the random set model to predict the volumetric capacitance of the electrodes and make a comparison with measurements.

State of the Art

Image Analysis and Segmentation

Reconstruction of highly porous three-dimensional pore spaces from FIB-SEM imaging is still in general an unsolved problem. Different methods have been tried, such as thresholding approaches as in [START_REF] Nanjundappa | A customized framework for 3-d morphological characterization of microporous layers[END_REF] or surface evolution as in [START_REF] Jørgensen | A framework for automatic segmentation in three dimensions of microstructural tomography data[END_REF]. Other methods include the ones shown in [START_REF] Salzer | A two-stage approach to the segmentation of FIB-SEM images of highly porous materials[END_REF], such as threshold backpropagation or valley detection. In this study we use the method presented in [START_REF] Prill | Morphological segmentation of FIB-SEM data of highly porous media[END_REF], based on mathematical morphology. Since it has been shown, that even accurately segmented FIB-SEM data can lead to false transport properties, we combine the automatic segmentation with stochastic modeling, as in [START_REF] Hutzenlaub | How coarsening the 3d reconstruction of a porous material influences diffusivity and conductivity values[END_REF]. Hence, we define a parametric stochastic model, which visually resembles the investigated materials. This allows for a better estimation of the transport propertied as well as a systematic alteration of the model parameters and the optimization of the model with respect to performance measures, such as capacitance and diffusion resistance.

Physical Modeling

Electrochemical double-layer capacitors are devices for the storage of electrical energy by means of an electrical double-layer forming on the surface of the electrode. Physically, energy is stored when charge carriers adsorb to the surface of the electrode and induce a strong electric field across the boundary. This constitutes the so called electrical double layer. The earliest theory of the formation of the double layer was developed by Helmholz in 1879. Later the theory was amended by Stern and then Guy and Chapman. An historical overview can be found in [START_REF] Conway | Electrochemical Supercapacitors: Scientific Fundamentals and Technological Applications[END_REF]. Modern approaches for modeling the double layer include the theories by Bazant, starting with [START_REF] Bazant | Diffuse-charge dynamics in electrochemical systems[END_REF].

As the capacitance of the electrode increases with the specific surface area of the electrode, micro-or nanoporous electrodes lend themselves as electrode materials. Different models for porous electrodes for EDLCs have been developed. The simplest ones are based on equivalent circuits, where an overview can be found in [START_REF] Barsoukov | Impedance Spectroscopy: Theory, Experiment, and Applications[END_REF]. A more detailed treatment of the porous electrode has been given in macrohomogeneous models pioneered by de [START_REF] De Levie | On porous electrodes in electrolyte solutions: I. capacitance effects[END_REF], and extended e.g. in [START_REF] Paasch | Theory of the electrochemical impedance of macrohomogeneous porous electrodes[END_REF] and [START_REF] Roßberg | The influence of porosity and the nature of the charge storage capacitance on the impedance behaviour of electropolymerized polyaniline films[END_REF]. A spatially resolving model is given, e.g., in [START_REF] Torquato | Physical interpretation of cyclic voltammetry for measuring electric double layer capacitances[END_REF].

In the present paper we will restrict ourselves to the most basic treatment of the double layer as an areal capacitance and focus on the geometrical aspects of the porous electrode. To this end we are using a stochastic model to represent the morphology of the pore space. The model is based on a Boolean Model of spheres combined with a convolution to achieve a smoothing effect on the geometry.

Materials and Imaging

The starting point for the study in this paper is segmented FIB-SEM images of two samples of porous carbon electrodes.

Materials

In this work, there are two samples of different electrode materials under investigation, denoted S12 and S14. Both have been synthesized by the same production process, yet with a slightly different composition, leading to a difference in pore sizes. Both materials have a hierarchical pore space, with mesopores, with a typical pore size of 2nm -50nm and micropores, with a typical pore size smaller than 2nm.

The samples have been produced by carbonization of resorcinol-formaldehyde (RF) gels at 800 • C in the presence of a cationic surfactant, cetyltrimethylammonium bromide (CTAB), used as a pore stabilizer. As the RF nanoparticles are sintered together during carbonization process, a pore stabilizer is added in order to avoid the collapse of the pores, leading to the formation of a well-connected mesopore network [START_REF] Balach | A direct and quantitative three-dimensional reconstruction of the internal structure of disordered mesoporous carbon with tailored pore size[END_REF]. The samples S12 and S14 have been produced using a CTAB to resorcinol molar ratio of, respectively, 0.12 and 0.14, leading to a difference in the morphology of the mesopore network. Based on [START_REF] Balach | A direct and quantitative three-dimensional reconstruction of the internal structure of disordered mesoporous carbon with tailored pore size[END_REF], in this study it is assumed that the micropore morphology is the same for both samples.

Imaging and Segmentation

The two samples were imaged by FIB-SEM Nanotomography [START_REF] Balach | A direct and quantitative image of the internal nanostructure of nonordered porous monolithic carbon using fib nanotomography[END_REF], using the secondary electron signal. The resolution of the SEM image stack of S12 was 3.57nm ×3.57nm and the slicing with the focused ion beam was carried out with 10 nm thickness. Additionally, the tilting of the sample by 52 • leads to a larger resolution in y-direction. This leads to a resulting voxel size in the threedimensional image of 3.57nm ×3.62nm ×10nm. A subwindow of the dimension 411 × 311 × 181 voxels was cropped from the original image stack for the analysis. To improve the segmentation result, the image was magnified by a factor of 1.5 yielding a data set of dimension 616 × 496 × 271 voxels with a voxel size of 2.38nm ×2.41nm ×6.67nm.

The sample S14 was imaged with a lateral resolution of 5nm ×5nm, with a 10nm slicing thickness, leading to a voxel size of 5nm ×6.27nm ×10nm. Analogously to S12, a subwindow of dimension 611 × 293 × 123voxels was cropped from the image stack. Fig. 1 shows two slice views of the cropped FIB-SEM image stack of S12 (on the left) and S14 (right).

After imaging and cropping, the images stacks were segmented using the morphological algorithm presented in [START_REF] Prill | Morphological segmentation of FIB-SEM data of highly porous media[END_REF], which uses the shading effect present in the SEM images, as the result of the large depth of field of the instrument. The segmentation parameters were chosen manually to optimize the visual impression of the segmentation. Slice views through the segmented images are shown in Fig. 2. The morphology of the porous network shows two tortuous components of a bipercolating medium, requiring the choice of adequate random textures models to generate a faithful description.

Finally, the images were scaled to yield an isotropic voxel edge length of 2.38 nm for S12 and 5 nm for S14.

Modeling the Morphology Nanoporous Carbon Electrodes

Since it has been shown, that computations of transport properties based on segmented FIB-SEM images can lead to large errors in the estimated properties resulting from the lack of spatial resolution in the z direction [START_REF] Hutzenlaub | How coarsening the 3d reconstruction of a porous material influences diffusivity and conductivity values[END_REF], it was decided to use a model of random set for the generation of 3D computational domains used as input data for estimation of the transport properties. This allows for correcting errors in the segmentation by, e.g., modeling an isotropic structure, and to explore systematic modifications of the microstructure. Starting from the segmented images, morphological measurements were made to estimate some probabilistic properties of the samples. This is the first step to generate representative domains matching the material morphology as best as possible. The random set model depends on a set of free parameters, which have to be determined by model fitting. Since both samples are fitted to the same model, this results in two sets of parameters representing the respective microstructures of the samples. It opens the possibility to make a systematic study of the impact of the morphology of the mesopores on the transport properties of electrodes.

Principle of Random Set Modeling

A random set A is completely known and identified from a functional, its Choquet capacity, defined on compact sets K in [START_REF] Matheron | Elements pour une théorie des milieux poreux[END_REF], [START_REF] Matheron | Elements pour une théorie des milieux poreux[END_REF][START_REF] Serra | Image Analysis and Mathematical Morphology[END_REF][START_REF] Jeulin | Random texture models for material structures[END_REF] defined by

T (K) = P (K ∩ A = ∅)
where P is a probability. If we denote K x the compact set K after translation to point x, we get

T (K x ) = P (K x ∩ A = ∅) = P {x ∈ A ⊕ Ǩ}
where A ⊕ Ǩ is the result of the dilation of set A by K:

A ⊕ Ǩ = {x, K x ∩ A = ∅} = ∪ y∈K A -y = ∪ x∈A,y∈K {x -y}
Similarly, we can define the erosion by K, A Ǩ, from

A Ǩ = {x, K x ⊂ A} = ∩ y∈K A -y = (A c ⊕ Ǩ) c
where A c is the complementary set of A. In principle, all types of compact sets should be used for a complete characterization of A. In what follows, we will use pair of points (K = {x, x + h}) and approximation of spheres in 3D by rhombo-cuboctaedra.

In addition to the Choquet capacity, size distributions of a random set can be accessed through the use of the opening transformation by convex compact sets (like a sphere with radius r, B(r)). This transformation starts from an erosion of the set by B(r), and is followed by a dilation by B(r):

γ B(r) (A) = (A B(r)) ⊕ B(r)
Alternatively, the size distribution of A c is accessed from the closing transformation of A by convex sets. For spheres B(r) the transformation starts from a dilation of A by B(r) and is followed by an erosion by B(r):

φ B(r) (A) = (A ⊕ B(r)) B(r)
It turns out that opening operations with spheres of increasing radius r progressively removes details of the random set A until its complete suppression. Voxels disappeared for size r correspond to details in A with size less than r, so that a cumulative size distribution can be easily obtained from the estimation of the probability for a point x to belong to γ B(r) (A), obtained by its volume fraction. Similarly, closing operations by spheres of increasing radius progressively fill the space, and removes details of A c , giving access to a cumulative size distribution of A c , estimated from the measurement of the volume fraction of φ B(r) (A). Applied to the segmented binary images of the nanoporous carbon electrodes of this study, both size distributions of the carbon phase and of the mesopores are estimated in 3D.

A typical size of A, namely a median radius r M , is obtained for

P {x ∈ γ B(r M ) (A)}/P {x ∈ A} = 0.5
. Similarly, a typical size of A c , is obtained for

P {x ∈ φ B(r M ) (A)} -P {x ∈ A}/(1 -P {x ∈ A}) = 0.5 .
The same operations and measurements are implemented on the 3D images of materials and of 3D realizations of parametric random sets in an iterative process aiming to minimize some distance between real and virtual specimens, as illustrated later.

Random Set Model

As a starting point for the modeling we use a Boolean model [START_REF] Matheron | Elements pour une théorie des milieux poreux[END_REF][START_REF] Matheron | Elements pour une théorie des milieux poreux[END_REF][START_REF] Serra | Image Analysis and Mathematical Morphology[END_REF][START_REF] Jeulin | Random texture models for material structures[END_REF]. It is built in two steps: a Poisson point process with intensity parameter θ (average number of points per unit volume) generates random germs in the 3D space; independent realizations of a random primary grains A are located on the Poisson points x k . The random set A is obtained by the union of A x k :

A = k A x k (1)
For this model the Choquet capacity is known in closed form. In the 3D space, we have, denoting V the mathematical expectation of the volume:

T (K) = 1 -exp(-θV (A ⊕ Ǩ) (2)
In the present case, we can use a monodisperse Boolean Model of Spheres (BSM) for the mesopores, the primary grain being a sphere of radius r, B(r). It turns out that for a volume fraction of spheres larger than 0.3 and less than 0.95, this model is bipercolating [START_REF] Jeulin | Multi-scale simulation of random spheres aggregatesapplication to nanocomposites[END_REF], which is required for the present materials.

A slice view through a realization of the Boolean model used for the modeling is shown in Fig. 3a.

Yet, clearly the segmented sample images do not resemble the morphology of the Boolean model, as the union of two spheres produces corners, which are not present in the segmented images. Hence, it was decided to modify the BSM by smoothing its boundary. This is done by using the indicator function of the BSM defined as

1 A (x) = 1 if x ∈ A 0 if x ∈ A c . (3)
To generate realizations of the modified model, the indicator function is convolved by a centered Gaussian kernel with width σ and integral 1 over the 3D space,

N (0, σ)(x) = 1 √ 2πσ 2 e -x 2 2σ 2 . ( 4
)
This yields a random function Z(x), with values in the interval [0, 1]: Z(x) = 1 A * N (0, σ).

(5)

A realization of Z(x) is shown on Fig. 3b. To obtain the indicator function of a random set A 1 , a threshold C is applied, as shown on Fig. 3c:

A 1 = {x, Z(x) ≥ C} (6)
The "convolution-thresholding" approach [START_REF] Grzhibovskis | A convolution-thresholding approximation of generalized curvature flows[END_REF] gives a good approximation of the evolution of the boundaries of the initial random set A by motion with a local velocity proportional to the local mean curvature in the present case. It has the effect to smoothen the irregularities of the boundary, as would be the case for some reaction-diffusion process.

In the present case, this sequence yields a random set depending on four parameters θ, r, σ and C.

The intensity θ can vary between zero and infinity but shows very small changes, when the volume fraction is close to zero or to unity. This hampers the optimization involved in the identification of parameters, which is based in a gradient descent. Hence, the parameter θ is replaced by the volume fraction of the BSM p using the relation derived from Eq. 2:

p = 1 -exp(-θ 4 3 πr 3 ). ( 7
)
The parameter p varies only between zero and unity and is hence more stable during the identification. To further simplify the fitting, the threshold C is chosen such that the volume fraction of the set A 1 is equal to the volume fraction of the BSM on each realization. This can be achieved by analyzing the histogram of realizations of the random function Z(x). This eliminates C as a free parameter, the parameter p being directly estimated from the volume fraction of segmented images.

Fitting Free Parameters to Images

To estimate the free model parameters giving the best match to the segmented images, a distance was defined, in order to measure the morphological similarity between model realizations and the segmented images. The distance Y is the sum of the squared difference of three morphological measurements [START_REF] Serra | Image Analysis and Mathematical Morphology[END_REF] performed on the segmented images and on realizations of the model. These morphological measurements are the set covariance, the opening curve and the closing curve.

Set Covariance

The set covariance is given by the probability of two points, separated by a distance h, to be included in the random set A:

Cov(h) = P (x ∈ A, x + h ∈ A) (8)
For a stationary random set with finite scale, the covariance reaches its sill (Cov(0)) 2 for a finite length a, called its range. For the studied materials, the covariance is invariant by rotation of vector h, so that their morphology is isotropic in 3D.

Opening curve

The opening curve Γ (r) is given by the volume fraction V V of the random set after a morphological opening by the spheres B(r)

Γ (r) = V V [γ B(r) (A)] (9)

Closing Curve

Analogously, the closing curve Φ(r) is defined by the volume fraction of the random set after closing by the spheres B(r)

Φ(r) = V V [φ B(r) (A)] (10)
For the identification, these functions were sampled at several distances h i and radii r i yielding the distance Y , defined as (the subscripts d and m stating for experimental data and for simulated model)

Y = h i [Cov d (h i ) -Cov m (h i )] 2 + l i [Γ d (r i ) -Γ m (r i )] 2 + l i [Φ d (r i ) -Φ m (r i )] 2 . (11)
The covariance was sampled on 37 equidistant points, spanning the interval h i = [0, 76] voxels. Each of the other two curves were sampled on 15 points covering the radii [0, 30] voxels. This gives the covariance roughly twice the weight of the other data.

It is assumed that the best parameter fit is reached for a minimum of the distance Y . Since, the measurements on the model are made on realizations, Y is prone to noise, and therefore a stochastic minimization algorithm was used, the SPSA algorithm in [START_REF] Spall | Multivariate stochastic approximation using a simultaneous perturbation gradient approximation[END_REF] in the present study.

The algorithm starts with a given set of parameters

θ 0 = (r 0 , σ 0 ). ( 12
)
Then a random sequence of computed according to the update rule

θ k+1 = θ k -a k g k (θ k ), (13)
with the stochastic gradient

g k (θ k ) =     Y (θ k +c k ∆ k )-Y (θ k -c k ∆ k ) 2c k ∆ k,1 . . . Y (θ k +c k ∆ k )-Y (θ k -c k ∆ k ) 2c k ∆ k,i     . ( 14
)
The sequences (a k ) and (c k ) control the convergence of the sequence, while the sequence (∆ k ) leads to a perturbation in parameter space, with the components of (∆ k ) following a Bernoulli distribution, i.e., P (∆ k,i = ±1) = 0.5.

In the present study, the sequences

a k = a (k + 1) α , a = 0.1, α = 0.101 c k = c (k + 1) γ , c = 0.03, γ = 0.602 (15)
have been used. This, leads to a slow but steady convergence to the minimum of the objective function Y (θ).

Experimental Results

The main 3D morphological properties are summarized in Table 1, clearly showing the differences in the scales of the two mesostructures.

Optimal Model Fit

For the model fitting, the stochastic optimization method presented in [START_REF] Spall | Multivariate stochastic approximation using a simultaneous perturbation gradient approximation[END_REF] was performed with the distance given in Eq. 11. The realizations of the modified BSM described in Sect. 4.2 for the fitting were generated on a voxel grid with dimension 256 × 256 × 256. The stochastic optimization converged to a minimum for each of the samples, leading to the respective best fit parameters. Overall, the model fitting yielded the following results.

Fitting S12

Optimal parameters for a model representing the sample S12 are shown in Table 2.

A comparison between the segmented image and a microstructure realization is shown in Fig. 4. As can be seen, the model reproduces the shape of the pores and the solid phase quite well. A slight difference can be observed in that the boundary of the model is slightly coarser, yet this might be due to noise induced by errors in the segmentation.

To quantify the similarity, the curves for the morphological measurements on the segmentation and on the model realizations are shown in Figs 567. As shown in Fig. 5 the covariance of the model realization is in good agreement with the one measured on the segmented image. Also, the opening curves in Fig. 6 are in good agreement. Yet, a small deviation can be observed in that the opening curve on the model shows a slightly lower probability to find larger details in the solid phase. The closing curves of the model and the segmentation (Fig. 7) virtually coincide, meaning that the size distribution of pores is recovered in the simulations of the model.

Fitting S14

Applying the same fitting procedure to the micrograph of sample S14, results in the parameters given in Tab. 3. A visual comparison of model and segmented image is shown in Fig. 8. Cov(h)

Covariance measured in x-direction The fitted model to S14 shows even better agreement in the measured characteristics. The covariances plotted in Fig. 9 and the opening curves plotted in Fig. 10 coincide almost perfectly, while a small deviation can be observed in the closing curves (Fig. 11), this time the model tending to less larger pore sizes.

Overall, as previously mentioned, the microstructure of S14 is much coarser than the S12s, as the fitted radii of the initial spheres as well as the width of the filter mask is about a factor 2 -3 times larger. Although the volume fraction is larger in S12, the amount of coarsening necessary to reproduce the microstructure for both samples is similar. This is indicated by the fact that the ratio between the initial radii and the width of the filter mask is about σ r = 1.4, for both materials. Finally, the surface area of both samples has been measured on the segmented images, as well as on the modeled microstructures. This is an important feature since the surface area of the mesopores plays an important role in capacitance of the materials. The resulting values are shown in Table 4. As can be seen, the values differ by up to 18% for the sample S14. This discrepancy is attributed to errors in the segmentation, since it is known that the segmentation algorithm does not reconstruct the surface of objects perfectly. 5 Estimating the Effective Transport Properties

Method of Numerical Homogenization

In order to predict the transport properties of the different mesopore structures, the effective conductivity and the effective diffusivity of the samples were estimated by numerical homogenization [START_REF] Papanicolau | Asymptotic Analysis for Periodic Structures[END_REF][START_REF] Sánchez-Palencia | Non-homogeneous media and vibration theory[END_REF], starting from model realizations of the modified BSM described in Sect. 4.2. In this study, the porous material is modeled as a two-phase heterogeneous material as shown in Fig. 12. The main difficulty concerns the contrast between the properties of components, which is infinite in the present case. Therefore, other estimation techniques like bounds derived from variational techniques can only provide an upper bound, the lower bound being equal to zero.

It is assumed that the local material properties are given by a constant conductivity σ 1 in the solid phase and a zero conductivity in the pore space. To estimate the effective conductivity of the material, we have to compute the electric field E(x) deriving from the potential Φ(x) in realizations of the medium, by solving the problem

∇ • [σ(x)∇Φ(x)] =0, ( 16
) ∇ × E(x) =0, (17)
with periodic boundary conditions on the domain boundaries and the condition E(x) = E, E(x) meaning the space average of the electric field. The problem is numerically solved using the "accelerated" Fourier scheme [START_REF] Eyre | A fast numerical scheme for computing the response of composites using grid refinement[END_REF]. It is an efficient method, which directly applies to images and does not require any meshing of the microstructure. The method is based on rewriting Eq. 16 as the implicit integral equation

E(x) = E -Γ (x )P (x -x )dx , ( 18
)
where Γ is the (second-rank, periodic) Green operator associated to the homogeneous conductivity σ 0 [START_REF] Duffy | Green's Functions with Applications[END_REF] and P (x) = J(x) -σ 0 E(x) is the polarization field associated to σ 0 .

In the accelerated FFT method, ( 18) is solved by a Neumann series, which is computed by an iterative fixed-point algorithm. Using explicit formula for the Green operator, the convolution product in the right-hand-side of ( 18) is determined in the Fourier domain. The use of fast Fourier transforms greatly reduces computation times. Additionally, for the Green operator, use is made of the "finitedifference" discretization proposed in [START_REF] Willot | Fourier-based schemes with modified green operator for computing the electrical response of heterogeneous media with accurate local fields[END_REF], which improves the convergence rate and also leads to more precise local fields. In FFT algorithms, the convergence rate generally depends on the reference conductivity σ 0 . In the present method, its optimal value is unknown. Based on numerical experiments, we set σ 0 σ 1 = 0.36 for all computations. We stress that this value is not necessarily optimal.

The effective conductivity σ eff is estimated using:

σ eff E i (x) = σ(x)E i (x) = J i (x) , ( 19
)
where a macroscopic potential gradient ∇φ(x) is applied along the direction e i . The diffusivity is estimated by the same approach, since we consider a steady state for the diffusion, so that the time derivative of the concentration is equal to zero. For this problem, the diffusivity of the solid phase is equal to zero, while the diffusivity of pores, filled with some electrolyte, is set to one. Therefore, the conductivity problem and the diffusivity problems are solved on the same mesostructure, after exchanging the roles of pores and of the solid phase.

Predicting Effective Transport Properties

Bounds of effective properties

As a first estimate for the effective transport properties, the Wiener, Hashin-Shtrikman, and third-order upper bounds are computed for the effective conductivity as well as the effective diffusivity. Since one component has a vanishing transport property, the lower Wiener, Hashin-Shtrikman and third-order bounds vanish.

The upper Wiener bound is valid for any microstructure of a heterogeneous material, while the upper Hashin-Shtrikman bound are valid for isotropic microstructures. Both depend only on the volume fractions of the phases p 1 and p 2 and the bulk properties of the respective phases, σ 1 and σ 2 . The upper Wiener bounds for a two-phase material with vanishing property for σ 2 are given by

σ u W = p 1 σ 1 (20)
For a two-phase material, the upper Hashin-Shtrikman bounds are given by:

σ u HS = σ 1 1 + p 2 1 + p 2 3 . ( 21
)
The third-order bounds [START_REF] Beran | Use of the variational approach to determine bounds for the effective permittivity in random media[END_REF]) make use of the 3-points probability functions P {x ∈ A, x + h 1 ∈ A, x + h 2 } and are generally model-dependent. For practical applications to two-components media with an isotropic geometry, they depend separately on a function calculated by some integral of the 3-points probability and on the property of each component [START_REF] Milton | Bounds on the elastic and transport properties of two-component composites[END_REF]. This probability is not known in a closed form for the present "convolution-thresholding" model, but it can be approximated by the corresponding function for the BSM model:

P {x ∈ A,x + h 1 ∈ A, x + h 2 } = exp (-θV (B(r) ∪ B(r) -h 1 ∪ B(r) -h 2 ))
This expression was used to estimate the function involved in the calculation of the third-order upper and lower bounds σ u BSM and σ l BSM of the BSM [START_REF] Torquato | Physical interpretation of cyclic voltammetry for measuring electric double layer capacitances[END_REF], which depends linearly on the volume fraction of spheres p. Many other models of random media show a similar behavior [START_REF] Jeulin | Random structures in physics[END_REF].

For the upper Wiener, Hashin Shtrikman, and third order bounds for the effective conductivity and diffusivity of the modeled microstructures are shown in Tables 5 and6, respectively.

Effective Conductivity

The effective conductivities were predicted by solving the problem in Eqs. 16-19 with the software morphhom [START_REF] Willot | Fourier-based schemes with modified green operator for computing the electrical response of heterogeneous media with accurate local fields[END_REF]. Hence, the conductivities are computed for the electrostatic case (i.e., at frequency ω = 0). This is justified, by the low frequency the capacitances being measured in the EIS measurement.

For the determination of the effective conductivities, 20 model realizations with different dimensions between 128 3 -512 3 voxels were generated for each of the parameter sets in Tables 2 and3. Then, the effective conductivities were computed for each realization. From the statistical dispersion of conductivities, it is possible to estimate the representativeness of the computed values with respect to variations in the model realizations. A study on the statistical representative volume element (as defined in [START_REF] Kanit | Determination of the size of the representative volume element for random composites: statistical and numerical approach[END_REF]) for S12's model is shown in Fig. 13. To characterize the statistical dispersion of the conductivity, the mean and the empirical standard deviation were computed for the 20 realizations of the model for each edge length of the volume. Fig. 13 shows the mean and twice the standard deviation of the effective conductivity as error bars. Hence, the error bars correspond to the 95% confidence interval of the individual realizations. As can be seen, the realizations with edge length greater than 256 voxel, corresponding to a physical volume with edge length of about l ≈ 0.6µm, show no bias (there is no edge effect for this size) and have a relative precision of less than ±2%.

When combining the 20 realizations with the largest edge length of 512 voxels, the confidence interval of the effective conductivity for the model for S12 is estimated as

σ eff,S12 σ bulk = 0.4976 ± 0.0034 (22)
where the errors indicate the ≈ 95% confidence interval of the average given by the sample mean ±2 SD √ n , where both quantities were estimated form the 20 realizations with edge length 512.

The corresponding study for the sample S14 is shown in Figure 14. Due to the larger structure size, the RVE for S14 is slightly larger, yet for an edge length for about 384 voxel corresponding to a physical volume with edge length ≈ 1, 9µm, the edge effects are negligible and the relative precision on 20 realizations is about ±2%. The average effective conductivity for S14 estimated from a volume of 512 3 voxels is

σ eff,S14 σ bulk = 0.308 ± 0.0044. ( 23
)
With the error again being the ≈ 95% confidence interval of the average.

The estimated conductivities of the two samples are much lower than their upper Wiener and Hashin-Shtrikman bounds (Table 5), which would provide poor estimates. Third-order bounds provide better estimates, but still with an overestimation of the conductivity.

Effective Diffusivity

The other important transport property considered in this work is the effective diffusivity in the pore space. Since the charge carriers in the electrolyte are transported by diffusion in the electrolyte, which in turn is filling the pore space, the effective diffusivity has a major impact on the charging speed of the EDLC electrode. The effective diffusivity is estimated using the same procedure as for the effective conductivity. In the model, it was assumed, that the solid phase is impregnable for charge carriers and hence the diffusivity of the solid phase vanishes. The diffusivity of the pore space was set to one. RVE studies for the effective diffusivity for both samples are shown in Fig. 15 and 16. Both studies show, that for edge lengths of the computational domain of more than 256 Voxels, edge effects essentially vanish and the relative precision reduces to less than ±1%. This corresponds to a representative volume with an edge length of about l ≈ 0.6µm for S12 and l ≈ 1.2µm for S14, respectively.

For the parameter set representing the model for S12, the effective diffusivity was estimated from 20 realization with edge length of 512 Voxels to

D eff,S12 D bulk = 0.097 ± 0.0048. ( 24
)
For the parameter set representing S14 the effective diffusivity was estimated to

D eff,S14 D bulk = 0.24 ± 0.0042. ( 25
)
In both cases, the errors represent the ≈ 95% confidence interval of the average, computed as for the effective conductivity.

Comparison to the bounds in Table 6 shows that the simulated values are closest to the BSM bounds.

Additionally, the diffusivities of the samples have been measured by electroimpedance spectroscopy resulting in the values in Table 7. A comparison shows that the qualitative behavior of the two samples is reproduced, yet quantitatively the level of agreement is low. RVE study for the effective diffusivity for S14

Parametric Study

Finally, a parametric study was performed, to assess the impact of the different model parameters on the effective conductivity and diffusivity of the resulting microstructures. To this end, the realizations of the model were generated with pore volume fractions ranging from 10% to 90%. Due to the linear nature of the conductivity problem, the absolute size of the radii and filter masks do not affect the effective properties. Thus, we limited the study to the estimation of the influence of the ratio between the radii of the spheres of the Boolean model and the size of the filter mask. The resulting effective conductivity and diffusivity are shown in Fig. 17 and18.

As can be seen, the ratio between the initial sphere radius and size of the filter mask has only a minute influence on the effective conductivity. Hence, it should in most cases be sufficient, to fit a Boolean model of spheres to the microstructure and estimate the effective conductivity and diffusivity on its realizations.

Comparison with Experiment

Experimentally, the gravimetric capacitance for both materials were measured using electroimpedance-spectroscopy. The corresponding capacitance values for both materials are:

C G,S12 = 157 F g , (26)
C G,S14 = 113 F g . ( 27
)
To compare the measurements with the calculated values, we calculate the ratio of the volumetric capacitances of the two samples through

C V,S12 = C G,S12 p S12 ρ, (28)
C V,S14 = C G,S14 p S14 ρ, (29)
where p S12 and p S14 are the volume fractions and ρ is the density of the microporous phase of the two samples in g cm 3 . Since, the density of the microporous phase (ρ) is not known, only the ratio of the volumetric capacitances of the two materials can be compared. The ratio of the measured volumetric capacitances is given by

C V,S12 C V,S14 = C G,S12 p S12 ρ C G,S12 p S14 ρ = 1.79, (30)
where the density ρ cancels out. By comparison, the ratio between the capacitances of the modeled structures is p S12 p S14 = 1.288. If the charge was stored only on the micropore surface, these ratios should be the same. Possible reasons for this discrepancy are that a double layer is forming on the mesopore surface, which has not been considered in the simulation. Also, the bulk volumetric capacitances of the solid phase might differ between the two samples, due to the possible presence of micropores with different volume fractions, unresolved by the FIB-SEM images. Since the sample S12 has a larger measured capacitance relative to S14 and also a larger surface area of the mesopores, which would explain a larger discrepancy if a double layer would be formed on the mesopore surface, the first possibility is further explored in this study.

8 Prediction of the Volumetric Capacitance Including a Double Layer on the Mesopore Surface It is assumed, that the total volumetric capacitance is a linear combination of the volumetric capacitance and the areal capacitance of the mesopore double layer multiplied by the specific surface area of the mesopores

C V,S12 =p S12 C V + a S12 C A , C V,S14 =p S14 C V + a S14 C A . (31)
The specific surface areas were measured on the segmented images as well as on the modeled microstructures. The resulting values are given in table 8.

Since the volumetric capacitances are known, only up to a factor ρ, Eqs. 31, 28 and 29 can be solved for the ratio of bulk volumetric capacitance coming from the micropores and volumetric surface capacitance coming from the mesopores. Using the measured surface area of the modeled microstructures, this leads to the following areal capacitance and volumetric capacitances depending on ρ

C V =99.3ρ F g , C A =1.12 * 10 -4 ρ F cm g . (32)
For the sample S12 this yields the following contributions to the capacitance:

C V,S12 = C V p S12 = 66.6ρ, (33)
C A V,S12 = C A a S12 = 38.6ρ. ( 34
)
Hence, in the sample S12 with a total capacitance of C S 12 = 105.9ρ F cm 3 , 36% of the total volumetric capacitance arises from the double layer on the mesopore surface and the remaining 66.6% from the micropores. For the sample S14 the corresponding relations are:

C V,S12 = C V p S12 = 51.7ρ, (35)
C A V,S12 = C A a S12 = 7.09ρ. (36)
Hence, in the sample S14, 12% of the capacitance of C S 14 = 58.76ρ F cm 3 comes from the mesopores and 87.9% comes from the microspores.

In order to estimate the optimal properties of the mesopore morphology, a parameter study was done to estimate the dependence of specific surface area on the model parameters. For the parameter study ten realizations with dimension 256 × 256 × 256 voxels were generated of the model for the same parameter values as in Sect. 6, i.e. with r + σ = 30 voxels. Then the surface area has been rescaled for a model with unit spheres of radius one, i.e. the expected surface area within a unit cube. The resulting specific surface areas are shown in Fig. 19.

For physical structures, the specific surface area for the radius r can be derived by the scaling law:

A(r) = A(r = 1) r (37)
Hence, for a sphere radius of 20nm the maximal specific surface area is

S A ≈ 1 voxel 2 voxel 3
20 * 10 -9 m = 0.05 * 10 9 m 2 m 3 = 5 * 10 5 cm 2 cm 3 (38)

To improve the microstructural properties using the model, we use Eqs. 31, 32, 37 and the computed data in Fig. 19 for a ratio of r σ = 1.5, which is close to the fitted value for both microstructures of r σ ≈ 1.4. We then compute the volumetric capacitance depending on d = r + σ between 10 nm and 50nm and the volume fraction p between 0.1 and 0.9. The resulting capacitances are shown in Fig. 20.

As can be seen, the largest capacitance is reached for the smallest radius, i.e., the smallest structure size of the model. This is due to the scaling in Eq. 37. For larger structures, the specific surface area is small, meaning that the volumetric term is dominant. Hence in the large structure regime, the maximal capacitance is computed for the largest volume fraction p. Yet, for smaller structures, the surface term dominates. This means, that the maximal capacitance is reached for a maximal surface area.

As shown in Fig. 19, the maximal surface area is reached for a volume fraction of p = 0.6 and minimal to no smoothing. For these parameters, the specific surface area peaks at around a(p = 0.6, r σ → ∞) ≈ 1.

Conclusion

In this paper, nanoporous materials for EDLC electrodes are investigated using FIB-SEM Nanotomography imaging, combined with automatic segmentation and random set modeling. It could be shown that transport properties could be predicted combining the structural model with simulation using the FFT method. The results agree qualitatively with experimental measurements. Furthermore, it could be shown that the capacitance values can be predicted when combining the mentioned techniques modeling based on the geometrical characteristics of the model. In this case, the experimental values do coincide very well with the experimental results. Finally, by varying the parameters of the structural model, improvements to the microstructure could be suggested. In conclusion, it was shown that by combining techniques of microimaging, image processing, stochastic modeling and simulations, material properties of EDLC electrodes could be predicted theoretically.

Outlook

The research presented in this paper focuses on the context of energy storage in EDLC's but can similarly be applied to lithium ion batteries or other nanoporous materials. In both cases, batteries as well as EDLC's the electrodes comprise mostly of a micro-or nanoporous carbon structure filled with an electrolyte. Furthermore, the techniques demonstrated here can be amended and improved individually. As the FIB technology improves, the image quality of finer structures like the ones investigated here can be improved, yielding a better segmentation. Also, since automatic segmentation for FIB-SEM images is a relatively new field, improved algorithms can improve the accuracy of the segmentation. Concerning the structural modeling, replacing the Boolean model, which serves as a basis for the modified model by a cherry-pit model, can improve the modeling of the sphere packing. Finally, improved physical simulations as mentioned in the introduction can be used to further enhance the prediction of the capacitance.

Fig. 1

 1 Fig. 1 Slice views of the FIB-SEM image stacks for the samples S12 (left) and S14 (right).

Fig. 2

 2 Fig. 2 Slice views of the segmentations of the images of samples S12 (left) and S14 (right).

Fig. 3

 3 Fig. 3 Realizations of a Boolean model of spheres (a), a smoothed indicator function Z(x) of the BSM (b) and a random set A 1 obtained by thresholding Z(x) (c).

Fig. 4

 4 Fig. 4 Comparison between the observed structure (left) and the fitted model (right) of dataset S12.

Fig. 5 Fig. 6

 56 Fig.5Comparison between the covariance of the fitted model (red) and the observed structure (green) of dataset S12.

Fig. 7

 7 Fig. 7 Comparison between the closing curve of the fitted model (red) and the observed structure (green) of dataset S12.

Fig. 8

 8 Fig. 8 Comparison between the observed structure (left) and the fitted model (right) of dataset S14.

Fig. 9 Fig. 10

 910 Fig.9Comparison between the covariance of the fitted model (red) and the observed structure (green) of dataset S14.

Fig. 11

 11 Fig. 11 Comparison between the closing curve of the fitted model (red) and the observed structure (green) of dataset S14.

Fig. 12 A

 12 Fig.12A sample domain on which the effective conductivity and diffusivity are computed.

Fig. 14

 14 Fig. 13 Study on the RVE of the conductivity the model with parameters for S12. The error bars corresponding to the 95% confidence interval for the individual realizations.

  Fig. 15 Study on the RVE of the effective diffusivity of the model with parameters for S12. The error bars corresponding to the 95% confidence interval for the individual realizations.

Fig. 16

 16 Fig. 16 Study on the RVE of the effective diffusivity of the model with parameters for S14. The error bars corresponding to the 95% confidence interval for the individual realizations.

Fig. 17

 17 Fig.17Parameter study of the effective conductivity of the model described in Section 4.2.

Fig. 18

 18 Fig.18Parameter study of the effective diffusivity of the model described in Section 4.2.

Fig. 19

 19 Fig.19Parameter study of the effective surface area of the model described in Section 4.2 with spheres of radius one.

Fig. 20

 20 Fig.20Parameter study of the total capacity of a structure of the model described in Sect. 4.2, with spheres of radius d, the volumetric and areal capacities from Eq. 32 and r σ = 1.5.

Table 1

 1 Basic morphological properties of the two samples: solid volume fraction (p); range of covariance (rc); median radius solid (med(rs)); median pore radius (med(rp))

	Sample	p	rc	med(rs)	med(rp)
	S12	67% 47.6 nm	23.8 nm	14.28 nm
	S14	52%	150 nm	50 nm	60 nm

Table 2

 2 Optimal parameters obtained by fitting to the observed structure of S12 Parameter Optimal value (voxel) Optimal value (nm)

	p	0.67	0.67
	r	7	16.7
	σ	10	23.8

Table 3

 3 Optimal parameters obtained by fitting to the observed structure of S14

	Parameter Optimal value (voxel) Optimal value (nm)
	p	0.52	0.52
	r	10	50
	σ	14

Table 4

 4 Specific surface area measurement on the segmented images and the modeled structures.

	Sample aseg [ cm 2 cm 3 ] a mod [ cm 2 cm 3 ]
	S12	3.1 × 10 5	3.28 × 10 5
	S14	7.15 × 10 4	6.02 × 10 4

Table 5

 5 Upper Wiener, Hashin-Shtrikman and third-order bounds of the conductivity of the two samples

	sample	σ u W	σ u HS	σ u BSM
	S12	0.67 0.629	0. 557
	S14	0.52 0.428 0.3914
	Table 6 Upper Wiener, Hashin-Shtrikman and third-order bounds of the diffusivity of the
	two samples
	sample	D u W	D u HS	D u BSM
	S12	0.33 0.137	0.117
	S14	0.48	0.37	0.244

Table 7

 7 Measured and simulated diffusivities of the two samples, assuming a bulk diffusivity of D bulk = 1.08 × 10 -6 cm 2 s .

	Sample	D Exp	D Sim
	S12 S14	1.07 × 10 -6 cm 2 s 4.89 × 10 -6 cm 2	0.097 cm 2 s

Table 8

 8 Specific surface of the samples measured on the segmented images aseg as well as on the modeled microstructures a mod .

	Sample	aseg	a mod
	S12 S14	3.1 × 10 5 cm 2 cm 3 7.1 × 10 4 cm 2 cm 3	3.18 × 10 5 cm 2 cm 3 6.02 × 10 4 cm 2 cm 3

Acknowledgements This work was partly funded by the German Federal Ministry of Education and Research, project AMiNa (03MS603D). The authors are also grateful for the support of DFDK-CDFA (grant CDFA-02-11).


Sébastien Salva
email: sebastien.salva@uca.fr

Loukmen Regainia
email: loukmen.regainia@uca.fr

Using Data Integration for Security Testing

Keywords: Security, Security patterns, Attack Defense Trees, Test case generation

The explosion of digitisation makes a plethora of security data publicly available for developers. These numerous (often complex) documents expose them to the difficulty of choosing the most appropriate solution for securing their applications. We propose in this paper a method based upon data acquisition and integration, which assists developers in the Threat modelling stage and in the security test case execution. The method firstly helps devise Attack Defense Trees by means of a data-store. These trees show attacks, steps and defenses given under the form of security patterns, which are re-usable solutions to design more secure applications. These trees are then used for the test case generation. The data-store integrates test case stubs, which make this generation easier and developers more efficient. We evaluate our approach on 24 participants and show encouraging results on the use of data integration in software engineering.

Introduction

Since a decade, it is well admitted that software security is essential and has to be considered through all the software life cycle. Many developers, researchers and organisations have hence made security their hobby-horse and brought several improvements with the proposal of numerous digitised security bases and documents. These take security into consideration at different stages of the software life cycle and are presented with different viewpoints, abstraction levels or contexts. This plethora of diverse documents makes difficult the choices of security solutions and afterwards their validations. Indeed, developers cannot be experts in any field and they clearly lack guidance for conceiving and implementing both secure software and tests.

This work focuses on this issue and studies the possibility of using publicly available security resources for helping developers devise more secure applications. We propose a method, which aims at assisting developers in the Threat modelling stage and in the test case generation and execution. More precisely, the contributions of this paper are highlighted in the following points:

we present a security data acquisition and integration method, which extracts data from various Web and publicly accessible sources to conceive a data-store storing relationships among attacks, security principles, security patterns and test case parts written with the Given When Then (GWT for short) template. The security pattern intuitively relates countermeasures to threats and attacks in a given context [START_REF] Schumacher | Security Engineering with Patterns: Origins, Theoretical Models, and New Applications[END_REF]. These security re-usable solutions often are presented textually or with UML schema and are characterised by a set of structural and behavioural properties; from the data-store, our method helps in the generation of Attack Defense Trees (ADTrees [START_REF] Kordy | Attack-defense trees[END_REF]) showing the attacker possibilities to compromise an application and the defenses that may be put in place to prevent attacks.

We have chosen the ADTree model because it offers the advantage of being easy to understand even for novices in security. These ADTrees are composed of defenses given under the form of security pattern combinations; -ADTrees serve here to the test case generation. These test cases aim to check whether an application is vulnerable against the attacks exposed in an ADTree and whether security pattern consequences are observed in the application behaviour. Pattern consequences are observable events resulting from the good contextualisation and implementation of the pattern in an application. From an ADTree, our method extracts attack scenarios, GWT test case stubs and related procedures composed of comments or blocks of code, which aim at guiding the developer in the test case completion. As ADTrees can be formalised with formal expressions, called ADTerms, we strictly define the test case generation and execution.

Besides, we concentrated our attention on quality criteria and on education while the design of this method. ADTrees are indeed constructed with concrete data extracted from the CAPEC base3 . ADTrees also express security pattern combinations with respect to several criteria, i.e., Unambiguity, Navigability and Comprehensibility, which are quality criteria proposed in [START_REF] Alvi | A Comparative Study of Software Security Pattern Classifications[END_REF], the last two respectively related to: the ability to direct a software designer among collaborative and related patterns; the ease to understand patterns by both a novice and expert developer. Test case stubs are also structured to ease Readability and Re-usability and to try to increase Effectiveness.

We have generated a data-store specialised to the context of Web applications (Web sites), which is composed of 215 CAPEC attacks, 26 security patterns and 669 test case parts. We employed this data-store to evaluate on 24 human subjects the benefits of using the notion of data acquisition and integration in the software life cycle. This evaluation shows encouraging results about Comprehensibility and Effectiveness.

The remainder of the paper is organised as follows: Section 2 presents some related work. The data integration step is shortly presented in Section 3. The next section shows how ADTrees are generated by means of the data-store. Section 5 addresses the test case generation and execution. We present our evaluation in Section 6 and finally conclude in Section 7.

A plethora of papers deals with model-based security testing. Due to lack of room, we only present some of them related to our work, which consider models not to describe the implementation behaviour but rather to express the attacker's goals or the vulnerability causes of the system [START_REF] Morais | Security protocol testing using attack trees[END_REF][START_REF] Marback | Security test generation using threat trees[END_REF][START_REF] Marback | A threat model-based approach to security testing[END_REF][START_REF] Shahmehri | An advanced approach for modeling and detecting software vulnerabilities[END_REF]. Some authors focused on trees (Attack trees, Security Activity Graphs, etc.), which express the treats or attacks or vulnerability causes that should be prevented in systems. From these models, test cases are then written to check whether attacks can be successfully executed. Morais et al. introduced a security testing approach specialised for protocols [START_REF] Morais | Security protocol testing using attack trees[END_REF]. Attack scenarios are extracted from an Attack tree and are converted to Attack patterns and UML specifications. From these, attack scripts are manually written and are completed with the injection of (network) faults. In [START_REF] Marback | Security test generation using threat trees[END_REF], data flow diagrams are converted into Attack trees from which sequences are extracted. These sequences are composed of events combined with parameters related to regular expressions allowing the generation of concrete values. These events are then replaced with blocks of code to produce test cases. In [START_REF] Marback | A threat model-based approach to security testing[END_REF], test cases are generated from Threat trees. The latter are completed with parameters associated to regular expressions. Security scenarios are extracted from the Threat trees and are manually converted to executable test scripts. Shahmehri et al. proposed a passive testing approach to detect vulnerabilities [START_REF] Shahmehri | An advanced approach for modeling and detecting software vulnerabilities[END_REF]. The undesired vulnerabilities are modelled with models called SGMs, which are specialised DAGs showing security goals, vulnerabilities and eventually mitigations. Detection conditions are then semi-automatically extracted and given to a monitoring tool, which returns test verdicts.

These works take Threat models as inputs, which are manually written. If these lack of details (parameters, attack steps, etc.), the final test cases will be too abstract as well. Furthermore, these methods do not give any recommendation on how to write tests and on how to structure them. Hence, developers lack guidance to write tests and to reuse them. This paper proposes a method, which semiautomatically infer Attack Defense Trees, composed of attacks steps, techniques and defenses. To ease the understanding and readability of the generated test case stubs, we associate in our data-store every attack and defense step to some test case sections, which are classified w.r.t. an application context and to an attack step or security pattern.

Few works tackled the testing of security patterns, which is another topic of this paper. Yoshizawa et al. introduced a method for testing whether behavioural and structural properties of patterns can be observed in the traces of instrumented applications [START_REF] Yoshizawa | Verifying implementation of security design patterns using a test template[END_REF]. Given a security pattern, two test templates (OCL expressions) are written, one to specify the pattern structure and another to describe its behaviour. Then, developers have to make templates concrete by writing Selenium scripts for experimenting the application. The latter returns traces on which the OCL expressions are verified. In contrast to the previous paper, our inferred ADTrees firstly help developers choose for every attack, the combinations of patterns that can be used as countermeasures. Then, our testing approach aims at testing the security pattern consequences. We do not check the structure of the patterns. Hence, our approach is complementary to the previous one.

We also proposed a semi-automatic data integration method in [START_REF] Regainia | A methodology of security pattern classification and of attack-defense tree generation[END_REF], in order to extract a pattern classification. We took inspiration from this paper to infer pattern combinations. In contrast, the notion of data-store, its architecture, the considered security properties and the test case generation and execution are new contributions.

3 Data Integration Figure 1(a) exposes the meta-model of the first part of the data-store used to integrate relationships among attacks, attack steps, techniques, security principles and security patterns. The entities of Figure 1(a) refer to these security properties. To increase the precision of the relations, we chose to decompose attacks into sub-attacks, and into attack steps. These steps are associated to countermeasures, allowing to prevent or counter attack steps. We also decompose security patterns into strong points, which are sub-properties expressing pattern key design features. Relying on a hierarchical organisation of security principles, the method maps countermeasure clusters to principles and strong points to principles. As countermeasures usually are detailed properties, we gather them into clusters (groups) to reach about the same abstraction levels as those of the security principles.

The meta-model of Figure 1(a) is extended with new entities and relations, which are required for the testing process. This extension is depicted in Figure 1(b). An attack step is also associated to a Test architecture and to one Application context. The context refers to an application family, e.g., Android applications, or Web sites. The "Test Architecture" entity refers to textual paragraphs explaining the points of observation and control, testers or tools required to execute the attack step on an application, which belongs to an Application context. Next, we map attack steps onto GWT test case sections. For readability and re-usability purposes, we chose to consider the "Given When Then" pattern to break up test cases into several parts:

the Given section aims at putting an application under test in a known state; the When section triggers some actions; the Then section is used to check whether the conditions of success of the test case are meet (assertions). In our context, we suppose that a Then section returns "P ass st " if an attack step st has been successfully executed on an implementation and "F ail st " otherwise.

Likewise, we map security pattern consequences onto Then sections to check whether the consequences of the pattern can be observed in the application traces. We assume that a Then section returns "F ail sp " if a consequence of the security pattern sp is not observed from the implementation. For instance, if an application is conceived with the "Input Guard" pattern, then unexpected inputs should bring the application to a quiescent state (no output) or outputs reflecting errors should be observed.

Each test case section is linked to one procedure stored in the Procedure table of the data-store, which implements the section. A Given, When or Then section can be reused with several attack steps or security patterns. With regard to the meta-model given in Figure 1(b), a GWT test case section (and procedure) is classified according to one application context and one attack step or pattern consequence.

In some specific application contexts, the procedures can be completed with comments or with blocks of code to ease the test case development. When the procedure content can be reused with any application in a precise context, we call it Generic procedure: Definition 1 (Generic procedure). Let C be an Application context. A Generic procedure is a block of code, related to a Given, When or Then test case section, that can be used with any application of C;

The data-store must only contain Generic procedures related to an application context. It worth mentioning that an empty procedure is generic.

Security data acquisition and integration

This section summarises a data integration example for the Web application (Web sites) context. We chose to focus on the CAPEC base to extract information about security attacks. The CAPEC base offers an open and rich catalogue of attacks in a comprehensive schema. We conceived a tool for data acquisition and extraction, based on text mining and on the ELT (Extraction, Load, Transform) tool Talend 4 . With it, we automatically scanned all the CAPEC base (Version 2.8) and collected 215 attacks, 209 steps, 448 techniques and 217 countermeasures, knowing that attacks can share steps, attack techniques and countermeasures. Among these, we gathered 75 attacks and 142 attack steps specialised for the context of Web sites.

As security patterns are described in an abstract manner with texts, we manually collected security patterns, their strong points and consequences from the catalogue given in [12]. We gathered 26 security patterns, 43 consequences and 36 strong points. We also integrated the inter-pattern relations given in [START_REF] Yskout | Do security patterns really help designers?[END_REF]. We organised 66 security principles found in the literature into a hierarchy composed of four levels, from the most abstract to the most concrete principles.

The data integration of the GWT test case sections was automatically performed. For a given attack, the CAPEC base provides two text sections called "Attack Prerequisites" and "Resources Required". We automatically scanned these paragraphs and completed 209 procedures including comments composed of the two previous paragraphs. Each procedure is associated to one Given test case section (one section for each attack step). For every step st, we added one When test case section and one procedure composed of comments listing the techniques related to st. Still in the CAPEC documents, the paragraphs "Indicators" and "Outcomes" provide some directives and conditions on an attack step realisation. In the same way as previously, we automatically scanned these paragraphs and, for every attack step, we completed the data-store with one Then section associated to one procedure, itself composed of the two previous paragraphs of the CAPEC base given as comments. In this way, we generated 627 GWT test case sections. For every security pattern consequence found in the data-store, a Then test case section and its related procedure are also automatically inserted. The procedure is composed of comments listing the consequence, which have to be observed from the application traces. In the context of Web applications, we observed that several procedures can be completed with blocks of code calling penetration testing tools. We completed 32 procedures, which cover 43 attack steps. We used the tools Selenium and ZAProxy5 , which is a penetration testing tool covering various Web vulnerabilities.

This data-store is available here6 .

Threat Modelling

Before presenting how our method assists developers in threat modelling, we recall some notions about the ADTree model.

Attack Defense Trees (ADTrees)

ADTrees are graphical representations of possible measures an attacker might take in order to attack a system and the defenses that a defender can employ to protect the system [START_REF] Kordy | Attack-defense trees[END_REF]. As illustrated in Figures 3(a) and 3(b), ADTrees have two different kinds of nodes: attack nodes (red circles) and defense nodes (green squares). A node can be refined with child nodes and can have one child of the opposite type (linked with a dashed line). Node refinements can be disjunctive (like in Figure 3(a)) or conjunctive. The former is recognisable by edges going from a node to its children. The latter is graphically distinguishable by connecting these edges with an arc. We extend these two refinements with the sequential conjunctive refinement of attack nodes, defined by the same authors in [START_REF] Jhawar | Attack trees with sequential conjunction[END_REF]. This operator expresses the execution order of child attack nodes. Graphically, a sequential conjunctive refinement is depicted by connecting the edges, going from a node to its children, with an arrow. For instance, the node "Attack Step" in Figure 3(b) is refined with a sequence of others steps. Alternatively, an ADTree T can be formulated with an algebraic expression called ADTerm and denoted ι(T ).

In short, the ADTerm syntax is composed of operators having types given as exponents in {o, p} with o modelling an opponent and p a proponent. ∨ s , ∧ s , -→ ∧ s , with s ∈ {o, p} respectively stand for the disjunctive refinement, the conjunctive refinement and the sequential conjunctive refinement of a node. A last operator c expresses counteractions (dashed lines in the graphical tree). The first stage of our method takes place in the Threat modelling phase of the software life cycle, which occurs while the requirement analysis. Threat modelling is a process consisting in identifying and describing the attacker's goals and capabilities, as well as identifying the potential threats of an application. Different methods can be followed, e.g., DREAD, or STRIDE. Our method starts to semi-automatically generate an ADTree by means of the data-store. The ADTree generation is illustrated in the the fourth first steps of Figure 2. We present them below.

Attack Defense Tree generation

Step 1: Initial ADTree design

The developer initially establishes a first ADTree T whose root node represents the attacker's goal, which may be refined with child nodes. The ADTree T describes attack combinations, which can be applied on the application. We here assume that T at least has leaves labelled by attacks kept in the data-store. Otherwise, a semantic alignment may be required to replace some labels by similar attack names.

Figure 3(a) depicts an ADTree example: the goal, given in the root node, is to inject malicious code into an application. This node is disjunctively refined with two children expressing two more concrete attacks, CAPEC-66: SQL Injection and CAPEC-244: Cross-Site Scripting via Encoded URI Schemes. Step 2: ADTree generation Usually, yhe ADTree T does not include enough details on how the attack is sequenced and on the defenses expressing how the application can be protected. Implementing a secure application and deriving test cases from this kind of tree remains a tedious task. The data-store can be used to augment T . For every node labelled with an attack Att, we automatically generate an ADTree denoted T (Att). The architecture of the data-store leads to the generation of ADTrees having the general form illustrated in Figure 3(b). The root of an ADTree T (Att) is labelled by Att. This node may have children expressing more concrete attacks and so forth. The most concrete attacks have step sequences (edges connected with an arrow). These steps are connected to techniques with a disjunctive refinement. The lowest attack steps in the ADTree are also linked to defense nodes, which may be the roots of sub-trees expressing security pattern combinations whose purpose is to counteract the attack step.

We implemented the ADTree generation with a tool, which takes attacks of the data-store and yields XML files. These can be edited with the tool ADTool [START_REF] Kordy | Attack-defense trees[END_REF]. For instance, Figure 4 depicts the ADTree of the attack CAPEC-66, which was exported from ADTool. Each lowest attack step has a defense node express-ing pattern combinations. Step 2.1, which identifies the possibilities to inject malicious code through the application inputs, requires more patterns than the other steps to filter these inputs. Some of them have relations: for instance "Application Firewall" can be replaced by "Input Guard" with "Output Guard". Step 3: Security pattern choice

The developer can now edit the ADTrees T (Att) to keep or remove attack steps w.r.t. the application context. He or she also has to choose the security patterns that have to be contextualised and implemented in the application. After this step, we assume that a defense node either is labelled by a security pattern (it does not have children) or has a conjunctive refinement of nodes labelled by security patterns. The lowest nodes labelled by attack steps, must be linked to a defense node.

As a result of the steps 2 and 3, the generated ADTrees have specific forms and are expressed with specific ADTerms: Proposition 1. An ADTree T (Att) achieved by the previous steps has an ADTerm ι(T (Att)) having one of these forms:

1. ∨ p (t 1 , . . . , t n ) with t i (1 ≤ i ≤ n) an ADTerm also having one of these forms: 2.

-→ ∧ p (t 1 , . . . , t n ) with t i (1 ≤ i ≤ n) an ADTerm having the form given in 2) or 3); 3. c p (st, sp), with st an ADTerm expressing an attack step and sp an ADTerm modelling a security pattern combination.

The first ADTerm expresses children nodes labelled by more concrete attacks.

The second one represents sequences of attack steps. The last expression is composed of an attack step st refined with techniques, which can be counteracted by a security pattern combination sp. In the remainder of the paper, we call the last expression Basic Attack Defence Step, shortened as BADStep. These shall be particularly useful to build GWT test case stubs: Step 4: Final ADTree generation In the initial ADTree T , each attack node labelled by Att is now automatically replaced with the ADTree T (Att). This can be done by substituting every term Att in the ADTerm ι(T ) by ι(T (Att)). We denote T f the resulting ADTree. It depicts a logical breakdown of the various options available to an adversary and the defences, materialised with security patterns, which have to be inserted into the application model.

In this step, we also extract from the data-store a description of the test architecture required to run the attacks on the application under test and to observe its reactions.

Test Suite Generation

The semantics of an ADTree can be defined in terms of attack-defense scenarios. In general terms, a scenario is a minimal combination of events leading to the root attack, minimal in the sense that, if any event is omitted from the attack scenario, then the root goal will not be achieved. The semantics of an Adtree T f , i.e. its scenario set, can be extracted from its ADTerm ι(T f ): Definition 3 (Attack scenarios). Let T f be an ADTree and ι(T f ) be its ADTerm. The set of Attack scenarios of T f , denoted SC(T f ) is the set of clauses of the disjunctive normal form of ι(T f ) over BADStep(T f ). An attack scenario s of SC(T f ) is an ADTerm over BADSteps. BADStep(s) denotes the set of BADSteps of s. We also denote SP (s) the security pattern set found in s: SP (s) = {sp | ∃b ∈ BADStep(s) : sp ∈ defense(b)}. By extension, BADStep(T f ) stands for the set of BADSteps found in ι(T f ); SP (T f ) is the security pattern set of ι(T f ), found in all its scenarios.

Step 5: Test suite generation

Let us consider a security scenario s ∈ SC(T f ). Given a BADStep b = c p (st, sp) ∈ BADStep(s), we generate the GWT test case T C(b), which aims at checking whether the application under test I is vulnerable to the attack step st and whether the consequences of the security patterns in defense(b) can be observed from I. T C(b) is assembled from the data-store by means of the following steps:

1. the data-store provides for st, with the relations testG, testW and testT , one Given section, one When section and one Then section, each related to one procedure. The Then section aims to assert whether the application is vulnerable to the attack step st; 2. the data-store provides from the security pattern set defense(b) a set of other Then sections, each related to procedures. These Then sections aim to check whether the security pattern consequences can be observed from the application behaviours; 3. all these sections are assembled to make up the GWT test case stub T C(b).

By applying these steps on all the scenarios of SC(T f ), we obtain the test suite T S with

T S = {T C(b) | b = c p (st, sp) ∈ BADStep(s) and s ∈ SC(T f )}.
We implemented these steps to yield GWT test case stubs compatible with the Cucumber framework7 , which supports a large number of languages. Figure 5 gives a test case stub example obtained with our tool from the first step of the attack CAPEC-66 depicted in Figure 4. The test case lists the Given When Then sections in a readable manner. Every section is associated to a Generic procedure stored into another file. The procedure related to the When section is given in Figure 6. The comments comes from the data-store and are extracted from the CAPEC base. This procedure includes a generic block of code; the "getSpider()" method relates to the call of the ZAProxy tool, which crawls a Web application to get its URLs. In this example, it only remains for the developer to complete the initial URL before testing whether the application can be explored.

Feature: CAPEC-66: SQL Injection #1. Explore Scenario: Step1.1 Survey application #The attacker first takes an inventory of the functionality exposed by the application. Given a new scanning session When spider the application Then the application is spidered #assertions for security pattern testing Then Output Guard security pattern is present Then Input Guard security pattern is present Fig. 5. A GWT test case example

Step 6: Test case stub completion Now, the developer has to complete the previous GWT test case stubs. We believe that the decomposition of the test case and its link to the ADTree T f (associations among steps, security patterns and procedures) make this step easier. In addition, the Generic procedures, composed of comments or blocks of code @When("^spider the application") public void theApplicationIsSpidered() { // Try one of the following techniques : // 1. Spider web sites for all available links // 2. Sniff network communications with application using a utility such as WireShark. getSpider().setMaxDepth [START_REF] Yoshizawa | Verifying implementation of security design patterns using a test template[END_REF]; url = "URL to be scanned"; try { spider(url); } catch (InterruptedException e) {e.printStackTrace();} waitForSpiderToComplete();} Fig. 6. The procedure related to the When section of Figure 5 should make him or her more effective in the test case writing.

Step 7: Test suite execution

Once the GWT test case stubs are completed, these can be executed on the application under test I. The test architecture allowing the experimentation of I is described in the report provided by Step 4.

After the execution of one test case T C(b) on I, denoted T C(b)||I, we obtain sets of verdict messages of the form "P ass st ", "F ail st " or "F ail sp ", resulting from its assertions (see Section 3. -{F ail st } (resp. {P ass st }) means I is (resp. does not appear to be) vulnerable to the attack step st and that the consequences of the security patterns are observed; -{P ass st , F ail sp1 , ..., F ail sp k } means I does not appear to be vulnerable to the attack step st but some consequences of the security patterns sp 1 , . . . , sp k are not observed; -{F ail st , F ail sp1 , ..., F ail sp k } means I is vulnerable to the attack step st and that some consequences of the security patterns sp 1 , . . . , sp k are not observed.

Evaluation

We empirically studied two scenarios on 24 participants to assess whether developers can take profit of our approach. The duration of each scenario was set at most to one hour and half. The participants are third to fourth year computer science undergraduate students, having good skills in the development and test of Web applications.

The participants were given the task of choosing security pattern combinations to prevent two attacks, CAPEC 244: Cross-Site Scripting via Encoded URI Schemes, and CAPEC 66: SQL Injection, on two deliberately vulnerable Web sites, RopeyTasks and The Bodgeit Store. We also asked the participants to write test cases with the tool Selenium in order to: show that both Web sites are vulnerable to the two attacks, show that the application behaviours do not include at least one consequence of the security pattern "Input Guard" and at least one consequence of "Output Guard".

In the first scenario, denoted Part 1, we supplied the CAPEC base, two concrete attack examples, the detailed steps showing how to manually perform them along with the expected outcomes and the security pattern catalogue given in [12]. In the second scenario, denoted Part 2, we also supplied the ADTrees of the two attacks (Figure 4 is one of them) along with the generated GWT test case stubs for each attack step. At the end of each scenario, the students were invited to fill in a form listing ten questions. Due to lack of room, we only present the questions and results concerning the test case generation:

-Q7: Was it easy to write test cases? -Q8: How long did you take for writing test cases? -Q9: How confident are you about your test cases? -Q10: Provide your test cases (or suites). These questions was asked in order to evaluate the following criteria:

-C1: Comprehensibility: does our method ease the test case development? -C2: Effectiveness: can the test cases detect defects? -C3: Efficiency: does our method help reduce the time needed to write tests?

Experiment results

We extracted the following results from the forms returned by the participants (available here8 ). We collected the answers of Question Q7, proposing this fourvalued scale: easy, fairly easy, difficult, very difficult. We finally analysed the test cases given by the participants and evaluated their correctness with regard to four aspects: 1&2: detection (with at least one test case) that both applications are vulnerable to the attacks CAPEC 66 and CAPEC 244; 3&4: detection that the application behaviours do not include the consequences of the patterns "Input Guard" and "Output Guard". As we considered this last aspect as difficult for students, we expected at least one Then test case section for every pattern. Figure 8 presents the number of participants who meet these aspects.

Result interpretation

C1 Comprehensibility: we chose this criteria to evaluate whether our method makes testing easier. Figure 7 shows that one quarter of the students found easier the test case writing with our test case stubs. After discussion, it turned out that the test case structure with the GWT template made test cases more readable and that the links between test case sections and Attack steps helped students understand what to develop. In the meantime, Question Q9 reveals that the confidence level of the participants about their test cases increases by 20,83%.

C2 Effectiveness: Figure 8 depicts the results about the test case correctness. The columns "SQLi" and "XSS" provide the number of test cases allowing to reveal that the attacks can be successfully executed on the applications. In Part 1, few participants developed complete test cases despite the detailed steps we provided (assertions were missing or incorrect in most of the test cases). The number of correct test cases strongly increases in Part 2 thanks to the comments the participants found in the procedures. The columns "Input Guard" and "Output Guard" give the number of Then sections (and procedures) allowing to show that the consequences of these security patterns are not observed from the application behaviours. This task was much more difficult for the students as they are not yet expert in security patterns. Hence, it is not surprising to see that only one student was able to write at least an assertion showing that the Input Guard consequences are not present. The number of correct Then sections rises to 14 (58,3%) in Part 2. With the pattern "Output Guard", the number of correct Then sections rises from 0 to 23 in Part 2. We can conclude that the test case correctness strongly increases with our approach.

C3 Efficiency: on average, the participants took 46 minutes for writing test cases from scratch and 60 minutes with the use of our method. The additional time spent in Part 2 can be explained when we alongside focus on Effectiveness and Comprehensibility. Indeed, after discussion with the participants, we deduced that they took more time to follow and analyse the ADTrees, to read the comments in procedures, etc. As a result, almost all the test cases are correct in Part 2 (more assertions, etc.).

Conclusion

We have presented a method taking advantage of data integration for guiding developers devise more secure applications from the Threat modelling stage to the testing one. The method generates ADTrees and test case stubs allowing to check whether an application is vulnerable to attacks and whether security pattern consequences are observable from the application behaviour. We conducted an evaluation of the method, which shows it makes the participants more effective on security testing. But, several issues remain open. For instance, our method does not take into consideration the size of the ADTrees in the Threat modelling stage. This is a strong limitation since large trees are usually unreadable, which contradicts the method purposes. The ADTree reduction could be a first solution on this problem. But, the literature does not yet provide a generic method for this kind of reduction. Besides the tree structure, the node meaning must be taken into account in the node aggregating process, which must preserve the semantics of the ADTree.

3. 1

 1 Fig. 1.

Fig. 2 .

 2 Fig. 2. Threat modelling and security testing steps

  Fig. 3.

Fig. 4 .

 4 Fig. 4. ADTree of the Attack CAPEC-66

Definition 2 (

 2 Basic Attack Defence Step (BADStep)). A BADStep b is an ADTerm of the form c p (st, sp), where st is an ADTerm modelling an attack step and sp an ADTerm of the form sp 1 or ∧ o (sp 1 , . . . , sp m ) modelling a security pattern conjunction. defense(b) = {sp 1 | b = c p (st, sp 1 )}∪{sp 1 , . . . , sp m | b = c p (st, ∧ o (sp 1 , . . . , sp m ))}

1 ) : Definition 4 (

 14 Test verdict sets). Let I be an application under test, b = c p (st, sp) ∈ BADStep(T f ) with defense(b) = {sp 1 , . . . , sp m }(m > 0) and T C(b) ∈ T S be a test case. The execution of T C(b) on I leads to a verdict set, denoted Verdict(T C(b)||I), which can be:

From these verdict

  messages, we define two first relations. The first relation vulnerable defines that an application I is vulnerable to a BADStep b if the message F ail st belongs to the verdict set Verdict(T C(b)||I). The relation unsat c defines that I does not satisfies the consequences of the pattern sp if the message F ail sp belongs to Verdict(T C(b)||I): Definition 5 (Test case verdicts). Let I be an implementation, T f be an ADTree, b = c p (st, sp) ∈ BADStep(T f ) and T C(b) ∈ T S be a test case. 1. I vulnerable b = true if {F ail st } ∈ Verdict(T C(b)||I), f alse otherwise; 2. I unsat c sp if sp ∈ defense(b) and F ail sp ∈ Verdict(T C(b)||I). We now define the relation effective on a scenario s ∈ SC(T f ), composed of the BADSteps b 1 , . . . , b n and on I to formally state whether s detects vulnerabilities on I. The relation effective is evaluated by substituting every BADStep term b i with the evaluation of I vulnerable b i . These relations help define and evaluate the final implementation relations. Definition 6 (Implementation relations). Let I be an implementation, T f be an ADTree, and s ∈ SC(T f ), with BADStep(s) = {b 1 , . . . , b n }. 1. σ : BADStep(s) → {true, f alse} is a substitution {b 1 → (I vulnerable b 1 ), . . . , b n → (I vulnerable b n )}; 2. s effective I, if the evaluation of the result sσ of applying σ to s returns true; 3. I vulnerable T f ⇔ def ∃s ∈ SC(T f ) s effective I; 4. I unsat c SP (T f ) ⇔ def ∃sp ∈ SP (T f ), I unsat c sp.

Figure 7

 7 depicts the distribution of the participant opinions.

Fig. 7 .

 7 Fig. 7. Response rates for Question Q7 We collected the time spent by the participants for writing test cases. The participants needed between 15 and 70 minutes in Part 1, while they took between 20 minutes and 86 minutes in Part 2. On average, they spent 46 minutes in Part 1 and 60 minutes in Part2. The levels of confidence of the participants is estimated with Question Q9. The possible answers were for both scenarios: very sure, sure, fairly sure, not sure.We finally analysed the test cases given by the participants and evaluated their correctness with regard to four aspects: 1&2: detection (with at least one test case) that both applications are vulnerable to the attacks CAPEC 66 and CAPEC 244; 3&4: detection that the application behaviours do not include the consequences of the patterns "Input Guard" and "Output Guard". As we considered this last aspect as difficult for students, we expected at least one Then test case section for every pattern. Figure8presents the number of participants who meet these aspects.

Fig. 8 .

 8 Fig. 8. Test case correctness (Question Q10)

https://capec.mitre.org/

https://talend.com/

https://www.owasp.org/index.php/OWASP Zed Attack Proxy Project

http://regainia.com/research/companion.html

https://cucumber.io/

http://regainia.com/research/companion.html