Abstract Most often, Digital Image Correlation (DIC) is used to analyze a sequence of images. Exploitation of an expected temporal regularity in the displacement fields can be used to enhance the performances of a DIC analysis, either in terms of spatial resolution, or in terms of uncertainty. A general theoretical framework is presented, tested on artificial and experimental image series. ∗

In the ﬁeld of solid mechanics, since the early 1980’s [1, 2, 3], most of the procedures dealing with digital image correlation (DIC) of 2D pictures (or 3D volumes) are based upon the registration of a pair of pictures, a ﬁrst one corresponding to the reference conﬁguration and the second one to the deformed conﬁguration [4]. They consist in subdividing the region of interest in the reference picture into a set of independent zones of interest (ZOIs). The latter ones are small interrogation windows that are registered, and may overlap since there is no spatial constraint on neighboring windows. The advantage is that the analysis of each ZOI can be run independently of the other ones. The drawback is that the continuity of initially contiguous ZOIs is not obtained because of measurement uncertainties. This lack of continuity is one of the criteria used to estimate the quality of the registration [6]. This type of approach is referred to as ‘local’ in the sense that the registration is performed with small ZOIs with no information exchange with neighboring ZOIs. The same type of approach was developed independently in the context of ﬂuid mechanics (i.e., particle image velocimetry or PIV [5]) starting in the late 1970’s [7, 8, 9]. Very early on, regularization techniques [10, 11] were proposed [12, 13] since the gray level conservation, i.e., the underlying conservation principle to determine the optical ﬂow, is very diﬃcult to achieve 2

when using particles in a 3D ﬂow observed with a single camera. Conversely, in solid mechanics, this hypothesis is easier to achieve when, say, spraying B/W paint onto the observed surface. One other aspect is related to the measurement uncertainties that require special attention in terms of gray level interpolation to achieve resolution on the decipixel or even centipixel range [14, 15, 16] to be applicable to situations dealing with small strains. One way of prescribing a priori, say, continuity requirements on the displacement ﬁeld, is to resort to global approaches of DIC. In that case, the registration is carried out over the whole region of interest (ROI). This type of analysis gives more freedom to the user in the choice of the displacement ﬁeld in comparison to local approaches where the type of kinematic ﬁeld is generally of low degree [17]. However, the computation time is higher. For instance, Fourier expansions of the displacement ﬁeld are considered. They allow for a very robust measurement of displacement ﬂuctuations [18], provided the non-periodic contribution is accounted for by using additional kinematic degrees of freedom [19]. Another alternative is to consider shape functions of ﬁnite element procedures [20, 21]. This type of kinematics is very useful when comparisons are sought between experimental data and numerical simulations. In that case, the mesh and the shape functions can be identical during the measurement and identiﬁcation or validation steps [22]. Further, it was shown that the continuity requirement leads to lower displacement uncertainties when compared to FFT-based local approaches [23] or with bilinear interpolations [35]. Up to now, only spatial features were discussed in terms of regularization procedures. The temporal axis may also be considered. It was shown that, for 3

very specialized kinematics, it is possible to extract in a robust way steady state 1D and 2D velocity ﬁelds [25], a spatiotemporal velocity ﬁeld when only one spatial component is considered. In the last case, the combined use of space and time kinematic bases allows the user to decrease the spatial resolution by increasing the temporal resolution to achieve similar uncertainty levels [26]. This property can be used, say, when monitoring high speed experiments for which the number of pictures can be greater than 100, or even 1,000, but the spatial deﬁnition is of the order of 10 to 100 kpixels as opposed to quasi-static experiments for which cameras have a high deﬁnition (typically greater than 10 Mpixels). More generally, with nowadays cameras, it becomes classical to have access to more than 100, 1,000 or even more pictures. Consequently, alternative techniques to incremental approaches may be considered, as was proposed when analyzing nonlinear problems [27]. Broggiato et al. [28] proposed to use 5 consecutive pictures for global DIC and a restricted (parabolic) time change of the displacement ﬁeld. This multi-frame procedure leads to a smoother time evolution of strain rates. It is used to have a more precise evaluation of the strain rates in the central frame (i.e., local in time). In the following developments, a global approach in space and time is proposed. It consists in using diﬀerent space-time decompositions of the displacement ﬁeld and in solving globally the minimization of the correlation residuals. In Section 2, the formulation of the measurement problem is presented and its practical implementation discussed. The latter is based on a ‘dyadic’ spatiotemporal decomposition. Particular features of the spatiotemporal analysis are addressed in Section 3. Among them, the initialization 4

of the spatiotemporal analysis, the link with digital volume correlation. A resolution analysis shows the beneﬁt to be expected from a spatiotemporal approach. The analysis of a tensile test on an aluminum alloy sample illustrates the fact that time can be transformed into a loading parameter, which may not be evenly distributed over the temporal axis. Last, a high speed experiment on a composite material is analyzed in Section 4.

