Dynamic grain models via fast heuristics for diagram representations

ABSTRACT

The present paper introduces a mathematical model for analysing dynamic grain growth. In particular, we show how the characteristic measurements grain volumes, centroids, and central second-order moments at discrete moments in time can be quickly turned into a continuous description of the grain growth process in terms of geometric diagrams (which largely generalize the well-known Voronoi and Laguerre tessellations). We give a theoretical analysis of common optimization-free heuristics in terms of discriminant analysis and evaluate the computational behaviour of our algorithm on real-world data.

KEYWORDS:

• Power diagrams
• generalized balanced power diagrams
• tessellations
• grains
• polycrystals
• grain growth

1. Introduction

Grain growth is an important field of study in materials science as the resulting grain structures strongly influence the mechanical and physical properties of metals, ceramics, and other polycrystalline materials [1, 2]. Although, starting with works of Smith, von Neumann, and Mullins in the 1950s [3–5], a number of different models have been introduced over the years, understanding and controlling grain growth remains challenging, both in theory and practice. For instance, the recent empirical study [6] demonstrates that neither the classical Hillert model [7] nor the MacPherson–Srolovitz model [8] capture anisotropic grain growth adequately. Even more drastically, [9] concludes that ‘a new model for grain boundary migration is needed to predict microstructure evolution’; see also [10–13].

The present paper introduces and analyses a novel tool which, when applied to real data, provides an empirical grain growth representation that can be used to analyse physical models and simulate grain growth. More precisely, we develop a mathematical interpolation/extrapolation model for quickly turning grain information acquired at discrete moments in time into a continuous description of the grain growth process. It employs diagram representations of grain maps (which will be formally introduced in Section 2) based on measured (or simulated) estimates of the grain volumes, centroids, and, if available, central second-order moments. Such characteristics can be obtained, for instance, by utilizing tomographic techniques [1, 14] or from otherwise available grain maps.

The diagram representations considered in this paper are anisotropic generalizations of previously used (isotropic) Laguerre and Voronoi tessellations [15–17] and allow for curved boundaries and non-convex grains. While these representations are based on a phenomenological model whose physical foundation is still not fully understood, they have consistently been reported to capture the grain boundaries in polycrystalline materials quite accurately; see [18–25].

Several sophisticated optimization algorithms are available for computing or approximating such diagrams with high accuracy, [18, 20], but they require relatively large computational cost. This can be a bottleneck, particularly in dynamic 3D studies. We therefore focus in the present paper on extremely fast optimization-free heuristics for fitting diagrams to measured data. Methods for the related but generally different task of generating synthetic grain structures that follow a given distribution are discussed in [26–30].

First, we use discriminant analysis to interpret, for the first time in literature, a heuristic (subsequently referred to as H1) which was introduced in [15] for the isotropic case and extended in [25] to the anisotropic case. Our analysis leads naturally to an even simpler and faster heuristic H2, proposed in [19], which, again for the first time, can now be explained to perform equally well under certain conditions. Also, the comparison of H1 and H2 and, in particular, our interpretation of the so-called size parameters, suggests that diagram parameters may provide additional insight into the stage of the growth process, and this will be used later on in the paper. Our theoretical foundation of the heuristics is then empirically verified. In fact, we apply H1 and H2 to real-world 3D grain map data in order to compare their quality. Then, we introduce an interpolation/extrapolation model to obtain a continuous description of the grain growth process based on measurements at only a few moments in time. Finally, we evaluate the computational behaviour of our algorithm on a real-world data set yielding an empirical grain growth model. Our study indicates that the growth process converges towards particularly simple diagrams, which seem to be energetically favourable.

The paper is organized as follows. Section 2 introduces the relevant diagrams and formalizes the concept of grain scans. Sections 3 and 4 address optimization-free heuristics, both, in terms of a theoretical analysis and a practical evaluation on real-world data. Section 5 introduces the interpolation/extrapolation model which is capable of extending measurements at a few time steps to a continuous grain growth model very quickly. Further, the results of our empirical study on practical grain data are presented. Section 6 concludes with some final remarks.

2. Diagrams and grain scans

Adequate representations of grain maps can be valuable tools for understanding the physical principles of grain growth and the properties of the resulting material. In fact, the diagram representations studied in this paper require only a few parameters per grain and are thus relatively easy to process. As they exhibit geometric and combinatorial features of grain maps which are not readily available otherwise, such representations may contribute to predicting and, potentially, even controlling the forming of new materials.

Some or all parameters which specify a diagram can be estimated from certain tomographic measurements, and optimizing the remaining ones then constitutes an inverse problem; see [31]. Previous studies [18–25] indicate that, in particular, grain volumes, centroids, and central second-order moments govern decisive properties of grain maps.

Let us now formally introduce the appropriate class of diagrams. Let $d\in \mathbb{N}$ denote the dimension of space. Here, we are mainly interested in the 3 dimensional case. However, our approach works for general d, and, in fact, planar sections or surfaces of polycrystals are also relevant in practice. Further, let $k\in \mathbb{N}$ denote the number of different grains in the underlying grain map. An anisotropic power diagram is specified by the following data for each $i\in \left[k\right]:=\{1,\dots ,k\}$:

• a positive definite symmetric matrix $A_i\in {\mathbb{R}}^{d\times d}$,

• a site $s_i\in {\mathbb{R}}^d$, and a

• size ${\gamma }_i\in \mathbb{R}$.

For simpler reference, we collect these parameters in the two families $\mathcal{A}=\{A_1,\dots ,A_k\}$ and $\mathrm{\Gamma }=\{{\gamma }_1,\dots ,{\gamma }_k\}$, where repetitions of elements are permitted, and the set $S=\{s_1,\dots ,s_k\}$. The sites specify the positions of the grains while the matrices in $\mathcal{A}$ describe characteristics of their shapes. In fact, each grain is equipped with its own ellipsoidal norm $\Vert {\Vert }_{A_i}$ given through $A_i$ by $\Vert x{\Vert }_{A_i}=\sqrt{x^{\mathrm{\top }}{A}_{i}x},(x\in {\mathbb{R}}^d).$The sizes ${\gamma }_i$ are used to control the grain volume.

