Transcript Use of Rectangular and Triangular Elements for Nearly

Application of neBEM to solve MPGD electrostatics
Supratik Mukhopadhyay, Nayana Majumdar, Sudeb Bhattacharya
Saha Institute of Nuclear Physics, Kolkata, India
E-mail: [email protected]
Presented by
Rob Veenhof
CERN
Field Solver
BEM
Solve Poisson’s
equation
.(mP)  S
 Reduced
dimension
 Accurate for both
potential and its
gradient
x Complex numerics
x Numerical boundary
layer
x Numerical and
physical singularities
FEM / FDM
 Nearly arbitrary
geometry
Analytic
 Exact
 Simple interpretation
x Restricted
x 2D geometry
x Small set of
geometries
 Flexible
x Interpolation for nonnodal points
x Numerical differentiation
for field gradient
x Difficulty in unbounded
domains
Finite Element Blues
•
•
•
•
•
•
•
•
The complete 3D volume is discretized using
nodes distributed throughout the volume
The governing equation is satisfied at the
nodal points to solve for potential
The variation of potential from node to node is
determined by a basis function – usually a low
order polynomial
Fields are thus represented by even lower
order polynomial
Value of potential at an arbitrary location is
obtained by interpolating values from
surrounding nodes – inaccurate
At a non-nodal location, field values suffer
even more
Serious problems in near-field
Artificial truncation of far-field boundary is a
necessity for problems with open domain
Mesh used for micromegas FEM solution by P Cwetansky
BEM Advantages
•
•
Drift plane
•
•
Micromegas
•
•
Anode strip
Mesh used for micromegas neBEM solution
Solves for charge density distribution on each
of the boundary elements by adopting a Green
function approach
Once the charge density distribution is
estimated correctly, potential and flux can be
easily obtained at any arbitrary location –
essentially mesh-less
No precision-eroding interpolation and
extrapolation are ever necessary
Near-field can be seamlessly handled if charge
density distribution is appropriate
Far-field conditions are satisfied naturally for
open-domain problems
Discretization is easier due to reduction of
dimensionality of the problem. It is also easier
to discretize complex devices.
BEM Basics
Green’s identities
Boundary Integral Equations
Potential u at any point y in the domain V enclosed by a surface S is given by
u ( y )   U ( x, y )q ( x)dS ( x)   Q ( x, y )u ( x)dS ( x)   U ( x, y )b( x)dV ( x)
S
S
V
where y is in V, u is the potential function, q = u,n, the normal derivative of u on the
boundary, b(x) is the body source, y is the load point and x, the field point. U and Q are
fundamental solutions
U2D = (1/2) ln(r), U3D = 1 / (4r), Q = -(1/2r) r,n
 = 1 for 2D and 2 for 3D. Distance from y to x is r, ni denotes the components of the
outward normal vector of the boundary.
2D Case
3D Case
r=0
r  0, r≠ 0
ln(r)
1/r
Weak singularity
Nearly weak singularity
1/r
1/r2
Strong singularity
Nearly strong singularity
1/r2
1/r3
Hyper singularity
Nearly hyper-singularity
Solution of 3D Poisson's Equation
using BEM
•
•
Numerical implementation of boundary integral equations (BIE) based on Green’s
function by discretization of boundary.
Boundary elements endowed with distribution of sources, doublets, dipoles, vortices.
Electrostatics BIE
Potential at r




(
r
)

(
r
,r
)

(
r
)
d
S
G
Green’s function
1

G
(
r
,r
)

4

r

r
S
Charge density at r’
Influence
Coefficient
Matrix
discretization
A
{ρ} = [A]-1{Φ}
 - permittivity of medium
Accuracy depends critically on the
estimation of [A], in turn, the
integration of G, which involves
singularities when r → r'.
Most BEM solvers fail here.
Conventional BEM
 While computing the influences of the singularities, the singularities are
modelled by a sum of known basis functions with constant unknown
coefficients.
 The strengths of the singularities are solved depending upon the
