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 nelemsnelems 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

refined_by(refine)[source]

create refined space by refining dofs in existing one

elem_eval(funcs, ischeme, separate=False, title='evaluating')[source]

element-wise evaluation

elem_mean(funcs, geometry, ischeme, title='computing mean values')[source]

element-wise integration

grid_eval(funcs, geometry, C, title='grid-evaluating')[source]

evaluate grid points

build_graph(func)[source]

get matrix sparsity

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

projection(fun, onto, geometry, **kwargs)[source]

project and return as function

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

refinedfunc(dofaxis, refine, degree, title='refining')[source]

create refined space by refining dofs in existing one

refine(n)[source]

refine entire topology n times

get_simplices(maxrefine, title='getting simplices')[source]

Getting simplices

get_trimmededges(maxrefine, title='getting trimmededges')[source]

Getting trimmed edges

class nutils.topology.StructuredTopology(structure, periodic=())[source]

structured topology

make_periodic(periodic)[source]

add periodicity

linearfunc()[source]

linears

rectilinearfunc(gridvertices)[source]

rectilinear func

refined[source]

refine entire topology

trim(levelset, maxrefine, lscheme='bezier3', finestscheme='uniform2', evalrefine=0, title='trimming', log=<module 'nutils.log' from '/var/build/user_builds/nutils/checkouts/v1.0/nutils/log.pyc'>)[source]

trim element along levelset

neighbor(elem0, elem1)[source]

Neighbor detection, returns codimension of interface, -1 for non-neighboring elements.

class nutils.topology.IndexedTopology(topo, elements)[source]

trimmed topology

splinefunc(degree)[source]

create spline function space

class nutils.topology.UnstructuredTopology(elements, ndims, namedfuncs={})[source]

externally defined topology

splinefunc(degree)[source]

spline func

linearfunc()[source]

linear func

bubblefunc()[source]

linear func + bubble

class nutils.topology.HierarchicalTopology(basetopo, elements)[source]

collection of nested topology elments

class nutils.topology.ElemMap(mapping, ndims)[source]

dictionary-like element mapping

nutils.topology.glue(master, slave, geometry, tol=1e-10, verbose=False)[source]

Glue topologies along boundary group __glue__.

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.

class nutils.function.Evaluable(args, evalf)[source]

Base class

recurse_index(data, operations, cbuild)[source]

compile

compile()[source]
graphviz(title='graphviz')[source]

create function graph

stackstr(values=None)[source]

print stack

asciitree()[source]

string representation

exception nutils.function.EvaluationError(etype, evalue, evaluable, values)[source]

evaluation error

class nutils.function.Tuple(items)[source]

combine

static vartuple(*f)[source]

evaluate

class nutils.function.PointShape[source]

shape of integration points

static pointshape(points)[source]

evaluate

class nutils.function.Cascade(ndims, side=0)[source]

point cascade: list of (elem,points) tuples

static cascade(elem, points, ndims, side)[source]

evaluate

class nutils.function.ArrayFunc(evalf, args, shape, dtype=<type 'float'>)[source]

array function

vector(ndims)[source]

vectorize

dot(weights, axis=0)[source]

array contraction

find(elem, C)[source]

iteratively find x for f(x) = target, starting at x=start

normalized()[source]

normalize last axis

normal(ndims=-1)[source]
curvature(ndims=-1)[source]
swapaxes(n1, n2)[source]

swap axes

transpose(trans=None)[source]
grad(coords, ndims=0)[source]

gradient

laplace(coords, ndims=0)[source]

laplacian

add_T(axes=(-2, -1))[source]

add transposed

symgrad(coords, ndims=0)[source]

gradient

div(coords, ndims=0)[source]

gradient

dotnorm(coords, ndims=0)[source]

normal component

ngrad(coords, ndims=0)[source]

normal gradient

nsymgrad(coords, ndims=0)[source]

normal gradient

T[source]

transpose

class nutils.function.ElemArea(weights)[source]