2 2.1

Mathematical formulation General principle

Let us brieﬂy review the global DIC methodology, which is further extended to the spatiotemporal framework. Readers are referred to Refs. [21, 29] for further details. The starting point of the analysis is a series of images, f (x, t) where x is the vector position, of pixel coordinates (x, y), and t the time index. The ﬁrst image is acquired at time t0 . The basic assumption, called “gray level conservation,” is expressed as f (x + u(x, t), t) = f (x, t0 )

(1)

where u denotes the sought planar displacement ﬁeld. The strategy used to tackle this problem is to choose a space of trial displacement ﬁelds V, a basis of which is denoted by ψ(x), which is expected to contain a good approximation of the actual (i.e., to be determined) displacement ﬁeld u(x, t). Instantaneous DIC consists in writing v(x, t) =

∑ i

5

ai (t)ψi (x)

(2)

where ai (t) are the unknown amplitudes (or kinematic degrees of freedom) of the estimated displacement ﬁeld v. They are obtained as the argument that minimizes a global measure (over space for a given time) of the violation of gray level conservation u(x, t) = arg min ∥f (x + v(x, t), t) − f (x, t0 )∥2ROI v

(3)

evaluated over the whole Region of Interest (ROI). In this expression, the norm ∥•∥ROI is commonly chosen as its L2 variant, or the quadratic diﬀerence between the reference image f (x, t0 ) and the“corrected”one at time t summed over space. The proposed spatiotemporal analysis consists in generalizing the previous space regularization (2) to space and time. Although this methodology can straightforwardly be extended to arbitrary functions ψ(x, t), the following discussion is specialized to the particular case where space and time are separated, and hence the trial function will simply be the ‘dyadic’ product of space ψi (x) and time ϕj (t) decompositions, namely, trial displacement ﬁelds v are written as v(x, t) =

∑∑ i

aij ψi (x)ϕj (t)

(4)

j

As soon as the choice of a particular representation is made, the displacement ﬁeld is restricted to belong to a given space. However, this is common to all discretizations and not speciﬁc to the dyadic product form. As will be shown below, a dyadic form will be shown to be very convenient for the present implementation, and allows for faster execution speeds. Let us note that the instantaneous (i.e., image by image) analysis is a particular case of the above choice when ϕj (t) = δ(t − tj ), where tj stands for all instants of time 6

where images were captured in addition to the reference picture (t = t0 ). Amplitudes aij are to be determined from the same problem as that written in Equation (3) with the norm ∥ • ∥ROI×[t0 ,t1 ] deﬁned over the time interval [t0 , t1 ] corresponding to the series of pictures

2

( )

∑ ∑

αij ψi (x)ϕj (t), t − f (x, t0 ) aij = arg min f x + αij

i

j

(5)

ROI×[t0 ,t1 ]

where t1 is the ﬁnal time of analysis. The actual implementation of the resolution of the nonlinear problem (5) follows exactly the footsteps proposed for global DIC approaches [21]. However, the major steps are here recalled in order to highlight the beneﬁts to be gained from the speciﬁc choice of a dyadic partition between space and time decompositions.

≈ f (x, t) − f (x, t0 ) + v(x, t) · ∇f (x, t0 ) where the latter approximation (legitimate when the displacement correction is small) allows for the computation of the gradient term once for all iteration steps. Taking into account the speciﬁc form of the space-time decomposition (4), the following linear system in the unknown vector {a} that

The solution to this linear system provides a displacement ﬁeld v (1) that can be used to correct the image sequence f (1) (x, t) = f (x + v (1) , t). The latter is to be used in the computation of a corrected second member {b} to evaluate the correction of the displacement ﬁeld. This loop is iterated until convergence. The speciﬁcity of the space/time decomposition as a dyadic product is now exploited from the expression of matrix [M]. Let us introduce matrix [N] the components of which are the following integrals over the ROI ¨ [ψi (x) · ∇f (x, t0 )][ψk (x) · ∇f (x, t0 )] dx

Nik =

(10)

ROI

that is computed only from the reference image, then [M] is expressed as Mijkl = Nik Φjl with

ˆ

(11)

t1

Φjl =

ϕj (t)ϕl (t) dt t0

8

(12)

It is to be noted that [N] is precisely the matrix that would be used in a standard (i.e., instantaneous) global DIC analysis. Similarly the second member {c} reads ¨ [ψi (x) · ∇f (x, t0 )][f (x, t0 ) − f (x, t)] dx

ci (t) =

(13)

ROI

and thus

ˆ

t1

bij =

ci (t)ϕj (t) dt

(14)

t0