boundary conditions, modeled by shape functions.
Singularities are assumed to be concentrated at centroids of elements
(thus avoiding integration), except for cases such as self influence.
r = r’: Mathematical singularities
can be removed; Sufficient to
satisfy the boundary conditions at
centroids of the elements.
r -> r’: Difficulties in modeling physical
singularities
geometric singularity:
Closely spaced surfaces,
corners, edges
boundary condition
singularity: Dirichlet
and Neumann conditions
close-by
Numerical boundary layer
Present Approach
Analytic expressions of the integration yielding both potential and flux field at any arbitrary
location due to a uniform distribution of source on flat rectangular and triangular elements have
been derived using symbolic tools. Using these elements, surfaces of any 3D geometry can be
discretized without requiring the singularities to be concentrated only at certain specific nodes.
Restatement of the approximations
 Singularities distributed uniformly on the surface of boundary elements.
Nodal concentrations, rather than the integrations, are avoided altogether.
 Strength of the singularity changes from element to element (unchanged).
 Strengths of the singularities solved depending upon the boundary
conditions, modeled by the shape functions (unchanged)
ISLES library and neBEM Solver
Foundation expressions are analytic and valid for the complete physical domain
Contrast of approaches
nodal (conventional) versus distributed (neBEM)
Unrealistic representation:
•Near-field solutions grossly incorrect
•Aspect ratio of element sides cannot be high
•Size of elements cannot be varied sharply
•Edges, corners, closely packed surfaces cannot be
modeled easily
•Proximity of Dirichlet and Neumann condition
cannot be allowed
•Flourish of special formulations!
Realistic representation:
•Accurate solutions everywhere, including near-field
•Aspect ratio of element sides can be much larger
•Size of elements can be varied easily
•Edges, corner, closely packed surfaces can be modeled
easily
•Proximity of Dirichlet and Neumann conditions allowed
•Single formulation for many problems
Precision in flux computation
comparison with quadrature (nodal approach of ususal BEM)
zMax = 10.0
Quadrature with only the highest discretization
produces results comparable to ISLES
Quadrature with even the highest discretization
fails!
Engineering Analysis with Boundary Elements (EABE), Elsevier, Available online 3 August 2008
Precision in flux computation
Comparison with multipole expansions
Comparison of flux along a line parallel to the
Z axis passing through the barycenter
Comparison of flux along a diagonal passing
through the barycenter
The quadrupole results are still far from precise
EABE, online, 3 Aug 2008
Expected features of a Field Solver for MPGDs
GEM
Typical dimensions:
Electrodes (5 μm thick)
Insulator (50 μm thick)
Hole size D ~ 60 μm
Pitch p ~ 140 μm
Induction gap: 1.0 mm,
Transfer gap: 1.5 mm
Micromegas
Typical dimensions:
Mesh size: 50 μm
Micromesh sustained
by 50 μm pillars
• Variation of field over length scales of a
micron to a meter needs to be precisely
estimated
• Fields at arbitrary locations should be
available on demand
• Intricate geometrical features – essential to
use triangular elements, if needed
• Multiple dielectric devices
• Nearly degenerate surfaces
• Space charge effects can be very significant
• Dynamic charging processes may be
important
• It may be necessary to calculate field for the
same geometry, but with different electric
configuration, repeatedly
Micro-wire: a test case
•
•
•
•
•
Source: Adeva et al, USC-FP/99-01
Representative length variation
Closely spaced elements / surfaces
Multiple dielectric device
Complex geometry
Field variation very close to anode is
important
• FEM solutions available (Peter
Cwetansky:
http://consult.cern.ch/writeup/
garfield/examples/micropattern/
microwire/index.htm)
• Please note that for computation, a
similar, but not identical geometry
has been used
Discretization
In FEM, the complete
3D domain needs to be
discretized
Far-field in FEM is
truncated and, possibly,
supplied with
Neumann boundary
conditions to maintain
periodicity and assign
drift field strength
In neBEM, the farfield is realistically
represented by a drift
plane at a certain
voltage.
Close up
Elements of various sizes and aspect ratios have been used
Since we are interested mostly in the field variation around the anode, it has been
discretized using very small elements
The lines join the element centroids – it is not the true mesh and, hence, has apparent gaps
Influence coefficient matrix
•
•
•
•
 I11