element area

static elemarea(weights)[source]

evaluate

class nutils.function.ElemInt(func, weights)[source]

elementwise integration

static elemint(w, f, ndim)[source]

evaluate

class nutils.function.Align(func, axes, ndim)[source]

align axes

static align(arr, trans, ndim)[source]
class nutils.function.Get(func, axis, item)[source]

get

class nutils.function.Product(func, axis)[source]

product

class nutils.function.IWeights[source]

integration weights

static iweights(elem, weights)[source]

evaluate

class nutils.function.OrientationHack(side=0)[source]

orientation hack for 1d elements; VERY dirty

static orientation(elem, side)[source]

evaluate

class nutils.function.Transform(fromcascade, tocascade, side=0)[source]

transform

static transform(fromcascade, tocascade, side)[source]
class nutils.function.Function(cascade, stdmap, igrad, axis)[source]

function

static function(cascade, stdmap, igrad)[source]

evaluate

class nutils.function.Choose(level, choices)[source]

piecewise function

static choose(level, *choices)[source]
class nutils.function.Choose2D(coords, contour, fin, fout)[source]

piecewise function

static choose2d(xy, contour, fin, fout)[source]

evaluate

class nutils.function.Inverse(func)[source]

inverse

class nutils.function.DofMap(cascade, dofmap, axis)[source]

dof axis

static evalmap(cascade, dofmap)[source]

evaluate

class nutils.function.Concatenate(funcs, axis=0)[source]

concatenate

static concatenate(iax, *arrays)[source]

evaluate

class nutils.function.Interpolate(x, xp, fp, left=None, right=None)[source]

interpolate uniformly spaced data; stepwise for now

class nutils.function.Cross(func1, func2, axis)[source]

cross product

class nutils.function.Determinant(func)[source]

normal

class nutils.function.DofIndex(array, iax, index)[source]

element-based indexing

static dofindex(arr, *item)[source]

evaluate

class nutils.function.Multiply(func1, func2)[source]

multiply

cdef()[source]

generate C code

class nutils.function.Negative(func)[source]

negate

class nutils.function.Add(func1, func2)[source]

add

cdef()[source]

generate C code

class nutils.function.BlockAdd(func1, func2)[source]

block addition (used for DG)

class nutils.function.Dot(func1, func2, naxes)[source]

dot

cdef()[source]

generate C code

class nutils.function.Sum(func, axis)[source]

sum

class nutils.function.Debug(func)[source]

debug

static debug(arr)[source]
class nutils.function.TakeDiag(func)[source]

extract diagonal

class nutils.function.Take(func, indices, axis)[source]

generalization of numpy.take(), to accept lists, slices, arrays

class nutils.function.Power(func, power)[source]

power

class nutils.function.ElemFunc(domainelem, side=0)[source]

trivial func

static elemfunc(cascade, domainelem)[source]

evaluate

find(elem, C)[source]

find coordinates

class nutils.function.Pointwise(args, evalf, deriv)[source]

pointwise transformation

class nutils.function.Sign(func)[source]

sign

class nutils.function.Eig(func, symmetric=False, sort=False)[source]

Eig

class nutils.function.Zeros(shape)[source]

zero

static zeros(points, shape)[source]

prepend point axes

class nutils.function.Inflate(func, dofmap, length, axis)[source]

inflate

static inflate(array, indices, length, axis)[source]
class nutils.function.Diagonalize(func)[source]

diagonal matrix

class nutils.function.Repeat(func, length, axis)[source]

repeat singleton axis

class nutils.function.Const(func)[source]

pointwise transformation

static const(points, arr)[source]

prepend point axes

nutils.function.asarray(arg)[source]

convert to ArrayFunc or numpy.ndarray

nutils.function.insert(arg, n)[source]

insert axis

nutils.function.stack(args, axis=0)[source]

stack functions along new axis

nutils.function.chain(funcs)[source]
nutils.function.merge(funcs)[source]