is a simple weighted sum of the instantaneous second member by the chosen time functions. In practice, gradients ∇f are evaluated as derivatives of the gray level interpolation function. Image corrections also need an interpolation scheme since sub-pixel resolutions are sought; spline interpolations were chosen. Alternative choices can be made as in local DIC approaches [15]. Exact gray level conservation is in general not strictly obeyed, and this principle is often relaxed in local DIC formulations, say, by resorting to zeromean normalized cross-correlations that amount to allowing for an aﬃne transformation of gray level without impact on the objective functional to be minimized [4, 17]. A similar relaxation can also be implemented here. It consists in rescaling the gray levels in each element so that mean and variance of both images over each element are transformed to 0 and 1, respectively. This transformation is to be performed on the reference image (prior to the computation of matrix [M] and vector {b}), and on the corrected deformed images after correction. Another alternative is to consider a ﬁeld of contrast and brightness that would vary from picture to picture, and would be corrected during the initialization procedure [29]. None of these procedures were 9

used and no correction was considered herein. Although among the global DIC approaches [30, 18] a ﬁnite-element based development is the most common [20, 21], other variants exist. They may exploit analytical elastic solution (e.g., a Brazilian test [23]) or numerical solutions [31]. The discussion was not limited herein to ﬁnite element discretizations because this reduction has no consequence in the temporal regularization. From the observed decoupling between space and time contributions, diﬀerent routes may be considered:

• First, a simple time decomposition using linear shape functions spanning over 2n images allows for a sparse matrix [M]. In that case, the shape functions are written as ϕj (t) = φ((t − tj )/∆t)

If nt denotes the number of ϕk functions, it is observed that each iteration step involves nt solutions (with the same kernel [M]) of a problem similar to an instantaneous computation. The exploitation of the time regularity means in practice that nt is much less than the number of images, and hence the computation load is much less an eﬀort than the sequential evaluation of each image. Note however that the most computationally-intensive part of the algorithm is rather the image correction than the solution to the linear systems. If each degree of freedom in the time domain is apparently decoupled from the others as far as matrix [M] is concerned, the entire problem is still to be treated simultaneously because of the correction of each image by the predetermined displacement ﬁeld. The freedom of selecting appropriate shape functions is still quite large. Fourier series is an appealing example, yet it is to be stressed that they are not well ﬁtted to non-periodic evolutions, and hence if the latter basis is chosen, it may be convenient to add a linear ramp function to Fourier modes to account for diﬀerent initial and ﬁnal states. Yet another example is given by Legendre polynomials.

In the following analyses, the spatiotemporal steps will use a spatial discretization made of 4-noded spatial elements with bilinear interpolation (as in a standard Q4-DIC analysis [21]). The surface of each element is equal to ℓ2 . The time discretization will use linear shape functions spanning over 2n pictures (i.e., two time increments) so that a spatiotemporal volume is

11

deﬁned, as would be a 3D volume element, as containing ℓ2 n voxels.

3

Discussion

3.1

Initialization

It is important to comment on one diﬃculty of the above displacement approach. It is well known that a major diﬃculty of DIC algorithms is to cope with large displacements or large strains. This diﬃculty is here of major importance as the time span of the analysis may be large, and hence so will be displacement amplitudes over the entire time interval. It is thus of importance to design a way to circumvent this diﬃculty that allows for an automated procedure, valid for a large class of image series. The way to address it is to design a proper initialization of the displacement ﬁeld so that the remaining details to be resolved deal with moderate or small displacement amplitudes. Although initialization is a critical feature, it does not have to be very precise. A coarse instantaneous DIC analysis of a sampling of images is performed with few degrees of freedom, and this displacement ﬁeld, projected onto the space V × T , is used to initialize the analysis. The “coarseness” of this initialization step is an issue that cannot be discussed in general as it will depend very much on the experiment to be analyzed. If, for instance, the origin of the large displacements comes from a global (i.e., rigid body) translation of the ﬁeld of view over the time interval to be studied, a very quick correlation procedure based on a single translation may be appropriate (i.e., two degrees of freedom per analyzed

12

image couple). This will concern the cases discussed in Section 3.4 and 4. On the contrary, if strains are very large, a multiscale approach based on a low-pass ﬁltering of images is recommended [25]. However, in that case, it may be unnecessary to carry out the pyramidal reﬁnement procedure down to the ﬁnest details. Thus, although more challenging, this case would still represent a small additional cost as compared to the entire analysis.

3.2

Connection with Digital Volume Correlation

