- Email: [email protected]
Improved Hermite finite element smoothing method for full-field strain measurement over arbitrary region of interest in digital image correlation
Optics and Lasers in Engineering 50 (2012) 1662–1671
Contents lists available at SciVerse ScienceDirect
Optics and Lasers in Engineering journal homepage: www.elsevier.com/locate/optlaseng
Improved Hermite finite element smoothing method for full-field strain measurement over arbitrary region of interest in digital image correlation Jia-qing Zhao a, Pan Zeng a,n, Bin Pan b, Li-ping Lei a, Hong-fei Du a, Wen-bin He a, Yang Liu a, Yue-jie Xu a a b
Key Laboratory for Advanced Materials Processing Technology of MOE, Department of Mechanical Engineering, Tsinghua University, Beijing 100084, China Institute of Solid Mechanics, Beihang University, Beijing 100191, China
a r t i c l e i n f o
a b s t r a c t
Article history: Received 13 February 2012 Received in revised form 18 April 2012 Accepted 24 April 2012 Available online 30 June 2012
We investigate the finite element method and the Tikhonov regularization for smoothing noisy displacement field obtained from the digital image correlation (DIC). Among various finite elements, the C2 continuous Hermite finite element has been proven to be an effective high-order element which ensures both continuity and smoothness of the final strain field. However, the standard Hermite element is required to be rectangular, which cannot tackle the region of interest of arbitrary shape. To overcome this constraint, we propose the improved Hermite finite element smoothing method (IHFESM) to extend the generality of Hermite finite element that is able to measure the arbitrary region of interest. The IHFESM is implemented by eliminating the contribution to the global regularization matrix from the invalid region. The proposed method could be employed to smooth the noisy displacement field obtained by DIC, and reconstruct the reliable strain field by differentiation. The effectiveness of IHFESM is verified by two tension experiments of plates with three and four holes. & 2012 Elsevier Ltd. All rights reserved.
Keywords: Digital image correlation Displacement smooth Hermite finite element Tikhonov regularization Generalized cross-validation
1. Introduction Digital image correlation (DIC) [1,2] has become an effective and widely used non-contact full-field strain measurement method in the field of experimental mechanics. Ever since DIC was proposed in 1980s, research on measurement accuracy of the strain has always been a hot issue in this field [3–7]. Due to various sources of unavoidable noise [8] (such as shot noise and thermal noise) during the imaging and transmission process, the digital images of the test specimen are always contaminated to certain extend. In addition, the gray values stored in the images are digitized discrete values of continuous natural light intensities, thus the final correlation results by DIC are always full of noise. On one hand, the image noise could be alleviated by employing the high-quality equipment, such as the CCD/CMOS camera of higher signal-noise ratio and higher gray-level resolution [9], however this will increase the hardware cost greatly. On the other hand, many researchers have concentrated on exploring the effective displacement and the strain smoothing methods to filter noise in the measured data, and these methods could be classified into two main catalogs: The local smooth method: in the local smooth method, only small part of the dense displacement data located in one calculation window will be considered at one time, and the selected data
n
Correseponding author. Tel.: þ 86 10 6279 4261; fax: þ86 10 6277 3637. E-mail address: [email protected] (P. Zeng).
0143-8166/$ - see front matter & 2012 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.optlaseng.2012.04.008
will be fitted by the given basis function (polynomial function for instance) and the resultant coefficients will be used to compute the smoothed strain of the center point of the calculation window. Two typical methods fall into this catalog are pointwise least squares [5,10] and diffuse approximation [11]. The advantages of local smooth methods are: (a) low computation and programming complexity. Generally, these methods could be realized with few lines of programming codes, and they often run faster than global methods due to low computation complexity and low data density being considered at each time, and (b) high adaptivity to the shape of region of interest. When using these methods, one could simply omit the invalid data in the calculation window to tackle the case of nonrectangular region of interest, especially on the edge of the region. However, the performance of the local smooth methods heavily relies on the proper size selection of the calculation window [5] and the form of the basis function [11], and this deficiency is highlighted in the inhomogeneous deformation measurement. The global smooth method: for this method, all the noisy displacements are considered as a whole and they are fitted by one globally defined basis function or pieces of interconnecting basis functions. After the smoothing of displacements, the strains are obtained by differentiation. The thin-plate spline smooth method [12,13] and the finite element based smooth method [3,6,11,14] are the two widely used approaches. The finite element based smooth methods could be further classified into two types: (a) the finite element based least squares method [11,14] and (b) the finite element based Tikhonov regularization
J.-q. Zhao et al. / Optics and Lasers in Engineering 50 (2012) 1662–1671
[6,15]. Tikhonov regularization is more effective than the least squares method in that it considers both the ‘‘closeness’’ and ‘‘smoothness’’ of the best fitting surface, and the generalized cross-validation (GCV) [6] is often chosen to determine the regularization parameter. Here, we focus on the finite element based Tikhonov regularization method. In terms of the element type, the existing highorder elements, such as 8-noded isoparametric element [14] and o the triangular Cowper element with 18 % of freedoms (also denoted as Cowper-18 element) [15] are frequently used to tackle the arbitrary region of interest, but they are still unable to reconstruct a continuous and smooth strain field. For example, the strain fields obtained by using Cowper-18 element are continuous but not smooth. More recently, Sadati et al. [16] proposed the Tikhonov regularization based on another higherorder element, the conforming C2 continuous (second order continuous) Hermite finite with quintic Hermite polynomials, which successfully guarantees the continuity and smoothness of the smoothed data and its gradient simultaneously. They also showed that the C2 continuous Hermite elements coupled with the minimization of a norm of the data’s third derivative (denoted as C2-R3 in [16]) is superior to other combinations, such as C0-LSQ (C0 continuous element based least squares method), C0-R2, C1-R2 and C2-R2. However, due to the nature of Hermite interpolation, each C2 continuous Hermite finite is required to be rectangular, which causes difficulty in applying Hermite element to the measurement of the arbitrarily shaped region of interest (ROI) in DIC. In this paper, we propose an improved Hermite finite element smoothing method (IHFESM) to tackle the arbitrary ROI (nonrectangular ROI for example), which is common in DIC. IHFESM is based on the C2 continuous Hermite finite element and the Tikhonov regularization, in which the minimization of norm of the data’s third derivative used in [16] is adopted. The difference between IHFESM and the method in [16] is that, IHFESM tackles the arbitrary ROI by omitting the contribution from the invalid region when calculating the global regularization matrix. The proposed IHFESM is simple and effective in smoothing the noisy displacement field for inhomogeneous strain field calculation in DIC, and the effectiveness of IHFESM is verified by two plate tension experiments in Section 4. The remainder of this paper is organized as follows. In Section 2, the basic principle of DIC, C2 continuous Hermite finite element, Tikhonov regularization and generalized cross-validation are presented. Section 3 demonstrates the details of IHFESM. The performance of IHFESM is investigated in Section 4. Finally, the discussions and conclusions are in Section 5.
2. Basic principles 2.1. Digital image correlation DIC [1,2] is a well-established non-contact optical full-field deformation measurement technique, which computes the displacement and the strain of points of interest (POI) in the predefined region of interest (ROI) on the speckled surface (see Fig.1). To calculate all desired deformation variables, the gray values of pixels in the deformed subset are compared with that of pixels in the reference subset. To describe the comparison more precisely, the zero-mean normalized cross-correlation criteria (ZNCC) [17] is widely used, which is defined as PM PM 0 0 x ¼ M y ¼ M ½f ðx,yÞf m ½gðx ,y Þg m qP ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ffi CðpÞ ¼ 1 qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi PM PM PM M 2 2 0 0 x ¼ M y ¼ M ½f ðx,yÞf m x ¼ M y ¼ M ½gðx ,y Þg m
ð1Þ
1663
x
ROI
Step size f(x, y) (ux ,u y,
y POI
(x, y)
vx ,v
y)
g(x', y')
v u Reference subset
(x', y')
Deformed subset
Fig. 1. Definition of the region of interest (ROI), point of interest (POI), deformation parameters, reference and deformed subset in DIC.
where f(x,y) represents the gray value of point (x,y) in the reference 0 0 0 0 subset and g(x ,y ) is the gray value of corresponding point (x ,y ) in the deformed subset, all the coordinates are depicted in coordinate system anchored at the center of the reference subset with size PM PM (2Mþ1) (2Mþ1) pixels, f m ¼ 1=ð2M þ1Þ2 x ¼ M y ¼ M f ðx,yÞ P PM 0 0 gðx ,y Þ are the average gray and g m ¼ 1=ð2M þ 1Þ2 M x ¼ M y ¼ M values of all points in two subsets, and p is the deformation vector which describes the relationship between coordinate (x,y) and 0 0 (x ,y ). One merit of ZNCC criteria is that it is insensitive to the scale and offset changes in illumination lighting fluctuations, which often occurs in the practical experiment. Generally, when the subset is small enough, (x,y) is mapped to 0 0 (x ,y ) by first-order shape function as @u @u xþ y @x @y @v @v y0 ¼ y þv þ x þ y @x @y
x0 ¼ x þu þ
ð2Þ
where u, v are the two displacement components of the reference subset center on the x and y direction, @u/@x, @u/@y, @v/@xand @v/@y are the four displacement gradient components, thus p ¼[u, v, @u/ @x, @u/@y, @v/@x, @v/@y]T is the corresponding parameter vector. The first-order shape function allows translation, rotation, shear, normal strains and their combinations of subset, and it provides all the necessary deformation information for the common measurement. 0 0 The coordinate (x ,y ) often locates at the sub-pixel position in the deformed image. Therefore, interpolation is necessary to be introduced to reconstruct their gray values from neighboring pixels. In this work, the bi-cubic spline interpolation is used to 0 0 determine the gray values at (x ,y ) by gðx,yÞ ¼
3 3 X X
am,n xm ,yn
ð3Þ
m¼0n¼0
where the coefficients am,n could be calculated to guarantee C2 continuity by using the bi-cubic spline functions [1]. The six deformation parameters could be resolved by minimizing the nonlinear function C(p) in Eq. (1) using the Newton–Raphson (NR) [5] or quasi-Newton (qN) [18] method. The iteration equation of NR is p ¼ p0 rrCðp0 Þ1 rCðp0 Þ
ð4Þ
where p is the approximation value after iteration, p0 denotes an initial guess of solution that is ‘‘sufficiently close’’ to the real solution,
1664
J.-q. Zhao et al. / Optics and Lasers in Engineering 50 (2012) 1662–1671
rC(p) and rrC(p) are gradient and Hessian matrix of function C(p). In this study, NR is chosen as the sub-pixel registration algorithm. The displacements and displacement gradients calculated by the NR are usually contaminated by various noises [8]. Especially when the displacements and the strains are small, the noise might be typically of the same order of the real data. Therefore a smooth procedure is needed to extract the true displacement and strain from the noisy correlation results. In this work, we use the improved C2 continuous Hermite finite element to build a parameterized displacement field, and then employ the Tikhonov regularization to approximate the true noiseless field. In the following, the basic principles of the Hermite finite element, Tikhonov regularization and generalized cross-validation will be introduced.
2.2. Hermite finite element
The six quintic Hermite polynomial shape functions corresponding to the six degrees of freedom of 1D element are 0
2
uðxÞ ¼ Ne ðxÞqe
ð5Þ
e
e
where q is of size d 1, N (x) is a row vector of size 1 d, and d is the number of nodal degrees of freedom. In Eq. (5), qe, Ne(x) and d are determined by the specific form of the interpolation function. To ensure higher displacement continuity than traditional elements (such as 4-noded or 8-noded isoparametric element [14] and triangular Cowper-18 element [15]), we use the C2 continuous Hermite interpolation function [16] for Eq. (5), and this element guarantees the continuity and smoothness of both displacement field and strain field. The definition of node numbering, nodal degrees of freedom and local coordinates for one-dimensional (1D) and two-dimensional (2D) Hermite elements are shown in Fig. 2. For 1D C2 continuous Hermite element, three degrees of freedom, ui, @ui/ @x, @2ui/@x2, i¼1, 2 are associated with each node. The connecting elements have these three common degrees of freedom; therefore it ensures continuity of slope and higher derivatives up to the second order between the elements.
2
N 2 ¼ ðx þ 1Þ3 ð3x 9x þ8Þ=16
1 N1 2 N1
1 N2 2 N2
¼ ðx1Þ3 ð3x þ 5Þðx þ 1Þ=16, 3
2
¼ ðx1Þ ðx þ 1Þ =16,
¼ ðx þ 1Þ3 ð3x5Þðx1Þ=16 ¼ ðx þ 1Þ3 ðx1Þ2 =16 ð6Þ
where the subscript indicates the node number and the superscript denotes the order of the derivative with respect to the local element coordinate xA[ 1, 1]. To ensure the C2 continuity of 2D rectangular Hermite element, the nine degrees of freedom at each node are needed: u,@u/@x,@u/@Z,@2u/@x@Z,@2u/@x2,@2u/@Z2,@3u/ @x2@Z,@3u/@x@Z2,@4u/@x2@Z2. Then, all the thirty six biquintic Hermite shape functions are obtained as products of above 1D shape functions [16] j
In the finite element method (FEM) [19], given the vector of nodal degree of freedom values qe and the shape function Ne(x), the displacement u(x) of point at xAR2 inside of one element could be expressed as
0
N1 ¼ ðx1Þ3 ð3x þ9x þ8Þ=16,
k
ðx, ZÞ ¼ Ni0 ðxÞNi00 ðZÞ, Nðj,kÞ i
i ¼ 1,. . .,4,
i0 , i00 ¼ 1,2,
j,k ¼ 0,1,2
ð7Þ
ðx, ZÞ is associated with the degree of freedom (@j þ ku/ where N ðj,kÞ i 0 00 @xj@yk)i of node i of the rectangular element, i and i refer to the node number in the corresponding 1D element and they take values of 1 or 2 if the local coordinate of node i equals 1 or 1. For 2 2 instance, Nð0,2Þ 4 2 ðx, ZÞ is associated with @ u/@y and it is the product 0 of N1 ðxÞ N2 ðZÞ. In fact, Eq. (7) assumes that all the Hermite elements in the ROI are of the same size of 2 2, while for the non-uniform rectangular elements, the 2D shape functions are of form: j
k
k
Nðj,kÞ ðx, ZÞ ¼ Ni0 ðxÞNi00 ðZÞwj h , i
i ¼ 1,. . .,4,
i0 , i00 ¼ 1,2
ð8Þ
where w and h are the half width and half height of the element respectively. In the context of C2 continuous Hermite element, qe in Eq. (5) is of size 36 1, and Ne(x) is of size 1 36, whose components are determined by Eq. (8). It should be noted that, x in Ne(x) is expressed in the global physical coordinate system of ROI, and it should be transformed into (x,Z) in the local Hermite element coordinate system before using Eq. (8). The assembly of Eq. (5) over all E elements in the ROI gives uc ¼
E X
N e ðxÞqe ¼ Nq
ð9Þ
e¼1
u
u
u
u
u
u
where q is the global parameter vector of size D 1 and D is the degrees of freedom of all elements over whole ROI, N is a matrix of size n D and n is the number of points of interest in the ROI, N is also called the global stiffness matrix, uc is the calculated displacements at the location of corresponding point of interest, P and it is a column vector of size n 1. Here ðdÞ is a symbolic denotation of the assembly process in FEM [19]. 2.3. Tikhonov regularization The unknown parameter vector q in Eq. (9) could be obtained as the general least squares solution [11,14] by minimizing: 2
CðqÞ ¼ 99uc u992 n
ui
ui
ui
ui
ui
ui
ui
ui
ui
i
Fig. 2. Definition of node numbering, nodal degrees of freedom and local coordinates for 1D and 2D Hermite elements.
ð10Þ
where uAR is a column vector collecting all the noisy displace2 ment values and 99d992 refers to the sum over all squared components. In fact, Eq. (10) only depicts the ‘‘closeness’’ of smoothed displacement field uc to the measured field u. The outcome of above least squares depends strongly on the unavoidable noise in image acquisition and correlation in DIC. However, this problem could be conquered by using the Tikhonov regularization method which enforces an additional restriction of ‘‘smoothness’’ for the
J.-q. Zhao et al. / Optics and Lasers in Engineering 50 (2012) 1662–1671
solution [6,13,16]. The smoothness restriction is often interpreted as specific order of derivative of the smoothed surface f(x), and it is then added to Eq. (10) as the roughness penalty term by Z 2 2 CðqÞ ¼ 99uc u992 þ l 99f ðkÞ uc ðxÞ992 dO ð11Þ O
where l is the regularization parameter, O is the domain of all FEM elements, f(k) is the kth derivative of f, and f(k)uc(x) defines the smoothness within the element. There are various choices of parameter k in f(k)uc(x), and k¼2 (i.e. the second derivative) is often used in the relevant studies [6,13]. However, as pointed in [16], one’s priori information about the final solution would provide a useful guideline for the selection of k. More precisely, if one is concerned with the smoothness of reconstructed data, the second derivative is a proper choice; while if the smoothness of reconstructed data gradient is favored, the selection of the third derivative will be more reasonable. This criterion is also validated by the comprehensive investigations in [16]. In this work, we are more concerned about the strain (i.e. displacement gradient) measurement by DIC, therefore the third derivative, i.e. k¼3, is actually adopted. 2 ð3Þ The specific form of 99f uc ðxÞ992 is !2 !2 !2 !2 @3 u @3 u @3 u @3 u 2 ð3Þ 99f uc ðxÞ992 ¼ þ 3 þ 3 þ ð12Þ @x3 @x2 @y @x@y2 @y3 Substituting Eq. (5) into Eq. (12), the whole roughness penalty term could be written as Z Z E X 2 ð3Þ 99f uc ðxÞ992 dO ¼ ðqe ÞT ðN e,xxx ÞT N e,xxx þ 3ðN e,xxy ÞT N e,xxy O
e
O e¼1 þ 3ðN e,xyy ÞT N e,xyy þ ðN e,yyy ÞT N e,yyy dOqe
ð13Þ
N e,kkk ðk ¼ x,yÞ
Here denotes the row vector of third partial derivatives of the shape function Ne of element e with respect to the global physical coordinate x or y. Eq. (13) could be rewritten into a compact form as Z 2 ð3Þ 99f uc ðxÞðxÞ992 dO ¼ qT Sq ð14Þ O
where S is the D D global regularization matrix. Substituting Eq. (14) into Eq. (11), the solution of Eq. (11) will be qðlÞ ¼ ðN T N þ lSÞ1 N T u
ð15Þ
Then the smoothed displacement field and two displacement gradient fields could be obtained using us ¼ NqðlÞ,
usx ¼ N ,x qðlÞ,
usy ¼ N ,y qðlÞ
3. Improved Hermite finite element smoothing method (IHFESM) Above C2 continuous rectangular Hermite element based Tikhonov regularization is effective in smoothing noisy data, which ensures the smoothness of the measured data and the data gradient. However, in the practical measurement by DIC, the user-defined ROI are often of arbitrary shape, for example, the arbitrary polygon region or the rectangular region with an invalid circular region inside (here, we define the invalid region as the region where there is no POI or measured data inside). In such cases, there might be an invalid region within the rectangular Hermite element, and this will produce an incorrect global regularization matrix S in Eq. (14). As the example shown in Fig. 3, three Hermite elements are used to cover the nonrectangular ROI and each Heremite element contains an invalid region, which is discretized into triangular elements. For this general case, some specific considerations should be taken for the general Tikhonov regularization formulation in Eq. (11). The first term on the right side of Eq. (11) depicts the error between the measured and the approximated data, which is not affected by the existence of the invalid region, while the second term (i.e. the roughness penalty term) is closely related to the invalid region. Since the essence of term Roughness Penalty is to describe the roughness of reconstructed data field within the ROI, the roughness metric within invalid region should not be considered. From this regard, we propose the following improved Hermite finite element smoothing method (IHFESM) for arbitrary ROI. Let OROI, O and OInvalid be the whole region of ROI, Hermite elements and the invalid regions respectively, then we have
O ¼ OROI [ OInvalid
where N,x and N,y are the global matrices of the first partial derivatives of the shape function N with respect to the global physical coordinate x and y respectively.
ð18Þ
Thus, the specific Tikhonov regularization formulation for arbitrary ROI using C2 continuous Hermite element will be Z 2 2 ð3Þ CðqÞ ¼ 99uc u992 þ l 99f uc ðxÞ992 dO ð19Þ OROI
Considering Eq. (14), the global regularization matrix for invalid region satisfies Z Z E X 2 ð3Þ 99f uc ðxÞ992 dO ¼ ðqe ÞT ðN e,xxx ÞT N e,xxx OInvalid
Invalid
O e¼1 þ 3ðN e,xxy ÞT N e,xxy þ3ðN e,xyy ÞT N e,xyy þðN e,yyy ÞT N e,yyy dOqe
¼ qT SInvalid q ð20Þ
ROI
ð16Þ
1665
Invalid
N e,kkk ðk
, ¼ x,yÞ has the same meaning as that in and S ¼S S Eq. (13). Then, parameter vector in Eq. (15) turns to be qðlÞ ¼ ðN T N þ lS ROI Þ1 N T u
ð21Þ
where N is the global stiffness matrix in Eq. (9). The smoothed data could be obtained by Eq. (16).
2.4. Generalized cross-validation It is apparent that, parameter vector q in Eq. (15) is a function of regularization parameter l since matrix N and S are constant for give noisy data u and Hermite elements. The regularization parameter l controls the tradeoff between the ‘‘roughness’’ of the solution and the infidelity to the data [12]. Generalized crossvalidation (GCV) [6,12,13,16] is a well-known and effective method for determining the optimum regularization parameter lopt. In short, GCV takes the ‘‘leave-one-out’’ principle [16] and lopt minimizes the following error function V(l): 2
n99NqðlÞu992 VðlÞ ¼ 2 tr½I-AðlÞ where A(l) ¼N(NTNþ lS) 1NT.
ð17Þ
Invalid region
Region of interest
Hermite element
Fig. 3. Relationship between the region of interest, invalid region (discretized into triangular elements) and Hermite elements in IHFESM.
J.-q. Zhao et al. / Optics and Lasers in Engineering 50 (2012) 1662–1671
v = 20mm
y
4×ϕ 80 x
160
Detailed description of the proposed IHFESM is as follows: Step 1: For the given ROI, measured data u (measured displacements x or y direction) and number of Hermite elements, generate corresponding rectangular Hermite elements. Step 2: For each Hermite element, discretize the invalid region inside into triangular elements. Step 3: Calculate S in Eq. (14) by the Gaussian quadrature with 4 4 points for rectangular Hermite elements, and calculate SInvalid in Eq. (20) by the 7-point Hammer quadrature [19] for triangular elements generated in Step 2, thus SROI ¼S SInvalid; calculate the matrix N in Eq. (9). Step 4: Find the optimum regularization parameter lopt in Eq. (17) and compute the smoothed displacement and the strain fields with Eq. (16). Compared with method in [16], the proposed IHFESM makes additional consideration for the computation of the new global regularization matrix SROI by removing the part SInvalid from the original regularization matrix S. With the help of this improvement, we can not only take the advantage of the high continuity of Hermite shape function to guarantee the smoothness of displacement and strain fields simultaneously, but also explore wider range of practical measurement problems with the effective Tikhonov regularization. Meanwhile, IHFESM does not bring in any addition computation burden beside the calculation of SInvalid. Thus, the proposed IHFESM can tackle arbitrarily shaped ROI efficiently and it is especially applicable to reduce the noise in the measured displacement and the strain fields by DIC. As a contrast, we call the smooth method based on Eq. (11) as the Hermite finite element smoothing method (HFESM), which will be compared with IHFESM in the following section.
480
1666
160
480 Fig. 4. Specimen geometry used in simulation.
4. Experiments and discussion To investigate the performance of the proposed IHFESM, two experiments were carried out in this section. In the first experiment, the tension of square plate with four holes was investigated, and the proposed IHFESM was used to smooth the displacement data from the FEM simulation, and comparison was also made between IHFESM and HFESM. Then, in the second experiment, IHFESM was used to smooth the noisy displacements field obtained in a practical tension experiment of one plate with three holes. 4.1. Tension of square plate with four holes The 480 480 mm2 plate with four holes under tension shown in Fig. 4 was analyzed with the finite element method. The upper boundary of the plate was subjected to the displacement v ¼20 mm on y direction, and displacement v of lower boundary was set to zero and the displacement u of middle point on the lower boundary was totally fixed. The plate was of 1 mm thickness, and the elastic modulus for the plate material was 200 GPa, and Poisson’s ration was 0.3. The whole model was meshed with 39,361 4-noded quadrilateral elements of 40,089 nodes. The displacement and the strain fields obtained by FEM were treated as the theoretical ones for validating the reconstructed fields obtained by the proposed method. After the FEM simulation, the u (or v) displacements of 20,045 nodes were randomly selected from all the nodes and added with normally distributed noise N(0, s), where s was the standard deviation. In this experiment, s was set to 0.00734 for u field and 0.04976 for v field, which was 0.5% of the average of absolute u and v values respectively. Then, the proposed IHFESM was employed to smooth two noisy displacement fields, and 8 8 uniform Hermite elements (with 9 9 9 ¼729 degree of freedom) and 1108
Fig. 5. 8 8 Uniform Hermite elements and triangular elements used by IHFESM for the first experiment.
triangular elements were used (see Fig. 5). For comparison, the Hermite finite element smoothing method (HFESM) that does not consider the invalid region for arbitrary ROI was also utilized to smooth the same noisy displacements, and HFESM uses the same Hermite elements shown in Fig. 5. Fig. 6 shows the contour plots of the theoretical and noisy u fields, and fields smoothed by HFESM and IHFESM. The similar results for v field are shown in Fig. 7. It reveals that both HFESM and IHFESM could dramatically decrease the noise and produce smooth displacement fields. In terms of regularization error, the average difference between the fields by IHFESM and the theoretical ones are 0.010994 mm for u field and 0.022920 mm for v field, while for HFESM, the errors are 0.011576 mm for u field and 0.025005 mm for v field, therefore IHFESM yields smaller error than HFESM. For IHFESM, the optimum regularization parameter is 1134.47 for u field and 7549.16 for v field; for IHFESM, the values are 416.41 for u field and 1828.09 for v field. The reconstructed strain fields for Ex, Ey and gxy are shown in Figs. 8, 9 and 10 respectively, and the regularization errors by HFESM and IHFESM for these three components are listed in Table 1. As Table 1 indicates the regularization errors for all five components by IHFESM are smaller than that by HFESM. Compared with the contour plots by HFESM (see Figs. 8(b), 9(b) and 10(b)), the strain distributions obtained by the proposed IHFESM (in Figs. 8(c), 9(c) and 10(c)) are in a better agreement with the theoretical distributions, especially for the strain field Ey. The detailed verification in Figs. 11 and 12 reveals the same conclusion. Fig. 11 demonstrates the reconstructed strain Ey (mapped onto one horizontal line through the center of the specimen
J.-q. Zhao et al. / Optics and Lasers in Engineering 50 (2012) 1662–1671
Fig. 6. u (mm) displacement: theoretical (a), noisy (b) and smoothed by HFESM (c) and IHFESM (d).
Fig. 7. v (mm) displacement: theoretical (a), noisy (b), smoothed by HFESM (c) and IHFESM (d).
Fig. 8. Theoretical Ex (e) field (a), and reconstructed by HFESM (b) and IHFESM (c).
Fig. 9. Theoretical Ey (e) field (a), and reconstructed by HFESM (b) and IHFESM (c).
Fig. 10. Theoretical gxy (e) field (a), and reconstructed by HFESM (b) and IHFESM (c).
1667
1668
J.-q. Zhao et al. / Optics and Lasers in Engineering 50 (2012) 1662–1671
Table 1 Average difference between the results obtained by two methods and the theoretical ones.
HFESM IHFESM
u/mm
v/mm
Ex/e
Ey/e
gxy/e
0.011576 0.010994
0.025005 0.022920
0.001743 0.001598
0.003599 0.003395
0.003072 0.002665
200
0.05
100 0
0.04
-200
0.03
-200
0
200
Theoretical Ey
0.02
Ey by IHFESM Ey by HFESM
0.01 0
100
200 300 400 Distance on line /mm
500
Fig. 11. Comparison between the reconstructed component Ey (e) (mapped onto the horizontal line, whose location is shown in the upper-right sub-image) by HFESM and IHFESM.
The aluminum alloy plate with three holes shown in Fig. 13 under tension was investigated in this section. The tension was performed by applying displacement load on the upper edge. Before experiment, the plate was painted with random speckle pattern, and the deformation was recorded by a CCD camera (1392 1040 pixels 8 bits). After experiment, a pair of reference and deformed images shown in Fig. 14 were selected and analyzed by DIC with the Newton–Raphson (NR) iteration. The full-field displacement and strain of the region of interest of size 50 50 mm2 (shown in Fig. 13) were computed. The subset size used was of 35 35 pixels and the distance between two neighboring subsets was 5 pixels, and there were a total of 13,892 points of interest in the ROI. For full-field measurement, the initial guess of deformation of the first point for NR was calculated by the DE–TCE algorithm in [18]. Readers are referred to [18] for details of the DE–TCE algorithm. The final correlation result of the first point was passed to the next point of interest as initial guess. The convergence conditions for NR was set to ensure that variations in displacements were no larger than 1E-4 pixels
x Theoretical E
0.05
y
E by IHFESM y
Ey by HFESM
0.048
200
F
Line
100
y
35 50
10
0.052
0 -100 -200 -200
3×ϕ 8 0
200
0.044 0.042
10 10
0.046
50
Ey /ε
4.2. Tension of plate with three holes
25
Ey /ε
-100
fields and reconstructing inhomogeneous strain fields. Meanwhile, due to the improvement in computation of new global regularization matrix SROI, IHFESM is capable of tackling the arbitrary region of interest in DIC.
20
0.04
45° 0
50
100 150 200 Distance on line /mm
250
300
10
0.038
35
Fig. 12. Comparison between the reconstructed component Ey (e) (mapped onto the vertical line, whose location is shown in the upper-right sub-image) by HFESM and IHFESM.
Fig. 13. Specimen geometry used in experiment.
displayed in the upper-right sub-image) by two methods. It is apparent that, both HFESM and IHFESM are able to restore a smooth and reasonable strain field on the whole, while IHFESM provides more accurate result at the two valleys of distribution. Fig. 12 demonstrates the result of strain Ey mapped on the vertical line near the edge of ROI, and comparison shows that the strain obtained by IHFESM has smaller variance than that by HFESM, even though both of the results fluctuate around the real value. It should be noted that, even if the HFESM is employed to smooth the displacement fields in the current experiment, as stated in Section 3, HFESM is not recommended to be used to smooth the arbitrary region of interest, and this is also proven by the above experiment results. The comparison between IHFESM and HFESM is conducted in order to show the difference between these two methods. In conclusion, the proposed IHFESM produces smaller regularization error than HFESM in smoothing complex displacement
Fig. 14. Reference image (left) and deformed image (right) recorded in tension experiment.
F
J.-q. Zhao et al. / Optics and Lasers in Engineering 50 (2012) 1662–1671
Fig. 15. Comparison between the displacement field u by DIC (left) and u field smoothed by IHFESM (right).
Fig. 16. Comparison between the displacement field v by DIC (left) and v field smoothed by IHFESM (right).
Fig. 17. Comparison between the displacement field Ex by DIC (left) and Ex field smoothed by IHFESM (right).
Fig. 18. Comparison between the displacement field gxy by DIC (left) and gxy field smoothed by IHFESM (right).
1669
1670
J.-q. Zhao et al. / Optics and Lasers in Engineering 50 (2012) 1662–1671
Fig. 19. Comparison between the displacement field Ey by DIC (left) and Ey field smoothed by IHFESM (right).
and variations in displacement gradients were no larger than 0.5E-6. After full-field analysis by DIC with NR iteration, the obtained noisy u and v displacement field were smoothed by the proposed IHFESM method. For IHFESM, the 8 8 uniform Hermite elements 1886 triangular elements were used, and the optimum regularization parameter obtained was 490.92 for u field and 150.50 for v field. Figs. 15 and 16 show the contour plots of the u and v fields smoothed by IHFESM, and the comparison of the three strain components, i.e. Ex, gxy and Ey, are demonstrated in Figs. 17, 18 and 19 respectively. Results show that the displacement and the strain fields directly calculated by NR iteration are contaminated by a high level of noise, while after being smoothed by IHFESM, the noises in all five components are dramatically reduced, and this could be clearly seen in 3D profile of Ey field depicted in Fig. 19. From Fig. 19 we find that, the strains near the two sides of the two lower holes are larger than that near the one upper hole, and this is reasonable because the holes harms the strength of plate, and strain concentration is serious near region with two holes in the current situation.
5. Conclusions In this work, the improved C2 continuous Hermite finite element based displacement smoothing method, i.e. IHFESM, was proposed to tackle the measurement of the arbitrary region of interest in DIC. IHFESM utilizes the well-known Tikhonov regularization to characterize the smooth problem and employs the generalized cross-validation to determine the optimum regularization parameter. The C2 continuous Hermite element is an effective high-order conforming element, which is capable of ensuring both continuity and smoothness of the final strain field. However, this element needs to be rectangular, which causes difficulty in handling the arbitrary region of interest. To address this problem, IHFESM makes an improvement in the computation of the global regularization matrix by ignoring the contribution to the matrix from the invalid region. Compared with the existing Hermite finite element smoothing method (HFESM), IHFESM involves little additional computation burden, while as revealed in the simulation experiment of the plate with four holes, IHFESM produces smaller regularization error than HFESM. The tension
experiment of the plate with three holes demonstrates that, IHFESM could remarkably reduce the noises in the displacement and the strain fields which are obtained directly by the NR iteration. Therefore, the proposed IHFESM is an effective method in smoothing complex displacement field and reconstructing inhomogeneous strain field for the arbitrary region of interest in DIC.
Acknowledgments The work is supported by the National Natural Science Foundation of China under Grant 50635050, 51175293 and 10972114. References [1] Sutton MA, Orteu JJ, Schreier HW. Image correlation for shape, motion and deformation measurements. Springer; 2009. [2] Pan B. Recent Prog Digital Image Correl 2010:1–13, http://dx.doi.org/ 10.1007/s11340-010-9418-3. [3] Sutton MA, Turner JL, Bruck HA, Chae TA. Full-field representation of discretely sampled surface deformation for displacement and strain analysis. Exp Mech 1991;31(2):168–77. [4] Schreier HW, Braasch JR, Sutton MA. Systematic errors in digital image correlation caused by intensity interpolation. Opt Eng 2000;39(11):2915–21, http://dx.doi.org/10.1117/1.1314593. [5] Pan B, Asundi A, Xie HM, Gao JX. Digital image correlation using iterative least squares and pointwise least squares for displacement field and strain field measurements. Opt Laser Eng 2009;47(7-8):865–74, http://dx.doi.org/ 10.1016/j.optlaseng.2008.10.014. [6] Meng LB, Jin GC, Yao XF. Application of iteration and finite element smoothing technique for displacement and strain measurement of digital speckle correlation. Opt Laser Eng 2007;45(1):57–63. [7] Schreier HW, Sutton MA. Systematic errors in digital image correlation due to undermatched subset shape functions. Exp Mech 2002;42(3):303–10, http:// dx.doi.org/10.1007/BF02410987. [8] Davis CQ, Freeman DM. Statistics of subpixel registration algorithms based on spatiotemporal gradients or block matching. Opt Eng 1998;37(4):1290–8, http://dx.doi.org/10.1117/1.601966. [9] Reu PL. Experimental and Numerical Methods for Exact Subpixel Shifting. Exp Mech 2010:1–10. [10] Pan B, Xie H, Guo Z, Hua T. Full-field strain measurement using a twodimensional Savitzky-Golay digital differentiator in digital image correlation. Opt Eng 2007;46(3):033601, http://dx.doi.org/10.1117/1.2714926. [11] Avril S, Feissel P, Pierron F, Villon P. Comparison of two approaches for differentiating full-field data in solid mechanics. Meas Sci Technol 2010;21:015703. [12] Craven P, Wahba G. Smoothing noisy data with spline functions. Numer Math 1978;31(4):377–403. [13] Wang CC, Deng J, Ateshian GA, Hung CT. An Automated approach for direct measurement of two-dimensional strain distributions within articular
J.-q. Zhao et al. / Optics and Lasers in Engineering 50 (2012) 1662–1671
cartilage under unconfined compression. J Biomech Eng 2002;124(5): 557–67, http://dx.doi.org/10.1115/1.1503795. [14] Yoneyama S. Smoothing measured displacements and computing strains utilising finite element method. Strain 2011;47(s2):258–66, http://dx.doi.org/ 10.1111/j.1475-1305.2010.00765.x. [15] Segalman DJ, Woyak DB, Rowlands RE. Smooth spline-like finite-element differentiation of full-field experimental data over arbitrary geometry. Exp Mech 1979;19(12):429–37. [16] Sadati M, Luap C, Kroger M, Gusev AA, Ottinger HC. Smooth full field reconstruction of velocity and its gradients from noisy scattered velocimetry data in a cross-slot flow. J Rheol 2011;55(2):353–77.
1671
[17] Pan B, Xie H, Wang Z. Equivalence of digital image correlation criteria for pattern matching. Appl Opt 2010;49(28):5501–9, http://dx.doi.org/10.1364/ AO.49.005501. [18] Zhao J, Zeng P, Lei L, Ma Y. Initial guess by improved population-based intelligent algorithms for large inter-frame deformation measurement using digital image correlation. Opt Laser Eng 2012;50(3):473–90, http://dx.doi.org/ 10.1016/j.optlaseng.2011.10.005. [19] Szilard R. Theories and applications of plate analysis: classical, numerical, and engineering methods. Wiley; 2004.
Contents lists available at SciVerse ScienceDirect
Optics and Lasers in Engineering journal homepage: www.elsevier.com/locate/optlaseng
Improved Hermite finite element smoothing method for full-field strain measurement over arbitrary region of interest in digital image correlation Jia-qing Zhao a, Pan Zeng a,n, Bin Pan b, Li-ping Lei a, Hong-fei Du a, Wen-bin He a, Yang Liu a, Yue-jie Xu a a b
Key Laboratory for Advanced Materials Processing Technology of MOE, Department of Mechanical Engineering, Tsinghua University, Beijing 100084, China Institute of Solid Mechanics, Beihang University, Beijing 100191, China
a r t i c l e i n f o
a b s t r a c t
Article history: Received 13 February 2012 Received in revised form 18 April 2012 Accepted 24 April 2012 Available online 30 June 2012
We investigate the finite element method and the Tikhonov regularization for smoothing noisy displacement field obtained from the digital image correlation (DIC). Among various finite elements, the C2 continuous Hermite finite element has been proven to be an effective high-order element which ensures both continuity and smoothness of the final strain field. However, the standard Hermite element is required to be rectangular, which cannot tackle the region of interest of arbitrary shape. To overcome this constraint, we propose the improved Hermite finite element smoothing method (IHFESM) to extend the generality of Hermite finite element that is able to measure the arbitrary region of interest. The IHFESM is implemented by eliminating the contribution to the global regularization matrix from the invalid region. The proposed method could be employed to smooth the noisy displacement field obtained by DIC, and reconstruct the reliable strain field by differentiation. The effectiveness of IHFESM is verified by two tension experiments of plates with three and four holes. & 2012 Elsevier Ltd. All rights reserved.
Keywords: Digital image correlation Displacement smooth Hermite finite element Tikhonov regularization Generalized cross-validation
1. Introduction Digital image correlation (DIC) [1,2] has become an effective and widely used non-contact full-field strain measurement method in the field of experimental mechanics. Ever since DIC was proposed in 1980s, research on measurement accuracy of the strain has always been a hot issue in this field [3–7]. Due to various sources of unavoidable noise [8] (such as shot noise and thermal noise) during the imaging and transmission process, the digital images of the test specimen are always contaminated to certain extend. In addition, the gray values stored in the images are digitized discrete values of continuous natural light intensities, thus the final correlation results by DIC are always full of noise. On one hand, the image noise could be alleviated by employing the high-quality equipment, such as the CCD/CMOS camera of higher signal-noise ratio and higher gray-level resolution [9], however this will increase the hardware cost greatly. On the other hand, many researchers have concentrated on exploring the effective displacement and the strain smoothing methods to filter noise in the measured data, and these methods could be classified into two main catalogs: The local smooth method: in the local smooth method, only small part of the dense displacement data located in one calculation window will be considered at one time, and the selected data
n
Correseponding author. Tel.: þ 86 10 6279 4261; fax: þ86 10 6277 3637. E-mail address: [email protected] (P. Zeng).
0143-8166/$ - see front matter & 2012 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.optlaseng.2012.04.008
will be fitted by the given basis function (polynomial function for instance) and the resultant coefficients will be used to compute the smoothed strain of the center point of the calculation window. Two typical methods fall into this catalog are pointwise least squares [5,10] and diffuse approximation [11]. The advantages of local smooth methods are: (a) low computation and programming complexity. Generally, these methods could be realized with few lines of programming codes, and they often run faster than global methods due to low computation complexity and low data density being considered at each time, and (b) high adaptivity to the shape of region of interest. When using these methods, one could simply omit the invalid data in the calculation window to tackle the case of nonrectangular region of interest, especially on the edge of the region. However, the performance of the local smooth methods heavily relies on the proper size selection of the calculation window [5] and the form of the basis function [11], and this deficiency is highlighted in the inhomogeneous deformation measurement. The global smooth method: for this method, all the noisy displacements are considered as a whole and they are fitted by one globally defined basis function or pieces of interconnecting basis functions. After the smoothing of displacements, the strains are obtained by differentiation. The thin-plate spline smooth method [12,13] and the finite element based smooth method [3,6,11,14] are the two widely used approaches. The finite element based smooth methods could be further classified into two types: (a) the finite element based least squares method [11,14] and (b) the finite element based Tikhonov regularization
J.-q. Zhao et al. / Optics and Lasers in Engineering 50 (2012) 1662–1671
[6,15]. Tikhonov regularization is more effective than the least squares method in that it considers both the ‘‘closeness’’ and ‘‘smoothness’’ of the best fitting surface, and the generalized cross-validation (GCV) [6] is often chosen to determine the regularization parameter. Here, we focus on the finite element based Tikhonov regularization method. In terms of the element type, the existing highorder elements, such as 8-noded isoparametric element [14] and o the triangular Cowper element with 18 % of freedoms (also denoted as Cowper-18 element) [15] are frequently used to tackle the arbitrary region of interest, but they are still unable to reconstruct a continuous and smooth strain field. For example, the strain fields obtained by using Cowper-18 element are continuous but not smooth. More recently, Sadati et al. [16] proposed the Tikhonov regularization based on another higherorder element, the conforming C2 continuous (second order continuous) Hermite finite with quintic Hermite polynomials, which successfully guarantees the continuity and smoothness of the smoothed data and its gradient simultaneously. They also showed that the C2 continuous Hermite elements coupled with the minimization of a norm of the data’s third derivative (denoted as C2-R3 in [16]) is superior to other combinations, such as C0-LSQ (C0 continuous element based least squares method), C0-R2, C1-R2 and C2-R2. However, due to the nature of Hermite interpolation, each C2 continuous Hermite finite is required to be rectangular, which causes difficulty in applying Hermite element to the measurement of the arbitrarily shaped region of interest (ROI) in DIC. In this paper, we propose an improved Hermite finite element smoothing method (IHFESM) to tackle the arbitrary ROI (nonrectangular ROI for example), which is common in DIC. IHFESM is based on the C2 continuous Hermite finite element and the Tikhonov regularization, in which the minimization of norm of the data’s third derivative used in [16] is adopted. The difference between IHFESM and the method in [16] is that, IHFESM tackles the arbitrary ROI by omitting the contribution from the invalid region when calculating the global regularization matrix. The proposed IHFESM is simple and effective in smoothing the noisy displacement field for inhomogeneous strain field calculation in DIC, and the effectiveness of IHFESM is verified by two plate tension experiments in Section 4. The remainder of this paper is organized as follows. In Section 2, the basic principle of DIC, C2 continuous Hermite finite element, Tikhonov regularization and generalized cross-validation are presented. Section 3 demonstrates the details of IHFESM. The performance of IHFESM is investigated in Section 4. Finally, the discussions and conclusions are in Section 5.
2. Basic principles 2.1. Digital image correlation DIC [1,2] is a well-established non-contact optical full-field deformation measurement technique, which computes the displacement and the strain of points of interest (POI) in the predefined region of interest (ROI) on the speckled surface (see Fig.1). To calculate all desired deformation variables, the gray values of pixels in the deformed subset are compared with that of pixels in the reference subset. To describe the comparison more precisely, the zero-mean normalized cross-correlation criteria (ZNCC) [17] is widely used, which is defined as PM PM 0 0 x ¼ M y ¼ M ½f ðx,yÞf m ½gðx ,y Þg m qP ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ffi CðpÞ ¼ 1 qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi PM PM PM M 2 2 0 0 x ¼ M y ¼ M ½f ðx,yÞf m x ¼ M y ¼ M ½gðx ,y Þg m
ð1Þ
1663
x
ROI
Step size f(x, y) (ux ,u y,
y POI
(x, y)
vx ,v
y)
g(x', y')
v u Reference subset
(x', y')
Deformed subset
Fig. 1. Definition of the region of interest (ROI), point of interest (POI), deformation parameters, reference and deformed subset in DIC.
where f(x,y) represents the gray value of point (x,y) in the reference 0 0 0 0 subset and g(x ,y ) is the gray value of corresponding point (x ,y ) in the deformed subset, all the coordinates are depicted in coordinate system anchored at the center of the reference subset with size PM PM (2Mþ1) (2Mþ1) pixels, f m ¼ 1=ð2M þ1Þ2 x ¼ M y ¼ M f ðx,yÞ P PM 0 0 gðx ,y Þ are the average gray and g m ¼ 1=ð2M þ 1Þ2 M x ¼ M y ¼ M values of all points in two subsets, and p is the deformation vector which describes the relationship between coordinate (x,y) and 0 0 (x ,y ). One merit of ZNCC criteria is that it is insensitive to the scale and offset changes in illumination lighting fluctuations, which often occurs in the practical experiment. Generally, when the subset is small enough, (x,y) is mapped to 0 0 (x ,y ) by first-order shape function as @u @u xþ y @x @y @v @v y0 ¼ y þv þ x þ y @x @y
x0 ¼ x þu þ
ð2Þ
where u, v are the two displacement components of the reference subset center on the x and y direction, @u/@x, @u/@y, @v/@xand @v/@y are the four displacement gradient components, thus p ¼[u, v, @u/ @x, @u/@y, @v/@x, @v/@y]T is the corresponding parameter vector. The first-order shape function allows translation, rotation, shear, normal strains and their combinations of subset, and it provides all the necessary deformation information for the common measurement. 0 0 The coordinate (x ,y ) often locates at the sub-pixel position in the deformed image. Therefore, interpolation is necessary to be introduced to reconstruct their gray values from neighboring pixels. In this work, the bi-cubic spline interpolation is used to 0 0 determine the gray values at (x ,y ) by gðx,yÞ ¼
3 3 X X
am,n xm ,yn
ð3Þ
m¼0n¼0
where the coefficients am,n could be calculated to guarantee C2 continuity by using the bi-cubic spline functions [1]. The six deformation parameters could be resolved by minimizing the nonlinear function C(p) in Eq. (1) using the Newton–Raphson (NR) [5] or quasi-Newton (qN) [18] method. The iteration equation of NR is p ¼ p0 rrCðp0 Þ1 rCðp0 Þ
ð4Þ
where p is the approximation value after iteration, p0 denotes an initial guess of solution that is ‘‘sufficiently close’’ to the real solution,
1664
J.-q. Zhao et al. / Optics and Lasers in Engineering 50 (2012) 1662–1671
rC(p) and rrC(p) are gradient and Hessian matrix of function C(p). In this study, NR is chosen as the sub-pixel registration algorithm. The displacements and displacement gradients calculated by the NR are usually contaminated by various noises [8]. Especially when the displacements and the strains are small, the noise might be typically of the same order of the real data. Therefore a smooth procedure is needed to extract the true displacement and strain from the noisy correlation results. In this work, we use the improved C2 continuous Hermite finite element to build a parameterized displacement field, and then employ the Tikhonov regularization to approximate the true noiseless field. In the following, the basic principles of the Hermite finite element, Tikhonov regularization and generalized cross-validation will be introduced.
2.2. Hermite finite element
The six quintic Hermite polynomial shape functions corresponding to the six degrees of freedom of 1D element are 0
2
uðxÞ ¼ Ne ðxÞqe
ð5Þ
e
e
where q is of size d 1, N (x) is a row vector of size 1 d, and d is the number of nodal degrees of freedom. In Eq. (5), qe, Ne(x) and d are determined by the specific form of the interpolation function. To ensure higher displacement continuity than traditional elements (such as 4-noded or 8-noded isoparametric element [14] and triangular Cowper-18 element [15]), we use the C2 continuous Hermite interpolation function [16] for Eq. (5), and this element guarantees the continuity and smoothness of both displacement field and strain field. The definition of node numbering, nodal degrees of freedom and local coordinates for one-dimensional (1D) and two-dimensional (2D) Hermite elements are shown in Fig. 2. For 1D C2 continuous Hermite element, three degrees of freedom, ui, @ui/ @x, @2ui/@x2, i¼1, 2 are associated with each node. The connecting elements have these three common degrees of freedom; therefore it ensures continuity of slope and higher derivatives up to the second order between the elements.
2
N 2 ¼ ðx þ 1Þ3 ð3x 9x þ8Þ=16
1 N1 2 N1
1 N2 2 N2
¼ ðx1Þ3 ð3x þ 5Þðx þ 1Þ=16, 3
2
¼ ðx1Þ ðx þ 1Þ =16,
¼ ðx þ 1Þ3 ð3x5Þðx1Þ=16 ¼ ðx þ 1Þ3 ðx1Þ2 =16 ð6Þ
where the subscript indicates the node number and the superscript denotes the order of the derivative with respect to the local element coordinate xA[ 1, 1]. To ensure the C2 continuity of 2D rectangular Hermite element, the nine degrees of freedom at each node are needed: u,@u/@x,@u/@Z,@2u/@x@Z,@2u/@x2,@2u/@Z2,@3u/ @x2@Z,@3u/@x@Z2,@4u/@x2@Z2. Then, all the thirty six biquintic Hermite shape functions are obtained as products of above 1D shape functions [16] j
In the finite element method (FEM) [19], given the vector of nodal degree of freedom values qe and the shape function Ne(x), the displacement u(x) of point at xAR2 inside of one element could be expressed as
0
N1 ¼ ðx1Þ3 ð3x þ9x þ8Þ=16,
k
ðx, ZÞ ¼ Ni0 ðxÞNi00 ðZÞ, Nðj,kÞ i
i ¼ 1,. . .,4,
i0 , i00 ¼ 1,2,
j,k ¼ 0,1,2
ð7Þ
ðx, ZÞ is associated with the degree of freedom (@j þ ku/ where N ðj,kÞ i 0 00 @xj@yk)i of node i of the rectangular element, i and i refer to the node number in the corresponding 1D element and they take values of 1 or 2 if the local coordinate of node i equals 1 or 1. For 2 2 instance, Nð0,2Þ 4 2 ðx, ZÞ is associated with @ u/@y and it is the product 0 of N1 ðxÞ N2 ðZÞ. In fact, Eq. (7) assumes that all the Hermite elements in the ROI are of the same size of 2 2, while for the non-uniform rectangular elements, the 2D shape functions are of form: j
k
k
Nðj,kÞ ðx, ZÞ ¼ Ni0 ðxÞNi00 ðZÞwj h , i
i ¼ 1,. . .,4,
i0 , i00 ¼ 1,2
ð8Þ
where w and h are the half width and half height of the element respectively. In the context of C2 continuous Hermite element, qe in Eq. (5) is of size 36 1, and Ne(x) is of size 1 36, whose components are determined by Eq. (8). It should be noted that, x in Ne(x) is expressed in the global physical coordinate system of ROI, and it should be transformed into (x,Z) in the local Hermite element coordinate system before using Eq. (8). The assembly of Eq. (5) over all E elements in the ROI gives uc ¼
E X
N e ðxÞqe ¼ Nq
ð9Þ
e¼1
u
u
u
u
u
u
where q is the global parameter vector of size D 1 and D is the degrees of freedom of all elements over whole ROI, N is a matrix of size n D and n is the number of points of interest in the ROI, N is also called the global stiffness matrix, uc is the calculated displacements at the location of corresponding point of interest, P and it is a column vector of size n 1. Here ðdÞ is a symbolic denotation of the assembly process in FEM [19]. 2.3. Tikhonov regularization The unknown parameter vector q in Eq. (9) could be obtained as the general least squares solution [11,14] by minimizing: 2
CðqÞ ¼ 99uc u992 n
ui
ui
ui
ui
ui
ui
ui
ui
ui
i
Fig. 2. Definition of node numbering, nodal degrees of freedom and local coordinates for 1D and 2D Hermite elements.
ð10Þ
where uAR is a column vector collecting all the noisy displace2 ment values and 99d992 refers to the sum over all squared components. In fact, Eq. (10) only depicts the ‘‘closeness’’ of smoothed displacement field uc to the measured field u. The outcome of above least squares depends strongly on the unavoidable noise in image acquisition and correlation in DIC. However, this problem could be conquered by using the Tikhonov regularization method which enforces an additional restriction of ‘‘smoothness’’ for the
J.-q. Zhao et al. / Optics and Lasers in Engineering 50 (2012) 1662–1671
solution [6,13,16]. The smoothness restriction is often interpreted as specific order of derivative of the smoothed surface f(x), and it is then added to Eq. (10) as the roughness penalty term by Z 2 2 CðqÞ ¼ 99uc u992 þ l 99f ðkÞ uc ðxÞ992 dO ð11Þ O
where l is the regularization parameter, O is the domain of all FEM elements, f(k) is the kth derivative of f, and f(k)uc(x) defines the smoothness within the element. There are various choices of parameter k in f(k)uc(x), and k¼2 (i.e. the second derivative) is often used in the relevant studies [6,13]. However, as pointed in [16], one’s priori information about the final solution would provide a useful guideline for the selection of k. More precisely, if one is concerned with the smoothness of reconstructed data, the second derivative is a proper choice; while if the smoothness of reconstructed data gradient is favored, the selection of the third derivative will be more reasonable. This criterion is also validated by the comprehensive investigations in [16]. In this work, we are more concerned about the strain (i.e. displacement gradient) measurement by DIC, therefore the third derivative, i.e. k¼3, is actually adopted. 2 ð3Þ The specific form of 99f uc ðxÞ992 is !2 !2 !2 !2 @3 u @3 u @3 u @3 u 2 ð3Þ 99f uc ðxÞ992 ¼ þ 3 þ 3 þ ð12Þ @x3 @x2 @y @x@y2 @y3 Substituting Eq. (5) into Eq. (12), the whole roughness penalty term could be written as Z Z E X 2 ð3Þ 99f uc ðxÞ992 dO ¼ ðqe ÞT ðN e,xxx ÞT N e,xxx þ 3ðN e,xxy ÞT N e,xxy O
e
O e¼1 þ 3ðN e,xyy ÞT N e,xyy þ ðN e,yyy ÞT N e,yyy dOqe
ð13Þ
N e,kkk ðk ¼ x,yÞ
Here denotes the row vector of third partial derivatives of the shape function Ne of element e with respect to the global physical coordinate x or y. Eq. (13) could be rewritten into a compact form as Z 2 ð3Þ 99f uc ðxÞðxÞ992 dO ¼ qT Sq ð14Þ O
where S is the D D global regularization matrix. Substituting Eq. (14) into Eq. (11), the solution of Eq. (11) will be qðlÞ ¼ ðN T N þ lSÞ1 N T u
ð15Þ
Then the smoothed displacement field and two displacement gradient fields could be obtained using us ¼ NqðlÞ,
usx ¼ N ,x qðlÞ,
usy ¼ N ,y qðlÞ
3. Improved Hermite finite element smoothing method (IHFESM) Above C2 continuous rectangular Hermite element based Tikhonov regularization is effective in smoothing noisy data, which ensures the smoothness of the measured data and the data gradient. However, in the practical measurement by DIC, the user-defined ROI are often of arbitrary shape, for example, the arbitrary polygon region or the rectangular region with an invalid circular region inside (here, we define the invalid region as the region where there is no POI or measured data inside). In such cases, there might be an invalid region within the rectangular Hermite element, and this will produce an incorrect global regularization matrix S in Eq. (14). As the example shown in Fig. 3, three Hermite elements are used to cover the nonrectangular ROI and each Heremite element contains an invalid region, which is discretized into triangular elements. For this general case, some specific considerations should be taken for the general Tikhonov regularization formulation in Eq. (11). The first term on the right side of Eq. (11) depicts the error between the measured and the approximated data, which is not affected by the existence of the invalid region, while the second term (i.e. the roughness penalty term) is closely related to the invalid region. Since the essence of term Roughness Penalty is to describe the roughness of reconstructed data field within the ROI, the roughness metric within invalid region should not be considered. From this regard, we propose the following improved Hermite finite element smoothing method (IHFESM) for arbitrary ROI. Let OROI, O and OInvalid be the whole region of ROI, Hermite elements and the invalid regions respectively, then we have
O ¼ OROI [ OInvalid
where N,x and N,y are the global matrices of the first partial derivatives of the shape function N with respect to the global physical coordinate x and y respectively.
ð18Þ
Thus, the specific Tikhonov regularization formulation for arbitrary ROI using C2 continuous Hermite element will be Z 2 2 ð3Þ CðqÞ ¼ 99uc u992 þ l 99f uc ðxÞ992 dO ð19Þ OROI
Considering Eq. (14), the global regularization matrix for invalid region satisfies Z Z E X 2 ð3Þ 99f uc ðxÞ992 dO ¼ ðqe ÞT ðN e,xxx ÞT N e,xxx OInvalid
Invalid
O e¼1 þ 3ðN e,xxy ÞT N e,xxy þ3ðN e,xyy ÞT N e,xyy þðN e,yyy ÞT N e,yyy dOqe
¼ qT SInvalid q ð20Þ
ROI
ð16Þ
1665
Invalid
N e,kkk ðk
, ¼ x,yÞ has the same meaning as that in and S ¼S S Eq. (13). Then, parameter vector in Eq. (15) turns to be qðlÞ ¼ ðN T N þ lS ROI Þ1 N T u
ð21Þ
where N is the global stiffness matrix in Eq. (9). The smoothed data could be obtained by Eq. (16).
2.4. Generalized cross-validation It is apparent that, parameter vector q in Eq. (15) is a function of regularization parameter l since matrix N and S are constant for give noisy data u and Hermite elements. The regularization parameter l controls the tradeoff between the ‘‘roughness’’ of the solution and the infidelity to the data [12]. Generalized crossvalidation (GCV) [6,12,13,16] is a well-known and effective method for determining the optimum regularization parameter lopt. In short, GCV takes the ‘‘leave-one-out’’ principle [16] and lopt minimizes the following error function V(l): 2
n99NqðlÞu992 VðlÞ ¼ 2 tr½I-AðlÞ where A(l) ¼N(NTNþ lS) 1NT.
ð17Þ
Invalid region
Region of interest
Hermite element
Fig. 3. Relationship between the region of interest, invalid region (discretized into triangular elements) and Hermite elements in IHFESM.
J.-q. Zhao et al. / Optics and Lasers in Engineering 50 (2012) 1662–1671
v = 20mm
y
4×ϕ 80 x
160
Detailed description of the proposed IHFESM is as follows: Step 1: For the given ROI, measured data u (measured displacements x or y direction) and number of Hermite elements, generate corresponding rectangular Hermite elements. Step 2: For each Hermite element, discretize the invalid region inside into triangular elements. Step 3: Calculate S in Eq. (14) by the Gaussian quadrature with 4 4 points for rectangular Hermite elements, and calculate SInvalid in Eq. (20) by the 7-point Hammer quadrature [19] for triangular elements generated in Step 2, thus SROI ¼S SInvalid; calculate the matrix N in Eq. (9). Step 4: Find the optimum regularization parameter lopt in Eq. (17) and compute the smoothed displacement and the strain fields with Eq. (16). Compared with method in [16], the proposed IHFESM makes additional consideration for the computation of the new global regularization matrix SROI by removing the part SInvalid from the original regularization matrix S. With the help of this improvement, we can not only take the advantage of the high continuity of Hermite shape function to guarantee the smoothness of displacement and strain fields simultaneously, but also explore wider range of practical measurement problems with the effective Tikhonov regularization. Meanwhile, IHFESM does not bring in any addition computation burden beside the calculation of SInvalid. Thus, the proposed IHFESM can tackle arbitrarily shaped ROI efficiently and it is especially applicable to reduce the noise in the measured displacement and the strain fields by DIC. As a contrast, we call the smooth method based on Eq. (11) as the Hermite finite element smoothing method (HFESM), which will be compared with IHFESM in the following section.
480
1666
160
480 Fig. 4. Specimen geometry used in simulation.
4. Experiments and discussion To investigate the performance of the proposed IHFESM, two experiments were carried out in this section. In the first experiment, the tension of square plate with four holes was investigated, and the proposed IHFESM was used to smooth the displacement data from the FEM simulation, and comparison was also made between IHFESM and HFESM. Then, in the second experiment, IHFESM was used to smooth the noisy displacements field obtained in a practical tension experiment of one plate with three holes. 4.1. Tension of square plate with four holes The 480 480 mm2 plate with four holes under tension shown in Fig. 4 was analyzed with the finite element method. The upper boundary of the plate was subjected to the displacement v ¼20 mm on y direction, and displacement v of lower boundary was set to zero and the displacement u of middle point on the lower boundary was totally fixed. The plate was of 1 mm thickness, and the elastic modulus for the plate material was 200 GPa, and Poisson’s ration was 0.3. The whole model was meshed with 39,361 4-noded quadrilateral elements of 40,089 nodes. The displacement and the strain fields obtained by FEM were treated as the theoretical ones for validating the reconstructed fields obtained by the proposed method. After the FEM simulation, the u (or v) displacements of 20,045 nodes were randomly selected from all the nodes and added with normally distributed noise N(0, s), where s was the standard deviation. In this experiment, s was set to 0.00734 for u field and 0.04976 for v field, which was 0.5% of the average of absolute u and v values respectively. Then, the proposed IHFESM was employed to smooth two noisy displacement fields, and 8 8 uniform Hermite elements (with 9 9 9 ¼729 degree of freedom) and 1108
Fig. 5. 8 8 Uniform Hermite elements and triangular elements used by IHFESM for the first experiment.
triangular elements were used (see Fig. 5). For comparison, the Hermite finite element smoothing method (HFESM) that does not consider the invalid region for arbitrary ROI was also utilized to smooth the same noisy displacements, and HFESM uses the same Hermite elements shown in Fig. 5. Fig. 6 shows the contour plots of the theoretical and noisy u fields, and fields smoothed by HFESM and IHFESM. The similar results for v field are shown in Fig. 7. It reveals that both HFESM and IHFESM could dramatically decrease the noise and produce smooth displacement fields. In terms of regularization error, the average difference between the fields by IHFESM and the theoretical ones are 0.010994 mm for u field and 0.022920 mm for v field, while for HFESM, the errors are 0.011576 mm for u field and 0.025005 mm for v field, therefore IHFESM yields smaller error than HFESM. For IHFESM, the optimum regularization parameter is 1134.47 for u field and 7549.16 for v field; for IHFESM, the values are 416.41 for u field and 1828.09 for v field. The reconstructed strain fields for Ex, Ey and gxy are shown in Figs. 8, 9 and 10 respectively, and the regularization errors by HFESM and IHFESM for these three components are listed in Table 1. As Table 1 indicates the regularization errors for all five components by IHFESM are smaller than that by HFESM. Compared with the contour plots by HFESM (see Figs. 8(b), 9(b) and 10(b)), the strain distributions obtained by the proposed IHFESM (in Figs. 8(c), 9(c) and 10(c)) are in a better agreement with the theoretical distributions, especially for the strain field Ey. The detailed verification in Figs. 11 and 12 reveals the same conclusion. Fig. 11 demonstrates the reconstructed strain Ey (mapped onto one horizontal line through the center of the specimen
J.-q. Zhao et al. / Optics and Lasers in Engineering 50 (2012) 1662–1671
Fig. 6. u (mm) displacement: theoretical (a), noisy (b) and smoothed by HFESM (c) and IHFESM (d).
Fig. 7. v (mm) displacement: theoretical (a), noisy (b), smoothed by HFESM (c) and IHFESM (d).
Fig. 8. Theoretical Ex (e) field (a), and reconstructed by HFESM (b) and IHFESM (c).
Fig. 9. Theoretical Ey (e) field (a), and reconstructed by HFESM (b) and IHFESM (c).
Fig. 10. Theoretical gxy (e) field (a), and reconstructed by HFESM (b) and IHFESM (c).
1667
1668
J.-q. Zhao et al. / Optics and Lasers in Engineering 50 (2012) 1662–1671
Table 1 Average difference between the results obtained by two methods and the theoretical ones.
HFESM IHFESM
u/mm
v/mm
Ex/e
Ey/e
gxy/e
0.011576 0.010994
0.025005 0.022920
0.001743 0.001598
0.003599 0.003395
0.003072 0.002665
200
0.05
100 0
0.04
-200
0.03
-200
0
200
Theoretical Ey
0.02
Ey by IHFESM Ey by HFESM
0.01 0
100
200 300 400 Distance on line /mm
500
Fig. 11. Comparison between the reconstructed component Ey (e) (mapped onto the horizontal line, whose location is shown in the upper-right sub-image) by HFESM and IHFESM.
The aluminum alloy plate with three holes shown in Fig. 13 under tension was investigated in this section. The tension was performed by applying displacement load on the upper edge. Before experiment, the plate was painted with random speckle pattern, and the deformation was recorded by a CCD camera (1392 1040 pixels 8 bits). After experiment, a pair of reference and deformed images shown in Fig. 14 were selected and analyzed by DIC with the Newton–Raphson (NR) iteration. The full-field displacement and strain of the region of interest of size 50 50 mm2 (shown in Fig. 13) were computed. The subset size used was of 35 35 pixels and the distance between two neighboring subsets was 5 pixels, and there were a total of 13,892 points of interest in the ROI. For full-field measurement, the initial guess of deformation of the first point for NR was calculated by the DE–TCE algorithm in [18]. Readers are referred to [18] for details of the DE–TCE algorithm. The final correlation result of the first point was passed to the next point of interest as initial guess. The convergence conditions for NR was set to ensure that variations in displacements were no larger than 1E-4 pixels
x Theoretical E
0.05
y
E by IHFESM y
Ey by HFESM
0.048
200
F
Line
100
y
35 50
10
0.052
0 -100 -200 -200
3×ϕ 8 0
200
0.044 0.042
10 10
0.046
50
Ey /ε
4.2. Tension of plate with three holes
25
Ey /ε
-100
fields and reconstructing inhomogeneous strain fields. Meanwhile, due to the improvement in computation of new global regularization matrix SROI, IHFESM is capable of tackling the arbitrary region of interest in DIC.
20
0.04
45° 0
50
100 150 200 Distance on line /mm
250
300
10
0.038
35
Fig. 12. Comparison between the reconstructed component Ey (e) (mapped onto the vertical line, whose location is shown in the upper-right sub-image) by HFESM and IHFESM.
Fig. 13. Specimen geometry used in experiment.
displayed in the upper-right sub-image) by two methods. It is apparent that, both HFESM and IHFESM are able to restore a smooth and reasonable strain field on the whole, while IHFESM provides more accurate result at the two valleys of distribution. Fig. 12 demonstrates the result of strain Ey mapped on the vertical line near the edge of ROI, and comparison shows that the strain obtained by IHFESM has smaller variance than that by HFESM, even though both of the results fluctuate around the real value. It should be noted that, even if the HFESM is employed to smooth the displacement fields in the current experiment, as stated in Section 3, HFESM is not recommended to be used to smooth the arbitrary region of interest, and this is also proven by the above experiment results. The comparison between IHFESM and HFESM is conducted in order to show the difference between these two methods. In conclusion, the proposed IHFESM produces smaller regularization error than HFESM in smoothing complex displacement
Fig. 14. Reference image (left) and deformed image (right) recorded in tension experiment.
F
J.-q. Zhao et al. / Optics and Lasers in Engineering 50 (2012) 1662–1671
Fig. 15. Comparison between the displacement field u by DIC (left) and u field smoothed by IHFESM (right).
Fig. 16. Comparison between the displacement field v by DIC (left) and v field smoothed by IHFESM (right).
Fig. 17. Comparison between the displacement field Ex by DIC (left) and Ex field smoothed by IHFESM (right).
Fig. 18. Comparison between the displacement field gxy by DIC (left) and gxy field smoothed by IHFESM (right).
1669
1670
J.-q. Zhao et al. / Optics and Lasers in Engineering 50 (2012) 1662–1671
Fig. 19. Comparison between the displacement field Ey by DIC (left) and Ey field smoothed by IHFESM (right).
and variations in displacement gradients were no larger than 0.5E-6. After full-field analysis by DIC with NR iteration, the obtained noisy u and v displacement field were smoothed by the proposed IHFESM method. For IHFESM, the 8 8 uniform Hermite elements 1886 triangular elements were used, and the optimum regularization parameter obtained was 490.92 for u field and 150.50 for v field. Figs. 15 and 16 show the contour plots of the u and v fields smoothed by IHFESM, and the comparison of the three strain components, i.e. Ex, gxy and Ey, are demonstrated in Figs. 17, 18 and 19 respectively. Results show that the displacement and the strain fields directly calculated by NR iteration are contaminated by a high level of noise, while after being smoothed by IHFESM, the noises in all five components are dramatically reduced, and this could be clearly seen in 3D profile of Ey field depicted in Fig. 19. From Fig. 19 we find that, the strains near the two sides of the two lower holes are larger than that near the one upper hole, and this is reasonable because the holes harms the strength of plate, and strain concentration is serious near region with two holes in the current situation.
5. Conclusions In this work, the improved C2 continuous Hermite finite element based displacement smoothing method, i.e. IHFESM, was proposed to tackle the measurement of the arbitrary region of interest in DIC. IHFESM utilizes the well-known Tikhonov regularization to characterize the smooth problem and employs the generalized cross-validation to determine the optimum regularization parameter. The C2 continuous Hermite element is an effective high-order conforming element, which is capable of ensuring both continuity and smoothness of the final strain field. However, this element needs to be rectangular, which causes difficulty in handling the arbitrary region of interest. To address this problem, IHFESM makes an improvement in the computation of the global regularization matrix by ignoring the contribution to the matrix from the invalid region. Compared with the existing Hermite finite element smoothing method (HFESM), IHFESM involves little additional computation burden, while as revealed in the simulation experiment of the plate with four holes, IHFESM produces smaller regularization error than HFESM. The tension
experiment of the plate with three holes demonstrates that, IHFESM could remarkably reduce the noises in the displacement and the strain fields which are obtained directly by the NR iteration. Therefore, the proposed IHFESM is an effective method in smoothing complex displacement field and reconstructing inhomogeneous strain field for the arbitrary region of interest in DIC.
Acknowledgments The work is supported by the National Natural Science Foundation of China under Grant 50635050, 51175293 and 10972114. References [1] Sutton MA, Orteu JJ, Schreier HW. Image correlation for shape, motion and deformation measurements. Springer; 2009. [2] Pan B. Recent Prog Digital Image Correl 2010:1–13, http://dx.doi.org/ 10.1007/s11340-010-9418-3. [3] Sutton MA, Turner JL, Bruck HA, Chae TA. Full-field representation of discretely sampled surface deformation for displacement and strain analysis. Exp Mech 1991;31(2):168–77. [4] Schreier HW, Braasch JR, Sutton MA. Systematic errors in digital image correlation caused by intensity interpolation. Opt Eng 2000;39(11):2915–21, http://dx.doi.org/10.1117/1.1314593. [5] Pan B, Asundi A, Xie HM, Gao JX. Digital image correlation using iterative least squares and pointwise least squares for displacement field and strain field measurements. Opt Laser Eng 2009;47(7-8):865–74, http://dx.doi.org/ 10.1016/j.optlaseng.2008.10.014. [6] Meng LB, Jin GC, Yao XF. Application of iteration and finite element smoothing technique for displacement and strain measurement of digital speckle correlation. Opt Laser Eng 2007;45(1):57–63. [7] Schreier HW, Sutton MA. Systematic errors in digital image correlation due to undermatched subset shape functions. Exp Mech 2002;42(3):303–10, http:// dx.doi.org/10.1007/BF02410987. [8] Davis CQ, Freeman DM. Statistics of subpixel registration algorithms based on spatiotemporal gradients or block matching. Opt Eng 1998;37(4):1290–8, http://dx.doi.org/10.1117/1.601966. [9] Reu PL. Experimental and Numerical Methods for Exact Subpixel Shifting. Exp Mech 2010:1–10. [10] Pan B, Xie H, Guo Z, Hua T. Full-field strain measurement using a twodimensional Savitzky-Golay digital differentiator in digital image correlation. Opt Eng 2007;46(3):033601, http://dx.doi.org/10.1117/1.2714926. [11] Avril S, Feissel P, Pierron F, Villon P. Comparison of two approaches for differentiating full-field data in solid mechanics. Meas Sci Technol 2010;21:015703. [12] Craven P, Wahba G. Smoothing noisy data with spline functions. Numer Math 1978;31(4):377–403. [13] Wang CC, Deng J, Ateshian GA, Hung CT. An Automated approach for direct measurement of two-dimensional strain distributions within articular
J.-q. Zhao et al. / Optics and Lasers in Engineering 50 (2012) 1662–1671
cartilage under unconfined compression. J Biomech Eng 2002;124(5): 557–67, http://dx.doi.org/10.1115/1.1503795. [14] Yoneyama S. Smoothing measured displacements and computing strains utilising finite element method. Strain 2011;47(s2):258–66, http://dx.doi.org/ 10.1111/j.1475-1305.2010.00765.x. [15] Segalman DJ, Woyak DB, Rowlands RE. Smooth spline-like finite-element differentiation of full-field experimental data over arbitrary geometry. Exp Mech 1979;19(12):429–37. [16] Sadati M, Luap C, Kroger M, Gusev AA, Ottinger HC. Smooth full field reconstruction of velocity and its gradients from noisy scattered velocimetry data in a cross-slot flow. J Rheol 2011;55(2):353–77.
1671
[17] Pan B, Xie H, Wang Z. Equivalence of digital image correlation criteria for pattern matching. Appl Opt 2010;49(28):5501–9, http://dx.doi.org/10.1364/ AO.49.005501. [18] Zhao J, Zeng P, Lei L, Ma Y. Initial guess by improved population-based intelligent algorithms for large inter-frame deformation measurement using digital image correlation. Opt Laser Eng 2012;50(3):473–90, http://dx.doi.org/ 10.1016/j.optlaseng.2011.10.005. [19] Szilard R. Theories and applications of plate analysis: classical, numerical, and engineering methods. Wiley; 2004.