A highly accurate high-order validated method to solve 3D Laplace - - PowerPoint PPT Presentation

A highly accurate high-order validated method to solve 3D Laplace equation

M.L.Shashikant, Martin Berz and Kyoko Makino Michigan State University, E Lansing, MI ,USA December 19, 2004

Approach We present an approach based on high-order quadrature and a high-order finite element method to find a validated solution of the Laplace equation when derivatives of the solution are specified on the boundary ∆ψ

¡−

→ r

¢ = 0 in volume Ω ³

R3´

∇ψ

¡−

→ r

¢ = −

→ f

¡−

→ r

¢ on surface ∂Ω ³

R3´

Where do we want to use this approach? In accelerator/spectrometer magnets where the magnet manufacturer pro- vides only discrete field data in the volume of interest

MAGNEX: A large acceptance MAGNetic spectrometer for EXcyt beams, at the Labora- tori Nazionali del Sud - Catania (Italy). (Fringe fields, high aspect ratio, discrete data)

What do we expect from this method/tool ?

• Provide validated local expansion of the field ( ψ

¡−

→ r

¢ and ∂n

xiψ

¡−

→ r

¢)

• Highly accurate (work for case with high aspect ratio)
• Computationally inexpensive
• Provide information about the field quality and if possible reduce noise

in experimentally obtained field data

Note about Laplace Equation

• Existence and uniqueness of the solution for 3D case can be shown

using Green’s formula

• Integral kernels that provides interior fields in terms of the boundary

fields or source are smoothing Interior fields will be analytic even if the field/source on the surface data fails to be differentiable

• Analytic closed form solution can be found for few problems with

certain regular geometries where separation of variables method can be applied

Numerical methods to solve Laplace equation

• Finite Difference, Method of weighted residuals and Finite element

methods — Numerical solution as data set in the region of interest — Relatively low approximation order — Prohibitively large number of mesh points and careful meshing re- quired

• Boundary integral methods or Source based field models

— Field inside of a source free volume due to a real sources outside

• f it can be exactly replicated by a distribution of fictitious sources
• n its surface. Error due to discretization of the source falls off

rapidly as the field point moves away from the source. ∗ Image charge method · Choose planes/grids to place point charges (or Gaussian dist) · Solve a large least square fit problem to find the charges · Lot of guess work and computation time involved in getting the solution ∗ Methods using Helmholtz’ theorem · Helmholtz’ theorem is used to find electric or magnetic field directly from the surface field data

· In our approach we make use of the Taylor model frame work to implement this

Helmholtz’ theorem Any vector field − → B that vanishes at infinity can be written as the sum of two terms, one of which is irrotational and the other, solenoidal − → B ( x) = ∇ × At ( x) + ∇φn ( x) φn ( x) = 1 4π

Z

∂Ω

• n (

xs) · − → B ( xs) | x − xs| ds − 1 4π

Z

Ω

• ∇ · −

→ B ( xv) | x − xv| dV

• At (

x) = − 1 4π

Z

∂Ω

• n (

xs) × − → B ( xs) | x − xs| ds + 1 4π

Z

Ω

• ∇ × −

→ B ( xv) | x − xv| dV

For a source free volume we have, ∇ × − → B ( xv) = 0 and ∇ · − → B ( xv) = 0 Volume integral terms vanish, φn ( x) and At ( x) are completely determined from the normal and the tangential magnetic field data on surface ∂Ω φn ( x) = 1

4π

R

∂Ω

• n(

xs)·− → B( xs) | x− xs|

ds

• At (

x) = − 1

4π

R

∂Ω

• n(

xs)×− → B( xs) | x− xs|

ds

• B is Electric or Magnetic field

∂Ω is a surface which bounds volume Ω

• xs and

xv denote points on ∂Ω and within Ω

• ∇ denote the gradient with respect to

xv

• n is a unit normal vector pointing away from ∂Ω

Implementation using Taylor Models

• Split domain of integration ∂Ω in to smaller regions Γi
• Expand them to higher orders in surface variables

rs and the volume variables r — Expanded in rs about the center of each surface element — Expanded in r about the center of each volume element — Field is chosen to be constant over each surface element

Z xNx

x0

Z yNy

y0

g (x, y) dxdy =

ix=Nx−1,iy=Ny−1,kx=∞,ky=∞

X

ix=0,iy=0,kx=0,ky=0

h2kx+1

x

(2kx + 1)! · 22kx h2ky+1

y

(2ky + 1)! · 22ky g2kx,2ky

Ã

xix+1 + xix 2 , yiy+1 + yiy 2

!

We can obtain: Scalar potential φn ( r) if we choose g(x, y) = ns ·

• f(

rs) | r− rs|

Vector potential At ( r) if we choose g(x, y) = ns ×

• f(

rs) | r− rs|

• — Benefits

∗ The dependence on the surface variables are integrated over sur- face sub-cells Γi, which results in a highly accurate integration formula ∗ The dependence on the volume variables are retained, which leads to a high order finite element method ∗ By using sufficiently high order, high accuracy can be achieved with a small number of surface elements

• Depending on the accuracy of the computation needed, we choose

step sizes, order of expansion in r (x, y, z) and order of expansion in rs (xs, ys, zs)

Validated Integration in COSY

xiu

Z

xil

f ( x) dxi ∈

³

Pn,∂−1f

³

• x|xi=xiu−xi0

´

− Pn,∂−1f

³

• x|xi=xil−xi0

´

, In,∂−1f

´

This method has following advantages: * No need to derive quadrature formulas with weights, support points xi, and an explicit error formula * High order can be employed directly by just increasing the order of the Taylor model limited only by the computational resources * Rather large dimensions are amenable by just increasing the dimensionality

