Baker-Atlas Report (May 04)

FULLY AUTOMATIC GOAL-ORIENTED hp-ADAPTIVITY FOR ELLIPTIC PROBLEMS D. Pardo∗, L. Demkowicz∗ ∗

Institute for Computational Engineering and Sciences The University of Texas at Austin

Abstract In this document, we present a Finite Element (FE) self-adaptive hp goal-oriented algorithm for elliptic problems. The algorithm delivers a sequence of optimal hp-grids. This sequence of grids minimizes the error of a prescribed quantity of interest with respect to the problem size. The refinement strategy is based on a fully automatic energy norm based hp-adaptive algorithm. We illustrate the efficiency of the method with 2D numerical results.

Key words: hp-Finite Elements, automatic adaptivity, goal-oriented adaptivity. AMS subject classification: 65N30, 35L15.

Acknowledgment A. Bespalov, L. Tabarovski, and C. Torres-Verdin.

1

1

Introduction

During the last decades, different algorithms intended to generate optimal grids for solving relevant engineering problems have been designed and implemented. Among those algorithms, a Finite Element (FE) self-adaptive hp-refinement strategy has recently (2001) been developed in the Institute for Computational Engineering and Sciences (ICES) at The University of Texas at Austin. The strategy produces automatically a sequence of hp-meshes that delivers exponential convergence rates in terms of the energy norm error with respect to the number of unknowns (for problems with and without singularities), providing high accuracy approximations of solutions corresponding to a variety of engineering applications. Furthermore, the self-adaptive strategy is problem independent, and can be applied to FE discretizations of H 1 -, H(curl)-, and H(div)-spaces, as well as to nonlinear problems (see [8, 17] for details). But the energy norm is a quantity of very little relevance for most engineering applications, especially when a particular objective is pursued, as for example, to simulate the electromagnetic behavior of a petroleum engineering logging tool in a borehole environment. Thus, a self-adaptive strategy intended to approximate a specific feature of the solution is needed. This type of refinement strategies are called goal-oriented adaptive algorithms [10, 15], and are based on minimizing the error of a prescribed quantity of interest mathematically expressed in terms of a linear functional. In this paper, we formulate and study numerically a number of self-adaptive goal-oriented algorithms intended to solve elliptic problems. These algorithms are an extension of the fully automatic (energy norm based) hp-adaptive strategy described in [8, 17]. In a forthcoming paper, a goal-oriented algorithm for solving Helmholtz and Maxwell’s equations will be presented. In Section 2, we introduce a family of goal-oriented h and hp self-adaptive algorithms. Section 3 is devoted toward presenting the most relevant implementation details. Notice that we intend to use not only the main ideas of energy norm based adaptivity, but (possibly) also the same lines of code. Comparative 2D numerical results illustrating the efficiency of each of the self-adaptive algorithms are shown in Section 4, along with the corresponding final h or hp-grids. Finally, in Section 5 we summarize the main conclusions, and we outline the future lines of research.

2

2

Formulation

We are interested in solving a variational problem in the standard form:   Find u ∈ V

 b(u, v) = f (v)

(2.1)

∀v ∈ V .

Here

• V is a Hilbert space. • f ∈ V 0 is a linear and continuous functional on V . • b is a bilinear and symmetric form (assumed to be coercive and continuous on space V ). Thus, we can define an inner product on V as (u, v) := b(u, v), with the corresponding norm denoted by kuk. Given an hp-FE subspace Vhp ⊂ V , we discretize (2.1) as follows:   Find uhp ∈ Vhp

 b(u , v ) = f (v ) hp hp hp

(2.2)

∀vhp ∈ Vhp .

The objective of goal-oriented adaptivity is to construct an optimal hp-grid, in the sense that it minimizes the problem size needed to achieve a given tolerance error in a quantity of interest L ∈ V 0 . By linearity of L, we have: E = L(u) − L(uhp ) = L(u − uhp ) = L(e) ,

(2.3)