Given $\mathcal{A}$, S, and $\mathrm{\Gamma }$, the anisotropic power diagram (in the following often simply referred to as diagram) $\mathcal{P}=\mathcal{P}(\mathcal{A},S,\mathrm{\Gamma })=\{P_1,\dots ,P_k\},$is a decomposition of ${\mathbb{R}}^d$ into cells $P_i$ defined by $P_i=\{x\in {\mathbb{R}}^d:\Vert x-s_i{\Vert }_{A_i}^2+{\gamma }_i\le \Vert x-s_{\ell }{\Vert }_{A_{\ell }}^2+{\gamma }_{\ell },\mathrm{\forall }\ell \in [k]\setminus \{i\left\}\right\}.$In particular, two different cells do not share interior points, i.e. $int\left(P_i\right)\cap int\left(P_{\ell }\right)=\mathrm{\varnothing }$ for $i\ne \ell$. See [32, 33] for examples and further properties of different classes of diagrams.

Diagrams can be viewed as continuous (or resolution-independent) representations of polycrystalline samples. Conventionally, the samples are represented as voxel-based digital images at some resolution, called grain maps. More generally, given a set of points $X=\{x_1,\dots ,x_n\}\subseteq {\mathbb{R}}^d$ in the sample and labels, i.e. $y_1,\dots ,y_n\in \left[k\right]$, with $k\in \mathbb{N}$ denoting again the number of different grains in the sample, we will call the labeled data set $\left\{\right(x_j,y_j):j\in [n\left]\right\}$ a grain scan or polycrystal scan. In the following, we refer to it as $(X,Y)$ with the understanding, that the indices of the labels in the family $Y=\{y_1,\dots ,y_n\}$ correspond to the indices of the points in X. A grain scan therefore provides a partitioning of X into $G_i=\{x_j:(x_j,i)\in (X,Y\left)\right\}$, the (discretizations of the) grains.

We can measure the fit of the diagram in terms of the symmetric difference of $G_i$ and $X\cap P_i$ for all $i\in \left[k\right]$. (Points on the boundary of two or more cells are assigned according to some preference criterion or, conservatively, can be viewed as classification error, but are essentially irrelevant in our context and hence generally discarded.)

Usual grain maps are special grain scans where the sample is assumed to be in the range $[0,1{]}^d$ and points of X are identified with the points of $[0,1{]}^d$ in the grid $\rho {\mathbb{Z}}^d$ for some resolution $\rho \in (0,1)$. Note that, with ${\kappa }_i=\left|G_i\right|$, the term ${\rho }^d{\kappa }_i$ can be regarded as an approximation of the volume ${\nu }_i$ of the ith grain.

If a grain scan is available the computation of a best fitting diagram requires the optimization of all three characteristics $\mathcal{A}$, S, and $\mathrm{\Gamma }$. While this can indeed be done, see [20, 31], it is, however, too time-consuming, given current optimization technology, for dynamic 3D applications. Hence, we will in the following assume that $\mathcal{A}$ and S are available through measurements and $\mathrm{\Gamma }$ will be determined by a heuristic which may, additionally, use grain volume information, i.e. (approximations of) ${\nu }_1,\dots ,{\nu }_k$. As such heuristics abstain from any optimization, they are computationally fast but potentially result in a higher misclassification error.

3. Optimization-free heuristics

This section considers two heuristics that derive all characteristic parameters $\mathcal{A},S,\mathrm{\Gamma }$ for an anisotropic power diagram directly from given measurements. These heuristics, introduced in [15, 19, 25], avoid optimization routines altogether, are fast and, as we will see, already provide a reasonable classification accuracy. More specifically, we provide a new interpretation of the heuristics within the framework of discriminant analysis.

3.1. Choices for $\mathcal{A}$ and S

Discriminant analysis is a popular classification technique in statistics and provides a different perspective on the a priori choices for the characteristic parameters $(\mathcal{A},S,\mathrm{\Gamma }$); see e.g. [34]. Within its realm, the grain scan is regarded as a sample of the multivariate random variable $(\mathcal{X},\mathcal{Y})$ where $\mathcal{X}$ specifies a point in $[0,1{]}^d$ and $\mathcal{Y}$ gives the index $i\in \left[k\right]$ of the corresponding grain. Discriminant analysis estimates the distribution of the multivariate random variable from samples and assigns a point to the class with the highest probability p. Hence, the decision boundary between two classes $i,\ell$ consists of those points x that have an equal conditional probability of lying in either class, i.e. $\{x\in {\mathbb{R}}^d:p(\mathcal{Y}=i|\mathcal{X}=x)=p(\mathcal{Y}=\ell |\mathcal{X}=x\left)\right\}.$For its estimation, discriminant analysis applies the maximum likelihood principle under the assumption that the members from each class are normally distributed. As each normal distribution $\mathcal{N}({\mu }_i,{\mathrm{\Sigma }}_i)$ is fully determined by its d-dimensional mean vector ${\mu }_i$ and its $(d\times d)$-covariance matrix ${\mathrm{\Sigma }}_i$, it is these characteristics that need to be estimated. It is well known that the maximum likelihood principle results in the choices of the sites S as the centroids and the matrices $\mathcal{A}$ as the inverse of the covariance matrices of the sample, i.e. more precisely, $\begin{array}{rl}{\mu }_i& =c_i=\frac{1}{\left|G_i\right|}\sum _{x_j\in G_i}{x}_j,\\ A_i& ={\left(\mathrm{\Sigma }\right(G_i\left)\right)}^{-1}={\left(\frac{1}{\left|G_i\right|}\right(x_j-c_i:x_j\in G_i\left)\right(x_j-c_i:x_j\in G_i{)}^T)}^{-1}\end{array}$for $i\in \left[k\right]$. Note that the latter requires that the matrices $\mathrm{\Sigma }\left(G_i\right)$ have full rank d. This can be assumed here as, in the underlying model, grains that are degenerate to volume 0 do not matter. The above maximum likelihood setting is precisely the choice for $\mathcal{A}$ and S suggested in [18], and subsequently also used in [25]. Also, estimates for these parameters are available through tomographic measurement.

Let us point out that, in practice, discriminant analysis is often applied even when the underlying assumption of a normal distribution is only met approximately (as is the case here). This is justified as the procedure is generally reported to be quite robust; see, e.g. [35].

3.2. Choices for Γ

We will now turn to the third set of characteristics, i.e. the sizes $\mathrm{\Gamma }$. We show first that the choice of $\mathrm{\Gamma }$ in heuristic H1 can also be interpreted within the paradigm of discriminant analysis with a specific prior. Building on this interpretation, a different natural choice will lead to H2.

