Comparison of two techniques for implementing the proper orthogonal decomposition method in damage detection problems

Available online at www.sciencedirect.com

8CImNCIE~DIREOT e ELSEVIER

MATHEMATICAL

AND COMPUTER MODELLING

Mathematical and Computer Modelling 40 (2004) 553-571 www.elsevier.com/locate/ mc m

Comparison of Two Techniques for Implementing the Proper Orthogonal Decomposition M e t h o d in Damage Detection Problems M. L. JOYNER State University of West Georgia Carrollton, G A 30118, U.S.A. (Received and accepted October 2003) A b s t r a c t - - T h e viability of using the reduced-order proper orthogonal decomposition (POD) (also called principal component analysis) methods to reduce the total computational time to detect damages using eddy current nondestructive evaluation techniques was demonstrated in a previous paper [i]. In this paper, we concentrate on various alternatives used to form the reduced-order solution to the forward problem, still in the context of eddy current damage detection. In particular, we focus on two different algorithms, a POD/Galerkin technique and a POD/interpolation technique. The POD/Galerkin method is a popular choice in the implementation of the POD method; however, in this paper, we will point out some of the problems in the traditional implementation of the reducedorder POD/Galerkin method when used in conjunction with eddy current damage detection. We will also compare the POD/Galerkin method to the POD/interpolation method and argue that in certain circumstances, the POD/interpolation method may be a better choice. (~) 2004 Elsevier Ltd. All rights reserved. K e y w o r d s - - R e d u c e d - o r d e r methods, Proper orthogonal decomposition, Karhunen Loeve decomposition, Principal component analysis, Eddy current nondestructive evaluation.

1. I N T R O D U C T I O N D a m a g e d e t e c t i o n p r o b l e m s , t y p i c a l l y referred to as n o n d e s t r u c t i v e e v a l u a t i o n , a r e q u i t e n a t u r a l l y f o r m u l a t e d in t h e c o n t e x t of inverse p r o b l e m s . I n t h e case of e l e c t r o m a g n e t i c i n t e r r o g a t i o n t e c h n i q u e s , such as t h e e d d y c u r r e n t m e t h o d , t h e p r o b l e m involves M a x w e l l ' s e q u a t i o n s a t s o m e level. C o n s e q u e n t l y , one c a n a n t i c i p a t e c o m p u t a t i o n a l a l g o r i t h m s t h a t a r e t i m e a n d c o m p u t e r m e m o r y intensive. M o t i v a t e d b y goals o f online, r e a l - t i m e a l g o r i t h m s t o b e u s e d in p o r t a b l e t e s t i n g devices, one is led to r e d u c e d - o r d e r c o m p u t a t i o n a l ideas. I n a n e a r l i e r p a p e r [1] (see also [2]), we d e m o n s t r a t e d t h e v i a b i l i t y of using t h e p r o p e r o r t h o g o n a l d e c o m p o s i t i o n ( P O D ) (also c a l l e d p r i n c i p a l c o m p o n e n t a n a l y s i s ) m e t h o d s t o r e d u c e t h e t o t a l c o m p u t a t i o n a l t i m e w h e n u s i n g e i t h e r s i m u l a t e d o r e x p e r i m e n t a l d a t a . I n t h e p r e v i o u s p a p e r , we also briefly m e n t i o n e d This research was stimulated by discussions while the author was a visitor at the Statistical and Applied Mathematical Sciences Institute (SAMSI) in Research Triangle Park~ NC~ which is supported in part by the National Science Foundation under Grant DMS-0112069. 0895-7177/04/$ - see front matter (~) 2004 Elsevier Ltd. All rights reserved. doi: 10.1016/j.mcm.2003.10.050

Typeset by . A ~ - T E X

554

M . L . JOYNER

some issues concerning the formulation of the reduced-order POD approximation. We address these issues in this paper and offer some helpful revisions to the popular method. In this paper, we begin by giving a brief overview of the eddy current method and the corresponding model governing this problem. We then present an overview of the POD method in the context of inverse problems. The rest of the paper focuses on the accuracy of the reducedorder POD approximation using two different methods for calculating the coefficients in the POD approximation. We first discuss the typical POD/Galerkin method used in most POD approximations. We then formulate a revised POD/Galerkin method which performs better in the case of this eddy current damage detection problem. We compare the revised POD/Galerkin method to a POD/interpolation method for calculating the coefficients for the POD approximation. We conclude with a final comparison of methods and the pros and cons of each method in the context of the eddy current damage detection problem. 2. E D D Y

CURRENT

METHOD

AND

MODEL

Nondestructive evaluation (NDE) is the process of detecting damages in structures such as airplanes and pipelines without destroying the structure's future usefulness. A particular type of damage detection technique is cMled the eddy current method. Eddy currents are currents found in any conducting material which is subjected to a time-varying magnetic field. They are useful for NDE purposes, because if a flaw is present within a conducting material, the flow of the eddy currents will be disrupted in some manner. From this disruption, and its effect on the magnetic flux density, one can often discern information about the damage or defect within the material. A thorough analysis of eddy currents and their behavior can be found in [3]. In our implementation of the eddy current method, we examine the process of inducing eddy currents within a sample by placing a thin conducting sheet carrying a uniform current above the sample. The current within the sheet produces a magnetic field perpendicular to it that in turn produces eddy currents within the sample. The presence of a flaw within the sample causes a disruption in the flow of the eddy currents and this disruption manifests itself in the magnetic flux density which can be measured by a device placed above the conducting sheet. A schematic of the inspection process is shown in Figure 1. Ultimately, we want to be able to determine parameters describing a damage given the magnetic flux density detected by the measuring device. To do this, we first formulate a forward problem for the magnetic flux density given a specific damage within the sample. Although the sensor detects the magnetic flux density B, we formulate the forward problem in terms of the magnetic vector potential A (for computational purposes) from which we can find B according to the equation B=VxA.