where e = u − uhp denotes the error function. Defining residual rhp ∈ V 0 as rhp (v) = f (v) − b(uhp , v), we relate the quantity of interest error with the residual by seeking for the solution of the following problem (if solution exists):   Find w ∈ V

(2.4)

 L(e) = r (w) . hp

Or equivalently:

  Find w ∈ V

 b(e, w) = L(e) .

(2.5) 3

The above problem is necessarily solved if we find solution of the following dual problem:   Find w ∈ V

 b(v, w) = L(v)

(2.6)

∀v ∈ V .

Using the Lax-Milgram theorem we conclude that problem (2.6) has a unique solution in V , and therefore, (2.4) is solvable. Solution w is usually referred as the influence function. Discretizing (2.6) using, for example, Vhp ⊂ V , we obtain:   Find whp ∈ Vhp

 b(v , w ) = L(v ) hp hp hp

(2.7)

∀vhp ∈ Vhp .

Notice that due to the symmetry of bilinear form b and the use of same space Vhp for solving both (2.2) and (2.7), it is only necessary to factorize the system of linear equations once. Thus, the extra cost of solving (2.7) reduces to only one backward and forward substitution. By orthogonality of e with respect Vhp (in the b-inner product), we have b(e, vhp ) = 0 for all vhp ∈ Vhp . Defining ² = ²(vhp ) = w − vhp , we obtain : E = L(e) = b(e, ²) .

(2.8)

Once error E in the quantity of interest has been determined in terms of the bilinear form, we wish to obtain a sharp upper bound for | E | that utilizes only local and computable quantities. Then, a self-adaptive algorithm intended to minimize this bound will be defined. First, using a similar procedure to the one of [8], we approximate u and w by finer grid functions u h , p+1 , w h , p+1 , that have been obtained by solving iteratively the corresponding 2 2 linear systems of equations associated to the FE subspace V h , p+1 . In the remainder of this 2 paper, u and w will denote the fine grid solutions of the direct and dual problems (u = u h , p+1 , 2 and w = w h , p+1 respectively), and we will restrict ourselves to discrete FE spaces only. 2

Next, we bound the error in the quantity of interest as a sum of element contributions. Let bK be defined as bK (uK , vK ) = b(uK , vK ) for all uK , vK ∈ VK , where VK ⊂ V . Then: | E |=| b(e, ²) |≤

X

| bK (e, ²) | ,

(2.9)

K

where summation over K indicates summation over elements. Inner product and norm associated to bilinear form bK will be denoted as (·, ·)K and k · kK respectively.

4

Projection Based Interpolation Operator. At this point, we introduce the projectionbased interpolation operator Πhp : V −→ Vhp defined in [7], and used in [8, 17] for construction of the fully automatic energy norm based hp-adaptive algorithm. We define Php : V −→ Vhp as the b-projection, and we denote uhp = Php u. Then, (2.9) becomes: | E |≤

X

| bK (e, ²) |=

K

X

| bK (u − Πhp u, ²) + bK (Πhp u − Php u, ²) | .

(2.10)

K

Given an element K, it is expected for | bK (Πhp u − Php u, ²) | to be negligible compared to | bK (u − Πhp u, ²) |. Under this assumption, we conclude the following. | E |≤

X

| bK (u − Πhp u, ²) | .

(2.11)

K

In particular, for ² = w − Πhp w, | E |≤

X

| bK (u − Πhp u, w − Πhp w) | .

(2.12)

K

Applying the parallelogram law to the last equation, we obtain the next upper bound for | E |: | E |≤

1X |k e˜ + ²˜ k2K − k e˜ − ²˜ k2K | , 4 K

(2.13)

where e˜ = u − Πhp u, and ²˜ = w − Πhp w. A second upper bound can be obtained from Cauchy-Schwartz inequality as follows: | E |≤k e˜ kK k ²˜ kK .

(2.14)

Now, we define two families of self-adaptive goal oriented algorithms: one for h-FE spaces, and a second one for hp-FE spaces.

2.1

A Self-Adaptive Goal-Oriented h-Refinement Algorithm