To begin with, let us briefly recall the volume argument of [15] and (in its anisotropic generality) [25] that leads to the specific choice for $\mathrm{\Gamma }$ in heuristic H1.

Let $A\in {\mathbb{R}}^{d\times d}$ be positive definite and symmetric. As A can be diagonalized there exist an orthogonal matrix $U\in {\mathbb{R}}^{d\times d}$ and a diagonal matrix $D=\mathrm{diag}({\lambda }_1,\dots ,{\lambda }_d)\in {\mathbb{R}}^{d\times d}$ with ${\lambda }_1,\dots ,{\lambda }_d>0$ such that $A=UD{U}^T$. Then, using the abbreviation $D^{\frac{1}{2}}=\mathrm{diag}(\sqrt{{\lambda }_1},\dots ,\sqrt{{\lambda }_d})$, we have for the unit ball ${\mathbb{B}}_A^d$ and $\rho \in (0,\mathrm{\infty })$ $\begin{array}{rl}\rho {\mathbb{B}}_A^d& =\{x\in {\mathbb{R}}^d:\Vert x{\Vert }_A\le \rho \}=\{x\in {\mathbb{R}}^d:x^{T}Ax\le {\rho }^2\}\\ & =\{x\in {\mathbb{R}}^d:x^{T}UD{U}^{T}x\le {\rho }^2\}=\rho U{D}^{-\frac{1}{2}}{\mathbb{B}}_2^d.\end{array}$Hence, the volume $\mathrm{vol}\left(\rho {\mathbb{B}}_A^d\right)$ of this ellipsoid is given by $\mathrm{v}\mathrm{o}\mathrm{l}\left(\rho {\mathbb{B}}_A^d\right)={\rho }^d\mathrm{v}\mathrm{o}\mathrm{l}\left(U{D}^{-\frac{1}{2}}{\mathbb{B}}_2^d\right)={\rho }^{d}det(A{)}^{-\frac{1}{2}}\mathrm{v}\mathrm{o}\mathrm{l}\left({\mathbb{B}}_2^d\right).$Regarding the ellipsoid $c_i+\rho {\mathbb{B}}_{A_i}^d$ as an approximation of the ith grain, it is now natural to choose ${\gamma }_i$ according to the volume condition $\mathrm{v}\mathrm{o}\mathrm{l}\left(\{x\in {\mathbb{R}}^d:\Vert x{\Vert }_{A_i}^2+{\gamma }_i\le 0\}\right)={\nu }_i$involving the available approximation ${\nu }_i$ for the volume of the ith grain. A simple computation yields ${\gamma }_i=-{\rho }^2=-{\left(\frac{{\nu }_i}{\mathrm{v}\mathrm{o}\mathrm{l}\left({\mathbb{B}}_{A_i}^d\right)}\right)}^{\frac{2}{d}}=-{\left(\frac{{\nu }_i\sqrt{det\left(A_i\right)}}{\mathrm{v}\mathrm{o}\mathrm{l}\left({\mathbb{B}}_2^d\right)}\right)}^{\frac{2}{d}}$For d = 3, and with the specific setting of $A_i=\left(\mathrm{\Sigma }\right(G_i){)}^{-1}$, this yields the choice ${\gamma }_i=-{\left(\frac{3{\nu }_i}{4\pi \sqrt{det\left(\mathrm{\Sigma }\right(G_i\left)\right)}}\right)}^{\frac{2}{3}}$from [15, 25] (the heuristic H1). In order to interpret this choice via discriminant analysis, we determine which probability distribution, called prior, will lead to the same decision boundaries between the classes. So, let ${\pi }_1,\dots ,{\pi }_k\in [0,\mathrm{\infty }),N=\sum _{i=1}^k{\pi }_i\ne 0.$Then $({\pi }_1,\dots ,{\pi }_k)$ specifies the a priori ‘belief’ about the prior probability ${\pi }_i/N$ that a point (regardless of its position) belongs to cluster i, i.e. $\frac{{\pi }_i}{N}=p(\mathcal{Y}=i),(i\in [k\left]\right).$The computation of the decision boundary is then facilitated via Bayes' Theorem. As each class is assumed to be normally distributed, Bayes' Theorem involves the conditional density functions $f_{\mathcal{X}|\mathcal{Y}=i}\left(x\right)=\frac{det(A_i{)}^{\frac{1}{2}}}{(2\pi {)}^{\frac{d}{2}}}\cdot \exp (-\frac{1}{2}(x-c_i{)}^{\mathrm{\top }}{A}_i(x-c_i))$of the k estimated normal distributions $\mathcal{N}({\mu }_i,{\mathrm{\Sigma }}_i),$ i.e. with ${\mu }_i=c_i,{\mathrm{\Sigma }}_i=\mathrm{\Sigma }\left(G_i\right)=A_i^{-1}.$In this situation, Bayes' Theorem states that $p(\mathcal{Y}=i|\mathcal{X}=x)=\frac{f_{\mathcal{X}|\mathcal{Y}=i}\left(x\right){\pi }_i}{\sum _{\ell =1}^k{f}_{\mathcal{X}|\mathcal{Y}=\ell }\left(x\right){\pi }_{\ell }}.$Therefore, the maximum likelihood decision boundaries between two classes $i,\ell$ can be determined by solving the equation $f_{\mathcal{X}|\mathcal{Y}=i}\left(x\right)\cdot {\pi }_i=f_{\mathcal{X}|\mathcal{Y}=\ell }\left(x\right)\cdot {\pi }_{\ell }.$Taking the logarithm, using that $(x-c_i{)}^{\mathrm{\top }}{A}_i(x-c_i)=\Vert x-c_i{\Vert }_{A_i}^2$, and simplifying, this condition turns into the equation $\Vert x-c_i{\Vert }_{A_i}^2-\log (det(A_i\left)\right)-2\log \left({\pi }_i\right)=\Vert x-c_{\ell }{\Vert }_{A_{\ell }}^2-\log (det(A_{\ell }\left)\right)-2\log \left({\pi }_{\ell }\right).$Hence, these decision boundaries coincide with that of the diagram with the parameters $(\mathcal{A},S,\mathrm{\Gamma })$ if and only if ${\gamma }_i={\gamma }_i^{\left(1\right)}=-\log (det(A_i\left)\right)-2\log \left({\pi }_i\right)+\gamma ,(i\in [k\left]\right),$where $\gamma \in \mathbb{R}$ is an arbitrary constant. (Recall that the diagrams are invariant under such common additions to the sizes.) For the above choice of the sizes in H1, we hence obtain ${\pi }_i={\pi }_i^{\left(1\right)}=\frac{{\mathrm{e}}^{\frac{\gamma }{2}}}{\sqrt{det\left(A_i\right)}}\cdot \exp \left(\frac{1}{2}{\left(\frac{{\nu }_i\sqrt{det\left(A_i\right)}}{\mathrm{v}\mathrm{o}\mathrm{l}\left({\mathbb{B}}_2^d\right)}\right)}^{\frac{2}{d}}\right),(i\in [k\left]\right).$Since $\frac{1}{\sqrt{det\left(A_i\right)}}=\sqrt{det\left(\mathrm{\Sigma }\right(G_i\left)\right)}=\frac{1}{\mathrm{v}\mathrm{o}\mathrm{l}\left({\mathbb{B}}_2^d\right)}\cdot \mathrm{v}\mathrm{o}\mathrm{l}\left({\mathbb{B}}_{A_i}^d\right),$the first factor is proportional to the volume of the ellipsoid ${\mathbb{B}}_{A_i}^d$. The second factor is a ‘correction term’ which increases exponentially if the measured grain volume ${\rho }^d{\kappa }_i$ exceeds $\mathrm{vol}\left({\mathbb{B}}_{A_i}^d\right)$. The first term seems natural while the second is quite unusual.