It is to be noted that the above analysis could be considered as a degenerate case of global Digital Volume Correlation (DVC) [32], mapping the three dimensional space coordinate (x, y, z) onto (x, y, t), namely, the reference volume would here be f1 (x, t) = f (x, t0 ) (i.e., an extrusion in time of the ﬁrst image), and the deformed volume the entire stack of images f2 (x, t) = f (x, t). The 3D displacement ﬁeld has to be specialized to have a zero component along the time axis vt = 0. Hence, even if such a reduction may not provide the most eﬃcient implementation of the space-time approach, it may be obtained quite simply from a DVC code. When an unstructured mesh in space and time is appropriate, this reduction is an easy route to follow. The above approach aims at measuring the displacement ﬁeld from the reference image. Alternatively the same approach could be tuned to measure “velocity” ﬁelds. Velocity is here to be understood as the incremental displacement ﬁeld between two consecutive images, separated by a time interval ∆t. To carry on the analogy with a three-dimensional problem, the reference volume would be f1 (x, t) = f (x, t) (i.e., the original stack of images but the

13

last one) while the deformed volume would be f2 (x, t) = f (x, t + ∆t), (i.e., the original stack of images but the ﬁrst one). This approach was followed in Ref. [26] for the speciﬁc case of 1D space and 1D time discretizations.

3.3

Benefit from time regularization

One speciﬁc feature of DIC is that the uncertainty on the measured displacement varies with the number of kinematic degrees of freedom, and hence with element size for a ﬁnite-element based approach [21]. The lower the number of degrees of freedom, i.e., the larger the element size, the lower the uncertainty level. The rationale behind this observation is that the displacements are statistical quantities whose noise susceptibility is dictated by the available number of pixel per degree of freedom. Spatial regularization is a way to reduce the number of degrees of freedom by additional constraints [35]. Therefore, spatial regularization is expected to reduce the displacement uncertainty for a ﬁxed number of spatial degrees of freedom. Alternatively, for a ﬁxed uncertainty, this technique may oﬀer the means to enhance the spatial resolution, i.e., to have a ﬁner spatial mesh, if the supporting space-time element extends over longer time increments. In the limit of a zero displacement, this enhancement can be compared to what is called “super-resolution”, i.e., a way to achieve sub-pixel gray level deﬁnition from a collection of images of the same scene [33]. To illustrate this property, let us consider an artiﬁcial case for which the reference picture is shown in Figure 1. Its deﬁnition is 512 × 200 pixels, and its dynamic range is equal to 200 gray levels. A set of deformed pictures is

14

created by applying a uniform strain rate of 10−3 per picture in the horizontal direction, and 5 × 10−4 per picture in the vertical direction. A third degree gray level interpolation is used herein. The spatial mesh is composed of 16 × 16-pixel elements. Three image series are considered, namely, with 19, 49 and 99 deformed pictures in which an additional white Gaussian noise √ is added with a standard deviation 2σ = 1 gray level. The maximum strain in the longitudinal direction is equal to 2 %, 5 %, 10 % in the horizontal direction, and 1 %, 2.5 %, 5 % in the vertical direction. The corresponding maximum displacement is equal to 10, 25, 50 pixels in the horizontal direction, and 5, 12.5, 25 pixels in the vertical direction. Figure 2 shows the result in terms of mean squared error. In Figure 2(a) the results are compared when the deformed pictures are re-encoded as 8-bit integers and left as 64-bit ﬂoating-point numbers. There is a clear diﬀerence in terms of variance. However, the same trend is observed, namely, there is a quasi linear variation of σu2 with the number of temporal elements. This is also true for the second (Figure 2(b)) and third series (Figure 2(c)) for which ﬂoating-point numbers were used, irrespective of the maximum displacement and strain amplitudes. It is worth noting that the last point of each series corresponds to an instantaneous approach with the same spatial discretization as in the spatiotemporal approach. These results can be understood by performing the following resolution analysis [21]. In the present case, it is assumed that the measurement uncertainties are dominated by the eﬀect of acquisition in comparison with other uncertainties related, for instance, to gray level interpolations. It is assumed that the main cause for error is due to acquisition noise compared to gray 15

level interpolation. The noise-free reference being unknown, it is equivalent to considering a noise in the diﬀerence f (x, t0 ) − f (x, t). The picture diﬀerence has two terms η(x, t0 ) − η(x, t), namely, a ﬁrst one η(x, t0 ) of variance σ 2 , and a second one, which is dependent upon the current time t, η(x, t), with the same variance. The spatiotemporal components of η are chosen to be Gaussian and uncorrelated in space and time. Matrices [M] and [N] are assumed to be unaﬀected by this noise. On the contrary, vector {c} is modiﬁed by an increment δcm (t) = p

where p is the physical size of one pixel. This expression shows that the eﬀect of noise on the spatial part at ﬁxed time (i.e., τ = t ̸= t0 ) simply consists in considering a noise of variance 2σ 2 added to f (x, t) only [21]. √ This argument explains why a 2σ standard deviation was considered in the previous analysis. The increment of the second member vector {b} is given by δbij =