nsor

S a m pie

Material Flaw

Figure 1. 3-D schematic of eddy current inspection process.

Comparison of Two Techniques

555

Table 1. Description of parameters used in boundary value problem. Symbol

Given Parameter or Variable Parameter

Parameter Parameter Parameter Parameter Variable Variable NA

~T

A V¢ Jt

Description Magnetic permeability Conductivity Angular frequency Electric permittivity Total current in conducting sheet (cs) Magnetic vector potential Electric scalar potential Total current density = (a + iwe)(-iwA - ~7¢)

Using Maxwell's equations in phasor form (we use phasor form [4,5] to suppress the time dependence; see [6] for justification), constitutive laws, Ohm's law, and simplifying assumptions (e.g., uniformity in the direction of the current flow is assumed to reduce the problem into a twodimensional problem; see [6] for details), we can derive the following boundary value problem for a computational domain ~, with a description of the parameters and variables given in Table 1.

v ×

~

× A(~, y)

= (.(~, y) + i ~ ( x , y ) ) ( - i ~ A ( ~ , y) - r e ) ,

~,u e n,

(1) (2)

= fcs Jr..da = f3(x, y)+ and V¢ --- 0,

x, y E ~ \ cs

(3)

with A(x, -35) = 0 = A(x, 35), V A . nl(0,y) = 0 = V A . nl(50,y ). Given the assumption that there is uniformity in the direction of the current flow, A has the form A = (0, 0, A3), and hence, we will denote the vector A by its scalar nonzero component A, i.e., the A3 component of A. We can then reduce the above system into a single integro-differential equation by noting (see [6] for details) that the term V¢, the electric scalar potential, is piecewise constant across all regions. Therefore, V¢ can be written in terms of A in (2). Hence, we can combine all the equations as done in [6] to obtain the equivalent form of the boundary value problem

02A(x,y) 02A(x,y)) Ox2 + Oy2 + #(x, y)(a(x, y) + iwe(x, y))(iwA(z, y) + K(x, y)) -- O, where

(4)

V¢ = K(x, y) is defined by K(x, y) =

{

Its Acs(acu ÷ iwec~) 0,

iw ~s ~ A da, for (x,y) E cs,

(5)

for (x, y) E ~ \ cs.

The effects of the damage are inherently embedded in the conductivity term, a. The sample is considered to be an aluminum sample with conductivity aAl ~ 3.72 X 107 S m -1. If no damage exists, the conductivity within the sample is uniform throughout. However, if a damage exists, the damage will cause a discontinuity in the conductivity within the sample since the damage is assumed to be composed of air which has conductivity 0 S m -1. Hence, at the edges of the damage, there is a significant change in conductivity (a difference on the order of 107). This

556

M.L. JOYNER

discontinuity in the conductivity will effect the magnetic flux density according to the boundary value problem given by equation (4). As shown in [1], the solution to the boundary value problem depends continuously on the parameters which describe the damage, and thus, the degree to which the magnetic flux density is altered depends on the "size" of the damage. The task, however, is to determine the parameters describing the damage given the magnetic flux density detected by a portable sensing device, a typical inverse problem. The inverse problem is a computationMly intensive iterative process in which the boundary value problem must be solved possibly numerous times. Using standard finite element methods, the inverse problem would be extremely time consuming and not practical in experimental settings. Therefore, we developed a computationd method using reduced order Karhunen-Loeve or proper orthogonal decomposition (POD) techniques to reduce the computationM time significantly (see [1,2,6]).

3. O V E R V I E W

OF T H E P O D

METHOD

The POD technique is an attractive order reduction method, because basis elements are formed in an "optimal" way which span a data set consisting of either numerical simulations or experimental data. Since the POD basis is formed so that each basis element captures important aspects of the data set, only a small number of POD basis elements are needed in general to describe the solution [7]. Consequently, the POD method results in a decrease of computational time [1]. In this section, we give a generM overview of the POD method in the context of inverse problems. For further details on the general POD method, we refer the reader to [7-18] and the extensive list of references contained therein. The first step in forming the POD basis is to collect "snapshots" or solutions across time, space or a varied parameter. In our case, we let q be the vector parameter characterizing physical properties of the damage. In our computational efforts, we consider only rectangular damages, and thus, q determines the geometry of the damage including the length, thickness, depth, N, etc. of the damage. For an ensemble of damages {qJ}j=l, we obtain corresponding solutions, N. of the boundary value problem given by (4), for magnetic vector potentials which {A (qi)}j=1, we call our "snapshots". We seek basis elements of the form N, ¢' = Z Vi(j)A(q3)' (6) j----1

where coefficients Vi(j) are chosen such that each POD basis element ¢i, i = 1, 2 , . . . , Ns, maximizes 1 N, --N8~.= I(A(q3)' @i)L2(n'C)]2 subject to (@i, @~)L2(n,C) = II@iH2 = 1. Standard arguments guarantee that coefficients Vi(j) can be found by solving the eigenvalue problem cv