Within the realm of discriminant analysis it appears reasonable to simply set ${\gamma }_i=0$ for all clusters i. Then we obtain the priors ${\pi }_i={\pi }_i^{\left(2\right)}=\frac{{\mathrm{e}}^{\frac{\gamma }{2}}}{\sqrt{det\left(A_i\right)}},(i\in [k\left]\right),$which are proportional to the volumes of the ellipsoids corresponding to the covariance matrices of the grains. Thus, this simpler prior distribution encodes the assumption that the ellipsoid volumes already represent the grain volumes well. This justifies the choice ${\gamma }_1^{\left(2\right)}=\cdots ={\gamma }_k^{\left(2\right)}=0,$proposed in [19]. We will refer to the respective heuristic as H2.

In the next section we compare the behaviour of H1 and H2 on available real-world grain scans.

4. Comparison of H1 and H2

We compare the two heuristics H1 and H2 on a real-world 3D grain scan and indicate how their difference can be interpreted. This may be of particular interest for understanding the state of a grain growth process and will be taken up in Section 5.

4.1. Experimental evaluation

We compare the methods with a data set taken from [15] that has been obtained by a synchrotron micro-tomograph experiment conducted on a metastable beta titanium alloy (Ti β21S). The material was scanned with the volume of 240μm × 240μm × 420μm at a resolution corresponding to a voxel size of 0.7 μm. This scan resulted in a 3D voxel image of size $339\times 339\times 599$ composed of a total of 591 grains. All our computations are carried out in 3D, i.e. d = 3.

In the following, we assess the quality of fit by means of different measures. Let $\mathcal{P}=\{P_1,\dots ,P_k\}$ be a diagram (obtained by H1, H2, or, for that matter, any other algorithm). We set ${\hat{C}}_i=\{x_j\in X:x_j\in int(P_i\left)\right\}(i\in [k\left]\right),\hat{C}\left(\mathcal{P}\right)=\hat{C}=({\hat{C}}_1,\dots ,{\hat{C}}_k).$Note that $\hat{C}$ is not a clustering of X but only of $X\cap \bigcup _{i=1}^{k}int\left(P_i\right)$; this is indicated by writing $\hat{C}$ rather than C.

A natural measure how well $\mathcal{P}$ captures the grain scan $G=(G_1,\dots ,G_k)$ is the relative fit or accuracy ${\mathrm{\Phi }}_G$ of $\hat{C}$ defined by ${\mathrm{\Phi }}_G\left(\hat{C}\right)=\frac{1}{n}\left|\bigcup _{i=1}^k(G_i\cap {\hat{C}}_i)\right|.$Let us point out that, as they are not assigned within $\hat{C}$, all points on the cell boundaries count as misclassified.

In Table 1 we will also report on the behaviour with respect to the relative cluster weight error ${\mathrm{\Psi }}_G\left(\hat{C}\right)=\frac{1}{n}\sum _{i=1}^k|{\kappa }_i-|{\hat{C}}_i\left|\right|,$and the relative deviation of the centroids $c\left(G_i\right)$ and $c\left({\hat{C}}_i\right)$, and the covariance matrices $A_i^{-1}={\mathrm{\Sigma }}_i=\mathrm{C}\mathrm{o}\mathrm{v}\left(G_i\right)$ and $\mathrm{C}\mathrm{o}\mathrm{v}\left({\hat{C}}_i\right)$, respectively. Formally, the relative centroid error and the relative covariance error are defined as $\frac{1}{n}\sum _{i=1}^k{\kappa }_i\Vert c\left(G_i\right)-c\left({\hat{C}}_i\right){\Vert }_2\mathrm{a}\mathrm{n}\mathrm{d}\frac{1}{n}\sum _{i=1}^k{\kappa }_i\Vert {\mathrm{\Sigma }}_i-\mathrm{C}\mathrm{o}\mathrm{v}\left({\hat{C}}_i\right){\Vert }_2.$For the latter, we used the spectral norm. This matrix norm measures the largest singular value of the covariance differences and, thus, intuitively captures the difference in the largest singular vectors of the covariance matrices. Note that, as defined, the latter two measures depend on the original physical dimensions of the sample. We refrain from any additional normalization as we do not conduct inter-sample studies here.

Table 1. Comparison of H1 and H2.

Performance characteristics	H1	H2
Accuracy (in %)	92.98652	92.98658
Relative cluster weight error	0.02869	0.02869
Relative centroid error (in μm)	0.64346	0.64347
Relative covariance error (in μm${}^2$)	15.25372	15.25376
Correct neighbourhoods (in %)	42.30118	42.30118
Correct neighbourhoods up to 1 Error (in %)	80.54146	80.54146
Correct neighbourhoods up to 2 Errors (in %)	94.92386	94.92386
Running time (in s)	1.49641	0.00000