This adaptive algorithm iterates along the following steps. 1. Solve the direct and dual problems in a given h grid and the corresponding globally h refined grid (h/2), to obtain uh , wh , uh/2 , and wh/2 .

5

2. For each element K, compute the following quantity as an indicator of the error in goal. |k (uh/2 − uh ) + (wh/2 − wh ) k2K − k (uh/2 − uh ) − (wh/2 − wh ) k2K |

(2.15)

3. Refine elements that contribute with 66% of the maximum error, performing anisotropic refinements if the error is almost one-dimensional. Notice that for h-adaptive algorithms, only an error indicator (one number per element) is neccesary to decide which elements should be refined. Thus, there is no need for the use of the projection based interpolation operator. On the other side, for hp-adaptive algorithms the full error function (not just one number) is needed to decide between different h and p refinements. For these hp-algorithms, the projection based interpolation operator is utilized.

2.2

A Self-Adaptive Goal-Oriented hp-Refinement Algorithm

In this section, we define two goal-oriented hp self-adaptive algorithms that utilize the main ideas of the fully automatic (energy-norm based) hp-adaptive algorithm presented in [8, 17]. Thus, we start by recalling the main objective of the self-adaptive (energy norm based) hp-refinement strategy, which consists of solving the following maximization problem.  ˜  Find an optimal hp-grid in the following sense:    X  ˜ = arg max  hp   b K hp

2 2 K | u − ΠK hp u |1,K − | u − Πhp b u |1,K

∆N

,

(2.16)

where

• u = u h , p+1 is the fine grid solution, 2

c • ∆N > 0 is the increment in the number of unknowns from the grid hp to the grid hp, and

• | · |1,K is the H 1 -seminorm for element K.

6

Similarly, for goal-oriented hp-adaptivity, we propose the following two algorithms based on estimates (2.13) and (2.14).  ˜  Find an optimal hp-grid in the following sense:     2 K 2   X | | (u + w) − ΠK  hp (u + w) |1,K − | (u − w) − Πhp (u − w) |1,K |

˜ = arg max hp

       

b hp

K

−

(2.17) ∆N 2 K 2 | | (u + w) − ΠK b (u + w) |1,K − | (u − w) − Πhp b (u − w) |1,K | hp ∆N

 ˜  Find an optimal hp-grid in the following sense:     K   X | u − ΠK  hp u |1,K · | w − Πhp w |1,K

˜ = arg max hp

       

b hp

K

−

∆N K K | u − Πhp b u |1,K · | w − Πhp b w |1,K ∆N

(2.18) ,

where

• u = u h , p+1 and w = w h , p+1 are the fine grid solutions corresponding to the direct and 2 2 dual problems, c • ∆N > 0 is the increment in the number of unknowns from the grid hp to the grid hp, and

• | · |1,K is the H 1 -seminorm for element K. Implementation of both goal-oriented hp-adaptive algorithms is based on the optimization procedure used for energy norm hp-adaptivity [17].

7

3

Implementation Details

In this section, we briefly discuss the main implementation details needed to extend the fully automatic (energy norm based) hp-adaptive algorithm to a fully automatic goal-oriented hp-adaptive algorithm. 1. First, solution w of the dual problem in the fine grid is needed. This goal can be attained either by using a direct (frontal) solver or an iterative (two grid) solver (see [12]). 2. Then, we should treat both solutions as corresponding to two different partial differential equations (PDE’s). In the case of algorithm (2.17), we select functions u + w and u − w as the solutions of the system of two (undefined) PDE’s. In the case of algorithm (2.18), u and w will be the corresponding solutions. 3. Then, we redefine evaluation of the error. For algorithm (2.17), we replace the H 1 seminorm error evaluation of a two dimensional function (| u1 − Πhp u1 |21 + | u2 − Πhp u2 |21 ) by the following quantity: | | (u + w) − Πhp (u + w) |21 − | (u − w) − Πhp (u − w) |21 |. For algorithm (2.18), we replace the H 1 -seminorm error evaluation of a two dimensional function by | u − Πhp u |1 · | w − Πhp w |1 . 4. After these small modifications, the energy norm based adaptive algorithm may now be utilized as a fully automatic goal-oriented hp-adaptive algorithm.