Combines unchained funcs into one function object.

nutils.function.vectorize(args)[source]
nutils.function.expand(arg, shape)[source]
nutils.function.repeat(arg, length, axis)[source]
nutils.function.get(arg, iax, item)[source]

get item

nutils.function.align(arg, axes, ndim)[source]
nutils.function.bringforward(arg, axis)[source]

bring axis forward

nutils.function.elemint(arg, weights)[source]

elementwise integration

nutils.function.grad(arg, coords, ndims=0)[source]

local derivative

nutils.function.symgrad(arg, coords, ndims=0)[source]

symmetric gradient

nutils.function.div(arg, coords, ndims=0)[source]

gradient

nutils.function.sum(arg, axes=-1)[source]

sum over multiply axes

nutils.function.dot(arg1, arg2, axes)[source]

dot product

nutils.function.determinant(arg, axes=(-2, -1))[source]
nutils.function.inverse(arg, axes=(-2, -1))[source]
nutils.function.takediag(arg, ax1=-2, ax2=-1)[source]
nutils.function.localgradient(arg, ndims)[source]

local derivative

nutils.function.dotnorm(arg, coords, ndims=0)[source]

normal component

nutils.function.kronecker(arg, axis, length, pos)[source]
nutils.function.diagonalize(arg)[source]
nutils.function.concatenate(args, axis=0)[source]
nutils.function.transpose(arg, trans=None)[source]
nutils.function.product(arg, axis)[source]
nutils.function.choose(level, choices)[source]
nutils.function.cross(arg1, arg2, axis)[source]

cross product

nutils.function.outer(arg1, arg2=None, axis=0)[source]

outer product

nutils.function.pointwise(args, evalf, deriv)[source]

general pointwise operation

nutils.function.multiply(arg1, arg2)[source]
nutils.function.add(arg1, arg2)[source]
nutils.function.negative(arg)[source]

make negative

nutils.function.power(arg, n)[source]
nutils.function.sign(arg)[source]
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.swapaxes(arg, axes=(-2, -1))[source]

swap axes

nutils.function.opposite(arg)[source]

evaluate jump over interface

nutils.function.function(fmap, nmap, ndofs, ndims)[source]

create function on ndims-element

nutils.function.take(arg, index, axis)[source]

take index

nutils.function.inflate(arg, dofmap, length, axis)[source]
nutils.function.pointdata(topo, ischeme, func=None, shape=None, value=None)[source]

point data

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

nutils.function.iwscale(coords, ndims)[source]

integration weights scale

nutils.function.supp(funcsp, indices)[source]

find support of selection of basis functions

class nutils.function.CBuilder(cachedir='/tmp/nutils')[source]

cbuilder

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().

class nutils.debug.Frame(frame, lineno=None)[source]

frame info

class nutils.debug.Explorer(exc, frames, intro)[source]

traceback explorer

show_context()[source]

show traceback up to index

do_s(arg)[source]

Show source code of the currently focussed frame.

do_l(arg)[source]

List the stack and exception type

do_q(arg)[source]

Quit traceback exploror.

do_u(arg)[source]

Shift focus to the frame above the current one.

do_d(arg)[source]

Shift focus to the frame below the current one.

do_w(arg)[source]

Show overview of local variables.

do_p(arg)[source]

Print local of global variable, or function evaluation.

onecmd(text)[source]

wrap command handling to avoid a second death

do_pp(arg)[source]

Pretty-print local of global variable, or function evaluation.

completedefault(text, line, begidx, endidx)[source]

complete object names

nutils.debug.exception()[source]

constructor

nutils.debug.write_html(out, exc, frames)[source]

write exception info to html file

nutils.debug.breakpoint()[source]

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

generate_ischeme(elem, maxrefine)[source]

generate integration scheme

class nutils.element.Transformation(fromdim, todim)[source]

transform points

eval(points)[source]

evaluate

class nutils.element.SliceTransformation(fromdim, start=None, stop=None, step=None)[source]