We will also address combinatorial features. More precisely, we consider the percentage of grains with a correct neighbourhood, i.e. with all grain neighbours in the ground truth being also neighbours in the representation. Here, a grain is a neighbour of another grain if a voxel of the former grain is adjacent (face-connected, edge-connected or vertex-connected) to a voxel of the later grain. We also provide this percentage when one and, respectively, two errors in the neighbourhood are allowed, i.e. when there is either one additional or one missing voxel.

The results of our evaluation are shown in Table 1.

Let us point out that the displayed running times involve the heuristics but, for a fairer comparison, do not include the computation of the performance criteria. Hence, the entry 0 for the running time of H2 refers to the fact that all parameters defining $\mathcal{P}$ are given upfront. Also, while we exercised appropriate care for the implementations, we do not claim that our code is fully optimized.

In terms of fit, we observe that both heuristics perform extremely well, in fact, almost identical, on the given real-world data set. We will now offer an explanation for the latter behaviour.

4.2. Similarity and difference

While both heuristics, H1 and H2, employ the same matrices $\mathcal{A}$ and sites S, they generally differ in the sizes $\mathrm{\Gamma }$ of the resulting diagrams. Hence, the nearly identical performance depicted in Table 1 for all applied quality criteria may come as a surprise. As observed before, the priors underlying H1 and H2 denoted ${\pi }_i^{\left(1\right)}$ and ${\pi }_i^{\left(2\right)},$ as before, differ by a factor that depends exponentially on the deviation of $\mathrm{vol}\left({\mathbb{B}}_{A_i}^d\right)$ from the measured grain volume ${\nu }_i$.

Now, let ${\alpha }_i$ be a positive real number such that ${\nu }_i={\alpha }_i\cdot \mathrm{v}\mathrm{o}\mathrm{l}\left({\mathbb{B}}_{A_i}^d\right)={\alpha }_i\cdot \frac{\mathrm{v}\mathrm{o}\mathrm{l}\left({\mathbb{B}}_2^d\right)}{\sqrt{det\left(A_i\right)}}.$Then, for each $i\in \left[k\right]$, $\begin{array}{rl}{\pi }_i^{\left(1\right)}& =\frac{{\mathrm{e}}^{\frac{\gamma }{2}}}{\sqrt{det\left(A_i\right)}}\cdot \exp \left(\frac{1}{2}{\left(\frac{{\nu }_i\sqrt{det\left(A_i\right)}}{\mathrm{v}\mathrm{o}\mathrm{l}\left({\mathbb{B}}_2^d\right)}\right)}^{\frac{2}{d}}\right)=\frac{{\mathrm{e}}^{\frac{\gamma }{2}}}{\sqrt{det\left(A_i\right)}}\cdot \exp \left(\frac{1}{2}{\alpha }_i^{\frac{2}{d}}\right)\\ & ={\pi }_i^{\left(2\right)}\cdot \exp \left(\frac{1}{2}{\alpha }_i^{\frac{2}{d}}\right),\end{array}$and, expressed in terms of the sizes ${\gamma }_i^{\left(1\right)}$ used in H1, we obtain (1) ${\gamma }_i^{\left(1\right)}=-{\left(\frac{{\nu }_i}{\mathrm{v}\mathrm{o}\mathrm{l}\left({\mathbb{B}}_{A_i}^d\right)}\right)}^{\frac{2}{d}}=-{\alpha }_i^{\frac{2}{d}}.$(1) Now, suppose that all factors ${\alpha }_i$ are identical, i.e. ${\nu }_i=\alpha \cdot \mathrm{v}\mathrm{o}\mathrm{l}\left({\mathbb{B}}_{A_i}^d\right)$ for some positive real number α independent of i. Then, of course, ${\pi }_i^{\left(1\right)}={\pi }_i^{\left(2\right)}\cdot \exp \left(\frac{1}{2}{\alpha }^{\frac{2}{d}}\right),{\gamma }_i^{\left(1\right)}=-{\alpha }^{\frac{2}{d}},(i\in [k\left]\right).$Hence, the (normalized) priors for H1 and H2 coincide, which implies that the maximum likelihood decision boundaries for the two heuristics are the same. Further, since diagrams are invariant under the addition of a constant to each size, the diagrams produced by H1 and H2 are identical.

We suspect that grain structures for which the volumes of the unit balls ${\mathbb{B}}_{A_i}^d$ with respect to $A_i=\left(\mathrm{\Sigma }\right(G_i){)}^{-1}$ are nearly proportional to the grain volumes ${\nu }_i$ might turn out to be energetically favourable. According to the above reasoning, the former may be quantified by the variance of the powers ${\mathrm{\Gamma }}^{\frac{d}{2}}=({\gamma }_1^{\frac{d}{2}},\dots ,{\gamma }_k^{\frac{d}{2}})$ computed in H1, i.e. $\mathrm{v}\mathrm{a}\mathrm{r}\left({\mathrm{\Gamma }}^{\frac{d}{2}}\right)=E({\mathrm{\Gamma }}^{\frac{d}{2}}-E({\mathrm{\Gamma }}^{\frac{d}{2}}\left)\right)=\frac{1}{k}\sum _{i=1}^k{({\gamma }_i^{\frac{d}{2}}-\frac{1}{k}\sum _{i=1}^k{\gamma }_i^{\frac{d}{2}})}^2.$As an example, we will in Section 5 compute the evolution of $\mathrm{v}\mathrm{a}\mathrm{r}\left({\mathrm{\Gamma }}^{\frac{3}{2}}\right)$ at different moments in time for a real-world 3D grain growth process. In terms of H1 and H2, we might expect that if the two heuristics perform differently in practice, the current grain structure will not be energetically minimal, and further growth will result in an improved structure.

We close this section with a ‘constructed’ theoretical example which sheds some additional light on the behaviour of the two heuristics. Since the general principle of the construction is the same in higher dimension and can, hence, already be illustrated in 2D we restrict the following description to the case d = 2.

As depicted in the top row of Figure 1, the image consists of two grains $G_1$ and $G_2$; one (orange [light gray]) is the square $[0,\text{β}]\times [0,\text{β}]$ for different values of $\text{β}\in (0,1)$, the other (blue [gray]) is its complement in $[0,1]\times [0,1].$ Of course, one would not expect such samples to originate from real-world grain scans. The bottom row of Figure 1 shows the corresponding ellipses $c_i+{\mathbb{B}}_{A_i}^2$ with respect to $A_i=\left(\mathrm{\Sigma }\right(G_i){)}^{-1}$.

