Welcome to Nutils’s documentation!¶
Nutils: open source numerical utilities for Python, is a collaborative programming effort aimed at the creation of a modern, general purpose programming library for Finite Element applications and related computational methods. Identifying features are a heavily object oriented design, strict separation of topology and geometry, and CAS-like function arithmetic such as found in Maple and Mathematica. Primary design goals are:
- Readability. Finite element scripts built on top of Nutils should focus on work flow and maths, unobscured by Finite Element infrastructure.
- Flexibility. The Nutils are tools; they do not enforce a strict work flow. Missing components can be added locally without loosing interoperability.
- Compatibility. Exposed objects are of native python type or allow for easy conversion to leverage third party tools.
- Speed. Nutils are self-optimizing and support parallel computation. Typical scripting inefficiencies are discouraged by design.
For latest project news and developments visit the project website at nutils.org.
Contents¶
Introduction¶
To get one thing out of the way first, note that Nutils is not your classical Finite Element program. It does not have menus, no buttons to click, nothing to make a screenshot of. To get it to do anything some programming is going to be required.
That said, let’s see what Nutils can be instead.
Design¶
Nutils is a programming library, providing components that are rich enough to handle a wide range of problems by simply linking them together. This blurs the line between classical graphical user interfaces and a programming environment, both of which serve to offer some degree of mixing and matching of available components. The former has a lower entry bar, whereas the latter offers more flexibility, the possibility to extend the toolkit with custom algorithms, and the possibility to pull in third party modules. It is our strong belief that on the edge of science where Nutils strives to be a great degree of extensibility is adamant.
For those so inclined, one of the lesser interesting possibilities this gives is to write a dedicated, Nutils powered GUI application.
What Nutils specifically does not offer are problem specific components, such as, conceivably, a “crack growth” module or “solve navier stokes” function. As a primary design principle we aim for a Nutils application to be closely readable as a high level mathematical problem description; i.e. the weak form, domain, boundary conditions, time stepping of Newton iterations, etc. It is the supporting operations like integrating over a domain or taking gradients of compound functions that are being kept out of sight as much as possible.
Quick demo¶
As a small but representative demonstration of what is involved in setting up a problem in Nutils we solve the Laplace problem on a unit square, with zero Dirichlet conditions on the left and bottom boundaries, unit flux at the top and a natural boundary condition at the right. We begin by creating a structured nelems ⅹ nelems Finite Element mesh using the built-in generator:
verts = numpy.linspace( 0, 1, nelems+1 )
domain, geom = mesh.rectilinear( [verts,verts] )
Here domain is topology representing an interconnected set of elements, and geometry is a mapping from the topology onto ℝ², representing it placement in physical space. This strict separation of topological and geometric information is key design choice in Nutils.
Proceeding to specifying the problem, we create a second order spline basis funcsp which doubles as trial and test space (u resp. v). We build a matrix by integrating laplace = ∇v · ∇u over the domain, and a rhs vector by integrating v over the top boundary. The Dirichlet constraints are projected over the left and bottom boundaries to find constrained coefficients cons. Remaining coefficients are found by solving the system in lhs. Finally these are contracted with the basis to form our solution function:
funcsp = domain.splinefunc( degree=2 )
laplace = function.outer( funcsp.grad(geom) ).sum()
matrix = domain.integrate( laplace, geometry=geom, ischeme='gauss2' )
rhs = domain.boundary['top'].integrate( funcsp, geometry=geom, ischeme='gauss1' )
cons = domain.boundary['left,bottom'].project( 0, ischeme='gauss1', geometry=geom, onto=funcsp )
lhs = matrix.solve( rhs, constrain=cons, tol=1e-8, symmetric=True )
solution = funcsp.dot(lhs)
The solution function is a mapping from the topology onto ℝ. Sampling this together with the geometry generates arrays that we can use for plotting:
points, colors = domain.elem_eval( [ geom, solution ], ischeme='bezier4', separate=True )
with plot.PyPlot( 'solution', index=index ) as plt:
plt.mesh( points, colors, triangulate='bezier' )
plt.colorbar()
Library¶
The Nutils are separated in modules focussing on topics such as mesh generation, function manipulation, debugging, plotting, etc. They are designed to form relatively independent units, though some components such as output logging run through all. Others, such as topology and element, operate in tight connection, but are divided for reasons of scope and scale. A typical Nutils application uses methods from all modules, although, as seen above, very few modules require direct access for standard computations.
What follows is an automatically generated API reference.
Topology¶
The topology module defines the topology objects, notably the StructuredTopology and UnstructuredTopology. Maintaining strict separation of topological and geometrical information, the topology represents a set of elements and their interconnectivity, boundaries, refinements, subtopologies etc, but not their positioning in physical space. The dimension of the topology represents the dimension of its elements, not that of the the space they are embedded in.
The primary role of topologies is to form a domain for nutils.function objects, like the geometry function and function bases for analysis, as well as provide tools for their construction. It also offers methods for integration and sampling, thus providing a high level interface to operations otherwise written out in element loops. For lower level operations topologies can be used as nutils.element iterators.
- class nutils.topology.Topology(ndims)[source]¶
topology base class
- elem_mean(funcs, geometry, ischeme, title='computing mean values')[source]¶
element-wise integration
- integrate(funcs, ischeme, geometry=None, iweights=None, force_dense=False, title='integrating')[source]¶
- integrate_symm(funcs, ischeme, geometry=None, iweights=None, force_dense=False, title='integrating')[source]¶
integrate a symmetric integrand on a product domain
- project(fun, onto, geometry, tol=0, ischeme=None, title='projecting', droptol=1e-08, exact_boundaries=False, constrain=None, verify=None, maxiter=0, ptype='lsqr')[source]¶
L2 projection of function onto function space
- class nutils.topology.StructuredTopology(structure, periodic=())[source]¶
structured topology
- class nutils.topology.UnstructuredTopology(elements, ndims, namedfuncs={})[source]¶
externally defined topology
Function¶
The function module defines the Evaluable class and derived objects, commonly referred to as nutils functions. They represent mappings from a nutils.topology onto Python space. The notabe class of ArrayFunc objects map onto the space of Numpy arrays of predefined dimension and shape. Most functions used in nutils applicatons are of this latter type, including the geometry and function bases for analysis.
Nutils functions are essentially postponed python functions, stored in a tree structure of input/output dependencies. Many ArrayFunc objects have directly recognizable numpy equivalents, such as Sin or Inverse. By not evaluating directly but merely stacking operations, coplex operations can be defined prior to entering a quadrature loop, allowing for a higher lever style programming. It also allows for automatic differentiation and code optimization.
It is important to realize that nutils functions do not map for a physical xy-domain but from a topology, where a point is characterized by the combination of an element and its local coordinate. This is a natural fit for typical finite element operations such as quadrature. Evaluation from physical coordinates is possible only via inverting of the geometry function, which is a fundamentally expensive and currently unsupported operation.
- exception nutils.function.EvaluationError(etype, evalue, evaluable, values)[source]¶
evaluation error
- class nutils.function.Interpolate(x, xp, fp, left=None, right=None)[source]¶
interpolate uniformly spaced data; stepwise for now
- class nutils.function.Take(func, indices, axis)[source]¶
generalization of numpy.take(), to accept lists, slices, arrays
- nutils.function.eig( arg, axes [ symmetric ] )[source]¶
Compute the eigenvalues and vectors of a matrix. The eigenvalues and vectors are positioned on the last axes.
- tuple axes The axis on which the eigenvalues and vectors are calculated
- bool symmetric Is the matrix symmetric
- int sort Sort the eigenvalues and vectors (-1=descending, 0=unsorted, 1=ascending)
- nutils.function.fdapprox(func, w, dofs, delta=1e-05)[source]¶
Finite difference approximation of the variation of func in directions w around dofs. Input arguments: * func, the functional to differentiate * dofs, DOF vector of linearization point * w, the function space or a tuple of chained spaces * delta, finite difference step scaling of ||dofs||_inf
Core¶
The core module provides a collection of low level constructs that have no dependencies on other nutils modules. Primarily for internal use.
Debug¶
The debug module provides code inspection tools and the “traceback explorer” interactive shell environment. Access to these components is primarily via breakpoint() and an exception handler in nutils.util.run().
Element¶
The element module defines reference elements such as the QuadElement and TriangularElement, but also more exotic objects like the TrimmedElement. A set of (interconnected) elements together form a nutils.topology. Elements have edges and children (for refinement), which are in turn elements and map onto self by an affine transformation. They also have a well defined reference coordinate system, and provide pointsets for purposes of integration and sampling.
- class nutils.element.TrimmedIScheme(levelset, ischeme, maxrefine, finestscheme='uniform1', degree=3, retain=None)[source]¶
integration scheme for truncated elements
- class nutils.element.SliceTransformation(fromdim, start=None, stop=None, step=None)[source]¶
take slice
- class nutils.element.Element(ndims, vertices, index=None, parent=None, context=None, interface=None)[source]¶
Element base class.
Represents the topological shape.
- class nutils.element.ProductElement(elem1, elem2)[source]¶
element product
- orientation[source]¶
Neighborhood of elem1 and elem2 and transformations to get mutual overlap in right location. Returns 3-element tuple: * neighborhood, as given by Element.neighbor(), * transf1, required rotation of elem1 map: {0:0, 1:pi/2, 2:pi, 3:3*pi/2}, * transf2, required rotation of elem2 map (is indep of transf1 in UnstructuredTopology.
- class nutils.element.TrimmedElement(elem, levelset, maxrefine, lscheme, finestscheme, evalrefine, parent, vertices)[source]¶
trimmed element
- class nutils.element.QuadElement(ndims, vertices, index=None, parent=None, context=None, interface=None)[source]¶
quadrilateral element
- class nutils.element.TriangularElement(vertices, index=None, parent=None, context=None)[source]¶
Triangular element. Conventions: * reference elem: unit simplex {(x,y) | x>0, y>0, x+y<1} * vertex numbering: {(1,0):0, (0,1):1, (0,0):2} * edge numbering: {bottom:0, slanted:1, left:2} * edge local coords run counter-clockwise.
- class nutils.element.TetrahedronElement(vertices, index=None, parent=None, context=None)[source]¶
tetrahedron element
- class nutils.element.PolyLine(poly)[source]¶
polynomial on a line
- class nutils.element.PolyTriangle[source]¶
poly triangle (linear for now) conventions: dof numbering as vertices, see TriangularElement docstring.
Library¶
The library module provides a collection of application specific functions, that nevertheless have a wide enough range of applicability to be useful as generic building blocks.
Log¶
The log module provides print methods debug, info, user, warning, and error, in increasing order of priority. Output is sent to stdout as well as to an html formatted log file if so configured.
- class nutils.log.ProgressContextLog(text, iterable=None, target=None, showpct=True, depth=1)[source]¶
progress bar
- nutils.log.context¶
alias of StaticContextLog
- nutils.log.iterate¶
alias of ProgressContextLog
- nutils.log.progress¶
alias of ProgressContextLog
Matrix¶
The matrix module defines a number of 2D matrix objects, notably the SparseMatrix() and DenseMatrix(). Matrix objects support indexed addition, basic addition and subtraction operations, and provide a consistent insterface for solving linear systems. Matrices can be converted to numpy arrays via asarray.
- nutils.matrix.krylov(matvec, b, x0=None, tol=1e-05, restart=None, maxiter=0, precon=None, callback=None)[source]¶
solve linear system iteratively
restart=None: CG restart=integer: GMRES
Mesh¶
The mesh module provides mesh generators: methods that return a topology and an accompanying geometry function. Meshes can either be generated on the fly, e.g. rectilinear(), or read from external an externally prepared file, gmesh(), igatool(), and converted to nutils format. Note that no mesh writers are provided at this point; output is handled by the nutils.plot module.
Numeric¶
The numeric module provides methods that are lacking from the numpy module. An accompanying extension module _numeric.c should be compiled to benefit from extra performance, although a Python-only implementation is provided as fallback. A warning message is printed if the extension module is not found.
- nutils.numeric.align(arr, trans, ndim)[source]¶
create new array of ndim from arr with axes moved accordin to trans
- nutils.numeric.dot(A, B, axis=-1)[source]¶
Transform axis of A by contraction with first axis of B and inserting remaining axes. Note: with default axis=-1 this leads to multiplication of vectors and matrices following linear algebra conventions.
- nutils.numeric.eig(A, sort=False)[source]¶
Compute the eigenvalues and vectors of a hermitian matrix sort -1/0/1 -> descending / unsorted / ascending
Parallel¶
The parallel module provides tools aimed at parallel computing. At this point all parallel solutions use the fork system call and are supported on limited platforms, notably excluding Windows. On unsupported platforms parallel features will disable and a warning is printed.
- class nutils.parallel.AlternativeFork(nprocs)[source]¶
single master, multiple slave fork context, unwinds at exit
Util¶
The util module provides a collection of general purpose methods. Most importantly it provides the run() method which is the preferred entry point of a nutils application, taking care of command line parsing, output dir creation and initiation of a log file.
Plot¶
The plot module aims to provide a consistent interface to various plotting backends. At this point matplotlib and vtk are supported.
- class nutils.plot.PyPlot(name, imgtype=None, ndigits=3, index=None, **kwargs)[source]¶
matplotlib figure
- mesh(points, colors=None, edgecolors='k', edgewidth=None, triangulate='delaunay', setxylim=True, **kwargs)[source]¶
plot elemtwise mesh
- slope_triangle(x, y, fillcolor='0.9', edgecolor='k', xoffset=0, yoffset=0.1, slopefmt='{0:.1f}')[source]¶
Draw slope triangle for supplied y(x) - x, y: coordinates - xoffset, yoffset: distance graph & triangle (points) - fillcolor, edgecolor: triangle style - slopefmt: format string for slope number
- nutils.plot.writevtu(name, topo, coords, pointdata={}, celldata={}, ascii=False, superelements=False, maxrefine=3, ndigits=0, ischeme='gauss1', **kwargs)[source]¶
write vtu from coords function
- class nutils.plot.PylabAxis(ax, title)[source]¶
matplotlib axis augmented with nutils-specific functions
- add_mesh(coords, topology, deform=0, color=None, edgecolors='none', linewidth=1, xmargin=0, ymargin=0, aspect='equal', cbar='vertical', title=None, ischeme='gauss2', cscheme='contour3', clim=None, frame=True, colormap=None)[source]¶
plot mesh