8

4

Numerical Results

In this section, we present a comparative study via numerical experimentation of different adaptive strategies, when a particular feature of the solution is pursued. More precisely, we are interested in minimizing the error of the FE solution at a given point in a smooth part of the domain. Typically, special averaging functions (called mollifiers) are used to define a continuous linear functional L ∈ V 0 utilized as the quantity of interest (see [15]) to be minimized. Also, especial post-processing techniques may be used. Here, for simplicity, we define the following quantity of interest. 1 Z e(x)dx , L(e) = |Ωi | Ωi

(4.19)

where Ωi is a small subdomain containing target point xi at which approximation of the solution is desired. Two elliptic problems are used to perform this numerical study: an L-shape domain problem, and an orthotropic heat conduction problem in a thermal battery. For each model problem, we describe the geometry, governing equations, material coefficients, and boundary conditions. We also display the exact or approximate solution, and we briefly explain the relevance of each problem in this research. Then, results comparing different adaptive strategies are displayed.

4.1

L-Shape Domain Problem

In this problem, we want to solve Laplace equation (−∆u = 0) with Dirichlet Boundary Conditions (possibly, non-homogeneous) corresponding to the exact solution u = r 2/3 sin(2θ/3 + π/3), expressed in terms of cylindrical coordinates. Geometry ([−1, 0]×[0, 1]+[0, 1]×[0, 1]+ [0, 1]×[−1, 0]) of the computational domain along with the exact solution are shown in figure 1. This problem has a corner singularity located at the origin, and we are interested in the FE solution at the proximity of point (−1/2, 1/2). More precisely, we select functional L given by (4.19) with Ωi = [−5625, −0.4375] × [0.4375, 0.5625] as our quantity of interest. Figure 2 displays convergence history for both the pointwise error and the quantity of interest error, by using different self-adaptive algorithms. Figure 3 shows convergence history

9

Figure 1: Exact solution of the L-shape domain problem of the upper bound U of | E |, given by:  X  | u − uhp |21,K h-goal: U =     K  X  2   | | (u + w) − ΠK  hp-goal 1 (2.17): U = hp (u + w) |1,K K

2   − | (u − w) − ΠK  hp (u − w) |1,K |   X   K  | u − ΠK  hp u |1,K · | w − Πhp w |1,K  hp-goal 2 (2.18): U =

(4.20)

K

Final hp-grids for different self-adaptive algorithms are shown in figures 4, 5,6, and 7.

4.2

Orthotropic Heat Conduction Battery

We present now a Sandia1 benchmark problem, in which we solve the heat equation in a thermal battery with large and orthotropic jumps in the material coefficients (up to six orders of magnitude). The computational domain of dimensions [0, 8.4] × [0, 24], as well as an approximation to the unknown exact solution are displayed in figure 8. This problem is governed by the heat conduction equation ∇(K∇u) = f (k) , where K is given by: K=K 1

(k)

=

"

Kx(k) 0 0 Ky(k)

#

.

(4.21)

Sandia National Laboratories, USA

10

−3

−3

10

10 h−energy h−goal hp−energy hp−goal 1 hp−goal 2

−4

10

−5

−5

10

10

−6

−6

10

10

−7

−7

10

10

−8

−8

10

10

−9

−9

10

10

−10

−10

10

10

−11

−11

10

10

−12

−12

10

10

−13

10

h−energy h−goal hp−energy hp−goal 1 hp−goal 2

−4

10

−13

216

512

1000

1728

2744

4096

5832

10

8000

216

512

1000

1728

2744

4096

5832

8000

Figure 2: L-shape domain problem. Pointwise convergence (left) and convergence of quantity of interest L(e) (right) using different automatic refinement strategies. ’hp-goal 1’ and ’hpgoal 2’ correspond to hp-algorithms (2.17) and (2.18) respectively. For each material k, we define:

f (k)

  0,       1,

k=1 k=2 =  1, k=3   0, k=4     0, k=5

  25,       7,

Kx(k) =       