Figure 1. A family of grain scans involving a parameter β to illustrate differences between H1 and H2. Top, from left to right: Grain scans for $\text{β}=0.1$, $\text{β}=0.5$, and $\text{β}=0.9$. Grain $G_1$ (orange [light gray]) fills the square $[0,\text{β}]\times [0,\text{β}]$, $G_2$ (blue [gray]) fills the rest of the unit square. The bottom row shows the computed ellipses $c_i+{\mathbb{B}}_{A_i}^2$ with $A_i=\left(\mathrm{\Sigma }\right(G_i){)}^{-1}$ for i = 1, 2.

As argued before, the variance of the ratios ${\nu }_i/\mathrm{v}\mathrm{o}\mathrm{l}\left({\mathbb{B}}_{A_1}^2\right)$ is a relevant parameter. Figure 2 therefore provides plots of $\frac{{\text{β}}^2}{\mathrm{v}\mathrm{o}\mathrm{l}\left({\mathbb{B}}_{A_1}^2\right)}\mathrm{a}\mathrm{n}\mathrm{d}\frac{1-{\text{β}}^2}{\mathrm{v}\mathrm{o}\mathrm{l}\left({\mathbb{B}}_{A_2}^2\right)}$as a function of β.

Figure 2. Ratios ${\nu }_i/\mathrm{v}\mathrm{o}\mathrm{l}\left({\mathbb{B}}_{A_i}^2\right)$ of the areas of the grain and the corresponding ellipsoid for i = 1 (orange [light gray]) and i = 2 (blue [gray]) as a function of β.

Note that these ratios differ considerably for larger values of β. Hence, we would expect that the heuristics H1 and H2 behave quite differently in these cases. This is confirmed by the resulting diagrams depicted in Figure 3.

Figure 3. Diagrams for the family of grain scans from Figure 1 obtained via H1 (top row) and H2 (bottom row).

The behaviour for $\text{β}=0.9$ suggests that H1 tries to preserve the areas of the grains, while H2 seems to favour their shapes. If the relevant volume ratios are close, as is the case for small β, both heuristics achieve both aims and perform similarly.

Figure 4 shows the ratios of the grain volume ${\nu }_i$ and the areas of the diagram cells under the two models H1 and H2, respectively.

Figure 4. Ratio of the area ${\nu }_i$ of the original grain and the area of the corresponding diagram cell obtained by H1 (left) and H2 (right) for i = 1 (orange [light gray]) and i = 2 (blue [gray]) as a function of β.

5. Dynamic grain models

We now describe a simple interpolation/extrapolation model that is capable of quickly turning information acquired at a few moments in time into a continuous description of the grain growth process. To demonstrate the potential of the model, we then report on the results of an empirical evaluation on a real-world grain scan time series. Finally, we study for this empirical model the correlation between the diagram size parameters Γ and the disappearance of grains in the process.

5.1. The model

Since well-chosen diagrams have been demonstrated to yield high-quality representations of grain structures of many types of polycrystalline materials, it is natural to exploit their rather low-dimensional parameter space for deriving continuous models from data measured at discrete times. Recall that a diagram is fully characterized by the triple $(\mathcal{A},S,\mathrm{\Gamma })$. In order to simulate the dynamic evolution of grain structures over time, we can fit a polynomial (or any other desired type of function) for each diagram parameter extracted from a discrete-time series polycrystal data set. These functions then represent the continuous change in diagram parameters during the growth process. Understanding the evolution of the diagram parameters over time may help predict changes in the material during an annealing process and may provide further insights into the material characteristics.

Suppose, that 3D-grain scans (or at least appropriate characteristic measurements) are available at the discrete times $t_0,t_1,\dots ,t_m$. Then we compute, for each $t\in \{t_0,\dots ,t_m\}$, the parameters of an anisotropic power diagram ${\mathcal{P}}_t=\mathcal{P}({\mathcal{A}}_t,S_t,{\mathrm{\Gamma }}_t)$ which represent the grain scans. Subsequently, we select a class $\mathcal{F}$ of parameterized functions and determine their parameters that best fit the data ${\mathcal{A}}_t,S_t,{\mathrm{\Gamma }}_t,t\in \{t_0,\dots ,t_m\}.$ In our experimental setting we choose polynomials of a small degree to avoid overfitting, and we measure the fit in terms of least-squares deviation.

In 3D, the interpolation/extrapolation involves a class $\mathcal{F}$ consisting of 10 functions for each grain, three for the coordinates of the diagram sites, three for the rotation part of the grain's covariance matrix, three for the length of the principal axis of the covariance matrix, and one for the size parameter. Of course, in order to set up the computation of the functions, we need to ‘track’ each grain, i.e. we need to identify which grain corresponds to which at the different moments t in time.

During the grain growth process, grains may disappear. We assume, however, that no new grains emerge. In the experimental series that will be described in more detail below, about $40\text{\%}$ of all initially present grains disappeared throughout the process. Hence we distinguish between such non-persistent grains and persistent grains that exist throughout the observed time interval. Of, course, for non-persistent grains, the functions are fitted only as long as the grain ‘definitely exists’. More specifically, if a previously existing grain disappears from the measurement, we choose as last diagram parameters the covariance matrix from the previous time step and compute a size parameter such that its grain volume is essentially zero. If, on the other hand, for some moment in time our prediction generates a covariance matrix and a size parameter which result in a grain volume close to zero, we exclude the grain from further consideration.

The diagrams ${\mathcal{P}}_t$ can, in principle, be computed by any suitable method. In addition to the obvious criteria that any employed algorithm should be computationally fast and produce diagrams of appropriate accuracy, one has to pay attention to its behaviour concerning the invariances within the parameter space to avoid corruption of the interpolation/extrapolation. As we will be using the heuristic H1 described in Section 3, all parameters are, however, uniquely determined by the measurements, and we can hence, in the following, ignore this issue.

5.2. The data