take slice

class nutils.element.AffineTransformation(offset, transform)[source]

affine transformation

invtrans[source]

inverse transformation

det[source]

determinant

nest(other)[source]

merge transformations

get_transform()[source]

get transformation copy

invapply(coords)[source]

apply inverse transformation

class nutils.element.Element(ndims, vertices, index=None, parent=None, context=None, interface=None)[source]

Element base class.

Represents the topological shape.

neighbor(other)[source]

level of neighborhood; 0=self

eval(where)[source]

get points

zoom(elemset, points)[source]

zoom points

intersected(levelset, lscheme, evalrefine=0)[source]

check levelset intersection:

+1 for levelset > 0 everywhere -1 for levelset < 0 everywhere =0 for intersected element

trim(levelset, maxrefine, lscheme, finestscheme, evalrefine)[source]

trim element along levelset

get_simplices(maxrefine)[source]

divide in simple elements

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.

eval(where)[source]

get integration scheme

class nutils.element.TrimmedElement(elem, levelset, maxrefine, lscheme, finestscheme, evalrefine, parent, vertices)[source]

trimmed element

children[source]

all 1x refined elements

edge(iedge)[source]
get_simplices(maxrefine)[source]

divide in simple elements

class nutils.element.QuadElement(ndims, vertices, index=None, parent=None, context=None, interface=None)[source]

quadrilateral element

children_by(N)[source]

divide element by n

children[source]

all 1x refined elements

ribbons[source]

ribbons

edge(iedge)[source]
refine(n)[source]

refine n times

select_contained(points, eps=0)[source]

select points contained in 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.

children[source]

all 1x refined elements

edge(iedge)[source]
refined(n)[source]

refine

select_contained(points, eps=0)[source]

select points contained in element

class nutils.element.TetrahedronElement(vertices, index=None, parent=None, context=None)[source]

tetrahedron element

children[source]

all 1x refined elements

edge(iedge)[source]
refined(n)[source]

refine

select_contained(points, eps=0)[source]

select points contained in element

class nutils.element.StdElem[source]

stdelem base class

extract(extraction)[source]

apply extraction matrix

class nutils.element.PolyProduct[source]

multiply standard elements

class nutils.element.PolyLine(poly)[source]

polynomial on a line

classmethod bernstein_poly(degree)[source]

bernstein polynomial coefficients

classmethod spline_poly(p, n)[source]

spline polynomial coefficients

classmethod spline(degree, nelems, periodic=False, neumann=0, curvature=False)[source]

spline elements, any amount

extract(extraction)[source]

apply extraction

class nutils.element.PolyTriangle[source]

poly triangle (linear for now) conventions: dof numbering as vertices, see TriangularElement docstring.

class nutils.element.BubbleTriangle[source]

linear triangle + bubble function conventions: dof numbering as vertices (see TriangularElement docstring), then barycenter.

class nutils.element.ExtractionWrapper(stdelem, extraction)[source]

extraction wrapper

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.ContextLog(depth=1)[source]

base class

disable()[source]

disable this logger

class nutils.log.HtmlLog(fileobj, title, depth=1)[source]

html log

class nutils.log.HtmlStream(chunks, attr, html)[source]

html line stream

write(text)[source]

write to out and buffer for html

class nutils.log.ProgressContextLog(text, iterable=None, target=None, showpct=True, depth=1)[source]

progress bar

text[source]

get text

update(current)[source]

update progress

nutils.log.SimpleLog(chunks=('', ), attr=None)[source]

just write to stdout

class nutils.log.StaticContextLog(text, depth=1)[source]

simple text logger

nutils.log.context

alias of StaticContextLog

nutils.log.iterate

alias of ProgressContextLog

nutils.log.progress

alias of ProgressContextLog

nutils.log.setup_html

alias of HtmlLog

nutils.log.stack(msg, frames=None)[source]

print stack trace

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

nutils.matrix.parsecons(constrain, lconstrain, rconstrain, shape)[source]