I
 21
 I 31

 I 41
 0
I12
I13
I14
I 22
I 32
I 23
I 33
I 24
I 34
I 42
I 43
I 44
0
A3
A4
0  1  V 
   

0   2  V 
 1  3    0 
   
 1  4   0 
0 VF   0 
•
•
•
•
Consider a system of two conductors, each having
been discretized into two elements.
In order to generalize the situation, let us even
consider one of the conductors to be floating. Thus,
one of these conductors is at a known voltage, V.
The other conductor is at a floating voltage VF,
which is unknown.
Number the elements on the conductor with known
voltage to be 1, 2, and those on the floating
conductor to be 3, 4.
Denote charge densities by i, area by Ai, on each
element
Resulting system of equation is as shown – the last
equation reflecting the fact that the total charge on
a floating conducting object is zero
In the above system, Iij denotes the influence of the
jth element on the ith element.
Please note that if we have more than one floating
conductor, they cannot be assumed to be at the
same potential, and one column and one row as
shown above needs to be added for each floating
conductor.
Please note that the matrix is not sparse.
Evaluation of influence coefficients
( X , Y , Z ) 

2  ( X | Z | x | z )  ln Di , j  ( X | Z  xi | z j ) 
i
j
 D  (X | Z  x | z ) 

m
n 
 m,n

1 

2 

 R j  iI i 


  tanh1  R j  iI i
 tanh1 

iS
Y

j

 Di , j Z  z j 
 Di , j Z  z j






 Di , j  ( Z  z j ) 

FX ( X , Y , Z )  ln

D

(
Z

z
)
n 
 m,n
FY ( X , Y , Z ) 

 R j  iI i 

S j t anh1 
 Di , j Z  z j 

i



 Sign(Y )  
2
 R j  iI i

1 
 S j t anh 

 Di , j Z  z j
 Di , j  ( X  xi ) 

FZ ( X , Y , Z )  ln

D

(
X

x
)
m 
 m,n




  2Y

 

 
 
2 terms



C




4 log
terms
4+4 complex
tanh-1 terms
4+4
terms
2 terms
Similar, but more complex, expressions have been derived for trianglular elements
Charge density distribution
Anode
Kapton
Cathode
• The neBEM solves for the charge
density
• On the cathode, the density is
negative and with relatively less
variation
• On the kaptons, the density is of both
polarities and the variation is least
• On the anode, the density is positive
and the variation is very sharp, the
edges and corners having very large
charge densities
• It is at these edges and corners where
the FEM fails, in general
Intricate geometries
The Micro Wire Detector
Total electric field contours on the central plane
across cathode and anode
Variation of total electric field along an axis
passing through the mesh hole
 The MWD has an intricate design. In this case:
Drift plane 785μm from the anode strip at 1.11kV.
JINST, 2007, 2 P09006
Case studies
These are the possible design variation discussed in the earlier study, the conventional one being
termed as mesh, other one as segmented. The latter consumes 40% less copper!
Potential
•There is significant variation in the potential contours of these two possible designs (points
representing `mesh’ configuration, broken lines for `segmented’).
•The right figure shows the potential for the segmented design only.
Flux
•There is significant variation in the flux contours as well (points representing `mesh’
configuration, broken lines for `segmented’, once again).
•The right figure shows the potential for the segmented design only.
Comparison with FEM
•
•
•
•
•
The mesh and segmented configurations
have significant differences except at the
mid-zone of the detector unit (please note
that the Electric field is on a log scale)
Close to the anode, the mesh
configuration achieves a higher field
The drift field is more for the segmented
configuration
Since the far-field for the FEM
computations is treated quite differently
from the neBEM, values at these zones
differ considerably (once again, note that
the scale is logarithmic and the variation
is, in fact, less than what is apparent
Besides the far-field, the two results agree
quite closely
Comparison with FEM (near-field)
•
•
•
•
•
•
•
Field around a line just 1μm away from the
anode surface is considered here –sampling
for neBEM is as small as 0.1μm!
The mesh configuration has higher field
values throughout
Sharp rise in the field values is observed at
all the four edges
Smooth variation of field is observed on each
of the four surfaces
Field values are found to decrease sharply
once the points are beyond anode surfaces
FEM computation is clearly unable to
produce correct results near and at the edges
FEM, although better on the surfaces, still
falls behind neBEM in performance
Effect of discretization (near-field)
•
•
•
•
In the earlier computation, we had used
20 elements to represent the top surface
and 10 elements on the side surface. The
elements were made successively smaller
towards the edges
In order to study effect of using coarse
discretization, we also used larger
elements of fixed size – only 3 elements
each to represent both top and side
surfaces
Although there is significant difference
between the results, the overall trend is
represented well by the larger elements
It is important to note that there is no
jaggedness (at 0.1μm sampling) despite
the use of unreasonably large elements!
Multiple dielectric devices
Resistive plate chambers
Strip width: 3.0cm, Strip length: 50.0cm
Layer height: 2.0mm
Layer-3 permittivity (r): 7.75 (~glass)
Layer-2 (middle) permittivity (r): 1.000513
(~Argon)
Successful validation with Riegler et al.
Layer-4,5 height: 200µm (~PET)
Layer-6,7 height: 20µm (~Graphite)
Layer-4,5 permittivity (r): 3.0 (~PET)
Layer-6,7 permittivity (r): 12.0
(~Graphite)
NIM, A 595 (2008) 346-352
Conclusions
• We have outlined the approach neBEM takes in order to compute the field
properties in a given device
• The approach has been compared with the more standard FEM method in
fair detail
• The micro-wire detector has been analyzed since it has features typical of
many of the MPGDs, and has good FEM results available for it
• FEM and neBEM results have been compared for far- and near-field
• neBEM results have been found to be more accurate in critical regions
• Large variation of discretization has not deteriorated neBEM results to any
great extent
• RPC weighting field results have been presented in order to demonstrate
the accuracy of neBEM – comparisons against analytical results have
turned out to be encouraging
Plans for 2009
• Development of an interface to ROOT so that devices built using ROOT
can be directly imported to neBEM and solved for
• For simple shapes, thanks to Andrei, we have already been able to extract
the surfaces that can be exported to neBEM using a ROOT script or a
stand-alone C++ code
• For composite shapes, thanks to Timur, we have been able to get the
elements being used for the geometrical rendering in ROOT. We are trying
to understand this mechanism to be finally able to export these surfaces in
the neBEM format
• A working interface to garfield and the new detailed detector simulation
framework is expected to be complete by the middle of 2009
• We will also try to set up an experimental laboratory for the development
of MPGDs during this period
Plans for beyond 2009
o The problem of dynamic charging will be addressed.
o Particles on Surface (ParSue), a new model for space-charge simulation based on
this formulation has been proposed recently. This model needs to be explored and
integrated properly.
o Problems related to magnetostatics will be addressed.
o Optimization and introduction of adaptivity in the process of mesh generation will
be implemented
o Implementation of improved and more efficient matrix solution algorithms will be
carried out.
o Parallel computation can help the overall detailed simulation and efforts may be
made in this direction
o A toolkit version of the field-solver may be developed for use in other areas
governed by the Poisson’s equation.
Acknowledgements
• We thank the organizers for giving us a chance to present this
work
• We thank our Director, Prof. Bikash Sinha, for his support
• We thank Rob for encouraging us, suggesting improvements
and test cases for neBEM and, finally, for presenting this
material
• We thank you all for your kind attention
Looking forward to a very friendly and successful collaboration!