/[escript]/trunk/doc/user/finley.tex
ViewVC logotype

Annotation of /trunk/doc/user/finley.tex

Parent Directory Parent Directory | Revision Log Revision Log


Revision 3330 - (hide annotations)
Mon Nov 1 04:29:09 2010 UTC (8 years, 2 months ago) by caltinay
File MIME type: application/x-tex
File size: 27352 byte(s)
Finished user's guide chapters up to Appendix.

1 ksteube 1811
2     %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
3 ksteube 1316 %
4 jfenwick 2881 % Copyright (c) 2003-2010 by University of Queensland
5 ksteube 1811 % Earth Systems Science Computational Center (ESSCC)
6     % http://www.uq.edu.au/esscc
7 gross 625 %
8 ksteube 1811 % Primary Business: Queensland, Australia
9     % Licensed under the Open Software License version 3.0
10     % http://www.opensource.org/licenses/osl-3.0.php
11 gross 625 %
12 ksteube 1811 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
13 jgs 102
14 caltinay 3293 \chapter{The \finley Module}\label{CHAPTER ON FINLEY}
15 caltinay 3306 %\declaremodule{extension}{finley}
16     %\modulesynopsis{Solving linear, steady partial differential equations using finite elements}
17 jgs 102
18 caltinay 3330 {\it finley} is a library of C functions solving linear, steady partial
19     differential equations\index{partial differential equations} (PDEs) or systems
20     of PDEs using isoparametrical finite elements\index{FEM!isoparametrical}.
21     It supports unstructured 1D, 2D and 3D meshes.
22     The module \finley provides access to the library through the \LinearPDE class
23     of \escript supporting its full functionality.
24     {\it finley} is parallelized using the OpenMP\index{OpenMP} paradigm.
25 jgs 102
26 gross 993 \section{Formulation}
27 caltinay 3330 For a single PDE that has a solution with a single component the linear PDE is
28     defined in the following form:
29 gross 993 \begin{equation}\label{FINLEY.SINGLE.1}
30     \begin{array}{cl} &
31     \displaystyle{
32 jfenwick 3295 \int_{\Omega}
33     A_{jl} \cdot v_{,j}u_{,l}+ B_{j} \cdot v_{,j} u+ C_{l} \cdot v u_{,l}+D \cdot vu \; d\Omega } \\
34     + & \displaystyle{\int_{\Gamma} d \cdot vu \; d{\Gamma} }
35     + \displaystyle{\int_{\Gamma^{contact}} d^{contact} \cdot [v][u] \; d{\Gamma} } \\
36     = & \displaystyle{\int_{\Omega} X_{j} \cdot v_{,j}+ Y \cdot v \; d\Omega }\\
37     + & \displaystyle{\int_{\Gamma} y \cdot v \; d{\Gamma}} +
38     \displaystyle{\int_{\Gamma^{contact}} y^{contact}\cdot [v] \; d{\Gamma}} \\
39 gross 993 \end{array}
40     \end{equation}
41    
42     \section{Meshes}
43 gross 2793 \label{FINLEY MESHES}
44 jgs 102
45 caltinay 3330 \begin{figure}
46     \centerline{\includegraphics{FinleyMesh}}
47     \caption{Subdivision of an Ellipse into triangles order 1 (\finleyelement{Tri3})}
48     \label{FINLEY FIG 0}
49     \end{figure}
50    
51     To understand the usage of \finley one needs to have an understanding of how
52     the finite element meshes\index{FEM!mesh} are defined.
53     \fig{FINLEY FIG 0} shows an example of the subdivision of an ellipse into
54     so-called elements\index{FEM!elements}\index{element}.
55     In this case, triangles have been used but other forms of subdivisions can be
56     constructed, e.g. quadrilaterals or, in the three-dimensional case, into
57     tetrahedra and hexahedra. The idea of the finite element method is to
58     approximate the solution by a function which is a polynomial of a certain order
59     and is continuous across its boundary to neighbour elements.
60     In the example of \fig{FINLEY FIG 0} a linear polynomial is used on each
61     triangle. As one can see, the triangulation is quite a poor approximation of
62     the ellipse. It can be improved by introducing a midpoint on each element edge
63     then positioning those nodes located on an edge expected to describe the
64     boundary, onto the boundary.
65     In this case the triangle gets a curved edge which requires a parameterization
66     of the triangle using a quadratic polynomial.
67     For this case, the solution is also approximated by a piecewise quadratic
68     polynomial (which explains the name isoparametrical elements),
69     see \Ref{Zienc,NumHand} for more details.
70     \finley also supports macro elements\index{macro elements}.
71     For these elements a piecewise linear approximation is used on an element which
72     is further subdivided (in the case of \finley halved).
73     As such, these elements do not provide more than a further mesh refinement but
74     should be used in the case of incompressible flows, see \class{StokesProblemCartesian}.
75     For these problems a linear approximation of the pressure across the element is
76     used (use the \ReducedSolutionFS) while the refined element is used to
77     approximate velocity. So a macro element provides a continuous pressure
78     approximation together with a velocity approximation on a refined mesh.
79     This approach is necessary to make sure that the incompressible flow has a
80     unique solution.
81    
82 jgs 102 The union of all elements defines the domain of the PDE.
83 caltinay 3330 Each element is defined by the nodes used to describe its shape.
84     In \fig{FINLEY FIG 0} the element, which has type \finleyelement{Tri3}, with
85     element reference number $19$\index{element!reference number} is defined by the
86     nodes with reference numbers $9$, $11$ and $0$\index{node!reference number}.
87     Notice that the order is counterclockwise.
88     The coefficients of the PDE are evaluated at integration nodes with each
89     individual element.
90     For quadrilateral elements a Gauss quadrature scheme is used.
91     In the case of triangular elements a modified form is applied.
92     The boundary of the domain is also subdivided into elements\index{element!face}.
93     In \fig{FINLEY FIG 0} line elements with two nodes are used.
94     The elements are also defined by their describing nodes, e.g. the face element
95     with reference number $20$, which has type \finleyelement{Line2}, is defined by
96     the nodes with the reference numbers $11$ and $0$.
97     Again the order is crucial, if moving from the first to second node the domain
98     has to lie on the left hand side (in the case of a two-dimensional surface
99     element the domain has to lie on the left hand side when moving
100     counterclockwise). If the gradient on the surface of the domain is to be
101     calculated rich face elements need to be used. Rich elements on a face are
102     identical to interior elements but with a modified order of nodes such that the
103     'first' face of the element aligns with the surface of the domain.
104     In \fig{FINLEY FIG 0} elements of the type \finleyelement{Tri3Face} are used.
105     The face element reference number $20$ as a rich face element is defined by the
106     nodes with reference numbers $11$, $0$ and $9$.
107     Notice that the face element $20$ is identical to the interior element $19$
108     except that, in this case, the order of the node is different to align the first
109     edge of the triangle (which is the edge starting with the first node) with the
110     boundary of the domain.
111 jgs 102
112 caltinay 3330 Be aware that face elements and elements in the interior of the domain must
113     match, i.e. a face element must be the face of an interior element or, in case
114     of a rich face element, it must be identical to an interior element.
115     If no face elements are specified \finley implicitly assumes homogeneous
116     natural boundary conditions\index{natural boundary conditions!homogeneous},
117     i.e. \var{d}=$0$ and \var{y}=$0$, on the entire boundary of the domain.
118     For inhomogeneous natural boundary conditions\index{natural boundary conditions!inhomogeneous},
119 caltinay 3293 the boundary must be described by face elements.
120 jgs 102
121 caltinay 3330 \begin{figure}
122     \centerline{\includegraphics{FinleyContact}}
123     \caption{Mesh around a contact region (\finleyelement{Rec4})}
124     \label{FINLEY FIG 01}
125     \end{figure}
126 jgs 102
127 caltinay 3330 If discontinuities of the PDE solution are considered, contact
128     elements\index{element!contact}\index{contact conditions} are introduced to
129     describe the contact region $\Gamma^{contact}$ even if $d^{contact}$ and
130     $y^{contact}$ are zero.
131     \fig{FINLEY FIG 01} shows a simple example of a mesh of rectangular elements
132     around a contact region $\Gamma^{contact}$\index{element!contact}.
133     The contact region is described by the elements $4$, $3$ and $6$.
134     Their element type is \finleyelement{Line2_Contact}.
135     The nodes $9$, $12$, $6$ and $5$ define contact element $4$, where the
136     coordinates of nodes $12$ and $5$ and nodes $4$ and $6$ are identical, with the
137     idea that nodes $12$ and $9$ are located above and nodes $5$ and $6$ below the
138     contact region.
139     Again, the order of the nodes within an element is crucial.
140     There is also the option of using rich elements if the gradient is to be
141     calculated on the contact region. Similarly to the rich face elements these
142     are constructed from two interior elements by reordering the nodes such that
143     the 'first' face of the element above and the 'first' face of the element below
144     the contact regions line up. The rich version of element $4$ is of type
145     \finleyelement{Rec4Face_Contact} and is defined by the nodes $9$, $12$, $16$,
146     $18$, $6$, $5$, $0$ and $2$.
147     \tab{FINLEY TAB 1} shows the interior element types and the corresponding
148     element types to be used on the face and contacts.
149     \fig{FINLEY.FIG:1}, \fig{FINLEY.FIG:2} and \fig{FINLEY.FIG:4} show the ordering
150     of the nodes within an element.
151 gross 2748
152 jgs 102 \begin{table}
153 caltinay 3293 \centering
154     \begin{tabular}{l|llll}
155 caltinay 3330 \textbf{interior}&\textbf{face}&\textbf{rich face}&\textbf{contact}&\textbf{rich contact}\\
156 caltinay 3293 \hline
157     \finleyelement{Line2} & \finleyelement{Point1} & \finleyelement{Line2Face} & \finleyelement{Point1_Contact} & \finleyelement{Line2Face_Contact}\\
158     \finleyelement{Line3} & \finleyelement{Point1} & \finleyelement{Line3Face} & \finleyelement{Point1_Contact} & \finleyelement{Line3Face_Contact}\\
159     \finleyelement{Tri3} & \finleyelement{Line2} & \finleyelement{Tri3Face} & \finleyelement{Line2_Contact} & \finleyelement{Tri3Face_Contact}\\
160     \finleyelement{Tri6} & \finleyelement{Line3} & \finleyelement{Tri6Face} & \finleyelement{Line3_Contact} & \finleyelement{Tri6Face_Contact}\\
161     \finleyelement{Rec4} & \finleyelement{Line2} & \finleyelement{Rec4Face} & \finleyelement{Line2_Contact} & \finleyelement{Rec4Face_Contact}\\
162     \finleyelement{Rec8} & \finleyelement{Line3} & \finleyelement{Rec8Face} & \finleyelement{Line3_Contact} & \finleyelement{Rec8Face_Contact}\\
163     \finleyelement{Rec9} & \finleyelement{Line3} & \finleyelement{Rec9Face} & \finleyelement{Line3_Contact} & \finleyelement{Rec9Face_Contact}\\
164     \finleyelement{Tet4} & \finleyelement{Tri6} & \finleyelement{Tet4Face} & \finleyelement{Tri6_Contact} & \finleyelement{Tet4Face_Contact}\\
165     \finleyelement{Tet10} & \finleyelement{Tri9} & \finleyelement{Tet10Face} & \finleyelement{Tri9_Contact} & \finleyelement{Tet10Face_Contact}\\
166     \finleyelement{Hex8} & \finleyelement{Rec4} & \finleyelement{Hex8Face} & \finleyelement{Rec4_Contact} & \finleyelement{Hex8Face_Contact}\\
167     \finleyelement{Hex20} & \finleyelement{Rec8} & \finleyelement{Hex20Face} & \finleyelement{Rec8_Contact} & \finleyelement{Hex20Face_Contact}\\
168 caltinay 3330 \finleyelement{Hex27} & \finleyelement{Rec9} & N/A & N/A & N/A\\
169     \finleyelement{Hex27Macro} & \finleyelement{Rec9Macro} & N/A & N/A & N/A\\
170     \finleyelement{Tet10Macro} & \finleyelement{Tri6Macro} & N/A & N/A & N/A\\
171     \finleyelement{Rec9Macro} & \finleyelement{Line3Macro} & N/A & N/A & N/A\\
172     \finleyelement{Tri6Macro} & \finleyelement{Line3Macro} & N/A & N/A & N/A\\
173 caltinay 3293 \end{tabular}
174 caltinay 3330 \caption{Finley elements and corresponding elements to be used on domain faces
175     and contacts.
176     The rich types have to be used if the gradient of the function is to be
177     calculated on faces and contacts, respectively.}
178 jgs 102 \label{FINLEY TAB 1}
179     \end{table}
180    
181     The native \finley file format is defined as follows.
182 caltinay 3330 Each node \var{i} has \var{dim} spatial coordinates \var{Node[i]}, a reference
183     number \var{Node_ref[i]}, a degree of freedom \var{Node_DOF[i]} and a tag
184     \var{Node_tag[i]}.
185     In most cases \var{Node_DOF[i]}=\var{Node_ref[i]} however, for periodic
186     boundary conditions, \var{Node_DOF[i]} is chosen differently, see example below.
187     The tag can be used to mark nodes sharing the same properties.
188     Element \var{i} is defined by the \var{Element_numNodes} nodes
189     \var{Element_Nodes[i]} which is a list of node reference numbers.
190     The order of these is crucial. Each element has a reference number
191     \var{Element_ref[i]} and a tag \var{Element_tag[i]}.
192     The tag can be used to mark elements sharing the same properties.
193     For instance elements above a contact region are marked with tag $2$ and
194     elements below a contact region are marked with tag $1$.
195     \var{Element_Type} and \var{Element_Num} give the element type and the number
196     of elements in the mesh.
197     Analogue notations are used for face and contact elements.
198     The following \PYTHON script prints the mesh definition in the \finley file
199     format:
200 jgs 102 \begin{python}
201 caltinay 3330 print("%s\n"%mesh_name)
202     # node coordinates:
203     print("%dD-nodes %d\n"%(dim, numNodes))
204     for i in range(numNodes):
205     print("%d %d %d"%(Node_ref[i], Node_DOF[i], Node_tag[i]))
206     for j in range(dim): print(" %e"%Node[i][j])
207     print("\n")
208     # interior elements
209     print("%s %d\n"%(Element_Type, Element_Num))
210     for i in range(Element_Num):
211     print("%d %d"%(Element_ref[i], Element_tag[i]))
212     for j in range(Element_numNodes): print(" %d"%Element_Nodes[i][j])
213     print("\n")
214     # face elements
215     print("%s %d\n"%(FaceElement_Type, FaceElement_Num))
216     for i in range(FaceElement_Num):
217     print("%d %d"%(FaceElement_ref[i], FaceElement_tag[i]))
218     for j in range(FaceElement_numNodes): print(" %d"%FaceElement_Nodes[i][j])
219     print("\n")
220     # contact elements
221     print("%s %d\n"%(ContactElement_Type, ContactElement_Num))
222     for i in range(ContactElement_Num):
223     print("%d %d"%(ContactElement_ref[i], ContactElement_tag[i]))
224     for j in range(ContactElement_numNodes): print(" %d"%ContactElement_Nodes[i][j])
225     print("\n")
226     # point sources (not supported yet)
227     print("Point1 0")
228 jgs 102 \end{python}
229    
230     The following example of a mesh file defines the mesh shown in \fig{FINLEY FIG 01}:
231     \begin{verbatim}
232     Example 1
233     2D Nodes 16
234     0 0 0 0. 0.
235     2 2 0 0.33 0.
236     3 3 0 0.66 0.
237     7 4 0 1. 0.
238     5 5 0 0. 0.5
239     6 6 0 0.33 0.5
240     8 8 0 0.66 0.5
241     10 10 0 1.0 0.5
242     12 12 0 0. 0.5
243     9 9 0 0.33 0.5
244     13 13 0 0.66 0.5
245     15 15 0 1.0 0.5
246     16 16 0 0. 1.0
247     18 18 0 0.33 1.0
248     19 19 0 0.66 1.0
249     20 20 0 1.0 1.0
250     Rec4 6
251     0 1 0 2 6 5
252     1 1 2 3 8 6
253     2 1 3 7 10 8
254     5 2 12 9 18 16
255     7 2 13 19 18 9
256     10 2 20 19 13 15
257     Line2 0
258     Line2_Contact 3
259     4 0 9 12 6 5
260     3 0 13 9 8 6
261     6 0 15 13 10 8
262     Point1 0
263     \end{verbatim}
264     Notice that the order in which the nodes and elements are given is arbitrary.
265 caltinay 3330 In the case that rich contact elements are used the contact element section
266     gets the form
267 jgs 102 \begin{verbatim}
268     Rec4Face_Contact 3
269     4 0 9 12 16 18 6 5 0 2
270     3 0 13 9 18 19 8 6 2 3
271     6 0 15 13 19 20 10 8 3 7
272     \end{verbatim}
273 caltinay 3330 Periodic boundary conditions\index{boundary conditions!periodic} can be
274     introduced by altering \var{Node_DOF}.
275     It allows identification of nodes even if they have different physical locations.
276     For instance, to enforce periodic boundary conditions at the face $x_0=0$ and
277     $x_0=1$ one identifies the degrees of freedom for nodes $0$, $5$, $12$ and $16$
278     with the degrees of freedom for $7$, $10$, $15$ and $20$, respectively.
279     The node section of the \finley mesh now reads:
280 jgs 102 \begin{verbatim}
281     2D Nodes 16
282     0 0 0 0. 0.
283     2 2 0 0.33 0.
284     3 3 0 0.66 0.
285     7 0 0 1. 0.
286     5 5 0 0. 0.5
287     6 6 0 0.33 0.5
288     8 8 0 0.66 0.5
289     10 5 0 1.0 0.5
290     12 12 0 0. 0.5
291     9 9 0 0.33 0.5
292     13 13 0 0.66 0.5
293     15 12 0 1.0 0.5
294     16 16 0 0. 1.0
295     18 18 0 0.33 1.0
296     19 19 0 0.66 1.0
297     20 16 0 1.0 1.0
298     \end{verbatim}
299    
300 jfenwick 1955 \clearpage
301     \input{finleyelements}
302     \clearpage
303 jgs 102
304 caltinay 3330 \section{Macro Elements}
305     \label{SEC FINLEY MACRO}
306    
307 gross 2793 \begin{figure}[th]
308     \begin{center}
309 caltinay 3330 \includegraphics{FinleyMacroLeg}\\
310     \subfigure[Triangle]{\label{FINLEY MACRO TRI}\includegraphics{FinleyMacroTri}}\quad
311     \subfigure[Quadrilateral]{\label{FINLEY MACRO REC}\includegraphics{FinleyMacroRec}}
312 gross 2793 \end{center}
313 caltinay 3330 \caption{Macro elements in \finley}
314 gross 2793 \end{figure}
315    
316 caltinay 3330 \finley supports the usage of macro elements\index{macro elements} which can be
317     used to achieve LBB compliance when solving incompressible fluid flow problems.
318     LBB compliance is required to get a problem which has a unique solution for
319     pressure and velocity. For macro elements the pressure and velocity are
320     approximated by a polynomial of order 1 but the velocity approximation bases on
321     a refinement of the elements. The nodes of a triangle and quadrilateral element
322     are shown in Figures~\ref{FINLEY MACRO TRI} and~\ref{FINLEY MACRO REC},
323     respectively. In essence, the velocity uses the same nodes like a quadratic
324     polynomial approximation but replaces the quadratic polynomial by piecewise
325     linear polynomials. In fact, this is the way \finley defines the macro elements.
326     In particular \finley uses the same local ordering of the nodes for the macro
327     element as for the corresponding quadratic element. Another interpretation is
328     that one uses a linear approximation of the velocity together with a linear
329     approximation of the pressure but on elements created by combining elements to
330     macro elements. Notice that the macro elements still use quadratic
331     interpolation to represent the element and domain boundary.
332     However, if elements have linear boundaries a macro element approximation for
333     the velocity is equivalent to using a linear approximation on a mesh which is
334     created through a one-step global refinement.
335     Typically macro elements are only required to use when an incompressible fluid
336     flow problem is solved, e.g. the Stokes problem in \Sec{STOKES PROBLEM}.
337     Please see \Sec{FINLEY MESHES} for more details on the supported macro elements.
338 jgs 102
339 caltinay 3330 \section{Linear Solvers in \SolverOptions}
340 gross 2748
341 caltinay 3330 Table~\ref{TAB FINLEY SOLVER OPTIONS 1} and
342     Table~\ref{TAB FINLEY SOLVER OPTIONS 2} show the solvers and preconditioners
343     supported by \finley through the \PASO library.
344     Currently direct solvers are not supported under \MPI.
345     By default, \finley uses the iterative solvers \PCG for symmetric and \BiCGStab
346     for non-symmetric problems.
347     If the direct solver is selected, which can be useful when solving very
348     ill-posed equations, \finley uses the \MKL\footnote{If the stiffness matrix is
349     non-regular \MKL may return without a proper error code. If you observe
350     suspicious solutions when using \MKL, this may be caused by a non-invertible
351     operator.} solver package. If \MKL is not available \UMFPACK is used.
352     If \UMFPACK is not available a suitable iterative solver from \PASO is used.
353 gross 2748
354 gross 2558 \begin{table}
355 caltinay 3330 \centering
356 jfenwick 2651 {\scriptsize
357 gross 2558 \begin{tabular}{l||c|c|c|c|c|c|c|c}
358     \member{setSolverMethod} & \member{DIRECT}& \member{PCG} & \member{GMRES} & \member{TFQMR} & \member{MINRES} & \member{PRES20} & \member{BICGSTAB} & \member{LUMPING} \\
359     \hline
360     \hline
361     \member{setReordering} & $\checkmark$ & & & & & &\\
362     \hline \member{setRestart} & & & $\checkmark$ & & & $20$ & \\
363     \hline\member{setTruncation} & & & $\checkmark$ & & & $5$ & \\
364     \hline\member{setIterMax} & & $\checkmark$& $\checkmark$ & $\checkmark$& $\checkmark$& $\checkmark$ & $\checkmark$ \\
365     \hline\member{setTolerance} & & $\checkmark$& $\checkmark$ & $\checkmark$& $\checkmark$& $\checkmark$ & $\checkmark$ \\
366     \hline\member{setAbsoluteTolerance} & & $\checkmark$& $\checkmark$ & $\checkmark$& $\checkmark$& $\checkmark$ & $\checkmark$ \\
367 gross 2573 \hline\member{setReordering} & $\checkmark$ & & & & & & & \\
368 gross 2558 \end{tabular}
369     }
370 caltinay 3330 \caption{Solvers available for \finley and the \PASO package and the relevant
371     options in \class{SolverOptions}.
372     \MKL supports \member{MINIMUM_FILL_IN}\index{linear solver!minimum fill-in ordering}\index{minimum fill-in ordering}
373     and \member{NESTED_DISSECTION}\index{linear solver!nested dissection ordering}\index{nested dissection}
374 caltinay 3293 reordering.
375     Currently the \UMFPACK interface does not support any reordering.
376 caltinay 3330 \label{TAB FINLEY SOLVER OPTIONS 1}}
377 caltinay 3293 \end{table}
378 gross 2558
379     \begin{table}
380     {\scriptsize
381     \begin{tabular}{l||c|c|c|c|c|c|c|c}
382 caltinay 3330 \member{setPreconditioner}&
383     \member{NO_PRECONDITIONER}&
384     \member{AMG}&
385     \member{JACOBI}&
386 caltinay 3293 \member{GAUSS_SEIDEL}&
387 gross 2558 \member{REC_ILU}&
388 caltinay 3330 \member{RILU}&
389     \member{ILU0}&
390     \member{DIRECT}\\
391 gross 2558 \hline
392 caltinay 3330 status:&later&later&$\checkmark$&$\checkmark$&$\checkmark$&later&$\checkmark$&later\\
393 gross 2558 \hline
394 caltinay 3330 \hline
395     \member{setCoarsening}& &$\checkmark$& & & & & &\\
396     \hline
397     \member{setLevelMax}& &$\checkmark$& & & & & &\\
398     \hline
399     \member{setCoarseningThreshold}& &$\checkmark$& & & & & &\\
400     \hline
401     \member{setMinCoarseMatrixSize}& &$\checkmark$& & & & & &\\
402     \hline
403     \member{setNumSweeps}& & &$\checkmark$&$\checkmark$& & & &\\
404     \hline
405     \member{setNumPreSweeps}& &$\checkmark$& & & & & &\\
406     \hline
407     \member{setNumPostSweeps}& &$\checkmark$& & & & & &\\
408     \hline
409     \member{setInnerTolerance}& & & & & & & &\\
410     \hline
411     \member{setDropTolerance}& & & & & & & &\\
412     \hline
413     \member{setDropStorage}& & & & & & & &\\
414     \hline
415     \member{setRelaxationFactor}& & & & & &$\checkmark$& &\\
416     \hline
417     \member{adaptInnerTolerance}& & & & & & & &\\
418     \hline
419     \member{setInnerIterMax}& & & & & & & &\\
420 gross 2558 \end{tabular}
421     }
422 caltinay 3330 \caption{Preconditioners available for \finley and the \PASO package and the
423     relevant options in \class{SolverOptions}.
424     \label{TAB FINLEY SOLVER OPTIONS 2}}
425 gross 2558 \end{table}
426    
427 gross 2793 \section{Functions}
428 gross 2690 \begin{funcdesc}{ReadMesh}{fileName \optional{, \optional{integrationOrder=-1}, optimize=True}}
429 caltinay 3330 creates a \Domain object from the FEM mesh defined in file \var{fileName}.
430     The file must be in the \finley file format.
431     If \var{integrationOrder} is positive, a numerical integration scheme is chosen
432     which is accurate on each element up to a polynomial of degree
433     \var{integrationOrder}\index{integration order}.
434     Otherwise an appropriate integration order is chosen independently.
435     By default the labeling of mesh nodes and element distribution is optimized.
436     Set \var{optimize=False} to switch off relabeling and redistribution.
437 gross 2690 \end{funcdesc}
438    
439 gross 2748 \begin{funcdesc}{ReadGmsh}{fileName \optional{, \optional{integrationOrder=-1}, optimize=True\optional{, useMacroElements=False}}}
440 caltinay 3330 creates a \Domain object from the FEM mesh defined in file \var{fileName}.
441     The file must be in the \gmshextern file format.
442     If \var{integrationOrder} is positive, a numerical integration scheme is chosen
443     which is accurate on each element up to a polynomial of degree
444     \var{integrationOrder}\index{integration order}.
445     Otherwise an appropriate integration order is chosen independently.
446     By default the labeling of mesh nodes and element distribution is optimized.
447     Set \var{optimize=False} to switch off relabeling and redistribution.
448     If \var{useMacroElements} is set, second order elements are interpreted as
449     macro elements\index{macro elements}.
450     Currently \function{ReadGmsh} does not support \MPI.
451 jgs 102 \end{funcdesc}
452    
453 gross 2748 \begin{funcdesc}{MakeDomain}{design\optional{, integrationOrder=-1\optional{, optimizeLabeling=True\optional{, useMacroElements=False}}}}
454 caltinay 3330 creates a \finley \Domain from a \pycad \class{Design} object using \gmshextern.
455 gross 2748 The \class{Design} \var{design} defines the geometry.
456 caltinay 3330 If \var{integrationOrder} is positive, a numerical integration scheme is chosen
457     which is accurate on each element up to a polynomial of degree
458     \var{integrationOrder}\index{integration order}.
459     Otherwise an appropriate integration order is chosen independently.
460     Set \var{optimizeLabeling=False} to switch off relabeling and redistribution
461     (not recommended).
462     If \var{useMacroElements} is set, macro elements\index{macro elements} are used.
463     Currently \function{MakeDomain} does not support \MPI.
464 gross 2748 \end{funcdesc}
465 gross 2690
466 gross 2417 \begin{funcdesc}{load}{fileName}
467 caltinay 3330 recovers a \Domain object from a dump file \var{fileName} created by the
468     \function{dump} method of a \Domain object.
469 gross 2417 \end{funcdesc}
470    
471 jgs 102 \begin{funcdesc}{Rectangle}{n0,n1,order=1,l0=1.,l1=1., integrationOrder=-1, \\
472 jfenwick 3301 periodic0=\False, periodic1=\False, useElementsOnFace=\False, useMacroElements=\False,\\ optimize=\False}
473 caltinay 3330 generates a \Domain object representing a two-dimensional rectangle between
474     $(0,0)$ and $(l0,l1)$ with orthogonal edges.
475     The rectangle is filled with \var{n0} elements along the $x_0$-axis and
476 caltinay 3293 \var{n1} elements along the $x_1$-axis.
477 caltinay 3330 For \var{order}=1 and \var{order}=2, elements of type \finleyelement{Rec4} and
478 caltinay 3293 \finleyelement{Rec8} are used, respectively.
479 caltinay 3330 In the case of \var{useElementsOnFace}=\False, \finleyelement{Line2} and
480 caltinay 3293 \finleyelement{Line3} are used to subdivide the edges of the rectangle, respectively.
481 caltinay 3330 If \var{order}=-1, \finleyelement{Rec8Macro} and \finleyelement{Line3Macro}\index{macro elements}
482     are used. This option should be used when solving incompressible fluid flow
483     problems, e.g. \class{StokesProblemCartesian}.
484     In the case of \var{useElementsOnFace}=\True (this option should be used if
485     gradients are calculated on domain faces), \finleyelement{Rec4Face} and
486 caltinay 3293 \finleyelement{Rec8Face} are used on the edges, respectively.
487 caltinay 3330 If \var{integrationOrder} is positive, a numerical integration scheme is chosen
488     which is accurate on each element up to a polynomial of degree
489     \var{integrationOrder}\index{integration order}.
490     Otherwise an appropriate integration order is chosen independently.
491     If \var{periodic0}=\True, periodic boundary conditions\index{periodic boundary conditions}
492     along the $x_0$-direction are enforced.
493     That means for any solution of a PDE solved by \finley the values on the line
494     $x_0=0$ will be identical to the values on $x_0=\var{l0}$.
495     Correspondingly, \var{periodic1}=\True sets periodic boundary conditions in the
496     $x_1$-direction.
497     If \var{optimize}=\True mesh node relabeling will be attempted to reduce the
498     computation and also ParMETIS will be used to improve the mesh partition if
499     running on multiple CPUs with \MPI.
500 jgs 102 \end{funcdesc}
501    
502 jfenwick 3301 \begin{funcdesc}{Brick}{n0,n1,n2,order=1,l0=1.,l1=1.,l2=1., integrationOrder=-1,
503     periodic0=\False, periodic1=\False, \\ periodic2=\False, useElementsOnFace=\False, useMacroElements=\False, optimize=\False}
504 caltinay 3330 generates a \Domain object representing a three-dimensional brick between
505 jgs 102 $(0,0,0)$ and $(l0,l1,l2)$ with orthogonal faces. The brick is filled with
506 caltinay 3293 \var{n0} elements along the $x_0$-axis,
507     \var{n1} elements along the $x_1$-axis and
508     \var{n2} elements along the $x_2$-axis.
509 caltinay 3330 For \var{order}=1 and \var{order}=2, elements of type \finleyelement{Hex8} and
510 caltinay 3293 \finleyelement{Hex20} are used, respectively.
511 caltinay 3330 In the case of \var{useElementsOnFace}=\False, \finleyelement{Rec4} and
512 caltinay 3293 \finleyelement{Rec8} are used to subdivide the faces of the brick, respectively.
513 caltinay 3330 In the case of \var{useElementsOnFace}=\True (this option should be used if
514     gradients are calculated on domain faces), \finleyelement{Hex8Face} and
515 caltinay 3293 \finleyelement{Hex20Face} are used on the brick faces, respectively.
516 caltinay 3330 If \var{order}=-1, \finleyelement{Hex20Macro} and \finleyelement{Rec8Macro}\index{macro elements}
517     are used. This option should be used when solving incompressible fluid flow
518     problems, e.g. \class{StokesProblemCartesian}.
519     If \var{integrationOrder} is positive, a numerical integration scheme is chosen
520     which is accurate on each element up to a polynomial of degree
521     \var{integrationOrder}\index{integration order}.
522     Otherwise an appropriate integration order is chosen independently.
523     If \var{periodic0}=\True, periodic boundary conditions\index{periodic boundary conditions}
524     along the $x_0$-direction are enforced.
525     That means for any solution of a PDE solved by \finley the values on the plane
526     $x_0=0$ will be identical to the values on $x_0=\var{l0}$.
527     Correspondingly, \var{periodic1}=\True and \var{periodic2}=\True sets periodic
528     boundary conditions in the $x_1$-direction and $x_2$-direction, respectively.
529     If \var{optimize}=\True mesh node relabeling will be attempted to reduce the
530     computation and also ParMETIS will be used to improve the mesh partition if
531     running on multiple CPUs with \MPI.
532 jgs 102 \end{funcdesc}
533    
534     \begin{funcdesc}{GlueFaces}{meshList,safetyFactor=0.2,tolerance=1.e-13}
535 caltinay 3330 generates a new \Domain object from the list \var{meshList} of \finley meshes.
536     Nodes in face elements whose difference of coordinates is less than
537     \var{tolerance} times the diameter of the domain are merged.
538     The corresponding face elements are removed from the mesh.
539     %TODO: explain \var{safetyFactor} and show an example.
540 jgs 102 \end{funcdesc}
541    
542     \begin{funcdesc}{JoinFaces}{meshList,safetyFactor=0.2,tolerance=1.e-13}
543 caltinay 3330 generates a new \Domain object from the list \var{meshList} of \finley meshes.
544     Face elements whose node coordinates differ by less than \var{tolerance} times
545     the diameter of the domain are combined to form a contact element\index{element!contact}.
546 caltinay 3293 The corresponding face elements are removed from the mesh.
547 caltinay 3330 %TODO: explain \var{safetyFactor} and show an example.
548     \end{funcdesc}
549 jgs 102

Properties

Name Value
svn:eol-style native
svn:keywords Author Date Id Revision

  ViewVC Help
Powered by ViewVC 1.1.26