parse constraints

class nutils.matrix.Matrix((nrows, ncols))[source]

matrix base class

cond(constrain=None, lconstrain=None, rconstrain=None)[source]

condition number

res(x, b=0, constrain=None, lconstrain=None, rconstrain=None)[source]

residual

class nutils.matrix.DenseSubMatrix(data, indices)[source]

dense but non-contiguous data

class nutils.matrix.SparseMatrix(graph, ncols=None)[source]

sparse matrix

reshape((nrows, ncols))[source]

reshape matrix

clone()[source]

clone matrix

matvec(other)[source]

matrix-vector multiplication

T[source]

transpose

toarray()[source]

convert to numpy array

todense()[source]

convert to dense matrix

rowsupp(tol=0)[source]

return row indices with nonzero/non-small entries

solve(b=0, constrain=None, lconstrain=None, rconstrain=None, tol=0, x0=None, symmetric=False, maxiter=0, restart=999, title='solving system', callback=None, precon=None)[source]
class nutils.matrix.DenseMatrix(shape)[source]

matrix wrapper class

clone()[source]

clone matrix

addblock(rows, cols, vals)[source]

add matrix data

toarray()[source]

convert to numpy array

matvec(other)[source]

matrix-vector multiplication

T[source]

transpose

solve(b=0, constrain=None, lconstrain=None, rconstrain=None, title='solving system', **dummy)[source]

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.

nutils.mesh.rectilinear(vertices, periodic=(), name='rect')[source]

rectilinear mesh

nutils.mesh.revolve(topo, coords, nelems, degree=3, axis=0)[source]

revolve coordinates

nutils.mesh.gmesh(path, btags={}, name=None)[source]
nutils.mesh.triangulation(vertices, nvertices)[source]
nutils.mesh.igatool(path, name=None)[source]

igatool mesh

nutils.mesh.fromfunc(func, nelems, ndims, degree=1)[source]

piecewise

nutils.mesh.demo(xmin=0, xmax=1, ymin=0, ymax=1)[source]

demo triangulation of a rectangle

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.normdim(ndim, n)[source]

check bounds and make positive

nutils.numeric.align(arr, trans, ndim)[source]

create new array of ndim from arr with axes moved accordin to trans

nutils.numeric.get(arr, axis, item)[source]

take single item from array axis

nutils.numeric.expand(arr, *shape)[source]
nutils.numeric.linspace2d(start, stop, steps)[source]

linspace & meshgrid combined

nutils.numeric.contract(A, B, axis=-1)[source]
nutils.numeric.contract_fast(A, B, naxes)[source]

contract last n axes

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.fastrepeat(A, nrepeat, axis=-1)[source]

repeat axis by 0stride

nutils.numeric.fastmeshgrid(X, Y)[source]

mesh grid based on fastrepeat

nutils.numeric.meshgrid(*args)[source]

multi-dimensional meshgrid generalisation

nutils.numeric.appendaxes(A, shape)[source]

append axes by 0stride

nutils.numeric.inverse(A)[source]

linearized inverse

nutils.numeric.determinant(A)[source]
nutils.numeric.eig(A, sort=False)[source]

Compute the eigenvalues and vectors of a hermitian matrix sort -1/0/1 -> descending / unsorted / ascending

nutils.numeric.eigh(A, sort=False)[source]

Compute the eigenvalues and vectors of a hermitian matrix sort -1/0/1 -> descending / unsorted / ascending

nutils.numeric.reshape(A, *shape)[source]

more useful reshape

nutils.numeric.mean(A, weights=None, axis=-1)[source]

generalized mean

nutils.numeric.norm2(A, axis=-1)[source]

L2 norm over specified axis

nutils.numeric.normalize(A, axis=-1)[source]

devide by normal

nutils.numeric.cross(v1, v2, axis)[source]

cross product

nutils.numeric.stack(arrays, axis=0)[source]

powerful array stacker with singleton expansion