= ),v,

where the "covariance" matrix C is defined by [C]~j = ~ (A(q~), A(qj))L2(n,C ). Matrix C is a Hermitian positive semidefinite matrix, and thus, it possesses a complete set of orthogonal eigenvectors with corresponding nonnegative reM eigenvalues. We order the eigenvalues along with their corresponding eigenvectors such that the eigenvalues are in decreasing order, )h > ,~2 > ... > )~N. _>0.

Comparison of Two Techniques

557

We furthermore normalize the eigenvectors corresponding to the rule

¼.yj = Ns~" Then, the ith POD basis element is defined by (6) where Vi(j) represents the j t h component of the ith eigenvector of C. It can also be shown that tI ~ *Ji=l .1g~ are orthonormal in L2(fl, C) and span{(I)i}/N__*l = span{A(qj)};=" 1. Indeed, given any A(qj), we have N~

A(qj)

= ~ C~k(qj)~k, k---1

where ak(q/) = (A(q/), Ck)L2(n,C). To determine the reduced number N of POD basis elements required to accurately portray the N. we consider ensemble of "snapshots" {A (qJ)}j=l, N j--1 N,

'

(7)

j=l

which represents the percentage of "energy" in span{g(qj)}g=~ 1 that is captured in span{~j }N=I. The reduced basis consists of only the first N elements &~, i -- 1,... ,N, where N is chosen according to the percentage "energy" desired. We intuitively argue that the "energy" we are referring to is related to the total electrostatic energy. Simply stated, matrix C contains terms of the form

L AA da = L 'A'2 da, which can be written in terms of the electric field E according to

E = -iwA -

V¢. Therefore,

~ ,A,2 da = Cl ~ 'E + V¢,2 da, where C1 = 1/w 2. Since V¢ is piecewise constant (proved in [6]), the terms in matrix C are a perturbation of terms associated with electrostatic energy given by 1

f

WE = ~eo Jv IEI2dV. Thus, we conclude that when we snapshot on the magnetic vector potential, the ratio in (7) is a measure of the electrostatic energy stored across ~ (see [19]). Employing only the first N POD basis elements, we obtain the approximation AN(qj) for A(qj) such that N

A(%.) ,~ AN(q) =- ~-~ak(q~)¢k. k---1

To approximate AN(q) where q is a given parameter approximation formula to obtain

not in the

N, set {qJ}j=l, we must extend the

N

AN(q)

= Z ak(q)Ok.

(8)

k----1

Two possible ways of computing OLk(q) are by using a POD/Galerkin method or a POD/interpolation method. While the traditional way to find ak(q) is to use a POD/Galerkin method, there are many advantages in choosing a POD/interpolation method. We next examine both methods and compare the advantages and disadvantages of each.

558

M.L. JOYNER 4. POD/GALERKIN

APPROACH

The first approach we consider, the P O D / G a l e r k i n approach, is essentially an application of Galerkin's method to the integro-differential equation given by equation (4). T h e P O D / G a l e r k i n method uses the approximation given in (8) in the variational form of the integro-differential equation where the test functions are chosen to be the reduced-order P O D basis elements g The system then reduces to a linear system which we can solve for the coefficients (2k(q), k = 1,...,N. In our computational efforts reported on here, we follow the literature [20-24] and neglect the displacement current in the numerical implementation. Therefore, we will use the variational form (see [6] for details)

# Ox' "~x + \ , Oy -~y

+ (iwaA, ¢)

- A¢--~

s

s

= ~

Substituting (8) into (9) and letting ~ = ~l, l = 1 , . . . , N, we obtain the system