For large values of n, Φi = n, and Φij ∝ n, hence the ﬁrst term of Equation (25) is proportional to n0 while the second term is proportional to n−1 . ∑ −1 In addition, p Φjp Φp = 1 (independent of j apart at the ends). The √ ﬁrst matrix is essentially 1 for all j, l. Φ−1 = (1/n)P where P ≈ 3, jl ii jl Pi,i+1 ≈ −0.46, Pi,i+2 ≈ 0.12, . . . i.e., an alternating series whose amplitude dies oﬀ exponentially. Consequently, the ﬁrst matrix plays a major role. If the correlation length (i.e., the characteristic distance over which the normalized pair correlation function decays to 0) of analyzed the texture is signiﬁcantly larger than the pixel size, but not too large so that matrix [M] is invertible, the standard displacement resolution σu reads √ B pσ σu = √ A+ n ⟨(∇f )2 ⟩ℓ

(26)

where ⟨· · ·⟩ denotes the spatial average over the considered ROI, ℓ2 the number of pixels in the considered spatial element, n the number of pictures in each time increment, A ≈ 1 and B dimensionless constants dependent on the interpolation functions [21]. The result of Equation (26) when n = 1 is consistent with earlier results obtained with local [34, 35] and global [21, 35] 2D-DIC approaches. The square-root dependence with time (i.e., n) is a 17

direct consequence of the central limit theorem applied to image analysis for which averaging techniques are known to improve the signal to noise ratio [36]. It is worth noting that when the reference picture is noiseless, A vanishes and σu is then inversely proportional to the square root of n pσC σu = √ √ ⟨(∇f )2 ⟩ℓ n

(27)

where C is a dimensionless constant. A noiseless reference image may be more than a theoretical concept in the sense that prior to an experiment, it is possible to perform a large number of image acquisitions of the same static scene so that their average may provide such a reference image with a much smaller noise level than in the active part of the test where often time or material does not allow for such noise reduction strategies. In addition this procedure provides estimates of the noise level in the actual test conditions. This last result shows the beneﬁt to be expected by using a time regularization. For the same displacement uncertainty σu , it is possible to decrease the spatial resolution ℓ by considering more than one image pair to perform √ a correlation. When the product ℓ n is kept constant, the standard displacement uncertainty does not change. This is true provided the sought displacement ﬁeld belongs to the measurement basis. This is not always the case, as will be shown in Section 4 in which the reference picture cannot be assumed noiseless. Figure 3 summarizes the results observed for the three picture series when the root mean squared diﬀerence is plotted against the time increment n used in the temporal discretization. All three series collapse onto a single curve that is inversely proportional to the square root of 18

n, as expected from Equation (27). The highest level of the standard displacement resolution corresponds to an incremental analysis with the same spatial discretization. When the noiseless reference picture is unknown, the gain to be expected from the time regularization is less important (Equation 26) than previously. Figure 4 shows the results for the series of 19 pictures. The trend predicted by Equation 26) is found in the numerical results. The very small resolution levels are underestimated. Diﬀerent causes may explain this deviation. One of them is related to gray level interpolations. The ratio of the two constants A/B is found to be of the order of 1.4. Consequently, the gain of a spatiotemporal analysis with a large number of pictures per time increment is √ √ of the order of 35 % (i.e., ( A/B + 1 − 1)/ A/B + 1) in terms of standard displacement resolution when compared to an incremental approach.

3.4

Extension to time-like parameters

Up to now t was interpreted as a time parameter. However, other choices may be more appropriate. For instance, a specimen may be subjected to a mechanical test with one degree of freedom, characterized by a load, F , so that the solid remains in its linear elastic regime, and images are captured at random instants of time. In such a case, the time or image number sequence would be irrelevant. If F is used instead of t, a simple linear response is anticipated. The substitution is straightforward for the method, and an uneven sampling of diﬀerent loads does not involve any new diﬃculties. Such a change of variable may however open new doors, namely, the hyperelastic

19

behavior of an elastomer for instance could be tackled based on a polynomial expansion in F rather than simply linear. To illustrate this case, let us consider a tensile test on a 2024 aluminum alloy. A special setup, initially designed to test ceramics [37], was used to minimize spurious ﬂexure when the sample is tested in the elastic regime. A long distance microscope monitors a small part of one lateral surface (i.e., 4 mm2 ). An annular lighting was used to minimize reﬂectivity variations of the black and white paints that were sprayed onto the observed surface. A total of 22 pictures was taken when the applied load F varies between 218 and 2003 N. Figure 5 shows the reference picture (i.e., when F = 218 N). Its deﬁnition is 1008 × 1016 pixels with an 8-bit digitization. In the present case, the pre-correction is simply a rigid body translation for each picture, and the temporal interpolation is linear since the behavior of the material is assumed to be elastic. Consequently only one time increment is considered (i.e., n = 22). Spatially, a single Q4 element is considered to capture a very simple kinematics. The main displacement gradient is obtained by spatial diﬀerentiation of the shape functions. Two analyses are run. The ﬁrst one when the reference load is F = 218 N (i.e., the minimum load level), and a second one when F = 2003 N (i.e., the maximum load level). Figure 6 shows the convergence of the √ correlation residual (i.e., ∥f (x + v(x, t), t) − f (x, t0 )∥ROI×[t0 ,t1 ] /∆f , with ∆f = maxROI f (x, t0 )−minROI f (x, t0 )) and the maximum displacement correction as a function of the iteration number. The former tends very quickly to about 1.6 % of the dynamic range of the picture. This value provides an √ estimate of the acquisition noise 2σ = 0.016∆f (i.e., σ ≈ 3 gray levels). 20