k=1 k=2 5, k=3 0.2, k=4 0.05, k=5

       

Ky(k) =       

25, 0.8, 0.0001, 0.2, 0.05,

k=1 k=2 k=3 k=4 k=5

(4.22)

Ordering each of the four sides of the boundary clockwise, and starting with the left hand side boundary, we impose the following boundary conditions: K∇u · n = g (i) − α(i) u ,

(4.23)

where

α(i)

 0,    

i=1 1, i=2 = 2, i=3    3, i=4

g (i)

 0,    

i=1 3, i=2 = . 2, i=3    0, i=4

(4.24)

This problem has several singularities of different strength, and we are interested in the FE solution at the proximity of point (8.2, 0.4), located at the southeast part of the domain. 11

−3

10

h−goal hp−goal 1 hp−goal 2

−4

10

−5

10

−6

10

−7

10

−8

10

−9

10

−10

10

−11

10

−12

10

−13

10

216

512

1000

1728

2744

4096

5832

8000

Figure 3: L-shape domain problem. Convergence history of upper bound U of | E | corresponding to different goal-oriented algorithms. More precisely, we select functional L given by (4.19) with Ωi = [8.0, 8.4] × [0.0, 0.8] as our quantity of interest. Figure 9 shows convergence history of the pointwise error and the quantity of interest error for different self-adaptive algorithms. Figure 10 shows convergence history of the upper bound U of | E |, given by (4.20). Final hp-grids (and their zooms) for different self-adaptive algorithms are shown in figures 11, 12, 13, 14, and 15.

12

2Dhp90: A Fully automatic hp-adaptive Finite Element code

y z x Figure 4: L-shape domain problem. Final h-grid containing 5177 unknowns (degrees of freedom) obtained by using the fully automatic energy norm based h-adaptive algorithm.

13

2Dhp90: A Fully automatic hp-adaptive Finite Element code

y z x Figure 5: L-shape domain problem. Final h-grid containing 6691 unknowns (degrees of freedom) obtained by using the fully automatic goal-oriented h-adaptive algorithm.

14

2Dhp90: A Fully automatic hp-adaptive Finite Element code

y z x Figure 6: L-shape domain problem. Final hp-grid containing 7392 unknowns (degrees of freedom) obtained by using the fully automatic energy norm based hp-adaptive algorithm.

15

2Dhp90: A Fully automatic hp-adaptive Finite Element code

y z x Figure 7: L-shape domain problem. Final hp-grid containing 6815 unknowns (degrees of freedom) obtained by using the fully automatic goal-oriented hp-adaptive algorithm.

16

k=1

k=2

k=3 k=5

k=4

Figure 8: Geometry (left) and FE solution (right) of the orthotropic heat conduction problem in a thermal battery.

−1

−1

10

10 h−energy h−goal hp−energy hp−goal 1 hp−goal 2

−2

10

−2

10

−3

−3

10

10

−4

−4

10

10

−5

−5

10

10

−6

−6

10

10

−7

10

h−energy h−goal hp−energy hp−goal 1 hp−goal 2

−7

64

216

512

1000

1728

2744

4096

5832

10

8000

64

216

512

1000

1728

2744

4096

5832

8000

Figure 9: Orthotropic heat conduction problem. Pointwise convergence (left) and convergence of quantity of interest L(e) (right) using different automatic refinement strategies. ’hp-goal 1’ and ’hp-goal 2’ correspond to hp-algorithms (2.17) and (2.18) respectively.

17

1

10

h−goal hp−goal 1 hp−goal 2 0

10

−1

10

−2

10

−3

10

−4

10

−5

10

64

216

512

1000

1728

2744

4096

5832

8000

Figure 10: Orthotropic heat conduction problem. Convergence history of upper bound U of | E | corresponding to different goal-oriented algorithms.

18

0: A Fully automatic hp-adaptive Finite Eleme

Figure 11: Orthotropic heat conduction problem. Final hp-grid containing 7353 unknowns (degrees of freedom) obtained by using the fully automatic energy norm based hp-adaptive algorithm.

19