For our experimental study, we use the data set from [6, 36] (details of the data processing chain, from experiment to grain scan, can be found there and, in particular, in their supplementary materials sections). The data set, obtained by a diffraction contrast tomography (DCT) [37, 38] experiment at beamline ID11 at the European Synchrotron Radiation Facility (ESRF), comprises 15 3D-grain scans (voxelized grain maps) of a 99.9% pure polycrystalline iron sample taken at different moments in time of an annealing process. More precisely, the material was first cold-rolled and annealed for 30 min at ${700}^{\circ }\text{C}$ to fully recrystallize before it was first scanned (time $t_0).$ Subsequently, the sample was annealed at a constant temperature of ${800}^{\circ }\text{C}$ for 10 min leading to grain growth, cooled to stabilize, and scanned again (time $t_1)$. This was repeated another 13 times with annealing periods of 5 min (time $t_2,\dots ,t_{14}$). The effect of the growth process is already visible by the number of grains in the sample. While the first scan exhibits 1327 grains, the last one contains only 776, and the number of interior grains strictly decreases from 387 to 189; see [6, Table 1] for more details. The identification of grains over the 15 moments in time is already provided by in [6, 36].

For this data set, i.e. for each 3D-grain scan at $t\in \{t_0,t_1,\dots ,t_{14}\},$ we computed an anisotropic power diagram ${\mathcal{P}}_t=\mathcal{P}({\mathcal{A}}_t,S_t,{\mathrm{\Gamma }}_t)$ representing the respective grain scan. For the computations, we used H1 as we were also interested in featuring the temporal development of the corresponding size parameters in view of the relation of H1 and H2 elaborated in Section 3.

5.3. Evaluation

We apply a cross-validation to assess the quality of the model. To compute the accuracy at time $t_{\ell }\in \{t_0,\dots ,t_{14}\},$ we exclude the data for $t_{\ell }$ and fit the functions of $\mathcal{F}$ on the data for the remaining times $\{t_0,\dots ,t_{14}\}\setminus \left\{t_{\ell }\right\}.$ Then, by evaluating these functions at time $t_{\ell },$ we obtain the diagram parameters ${\mathcal{A}}_{\ell }^{\prime },S_{\ell }^{\prime },{\mathrm{\Gamma }}_{\ell }^{\prime },$ and hence the diagram ${\mathcal{P}}_{\ell }^{\prime }=\mathcal{P}({\mathcal{A}}_{\ell }^{\prime },S_{\ell }^{\prime },{\mathrm{\Gamma }}_{\ell }^{\prime }).$ Computing ${\mathrm{\Phi }}_G\left(\hat{C}\right({\mathcal{P}}_{\ell }^{\prime }\left)\right),$ the accuracy of ${\mathcal{P}}_{\ell }^{\prime }$ representing the 3D-grain scan G at $t_{\ell },$ gives us finally the accuracy (of the model) at time $t_{\ell }.$ In other words, we perform a leave-one-out cross validation using (in all combinations) one observation for validation and 14 observations for fitting.

The results for polynomials of degree 1, 2, and 3 are depicted in Figure 5 in yellow [light gray], orange [gray], and red [dark gray], respectively. Let us point out that the piecewise linear curves which connect the 15 data points for each colour are only drawn for a better perception of the trends over time.

Figure 5. Accuracy (within the cross-validation paradigm) of the interpolation/extrapolation model plotted for each $t_0,\dots ,t_{14}$ using polynomials of degree 1, 2, and 3, respectively for the fit.

We observe that the degree of the employed polynomials does not seem to have a strong impact on the fit. It would be interesting to link this observation to the ‘uniform setup’ in terms of the actions over time of the underlying grain growth process.

As an additional reference the blue [black] markers give an upper ‘orientation’ for the accuracy that one may expect from the model. In fact, the blue [black] markers depict the accuracy of H1 at each of the 15 time steps. The curves are dashed at both ends to signify that the model changes from interpolation to extrapolation at these extremes. Also, note that the earlier annealing times differ from the latter ones.

Figure 5 shows that the accuracy of our models is generally high. With quadratic polynomials, we obtain an initial accuracy of 74%, which increases over time to approximately 86% for later time steps. Note that this is not much below the accuracy obtained by H1 on the ground truth. The general trend that accuracy increases with time (and hence prediction becomes easier) for this data set is in perfect agreement of the findings in [6, 36], which show that over time the system approaches the self-similar regime.

Figure 6 shows a fixed 2D slice of the 3D polycrystal at various moments in time and the diagrams resulting from our model (in the cross-validation scheme described above).

Figure 6. Visual comparison of diagrams obtained from full data (available for $t_0,\dots ,t_{14}$) and the interpolation/extrapolation model that uses all data except for the respective moment in time. A fixed 2D slice through the 3D data set at times $t_0$, $t_5$, $t_9$ and $t_{13}$ (top row); the corresponding diagram representation obtained via H1 (2nd row) and, respectively, via the interpolation/extrapolation model that uses only data from the other moments in time and best-fitting quadratic polynomials (3rd row); symmetric difference between the images from the 2nd and 3rd row (bottom row).

A discussion of the question how much the fit can be further improved by optimizing over all three sets of parameters that define the anisotropic power diagrams can be found in [31].

5.4. Size variation

In the remaining part of this section, we study how the growth process affects the size parameters ${\mathrm{\Gamma }}_t$ (recall that they need not correlate with grain sizes and they should not, at this point, be confused with often similarly denoted interfacial energy terms; we have provided an interpretation in terms of scaling factors in Equation (1(1) ${\gamma }_i^{\left(1\right)}=-{\left(\frac{{\nu }_i}{\mathrm{v}\mathrm{o}\mathrm{l}\left({\mathbb{B}}_{A_i}^d\right)}\right)}^{\frac{2}{d}}=-{\alpha }_i^{\frac{2}{d}}.$(1) )). As Figure 5 shows, the accuracy of the model increases with the duration of the process. This may be seen as an indication that the employed diagram structures are energetically favourable. Complementing the analysis of Section 4, we will now empirically study how the size parameters $\mathrm{\Gamma }$ of the diagrams change as the grains grow.

Let us begin by recalling the interpretation for the similarities and differences of H1 and H2. In fact, based on the analysis of Section 4 we suspected that for persistent grains the variances of the powers ${\mathrm{\Gamma }}_t^{\frac{3}{2}}$ of the sizes chosen by H1 should be significantly smaller than for non-persistent grains (as the latter might not be energetically favourable). As Tables 2 and 3 show, these variances are indeed all small and significantly smaller for the persistent grains. This means, on the one hand, that the scaling factors ${\alpha }_i$ from Equation (1(1) ${\gamma }_i^{\left(1\right)}=-{\left(\frac{{\nu }_i}{\mathrm{v}\mathrm{o}\mathrm{l}\left({\mathbb{B}}_{A_i}^d\right)}\right)}^{\frac{2}{d}}=-{\alpha }_i^{\frac{2}{d}}.$(1) ) do not vary much, and we might as well have employed the simpler heuristic H2. On the other hand, in terms of the ${\alpha }_i,$ non-persistent grains are less ‘regular’ than persistent grains.