nutils.numeric.bringforward(arg, axis)[source]

bring axis forward

nutils.numeric.diagonalize(arg)[source]

append axis, place last axis on diagonal of self and new

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.Fork(nprocs)[source]

nested fork context, unwinds at exit

class nutils.parallel.AlternativeFork(nprocs)[source]

single master, multiple slave fork context, unwinds at exit

nutils.parallel.fork(func, nice=19)[source]

fork and run (return value is lost)

nutils.parallel.shzeros(shape, dtype=<type 'float'>)[source]

create zero-initialized array in shared memory

nutils.parallel.pariter(iterable)[source]

iterate parallel

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.

class nutils.util.ImmutableArray[source]

immutable array

nutils.util.delaunay(points)[source]

delaunay triangulation

nutils.util.withrepr(f)[source]

add string representation to generated function

class nutils.util.Cache(*args)[source]

cache

nutils.util.getpath(pattern)[source]

create file in dumpdir

nutils.util.sum(seq)[source]

a better sum

nutils.util.product(seq)[source]

multiply items in sequence

nutils.util.clone(obj)[source]

clone object

nutils.util.iterate(context='iter', nmax=-1)[source]

iterate forever

class nutils.util.NanVec[source]

nan-initialized vector

where[source]

find non-nan items

class nutils.util.Clock(interval)[source]

simple interval timer

nutils.util.tensorial(args)[source]

create n-dimensional array containing tensorial combinations of n args

nutils.util.arraymap(f, dtype, *args)[source]

call f for sequence of arguments and cast to dtype

nutils.util.objmap(func, *arrays)[source]

map numpy arrays

nutils.util.fail(msg, *args)[source]

generate exception

class nutils.util.Locals[source]

local namespace as object

nutils.util.getkwargdefaults(func)[source]

helper for run

class nutils.util.Statm(rusage=None)[source]

memory statistics on systems that support it

nutils.util.run(*functions)[source]

call function specified on command line

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.BasePlot(name, ndigits=0, index=None)[source]

base class for plotting objects

class nutils.plot.PyPlot(name, imgtype=None, ndigits=3, index=None, **kwargs)[source]

matplotlib figure

save(name)[source]

save images

mesh(points, colors=None, edgecolors='k', edgewidth=None, triangulate='delaunay', setxylim=True, **kwargs)[source]

plot elemtwise mesh

polycol(verts, facecolors='none', **kwargs)[source]

add polycollection

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

slope_trend(x, y, lt='k-', xoffset=0.1, slopefmt='{0:.1f}')[source]

Draw slope triangle for supplied y(x) - x, y: coordinates - slopefmt: format string for slope number

rectangle(x0, w, h, fc='none', ec='none', **kwargs)[source]
griddata(xlim, ylim, data)[source]

plot griddata

cspy(A, **kwargs)[source]

Like pyplot.spy, but coloring acc to 10^log of absolute values, where [0, inf, nan] show up in blue.

class nutils.plot.DataFile(name, index=None, ndigits=0, mode='w')[source]

data file

class nutils.plot.VTKFile(name, index=None, ndigits=0, ascii=False)[source]

vtk file

unstructuredgrid(points, npars=None)[source]

add unstructured grid

celldataarray(name, data)[source]

add cell array

pointdataarray(name, data)[source]

add cell array

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.Pylab(title, name='graph{0:03x}')[source]

matplotlib figure

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

add_quiver(coords, topology, quiver, sample='uniform3', scale=None)[source]

quiver builder

add_graph(xfun, yfun, topology, sample='contour10', logx=False, logy=False, **kwargs)[source]

plot graph of function on 1d topology

add_convplot(x, y, drop=0.8, shift=1.1, slope=True, **kwargs)[source]

Convergence plot including slope triangle (below graph) for supplied y(x), drop = distance graph & triangle, shift = distance triangle & text.

nutils.plot.preview(coords, topology, cscheme='contour8')[source]

preview function

Indices and tables