0: A Fully automatic hp-adaptive Finite Eleme

Figure 12: Orthotropic heat conduction problem. Final hp-grid containing 6605 unknowns (degrees of freedom) obtained by using a fully automatic goal-oriented hp-adaptive algorithm.

20

2Dhp90: A Fully automatic hp-adaptive Finite Element code

2Dhp90: A Fu 2Dhp90 Fully y automatic u om hp-adaptive hp d p v F Finite n Element E m n code od

y z 2Dhp90 A Fu y u om

hp d p v F n

E m n

od

2Dhp90 A Fu y u om

hp d p v F n

E m n

od

y z

x 2Dhp90 A Fu y u om

hp d p v F n

E m n

od

2Dhp90 A Fu y u om

hp d p v F n

E m n

od

y

y x

x

y

y x

2Dhp90 A Fu y u om

hp d p v F n

E m n

x

x

od

2Dhp90 A Fu y u om

hp d p v F n

E m n

od

y

y x

x

Figure 13: Orthotropic heat conduction problem. Zooms toward a singularity located at point (18.8, 6.1) by factors of 1, 10, 100, and 10000 (from top to bottom) of the final h-grids obtained by using the fully automatic energy norm based (left) and goal-oriented based (right) h-adaptive algorithms. 21

2Dhp90: A Fully automatic hp-adaptive Finite Element code

2Dhp90: A Fu 2Dhp90 Fully y automatic u om hp-adaptive hp d p v F Finite n E Element m n code od

y z x 2Dhp90 A Fu y u om

hp d p v F n

E m n

y z x

od

2Dhp90 A Fu y u om

hp d p v F n

E m n

od

y

y x

2Dhp90 A Fu y u om

hp d p v F n

E m n

x

od

2Dhp90 A Fu y u om

hp d p v F n

E m n

od

y

y x

x

Figure 14: Orthotropic heat conduction problem. Zooms toward a singularity located at point (18.8, 6.1) by factors of 1 (top), 10 (middle), and 100 (bottom) of the final hp-grids obtained by using the fully automatic energy norm based (left) and goal-oriented based (right) hp-adaptive algorithms.

22

2Dhp90: A Fully automatic hp-adaptive Finite Element code

2Dhp90: A Fu 2Dhp90 Fully y automatic u om hp-adaptive hp d p v F Finite n E Element m n code od

y z x 2Dhp90 A Fu y u om

hp d p v F n

E m n

y z x

od

2Dhp90 A Fu y u om

hp d p v F n

E m n

od

y

y x

2Dhp90 A Fu y u om

hp d p v F n

E m n

x

od

2Dhp90 A Fu y u om

hp d p v F n

E m n

od

y

y x

x

Figure 15: Orthotropic heat conduction problem. Zooms toward a singularity located at point (18.8, 6.1) by factors of 1000 (top), 10000 (middle), and 100000 (bottom) of the final hp-grids obtained by using the fully automatic energy norm based (left) and goal-oriented based (right) hp-adaptive algorithms.

23

5

Conclusions

In this report, we have designed and implemented a set of goal-oriented self-adaptive algorithms. These algorithms produce automatically a sequence of optimal h or hp-grids, that minimizes the error of a prescribed quantity of interest with respect to the problem size. Preliminary numerical results indicate the superior convergence of these goal-oriented algorithms for elliptic problems, especially when hp-elements are used. A generalization of these algorithms for Helmholtz and Maxwell’s equations will be formulated and implemented in a forthcoming paper. For electromagnetic simulations of logging tools in the area of Petroleum Engineering, it is expected that it will be essential the use of goal-oriented adaptivity. Since fields decay exponentially as we move away from the source antenna/s (and thus, energy norm based adaptivity becomes meaningless) there is an increasing need for goal-oriented adaptive algorithms.

24