The displacement corrections have an exponential decay with the number of iterations. The correlation was stopped when the value was less than 10−5 pixel. It required ten iterations of the spatiotemporal procedure. These two results indicate that the measurements are trustworthy in both cases with about the same order of conﬁdence. All the results reported herein are obtained with a Core i7 (2.9 GHz) computer. A total analysis of the 22 pictures with a 900 × 900-pixel ROI, i.e., a spatiotemporal volume of the order of 18 Mvoxels required less than 130 s of CPU time. Because a modiﬁed Newton scheme is used for all the analyses (initialization and spatiotemporal), a total of 173 (fast) iterations was needed for the initialization step. This is due to the fact that large displacements occurred in both directions (Figure 7). Conversely, the spatiotemporal analysis required only ten (slow) iterations because the strain levels remained modest, i.e., at most of the order of 10−3 when expressed in absolute values. This was possible thanks to the initialization step. Figure 7 shows the change of the mean displacement as a function of the applied load for both directions. There is a gradual change of the transverse displacement uy . In the longitudinal direction, the displacement ux is more erratic. This is due to the fact that the long distance microscope was manually repositioned during the experiment. It is especially true at the beginning of the experiment when the sample repositions itself thanks to the elastic joint of the setup. This result shows the need for a rigid body precorrection, picture by picture, before performing the spatiotemporal analysis. Had this pre-correction step not been implemented, low value of the correlation residuals could not have been reached and the spatiotemporal analysis 21

would have failed. Figure 8 shows the displacement ﬁelds in the both directions for the ﬁrst picture (F = 261 N) when the reference load is F = 218 N. It is worth noting that the dynamic range of the measured displacements is very small. There is a small rotation component of the same order of magnitude as the longitudinal and transverse strains. This is made possible thanks to the spatiotemporal analysis performed herein. By computing the mean eigen strains, it is possible to get Poisson’s ratio ν of the analyzed material. Two values are obtained, namely, ν = 0.36 when the reference picture is for F = 218 N, and ν = 0.30 when the reference picture is for F = 2003 N. Consequently, an estimate of Poisson’s ratio is 0.33 ± 0.03 in very good agreement with strain gauge measurements [38].

4

Experimental image sequence

The last example treated herein corresponds to an image sequence of a tensile test on a T700 carbon ﬁber / M21 epoxy matrix composite with plies oriented at ±45◦ with respect to the loading direction. A high speed camera was used to monitor the experiment with a frame rate of 65,000 fps. Consequently, the picture deﬁnition is reduced to 224 × 449 pixels (Figure 9). This is typical of high speed experiments for which the full deﬁnition is not reachable when the frame rate increases. This is a case where the spatiotemporal analysis proposed herein is useful to have a small spatial resolution by enlarging the temporal resolution since 155 pictures are available before the ﬁnal failure of the composite. However, after picture no. 100 the damage state of the ma22

terial becomes very important (Figure 9(b)). The analysis will be restricted to deformed pictures no. 1 to 99 in the sequel. The ROI is here composed of 336 × 175 pixels.

4.1

Resolution analysis

A ﬁrst analysis indicates that from deformed picture no. 1 to 19 the displacements and strains remain vanishingly small. Consequently, a spatiotemporal analysis can be performed to assess the resolution of the technique for different spatial and temporal discretizations. Figure 10 shows the change of the standard displacement resolution as a function of the volume of the spatiotemporal elements. The resolution is estimated from the values of the nodal displacements. Contrary to the artiﬁcial case, the standard resolution saturates in the centipixel range. This is presumably due to the digitization of the pictures that leads to a lower bound of the displacement resolution. When the volume of the spatiotemporal elements is recast in terms of an equivalent length (i.e., the cube root of the volume), it is observed that the standard displacement resolution is less than 4 centipixels for element lengths greater than or equal to 5 pixels. The trend predicted by Equation (26) is obtained for the data before saturation. A good agreement is observed. The ratio of the interpolation coeﬃcients A and B deﬁned in Equation (26) is of the order of 0.2, thereby indicating a beneﬁt to be expected from the spatiotemporal analysis for large time increments of the order of 10 % (i.e., √ √ ( A/B + 1 − 1)/ A/B + 1) for small spatial elements, see Figure 10. This situation will precisely be encountered in the following analysis in which very

23