for (2 ---- [(21, (22,..., (2N]T, where

[K] k = [M]lk = ~

f, t lock ocz oak 7, do) , ; \ 0x a~k~Z da,

and

[b]t = Its ~l da. Next, we consider the accuracy of the POD approximation B N ( q ) = V x AN(q)

(11)

by comparing the approximation to the finite element simulation using Ansoft Maxwell 2D Field Simulator for a sample containing a damage of length 1 = 1.3 mm or equivalently q = l = 1.3 mm with a fixed thickness of 2 mm and depth of 9 mm. In the examples below, snapshots were taken of the magnetic vector potential for damages having lengths 1 = 0 m m to I = 4 m m in increments of 0.2 mm. Therefore, l = 1.3 mm is a length not included in the snapshots. T h e approximate energy captured in the POD approximation for values of N ranging from 1 to 4, according to (7), is given in Table 2. Table 2. Energy captured with N basis elements using snapshots of A while varying length. N 1 2 3 4

Energy Captured 99.999469% 99.999998% 99.999999% 99.999999%

Comparison of Two Techniques

559

However, as discussed previously, this energy is related to the electrostatic energy, and hence, not the energy stored in the magnetic field given by Wm = -~o

IBI 2 = " ~ o

B . f3 d V

as stated in [19]. Since we are comparing the reduced approximation for the magnetic flux density for various values of N, we can argue t h a t N should be chosen based on the total magnetic energy captured (or similar energy) as opposed to electrostatic energy. In other words, we should consider the energy captured if we were to snapshot on B2 (for example) instead of A. Table 3 gives the energy captured when taking snapshots of the B2 field. We note t h a t although the choice of N should intuitively be based on snapshots of B2, the approximation for B2 behaves in the same manner whether we snapshot on A first and then form B2 according to (11) or if we snapshot on B2 directly. Table 3. Energy captured with N basis elements using snapshots of B2 while varying length. N 1 2 3 4 5

Energy Captured 95.752844% 98.938760% 99.515414% 99.680555% 99.749219%

N 6 7 8 9 10

Energy Captured 99.789119% 99.822617% 99.853187% 99.871258% 99.888757%

Table 4. Relative error in the POD approximation using the traditional POD/Galerkin method depending on the value of N. N 1 2 3 4 5 6 7 8 9 10 11

Relative Error 59.97% 6.88% 5.99% 10.63% 11.00% 14.89% 12.10% 13.64% 15.89% 24.74% 24.41%

N 12 13 14 15 16 17 18 19 20 21

Relative Error 26.69% 21.74% 23.99% 65.00% 90.93% 120.77% 96.68% 80.10% 199.88% 199.05%

Figures 2-4 give plots comparing the reduced-order P O D approximation for various values of N related to the finite element simulation for length l = 1.3 mm. We first note t h a t N = 1 P O D basis element (Figure 2) does not approximate the finite element solution well at all; however, when we use N = 3 P O D basis elements (Figure 3), there is a considerable improvement in the approximation. Using only N -- 3 P O D basis elements, we were able to fairly accurately a p p r o x i m a t e the finite element solution which uses over 7000 finite elements in its approximation. However, after N = 3, as the value of N increases, the approximation worsens and induces more error as shown in Table 4. Initially, as stated in [1], the problem with the P O D / G a l e r k i n approximation was thought to be the conditioning of the linear system used to solve for the coefficients, given by equation (10). Indeed, the conditioning of the linear system progressively worsens as N increases; however, this is not the reason for the poor approximation. W h e n N is equal to 1, 2, 3, 4, or 8, as shown in Figures 2-4, the P O D approximation has the same shape as the finite element solution;

560

M . L . JOYNER P O D a p p r o x i m a t i o n w i t h N = 1 vs. Ansoft F E M simulation for d a m a g e w i t h length l = 1.3mm. x 10 -a 3.5~

I

a

I

•

I

•

I #

2.5

b.c m

T

't,

%

•

# # I l

• • •

2

v

%

E

%

1.5

0.5

00

5

10

15

20

25

1

315

30

L

40

45

50

x (in mm)

(a) P O D a p p r o x i m a t i o n w i t h N = 2 vs. Ansoft F E M simulation for d a m a g e w i t h length l = 1.3mm.

x 10 4

..--

3.5 •- - - -

Ansoft FEM Simulation POD Approx.

2.5 %

2

1.5 %%

0.5

0

0

5

10

15

20

25 X (in mm)

30

35

40

45

(b) Figure 2. Finite element simulations vs. P O D a p p r o x i m a t i o n s formed using t h e P O D / G a l e r k i n m e t h o d w i t h N -- 1 and N = 2 basis elements.

50

C o m p a r i s o n of T w o T e c h n i q u e s

561

P O D a p p r o x i m a t i o n w i t h N = 3 vs. A n s o f t F E M s i m u l a t i o n for d a m a g e w i t h l e n g t h l = 1.3mm. 3.5 X 10"-e

r

~

3

2.5

/

2

"%,

.~=_ 1.5 m

0.5

0

0

J.

L

5

10

15

20

25 x

30

35

40

45

50

(a) P O D a p p r o x i m a t i o n w i t h N = 4 vs. A n s o f t F E M s i m u l a t i o n for d a m a g e w i t h l e n g t h l ---- 1.3 m m . x 10 -s 3 . 5 ~

r

r

~

-

~

r

-

-

[- -

POD Approx.

2.5 qt

s

,,

f',<,,, 2

1.5

0.5

0

0

5

10

15

20

25 x (in mm)

30

35

40

45

(b) F i g u r e 3. F i n i t e e l e m e n t s i m u l a t i o n s vs. P O D a p p r o x i m a t i o n s f o r m e d u s i n g t h e POD/Galerkin method with N = 3 and N = 4 basis elements.

50

562

M . L . JOYNER P O D a p p r o x i m a t i o n with N -- 8 vs. Ansoft F E M simulation for d a m a g e w i t h length I --= 1.3 mm.

3.5 x 1-0"e

T

T

r

r

r

[- -

POD Approx.

p.1%%

2.5

2

1.5

t

"k,

r

-

0.5

0

0

5

10

15

20

25

30

35

40

45

50

x Ca) P O D a p p r o x i m a t i o n w i t h N = 21 vs. Ansoft F E M simulation for d a m a g e w i t h length 1 = 1.3 m m .

a.sxlo-'

2.5

2 b-A

m 1.5 E 1

0.5

00

5

10

15

20

25 X (in mm)

30

.L

L

35

40

45

(b) Figure 4. Finite element simulations vs. P O D a p p r o x i m a t i o n s formed using t h e P O D / G a l e r k i n m e t h o d w i t h N -- 8 a n d N = 21 basis elements.

50

C o m p a r i s o n of Two Techniques

563

however, the POD approximation overestimates the maximum amplitude of the plot. On the other hand, if one examines the case of N -- 21 more closely as shown in Figure 4, we notice that although there appears to be considerable "noise-like" behavior in the plot, both the shape and magnitude of the POD approximation, on average, is the same as that for the finite element solution throughout the entire graph. Hence, we are led to the idea that although most of the energy for the system is contained in the first few basis elements (as illustrated in Table 3), the energy and information contained in the last basis elements are essential to obtaining a POD approximation which approximates both the shape and magnitude of the solution. However, the noise-like behavior causes a large relative error as shown in Table 4, and therefore, would not be appropriate in an inverse problem. The cause of the "noise-like" behavior is due to the basis elements themselves. Figures 5 and 6 show the real and imaginary portions of all the basis elements, respectively. As N increases, the graphs of the POD basis elements have more irregular behavior, creating the "noise-like" behavior in the approximation. Therefore, to make the POD approximation more accurate, it makes sense to use all Ns = 21 basis elements in equation (10) and solve for all the coefficients ~k, k = 1 , 2 , . . . , N s = 21. However, to reduce the "noise-like" behavior, once the coefficients have been calculated, we will use only the first few of them in the approximation of the solution given by (8) N

A N(q) -- ~ O~k(q)Ok. k=l

Figure 7 is a plot of the POD approximation to the solution using all the basis elements to solve for the coefficients C~k, k = 1, 2 , . . . , Ns and then using only the first three basis elements (and

\/ (a) N ----1.

(b) N --- 2.

t (c) N -~ 3.

(d) N = 4.

(e) N -- 5.

(f) N = 6.

I

;L (h) N -- 8.

(g) N = 7.

I

(i) N --- 9.

(j) N ----10.

(k) N = 11.

(1) N = 12.

( m ) N = 13.

(n) N =

14.

I I

(o) N - - 1 5 .

(p) N = 1 6 .

(q) N = l T .

(r) N = 1 8 .

(s) N----19.

Figure 5. Real p a r t of P O D basis elements.

(t) N----20.

(u) N =

21.

564

M.L. JOYNER

/

F1

f

/l

vt

]

(a) N - - 1.

(b) N----2.

(c) N = 3 .

f

7

L/ (d) N - - 4 .

(e) N = 5 .

(f) N = 6 .

(g) N - - 7 .

/\ (h) N -- 8.

(i) N = 9.

(j) N

-- 10.

(k) N = 11.

(l) N :

12.

14.

( m ) N = 13.

(n) N =

(t) N = 20.

(u) N -- 21.

!

i (o) N = 15.

(p) N ----16.

(q) N ----17.

(r) N = 18.

(s) N ----19.

Figure 6. Imaginary part of POD basis elements. associated coefficients) in equation (8). approximation.

As shown, there is a remarkable improvement in the

Indeed, we can examine the reduction in error in the P O D approximation when we solve an m x m system for the first m coefficients and then use a reduced number, N , of coefficients in the approximation in (8). Table 5 gives the relative error between the P O D approximation using N = 3 basis elements in the formation of the solution when we vary the size of the linear system we solve, m, m -- 3 , . . . , Ns where m = N8 -- 21 is the graph given in Figure 7. There is not a whole lot of variation in the relative error until m is greater than 17. However, as suggested earlier, the best approximation is f o u n d when all the basis elements are used in solving for the coefficients, i.e., when we solve an N8 x N8 system for all Ns coefficients initially before reducing the number of basis elements used in the P O D approximation. Although N = 3 produced the best approximation with the traditional P O D / G a l e r k i n m e t h o d as noted earlier (see Table 4), it is of interest to see if the same result is true when we use the revised P O D / G a l e r k i n method, i.e., solving for all the coefficients first before choosing the reduced number for the approximation. The results are given in Table 6. Comparing Tables 4 and 6, it can be seen t h a t in almost every case, using all m = Ns of the basis elements to solve the linear system produced a more accurate result. However, as we increase the number of basis elements used in the P O D approximation, we start adding in the basis elements with the "noiselike" behavior (see Figures 5 and 6). When N is greater than six, this addition of basis elements results in a more rapid increase in the relative error between the finite element solution and P O D approximation. The best approximation is when N = 2 basis elements are used; however, there is not much variation in the relative error when using N -- 2, 3, or 4.

Comparison of Two Techniques

565

POD approximation with N -- 3; m = 21 vs. Ansoft FEM simulation for damage with length l ----1.3 mm. X 10-e 3.5 ~ " r

T

T

I

T

~

[- -

P O D Approx.

2.5

2

"~1.5

1

0.5

0

5

10

15

20

25 30 x (in ram)

35

40

45

50

Figure 7. Modified POD/Galerkin approach using N = 3 POD basis elements. Table 5. Relative error in the approximation depending on the size m of the linear system in (10). m

Relative Error

m

Relative Error

5.99%

13

5.36%

4

9.75%

14

5.14%

5

10.46%

15

5.37%

6

9.13%

16

3.84%

7

7.13%

17

1.58%

8

7.19%

18

1.34%

3

9

7.68%

19

1.44%

i0

5.62%

20

0.64%

11

5.16%

21

0.56%

12

7.49%

We have shown that we can obtain a better P O D approximation by using a revised P O D / G a l erkin approach instead of the standard P O D / G a l e r k i n approach. However, we need to address the impact of this alteration on the damage detection problem discussed in this paper. We noted earlier that we initially considered the reduced-order P O D methodology in order to create a method which could run in real-time, allowing the possibility of a portable scanning device. The traditional P O D / G a l e r k i n method, although extremely fast due to the small system which needed to be solved, did not produce accurate enough results for use in the inverse problem. On the other hand, the revised P O D / G a l e r k i n approach produced much more accurate results and would suffice for use in the inverse problem. However, the time required to solve the full linear system with 21 basis elements was still on the order of one to two minutes for each forward calculation. Although this is a reduction in time over the 5-7 minutes required to solve the

566

M. L. JOYNER

Table 6. Relative error in the approximation depending on N in approximation (8) using m = 21 for solving the linear system in (10).

6

1.12%

N 12 13 14 15 16 17

7

2.39%

18

86,16%

8

5.06% 7.09% 17.72% 18.33%

19 20 21

64.45% 200.20% 199.05%

N 1 2 3 4 5

9

10 II

Relative Error 52.06% 0.34% 0.56% 0.87% 1.34%

Relative Error 21.99% 21.03% 21.06% 46.23% 84.73% 116.78%

boundary value problem using approximately 7000 finite elements (the number when using the commercial software Ansoft 2D Field Simulator), the reduction in time is not large enough to produce real-time results of the inverse problem. In the next section, we will discuss the use of a P O D / i n t e r p o l a t i o n method to find the coefficients used in the P O D approximation which requires even less computational time than the P O D / G a l e r k i n method, and hence, is a viable alternative to the traditional P O D / G a l e r k i n approach.

5. P O D / I N T E R P O L A T I O N M E T H O D Unlike the P O D / G a l e r k i n method, the POD/interpolation method does not take into account the boundary value problem which A satisfies. It relies entirely on the values of the coefficients N~ upon which we snapshot. Using linear interpolation for the ak(q) for q in the set {qJ}j=l one-parameter case as considered in the previous section, i.e., for q = q (the scalar case), the coefficients c~k(q) are defined by , ,, q - q j (~k(q) ------~k(qj) + (ak(qj+l) - O l k ( q j ) ) q ~ - 1 --q , j -- j

(12)

where there exists j in { 1 , . . . , Ns - 1} such that qj < q < qj+l. Doing the same comparison as we did for the P O D / G a l e r k i n method, we snapshot on length and compare the P O D approximation using the POD/interpolation m e t h o d to the finite element solution for a length of I = 1.3 mm as depicted in Figures 8 and 9. W i t h N = 2 basis elements, there is still some visible error in the approximation, but using N = 3 P O D basis elements with POD/interpolation, the approximation appears to lie directly over the finite element simulation. Table 7 gives the relative error in the POD approximation using the P O D / i n t e r p o l a t i o n method. The best approximation is found when using N = 5, although the relative error is less than 1% for any N above four. We can compare the relative error when using the P O D / i n t e r p o l a t i o n m e t h o d (Table 7) to the error using the revised P O D / G a l e r k i n method (Table 6). First note t h a t the overall best approximation is found when using the revised P O D / G a l e r k i n method for approximating the coefficients with two basis elements in the approximation. However, the approximation is quite comparable to any of the POD/interpolation approximations with N larger than four. T h e more notable difference is that the error in the approximation when using the P O D / i n t e r p o l a t i o n method does not continually worsen as N increases as it does when using the revised P O D / G a l e r k i n method. At first glance, this would seem to contradict the previous argument we made, i.e., the argument t h a t adding in the extra basis elements creates a "noise-like" behavior causing higher relative error in the approximation. However, if one compares the coefficients found using the P O D / G a l e r k i n

C o m p a r i s o n of T w o T e c h n i q u e s

567

P O D a p p r o x i m a t i o n w i t h N = 1 vs. A n s o f t F E M s i m u l a t i o n for d a m a g e w i t h l e n g t h l = 1.3 m m .

3,5X107_~ ~ - ~ Ansoft FEM Simulation L" " POD Approx.

I

I I

I

I I

2"51J

$

mmu

2 tO

E

1.5

0.5

0 ~

0

5

10

15

20

25 x

30

35

40

45

50

(a) P O D a p p r o x i m a t i o n w i t h N = 2 vs. A n s o f t F E M s i m u l a t i o n for d a m a g e w i t h l e n g t h l = 1.3 ram.

3.5

x 104

1

~

T

~

T

5

10

15

20

25 X

2.5

2

1.5

0.5

0

30

35

40

45

(b) F i g u r e 8. F i n i t e e l e m e n t s i m u l a t i o n s vs. P O D a p p r o x i m a t i o n s f o r m e d u s i n g t h e P O D / i n t e r p o l a t i o n m e t h o d w i t h N ---- 1 a n d N = 2 b a s i s e l e m e n t s .

50

568

M . L . JoYNER P O D a p p r o x i m a t i o n w i t h N = 3 vs. A n s o f t F E M s i m u l a t i o n for d a m a g e w i t h l e n g t h l -----1 . 3 m m . 3.5 x 1 0 4 i -[- -

,-,,,~v,, r ~,v, ,.~,,,u,=,,v,, POD Approx.

2.5

2

1.5

0.5

0

0

5

10

15

20

25

30

35

40

45

50

X

(a) P O D a p p r o x i m a t i o n w i t h N = 4 vs. A n s o f t F E M s i m u l a t i o n for d a m a g e w i t h l e n g t h l = 1.3mm. 3.5

x 10 -s 1

T

I" "

POD Approx.

2.5

2

1,5

0.5

0

0

5

10

15

20

25

30

35

40

45

X

(b) F i g u r e 9. Finite e l e m e n t s i m u l a t i o n s vs. P O D a p p r o x i m a t i o n s f o r m e d u s i n ~ t h e POD/interpolation method with N = 3 and N = 4 basis elements.

50

Comparison of Two Techniques

569

Table 7. Relative error in the P O D approximation using t h e P O D / i n t e r p o l a t i o n method depending on the value of N. N

Relative Error

N

Relative Error

1

51.69%

12

0.72%

2

1.95%

13

0.51%

3

1.14%

14

0.52%

4

0.63%

15

0.56%

5

0.55%

16

0.57%

6

0.79%

17

O.7O%

7

0.79%

18

0.70%

8

0.72%

19

0.74%

9

0.72%

20

0.76%

i0

0.72%

21

0.91%

11

0,74%

Table 8. Comparison of P O D coefficients a k using P O D / G a l e r k i n method vs. P O D /interpolation method (order 10-13). k

P O D / G a l e r k i n Method

P O D / I n t e r p o l a t i o n Method

1

4 9 7 2 1 . 6 2 1 0 - 160.1739i

49572.9504 - 28.7329i

2

-66.3695 - 19.7766i

-64.6543 - 2.2641i

3

-2.4842 - 7.8178i

-3.5192 - 0.0383i

4

2.4126 + 2.0794i

1.4691 + 0.0230i

5

1.5355 - 2.9334i

-0.2024 - 0.0034i

6

0.7200 - 1.8792i

-0.1521 - 0.0034i

7

-2.4552 - 2.2961/

0.0508 + 0.0023i

8

-0.0316 ÷ 2.3680i

-0.0747 + 0.0007i

9

-1.1865 - 1.4410i

0.0082 - 0.0060i

10

-3.8374 - 1.7443i

-0.0181 + 0.0011i

11

-0.9876 - 0.8287i

-0.0890 - 0.0314i

12

-1.0773 - 3.0898i

-0.0389 -{- 0.0081/

13

-6.5288 + 0.6933i

0.0052 -I- 0.0257/ -0.0129 -{- 0.0512i

14

-0.1199 + 0.2753i

15

8.5887 + 1.9499/

16

5.0858 - 10.8746i

0.0842 -}- 0.0386i

17

9.2571 - 7.4390i

-0.0019 - 0.0025i

18

0.3887 - 2.6529i

-0.0289 + 0.0097i

19

-2.6303 + 2.5024i

0.0020 ÷ 0.0016i

20

-7.1746 - 7.0364i

0.0106 - 0.0110i

21

0.6766 + 0.3035i

0.0009 - 0,0017i

method to those found using the POD/interpolation

-0.0000

-

0.0059i

m e t h o d a s d o n e i n T a b l e 8, t h e r e is a n o t a b l e

d i f f e r e n c e i n t h e v a l u e o f t h e c o e f f i c i e n t s . T h e l a r g e s t c o e f f i c i e n t is ~1 s i n c e t h e f i r s t b a s i s e l e m e n t c o n t a i n s t h e m o s t e n e r g y o f t h e s y s t e m . A f t e r k = 1, t h e o r d e r o f t h e c o e f f i c i e n t d r o p s off s i g n i f i c a n t l y s i n c e less e n e r g y is b e i n g a d d e d i n t o t h e s y s t e m w i t h e a c h a d d i t i o n a l b a s i s e l e m e n t . the POD/interpolation

method

With

t h e o r d e r c o n t i n u e s t o d r o p off a s k c o n t i n u a l l y i n c r e a s e s . T h e r e -

fore, v e r y l i t t l e i n f o r m a t i o n is a d d e d i n f r o m t h e l a t t e r b a s i s e l e m e n t s , a n d h e n c e , less " n o i s e - l i k e " b e h a v i o r is b e i n g a d d e d t o t h e o v e r a l l a p p r o x i m a t i o n . POD/Galerkin

method

e l e m e n t is c o m p a r a b l e

However, the coefficients found using the

do not behave in the same manner. to contribution

Contribution

from the third basis

f r o m t h e 20 th b a s i s e l e m e n t . A s a c o n s e q u e n c e , t h e n o i s e -

like b e h a v i o r is b e i n g a d d e d i n a t a m o r e s u b s t a n t i a l l e v e l u s i n g t h e P O D / G a l e r k i n

method

than

570

M.L. JOYNER

when using the P O D / i n t e r p o l a t i o n method. The reason for this difference is not clear although one idea is t h a t we are solving a poorly conditioned linear system for coefficients on the order of 10 -13. Nonetheless, b o t h methods, the revised P O D / G a l e r k i n m e t h o d and P O D / i n t e r p o l a t i o n method, give a good approximation to the finite element solution if implemented as discussed above.

6. C O N C L U S I O N In this paper, we began by giving a brief overview of the eddy current m e t h o d and the corresponding model. Next, we presented an overview of the P O D method. We then discussed the reduced-order P O D methodology and focused on two different methods, the P O D / G a l e r k i n method and P O D / i n t e r p o l a t i o n method, for calculating the coefficients used in the reduced-order P O D approximation. We began by discussing the typical P O D / G a l e r k i n m e t h o d in which the weak form of the b o u n d a r y value problem was used to find the first N coefficients of the approximation. This was done by substituting the reduced-order approximation into the weak form and using the first N P O D basis elements as the basis elements in the weak form. An N × N linear system resulted which could be solved for the first N coefficients. However, using the P O D / G a l e r k i n m e t h o d as described resulted in significant relative error in the approximation when compared to the finite element solution. Although using all the basis elements to solve the system resulted in an extremely high relative error, the graph seemed to indicate t h a t the shape and m a g n i t u d e of the resulting approximation was on the same order as t h a t of the finite element solution. T h e "noise-like" behavior was believed to be responsible for the poor approximation. Therefore, we formulated a revised P O D / G a l e r k i n m e t h o d in which all 21 basis elements were used to solve for the coefficients and t h e n only the first few basis elements were used in the approximation. This resulted in an approximation with the same shape and m a g n i t u d e without the "noise-like" behavior. Hence, the error in the approximation was greatly reduced. However, using all 21 basis elements to find the coefficients could still be time consuming. Hence, we introduced the P O D / i n t e r p o l a t i o n method which was b o t h easier to implement and required less computational time. T h e P O D / i n t e r p o l a t i o n m e t h o d resulted in approximations comparable to the best approximation found when using the revised P O D / G a l e r k i n method. However, one distinct advantage of the P O D / i n t e r p o l a t i o n method is t h a t it does not rely on the equations describing the system. This can be very useful in some experimental applications in which d a t a is available but it is not easy to model the physical process corresponding to the data. In addition, for real-time inverse problem applications, such as the eddy current d a m a g e detection problem, the ease and computational savings of the P O D / i n t e r p o l a t i o n m e t h o d make this m e t h o d a wise choice. In either case, the reduced-order P O D m e t h o d is a an excellent reduced-order approximation to the finite element solution and since it is a reduced-order methodology, it will result in reduced computational time for the appropriate inverse problem.

REFERENCES 1. H.T. Banks, M.L. Joyner, B. Wincheski and W.P. Winfree, Real-time compuational algorithms for eddycurrent-based damage detection, Inverse Problems 18, 795-823, (2002). 2. H.T. Banks, M.L. Joyner, B. Wincheski and W.P. Winfree, Nondestructive evaluation using a reduced-order computational methodology, Inverse Problems 16, 929-945, (2000). 3. R.L. Stoll, The Analysis of Eddy Currents, Clarendon Press, Oxford, (1974). 4. Ansoft Corporation, Maxwell 2D Field Simulator--Technical Notes, (1995-1999). 5. D.K. Cheng, Field and Wave Electromagnetics, Second Edition, Addison-Wesley, Reading, MA, (1992). 6. M.L. Joyner, An application of a reduced order computational methodology for eddy current based nondestructive evaluation techniques, Ph.D. Thesis, North Carolina State University, (2001). 7. H.V. Ly and H.T. Tran, Proper orthogonal decomposition for flow calculations and optimal control in a horizontal CVD reactor, CRSC Tech. Rep. CRSC-TR98-13, North Carolina State University, (1998).

Comparison of Two Techniques

571

8. H.T. Banks, R.C. delRosario and R.C. Smith, Reduced order model feedback control design: Numerical implementation in a thin shell model, IEEE Trans. Auto. Control 45, 1312-1324, (July 2000). 9. G. Berkooz, Observations on the proper orthogonal decomposition, In Studies in Turbulence, pp. 229-247, Springer-Verlag, New York, (1992). 10. G. Berkooz, P. Holmes and J.L. Lumley, The proper orthogonal decomposition in the analysis of turbulent flows, Annual Review of Fluid Mechanics 25 (5), 539-575, (1993). 11. R.C. delRosario, Computational methods for feedback control in structural systems, Ph.D. Thesis, North Carolina State University, (1998). 12. K. Karhunen, Zur spektral theorie stochasticher prozesse, Ann. Acad. Sei. Fennicae 37 (A1), (1946). 13. M. Kirby, J.P. Boris and L. Sirovich, A proper orthogonal decomposition of a simulated supersonic shear layer, International Journal for Numerical Methods in Fluids 10, 411-428, (1990). 14. M. Kirby and L. Sirovich, Application of the Karhunen-Loeve procedure for the characterization of human faces, IEEE Transactions on Pattern Analysis and Machine Intelligence 12 (1), 103-108, (1990). 15. K. Kunisch and S. Volkwein, Control of Burgers' equation by a reduced-order approach using proper orthogonal decomposition, J. Optimization Theory and Applic. 102 (2), 345-371, (1999). 16. M. Loeve, Functions aleatoire de second ordre, Compte Rend. Aead. Sci. Paris, (1945). 17. J.L. Lumley, The structure of inhomogeneous turbulent flows, Atmospheric Turbulence and Radio Wave Propagation, 166-178, (1967). 18. J.L. Lumley, Stochastic Tools in Turbulence, Academic Press, New York, (1970). 19. R.S. Elliott, Electromagnetics: History, Theory, and Applications, IEEE Press, New York, (1993). 20. A. Konrad, Integrodifferential finite element formulation of two-dimensional steady-state skin effect problems, IEEE Transactions on Magnetics 18 (1), 284-292, (1981). 21. A. Konrad, The numerical solution of steady-state skin effect problems--An integrodifferential approach, IEEE Transactions on Magnetics 17 (1), 1148-1152, (1981). 22. A. Krawczyk and J.A. Tegopoulos, Numerical Modelling of Eddy Currents, Oxford University Press, Oxford, (1993). 23. E.E. Kriezis, T.D. Tsiboukis, S.M. Panas and J.A. Tegopoulos, Eddy currents: Theory and applications, Proceedings of the IEEE 80 (10), 1559-1589, (1992). 24. J. Weiss and Z.J. Csendes, A one-step finite element method for multiconductor skin effect problems, IEEE Transactions on Power Apparatus and Systems 101 (10), 3796-3803, (1982).