References [1] I. Babuˇska, B. Guo, ”Approximation properties of the hp-version of the finite element method”, Comp. Meth. Appl. Mech. Engg., 133 (1996), 319-346. [2] I. Babuˇska, T. Strouboulis, K. Copps, S.K. Gangaraj, C. Upadhyay, ”A-Posteriori Error Estimation for Finite Element and Generalized Finite Element Method”, TICAM Report, 98-01, 1998. [3] R. Becker, R. Rannacher, ”Weighted A-Posteriori Error Control in FE Method”, ENUMATH-95, Paris, September 1995. [4] R. Becker, R. Rannacher, ”A Feedback Approach to Error Control in Finite Elements Methods: Basic Analysis and Examples”, East-West J. Numer. Math 4, pp. 237-264, 1996. [5] L. Demkowicz, ”Fully Automatic hp-Adaptivity for Maxwell’s Equations”, TICAM Report 03-45, 2003. [6] L. Demkowicz, ”2D hp-Adaptive Finite Element Package (2Dhp90). Version 2.0”, TICAM Report 02-06, January 02. [7] L. Demkowicz, A. Buffa, ”H 1 , H(curl), and H(div) Conforming Projection-Based Interpolation in Three Dimensions: Quasi Optimal p-Interpolation Estimates”, ICES Report, 2004. [8] L. Demkowicz, W. Rachowicz, P. Devloo, ”A Fully Automatic hp-Adaptivity”, Journal of Scientific Computing, 17, 1-3, 127-155, 2002. [9] V. Heuveline, R. Rannacher, ”Duality-Based Adaptivity in the hp-Finite Element Method”, J. Numer. Math. 11 No. 2, pp. 95-103, 2003. [10] J.T. Oden, S. Prudhomme, ”Goal-Oriented Error Estimation and Adaptivity for the Finite Element Method”, Comp. Math. Appl. 41, pp. 735-756, 2001. [11] M. Paraschivoiu, A.T. Patera, ”A Hierarchichal Duality Approach to Bounds for the Outputs of Partial Differential Equations”, Comput. Methods Appl. Mech. Engrg. 158, pp. 398-407, 1998. [12] D. Pardo, ”Integration of hp-Adaptivity with a Two Grid Solver: Applications to Electromagnetics”, Ph.D. Dissertation, The University of Texas at Austin, 2004. 25

[13] D. Pardo, L. Demkowicz, ”Integration of hp-Adaptivity and Multigrid. I. A Two Grid Solver for hp-Finite Elements.”, TICAM Report 02-33, September 2002 [14] J. Peraire, A.T. Patera, ”Bounds for Linear-Functional Outputs of Coercive Partial Differential Equations: Local Indicators and Adaptive Refinement”, Advances in Adaptive Computational Methods in Mechanics(eds. P. Ladeveze, J.T. Oden), pp. 43-63, Amsterdam, 1998. [15] S. Prudhomme and J.T. Oden, ”On Goal-Oriented Error Estimation for Elliptic Problems: Application to the Control of Pointwise Errors”, Comput. Methods Appl. Mech. Engrg. 176, pp. 313-331, 1999. [16] W. Rachowicz, L. Demkowicz, ”An hp-Adaptive Finite Element Package for Electromagnetics. Part 1: Data structure and constrained approximation”, Computer Methods in Applied Mechanics and Engineering, 187 (1-2) (2000) pp. 307-335. [17] W. Rachowicz, D. Pardo, L.F. Demkowicz, ”Fully Automatic hp-Adaptivity in Three Dimensions”, ICES Report 04-22, 2004. [18] R. Rannacher, F.T. Stuttmeier, ”A Posteriori Error Control in Finite Element Methods Via Duality Techniques: Application to Perfect Plasticity”, Computational Mechanics 21, pp. 123-133, 1998. [19] P. Solin, L.F. Demkowicz, ”Goal-Oriented hp-Adaptivity for Elliptic Problems”, (TICAM Report 02-32) Comput. Methods Appl. Mech. Engrg., accepted for publication, 2002. [20] T. Walsh, L.F. Demkowicz, ”hp-Boundary Element Modeling of the External Human Auditory System – Goal-Oriented Adaptivity with Multiple Load Vectors”, Comput. Methods Appl. Mech. Engrg. 192, pp. 125-146, 2003.

26