• f the Taylor models, limited only by computational resources

An Analytical Example: the Bar Magnet x1 ≤ x ≤ x2, |y| ≥ y0, z1 ≤ z ≤ z2 As a reference problem we consider the magnetic field of rectangular iron bars with inner surfaces (y = ±y0) parallel to the mid-plane (y = 0)

From this bar magnet one can obtain analytic solution for the magnetic field B (x, y, z) of the form By (x, y, z) = B0 4π

X

i,j

(−1)i+j

⎡ ⎣arctan ⎛ ⎝ Xi · Zj

Y+ · R+

ij

⎞ ⎠ + arctan ⎛ ⎝ Xi · Zj

Y− · R−

ij

⎞ ⎠ ⎤ ⎦

Bx (x, y, z) = B0 4π

X

i,j

(−1)i+j

⎡ ⎣ln ⎛ ⎝Zj + R−

ij

Zj + R+

ij

⎞ ⎠ ⎤ ⎦

Bz (x, y, z) = B0 4π

X

i,j

(−1)i+j

⎡ ⎣ln ⎛ ⎝Xj + R−

ij

Xj + R+

ij

⎞ ⎠ ⎤ ⎦

where i, j = 1, 2 , Xi = x − xi, Y± = y0 ± y, Zi = z − zi and R± =

³

X2

i + Y 2 j + Z2 ±

´1

2

We note that only even order terms exist in the Taylor expansion of this field about the origin.

• 0.003
• 0.002
• 0.001

0.001 0.002 0.003 X

• 0.004
• 0.002

0.002 0.004 Y 2 4 6 8 10 12 BY

By component of the field on the mid-plane.

Results

• 1. To study the dependency of the Interval part of the potentials and

B field on the surface element length

• All of the volume is considered as just one volume element
• Examine contributions of each surface element towards the total

integral — Expansion is done at r = (.1, .1, .1) and — Plot of interval width VS surface element length for scalar po- tential — Plot of interval width VS surface element length for vector po- tential (x component)

• Plot of interval width VS Order for different surface element length

for x component of Magnetic field

• 22
• 20
• 18
• 16
• 14
• 12
• 10
• 8
• 6
• 4
• 2
• 13
• 12
• 11
• 10
• 9
• 8
• 7
• 6
• 5
• 4
• 3

LOG10(Interval Width) LOG2(Surface Element Length) Interval Width VS Surface Element Length for Scalar Potential Order 8 Order 7 Order 6 Order 5 Order 4 Order 3 Order 2

Figure 1: Integration over single surface element (for φ)

• 22
• 20
• 18
• 16
• 14
• 12
• 10
• 8
• 6
• 4
• 2
• 13
• 12
• 11
• 10
• 9
• 8
• 7
• 6
• 5
• 4
• 3

LOG10(Interval Width) LOG2(Surface Element Length) Interval Width VS Surface Element Length for X Vector Potential Order 8 Order 7 Order 6 Order 5 Order 4 Order 3 Order 2

Figure 2: Integration over single surface element (for Ax)

• 11
• 10
• 9
• 8
• 7
• 6
• 5
• 4

4 4.5 5 5.5 6 6.5 7 7.5 8 LOG(Interval Width) ORDER Interval Width VS Order for different stepsize Step size = 0.0166 Step size = 0.0192 Step size = 0.0147

Figure 3: Interval width VS Order (for different step size)

2 Study the dependency of the Polynomial part and Interval part of the B field on the volume element length

• The surface element length is locked at 1/128
• Plot of the error calculated for the polynomial part VS the volume

element length

• Plot of interval width VS volume element length for y component
• f Magnetic field

• 11
• 10
• 9
• 8
• 7
• 6
• 5
• 4
• 3
• 2
• 1
• 7.5
• 7
• 6.5
• 6
• 5.5
• 5
• 4.5
• 4
• 3.5
• 3

LOG(Interval Width) LOG2(Length of Volume Element) Interval Width VS Length of Volume Element for Scalar Potential Order 8 Order 7 Order 6

Figure 4: Interval width VS Volume element length (for φ)

• 11
• 10
• 9
• 8
• 7
• 6
• 5
• 4
• 3
• 2
• 1
• 7.5
• 7
• 6.5
• 6
• 5.5
• 5
• 4.5
• 4
• 3.5
• 3

LOG(Interval Width) LOG2(Length of Volume Element) Interval Width VS Length of Volume Element for X Component of Vector Potential Order 8 Order 7 Order 6

Figure 5: Interval width VS Volume element length (for Ax)

• 12
• 11
• 10
• 9
• 8
• 7
• 6
• 5
• 4
• 3

0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 LOG(Average error) Volume element length Error VS Volume element length Order 3 Order 5 Order 7 Order 9

Figure 6: Error VS Volume element length for polynomial part (for By)

• 10
• 9.5
• 9
• 8.5
• 8
• 7.5
• 7
• 6.5
• 6
• 5.5
• 5
• 1.95
• 1.9
• 1.85
• 1.8
• 1.75
• 1.7
• 1.65
• 1.6

LOG(Interval Width) LOG(Surface Element Length) Interval Width VS Surface Element Length for Y Component of Magnetic Field Order 8 Order 7 Order 6

Figure 7: Interval width VS Volume element length for By

Summary

• Helmholtz’ theorem implemented using the Taylor Model tools provide

a promising approach to find local expansion of the field in the volume

• f interest
• Accuracy achieved is very high compared to conventional numerical

field solvers

• Provides a good way to check the field quality
• This method can be extended for PDE’s