small spatial elements are considered.

4.2

Analysis of the experiment

The experiment reported herein is deemed diﬃcult because cracks appear during the analyzed sequence (Figure 9(b)); the hypothesis of displacement continuity in the spatial and temporal domains will be violated at the end of the sequence. Consequently, the spatial and temporal discretizations will be adapted to capture these events. A set of analyses is run with diﬀerent spatial and temporal discretizations. Figure 11 shows the change of the correlation with the volume of the spatiotemporal elements. It is observed that the ﬁnest spatial discretization is needed to get the smallest correlation residual. It corresponds to spatial elements of equivalent length less than 5.4 pixels. This length can be used since the temporal axis is composed of 5 increments of 18 pictures each. In the present case, the temporal discretization does not have to be very ﬁne. The correlation residual level is less than 3 %, which is higher than those observed in the resolution analysis (i.e., less than 1 %), but still low enough to study the measured results. This discretization is chosen for the remainder of the section. With these parameters, the spatiotemporal volume of interest is of the order of 5.3 Mvoxels. In the present case, the initialization step was very fast and required only a few (fast) iterations for the 90 analyzed pictures. On the other hand, 20 spatiotemporal iterations were needed because of the complexity of the analyzed kinematics due to the presence of cracks. The total CPU time was less than 70 s for the whole analysis. Further, the

24

54 analyses with various spatiotemporal discretizations shown in Figure 11 were obtained in less than two hours. This corresponds to a cumulative spatiotemporal volume of 286 Mvoxels. Figure 12(a) shows the change of the mean displacement in both directions. The longitudinal component has a steady change with the picture number. This is not the case for the transverse displacement that starts to ﬂuctuate, and becomes steadier from picture no. 60 on. This regime is due to the loading conditions prescribed by the high speed tensile testing machine. Again, the picture by picture displacement pre-correction is needed to enable for the spatiotemporal analysis. In Figure 12(b), the mean strain components are plotted. The shear strains remain small compared to the two normal components throughout the analysis. The mean normal strain components have a smooth history, even though rigid body motions ﬂuctuate. There are two regimes, namely a ﬁrst one prior to picture no. 50, presumably the elastic domain (e.g., when the ﬁrst 11 pictures are used, the mean ratio of the transverse to longitudinal strains is equal to -0.78 which is in complete agreement with the values measured by strain gauges), and a second one during which the strain rate increases. This trend is due to nonlinear phenomena such as multiple cracking (Figure 9(b)) as will be proven in the sequel. Since the strain gauges are not mounted on the same face, the fact that the transverse signal departs from the value evaluated from the spatiotemporal approach may be caused by the heterogeneity of the damage ﬁeld. From this plot, it is deduced that the mean longitudinal strain rate of the experiment is 43 s−1 . Further, the apparent Poisson’s ratio has small ﬂuctuations around a mean value of 0.82 during the whole experiment. These 25

values are in very good agreement with known levels observed for this type of material [39]. The global residuals were observed to be of the order of 3 % for the analyzed sequence. However, the mean value per picture is not constant. Figure 13 shows its change with the picture number. The ﬁrst part corresponds to displacements that are vanishingly small. Part of these pictures were analyzed when the resolution was assessed. There is a second regime during which the mean correlation residual is of the order of 2.4 %. It corresponds to a situation where no cracks are observed. The reason for the increase of the mean correlation residual in the third regime (i.e., after picture no. 60) is due to cracking. This last conclusion can be drawn from the analysis of the normalized residual ﬁelds shown in Figure 14. Cracks are clearly visible for the residual ﬁeld corresponding to picture no. 53. Multiple cracking zones can be observed in the residual of picture no. 71. Two fully damaged zones can be seen in the residual ﬁeld of picture no. 99. This analysis allows one to understand the gradual increase of the mean correlation per picture shown in Figure 13. For the same pictures, the two components of the displacement ﬁeld are shown in Figure 15. Multiple cracks are observed in the maps corresponding to picture no. 53. Their number increases with the strain level and a lot of them are present in the maps corresponding to picture no. 99. Their exact numbering is not addressed herein. The fact that numerous cracks appear explains why the mean correlation residuals gradually increase (Figure 13). This is also the reason for requiring a ﬁne spatial discretization.

26

5

Conclusions