Table 2. Variance of ${\mathrm{\Gamma }}_{t_{\ell }}^{\frac{3}{2}}$ for $\ell \in \{0,1,\dots ,14\}$ of persistent grains.

	$t_{\ell +0}$	$t_{\ell +1}$	$t_{\ell +2}$	$t_{\ell +3}$	$t_{\ell +4}$
$\ell =0$	$6.409\cdot {10}^{-8}$	$6.214\cdot {10}^{-8}$	$5.130\cdot {10}^{-8}$	$3.016\cdot {10}^{-6}$	$2.544\cdot {10}^{-7}$
$\ell =5$	$6.719\cdot {10}^{-6}$	$5.839\cdot {10}^{-6}$	$1.057\cdot {10}^{-5}$	$7.025\cdot {10}^{-4}$	$9.815\cdot {10}^{-7}$
$\ell =10$	$2.286\cdot {10}^{-4}$	$7.425\cdot {10}^{-5}$	$5.749\cdot {10}^{-5}$	$1.185\cdot {10}^{-3}$	$1.823\cdot {10}^{-5}$

Table 3. Variance of ${\mathrm{\Gamma }}_{t_{\ell }}^{\frac{3}{2}}$ for $\ell \in \{0,1,\dots ,14\}$ of non-persistent grains.

	$t_{\ell +0}$	$t_{\ell +1}$	$t_{\ell +2}$	$t_{\ell +3}$	$t_{\ell +4}$
$\ell =0$	$7.570\cdot {10}^{-3}$	$1.040\cdot {10}^{-2}$	$4.915\cdot {10}^{-3}$	$1.974\cdot {10}^{-2}$	$2.993\cdot {10}^{-3}$
$\ell =5$	$1.278\cdot {10}^{-2}$	$1.834\cdot {10}^{-2}$	$2.790\cdot {10}^{-2}$	$4.158\cdot {10}^{-3}$	$8.959\cdot {10}^{-3}$
$\ell =10$	$2.608\cdot {10}^{-2}$	$2.804\cdot {10}^{-2}$	$2.955\cdot {10}^{-2}$	$3.541\cdot {10}^{-3}$	−

A closer look at Table 2 shows that, perhaps unexpectedly, the variance of ${\mathrm{\Gamma }}_t^{\frac{3}{2}}$ does not decrease monotonously in time, not even for persistent grains. One might have anticipated that over time the behaviour of the volumes of the unit balls ${\mathbb{B}}_{A_i\left(t\right)}^3$ with respect to $A_i\left(t\right)=\left(\mathrm{\Sigma }\right(G_i\left(t\right)){)}^{-1}$ become increasingly similar to that of the grain volumes ${\nu }_i\left(t\right)$. The actual behaviour might be explained by grain disappearance with neighbouring grains filling the newly available space such that it might require additional energy to reach a more regular mass distribution.

5.5. Grain disappearance

Next, we use a standard logistic regression model for predicting the probability $p_i\left(t\right)$ of the binary response variable which signifies whether grain i disappears (response value 1) at the next observed moment t + 1 or not (response value 0). As predictor variables we use the grain's volume and its size parameter as computed by H1. Hence, we estimate three parameters, factors ${\text{β}}_1$ and ${\text{β}}_2$ for each of the two variables and one parameter ${\text{β}}_0$ for the intercept, i.e. we use the regression model $p_i\left(t\right)=\frac{1}{1+{\mathrm{e}}^{-({\text{β}}_0+{\text{β}}_1{\nu }_i(t)+{\text{β}}_2{\gamma }_i(t\left)\right)}}.$Note that the parameters ${\text{β}}_0,{\text{β}}_1,$ and ${\text{β}}_2$ are independent of t. Hence, the model is identical for all time steps in the sense that it depends only on the grain characteristic. Our null hypothesis is that neither the size nor volume parameter correlate with grain disappearance. As a significance level we choose 0.05.

Table 4 depicts the computed coefficients and p-values. Let us, in addition, point out that the McFadden's Pseudo $R^2$ value is 0.322, indicating a good model fit. As it turns out, we can reject the null hypothesis at a significance level of 0.05; the small p-value for ${\text{β}}_2$ suggests that the diagram size parameter ${\mathrm{\Gamma }}_t$ is statistically significant in this model for predicting grain disappearance.

Table 4. Coefficients estimated by maximum likelihood estimation and their p-values for the described logistic regression model.

	Coefficient	p-value
${\text{β}}_0$ (intercept)	−4.0537	<0.0001
${\text{β}}_1$ (volume)	−0.0010	<0.0001
${\text{β}}_2$ (size parameter)	−1.0924	<0.0001

We believe it is promising to investigate this further, not least because it agrees with recent experimental evidence from [39, 40], which demonstrates that the volume parameter alone is not a good predictor of grain growth (or shrinkage).

6. Final remarks

As we have seen, measurements at discrete times of a grain growth process can be turned into a continuous description of the process using diagram representations. Our findings about the increasing fit of anisotropic diagrams, the variance of ${\mathrm{\Gamma }}_t^{\frac{3}{2}},$ and the relevance of the ${\mathrm{\Gamma }}_t^{\frac{3}{2}}$ for predicting grain disappearance indicate that, in the course of the growth process, the grain scans tend towards particularly ‘simple’ diagrams with essentially uniform size parameters.

Since all computations are fast, and can, if necessary, generally even be speeded up by using H2 rather than H1, we believe that the described model has the potential of becoming a practical tool for simulation and further studies of grain growth in practice. A next step is to evaluate the model more extensively on various types of grain scan data. It would be desirable if a consolidating foundation of the model can be achieved based on first principles.

Finally, let us point out that the model allows for some flexibility to incorporate additional insight. It is possible, for instance, to incorporate available characteristics of grain growth into the curve fitting process.

Acknowledgments

We thank Henning Friis Poulsen for fruitful discussions and valuable comments on an earlier version of the present paper. We are also grateful to Henning Friis Poulsen and Allan Lyckegaard for providing the real-world data set used in this paper.

Disclosure statement

No potential conflict of interest was reported by the author(s).