A spatiotemporal analysis was proposed to analyze a picture series with a global digital image correlation technique. By extending to the temporal axis the spatial regularization provided by a global approach (e.g., Q4-DIC [21]), it is shown that a series of pictures can be correlated by directly considering the whole sequence as opposed to an incremental approach generally performed. It is worth remembering that the spatiotemporal basis may be completed by a basis of rigid body motions that vary from picture to picture. This is important when analyzing experiments during which the ﬁrst part of the loading leads to a repositioning of the sample in the testing device. This was the case of the two experiments reported herein. A resolution analysis shows that the standard displacement uncertainty is inversely proportional to the square root of the spatiotemporal element size. This result shows the beneﬁt to be expected from this type of approach. When the picture deﬁnition is small, as often happens when high speed cameras are used, it is possible to have a small spatial resolution (i.e., the element size in the spatial domain) by enlarging the temporal resolution (i.e., the time increment) to achieve the same overall measurement uncertainty. This result was validated against an artiﬁcial and a real case. The previous result can also be used to analyze the elastic response of a material for which a linear relationship between the measured displacement and the applied load can be enforced during the measurement stage. The time axis then corresponds to the load axis with no change in the algorithm, even when the load data are not evenly distributed along the temporal axis.

27

It was possible to extract the Poisson’s ratio of the tested aluminum alloy whose mean value was in very good agreement with a value obtained by using strain gauge data. The analysis of a dynamic tensile test on a composite material was also performed. Even though cracks appear during the applied load, the spatiotemporal analysis could be performed. The level of correlation residuals was low enough to have conﬁdence in the measurement results. The convergence is slower because of discontinuities that appear at the middle of the picture sequence. The spatial and temporal discretizations had to be adapted to capture the main displacement and strain features. This was possible by analyzing the correlation residuals. The maximum mean longitudinal strain was found to be of the order of 5 % for 90 analyzed pictures. The inherent ﬂexibility of global spatial approaches to DIC is to be expected from the spatiotemporal approach as well. In the present paper, bilinear spatial and linear temporal shape functions were considered. This is not a limitation. If other discretizations are to be considered, it is possible to implement them. This is in particular true for the temporal discretization that may tailored to the analyzed case as was performed with a uniform discretization in the analysis of the dynamic tensile test of the composite material. It is a straightforward perspective to implement automatic mesh generation along the temporal axis. The ﬁnite-element based approach developed herein allows for a direct coupling with FE simulations. The measured displacements on the boundary of the spatiotemporal ROI could be prescribed to an explicit ﬁnite element analysis of the dynamic experiment. In that case, higher degree time inter28

polations may be considered. An alternative route is to measure the velocity ﬁeld as was initially proposed [26] to analyze a high speed movie.

Acknowledgments This work was funded by Agence Nationale de la Recherche under the grant ANR-2006-MAPR-0022-01 (VULCOMP Phase 1 Project).

Reference picture (a) and 99th deformed picture (b) of the tensile test on a composite material. Cracks are clearly visible on the lower edge of the sample in its deformed conﬁguration.

10

47

Resolution analysis of the tensile test on a composite material. Log-log plot of the standard displacement resolution as a function of the volume ℓ2 n of the spatiotemporal elements. The interpolation (solid line) described by Equation (26) is used only for the two series that do not saturate. . . . . . . . . . . 48

Figure 2: Resolution analysis for diﬀerent numbers of deformed pictures. 40 of the number of temporal elements. Displacement variance σu2 as a function A quasi linear trend is observed for all analyzed cases.

−2

σu (pixel)

10

data √ B/ n

−3

10

−4

10

0

10

1

10 time increment, n (number of pictures)

2

10

Figure 3: Standard Displacement resolution σu as a function of the time increment n for the three series analyzed in Figure 2. The solid line corresponds to the trend predicted by Equation (27).

41

−1

standard displacement resolution (pixel)

10

Data p α/ℓ2 + β/ℓ2 n

−2

10

−3

10

−4

10

2

3

4

5

10 10 10 10 10 volume of spatiotemporal element (voxels)

6

Figure 4: Standard Displacement resolution σu as a function of the spatiotemporal volume ℓ2 n for a series 19 pictures. The solid line corresponds to the trend predicted by Equation (26).

Figure 6: Convergence study. Correlation residual (a) and maximum displacement correction (b) as functions of the iteration number. The iterations were stopped when the 10−5 -pixel limit was reached by the maximum displacement correction.

Figure 9: Reference picture (a) and 99th deformed picture (b) of the tensile test on a composite material. Cracks are clearly visible on the lower edge of the sample in its deformed conﬁguration.

47

standard displacement resolution (pixel)

0.05

Data p α/ℓ2 + β/ℓ2 n

0.04

0.03

0.02

0.01 1 2 3 4 5 6 10 10 10 10 10 10 volume of spatiotemporal element (voxels) Figure 10: Resolution analysis of the tensile test on a composite material. Log-log plot of the standard displacement resolution as a function of the volume ℓ2 n of the spatiotemporal elements. The interpolation (solid line) described by Equation (26) is used only for the two series that do not saturate.

48

0.05

correlation residual

0.045

0.04

0.035

0.03

0.025 2 10

4

6

8

10 10 10 volume of spatiotemporal element (voxels)

Figure 11: Correlation residual as a function of the volume ℓ2 n of spatiotemporal elements when 90 pictures are analyzed. A minimum value is observed when n = 5.