# BEAST.jl documentation

BEAST provides a number of types modelling concepts and a number of algorithms for the efficient and simple implementation of boundary and finite element solvers. It provides full implementations of these concepts for the LU based solution of boundary integral equations for the Maxwell and Helmholtz systems.

Because Julia only compiles code at execution time, users of this library can hook into the code provided in this package at any level. In the extreme case it suffices to provide overwrites of the `assemble`

functions. In that case, only the LU solution will be performed by the code here.

At the other end it suffices that users only supply integration kernels that act on the element-element interaction level. This package will manage all required steps for matrix assembly.

For the Helmholtz 2D and Maxwell 3D systems, complete implementations are supplied. These models will be discussed in detail to give a more concrete idea of the APIs provides and how to extend them.

Central to the solution of boundary integral equations is the assembly of the system matrix. The system matrix is fully determined by specifying a kernel G, a set of trial functions, and a set of test functions.

## Basis

Sets of both trial and testing functions are implemented by models following the basis concept. The term basis is somewhat misleading as it is nowhere required nor enforced that these functions are linearly independent. Models implementing the Basis concept need to comply to the following semantics.

`numfunctions(basis)`

: number of functions in the Basis.`coordtype(basis)`

: type of (the components of) the values taken on by the functions in the Basis.`scalartype(d)`

: the scalar field underlying the vector space the basis functions take value in.`refspace(basis)`

: returns the ReferenceSpace of local shape functions on which the Basis is built.`assemblydata(basis)`

:`assemblydata`

returns an iterable collection`elements`

of geometric elements and a look table`ad`

for use in assembly of interaction matrices. In particular, for an index`element_idx`

into`elements`

and an index`local_shape_idx`

in basis of local shape functions`refspace(basis)`

,`ad[element_idx, local_shape_idx]`

returns the iterable collection of`(global_idx, weight)`

tuples such that the local shape function at`local_shape_idx`

defined on the element at`element_idx`

contributes to the basis function at`global_idx`

with a weight of`weight`

.`geometry(basis)`

: returns an iterable collection of Elements. The order in which these Elements are encountered corresponds to the indices used in the assembly data structure.

## Reference Space

The *reference space* concept defines an API for working with spaces of local shape functions. The main role of objects implementing this concept is to allow specialization of the functions that depend on the precise reference space used.

The functions that depend on the type and value of arguments modeling *reference space* are:

`numfunctions(refspace)`

: returns the number of shape functions on each element.

## Kernel

A kernel is a fairly simple concept that mainly exists as part of the definition of a Discrete Operator. A kernel should obey the following semantics:

In many function definitions the kernel object is referenced by `operator`

or something similar. This is a misleading name as an operator definition should always be accompanied by the domain and range space.

## Discrete Operator

Informally speaking, a Discrete Operator is a concept that allows for the computation of an interaction matrix. It is a kernel together with a test and trial basis. A Discrete Operator can be passed to `assemble`

and friends to compute its matrix representation.

A discrete operator is a triple `(kernel, test_basis, trial_basis)`

, where `kernel`

is a Kernel, and `test_basis`

and `trial_basis`

are Bases. In addition, the following expressions should be implemented and behave according to the correct semantics:

`quaddata(operator,test_refspace,trial_refspace,test_elements,trial_elements)`

: create the data required for the computation of element-element interactions during assembly of discrete operator matrices.`quadrule(operator,test_refspace,trial_refspace,p,test_element,q_trial_element,qd)`

: returns an integration strategy object that will be passed to`momintegrals!`

to select an integration strategy. This rule can depend on the test/trial reference spaces and interacting elements. The indices`p`

and`q`

refer to the position of the interacting elements in the enumeration defined by`geometry(basis)`

and allow for fast retrieval of any element specific data stored in the quadrature data object`qd`

.`momintegrals!(operator,test_refspace,trial_refspace,test_element,trial_element,zlocal,qr)`

: this function computes the local interaction matrix between the set of local test and trial shape functions and a specific pair of elements. The target matrix`zlocal`

is provided as an argument to minimise memory allocations over subsequent calls.`qr`

is an object returned by`quadrule`

and contains all static and dynamic data defining the integration strategy used.

In the context of fast methods such as the Fast Multipole Method other algorithms on Discrete Operators will typically be defined to compute matrix vector products. These algorithms do not explicitly compute and store the interaction matrix (this would lead to unacceptable computational and memory complexity).

`BEAST.elements`

— Functionelements(geo)

Create an iterable collection of the elements stored in `geo`

. The order in which this collection produces the elements determines the index used for lookup in the data structures returned by `assemblydata`

and `quaddata`

.

`BEAST.numfunctions`

— Function`numfunctions(r)`

Return the number of functions in a `Space`

or `RefSpace`

.

`CompScienceMeshes.coordtype`

— Function`coordtype(mesh)`

Returns `eltype(vertextype(mesh))`

`coordtype(simplex)`

Return coordinate type used by simplex.

`BEAST.scalartype`

— Function`scalartype(x)`

The scalar field over which the values of a global or local basis function, or an operator are defined. This should always be a scalar type, even if the basis or operator takes on values in a vector or tensor space. This data type is used to determine the `eltype`

of assembled discrete operators.

`BEAST.assemblydata`

— Function`charts, admap = assemblydata(basis)`

Given a Basis this function returns a data structure containing the information required for matrix assemble. More precise the following expressions are valid for the returned object `ad`

:

```
ad[c,r,i].globalindex
ad[c,r,i].coefficient
```

Here, `c`

and `r`

are indices in the iterable set of geometric elements and the set of local shape functions on each element. `i`

ranges from 1 to the maximum number of basis functions local shape function `r`

on element `r`

contributes to.

For a triplet `(c,r,i)`

, `globalindex`

is the index in the Basis of the `i`

-th basis function that has a contribution from local shape function `r`

on element `r`

. `coefficient`

is the coefficient of that contribution in the linear combination defining that basis function in terms of local shape function.

*Note*: the indices `c`

correspond to the position of the corresponding element whilst iterating over `geometry(basis)`

.

`BEAST.geometry`

— Function`geometry(basis)`

Returns an iterable collection of geometric elements on which the functions in `basis`

are defined. The order the elements are encountered needs correspond to the element indices used in the data structure returned by `assemblydata`

.

`BEAST.refspace`

— Function`refspace(basis)`

Returns the ReferenceSpace of local shape functions on which the basis is built.

`BEAST.quaddata`

— Function`quaddata(operator, test_refspace, trial_refspace, test_elements, trial_elements)`

Returns an object cashing data required for the computation of boundary element interactions. It is up to the client programmer to decide what (if any) data is cached. For double numberical quadrature, storing the integration points for example can significantly speed up matrix assembly.

`operator`

is an integration kernel.`test_refspace`

and`trial_refspace`

are reference space objects.`quadata`

is typically overloaded on the type of these local spaces of shape functions. (See the implementation in `maxwell.jl`

for an example).

`test_elements`

and`trial_elements`

are iterable collections of the geometric

elements on which the finite element space are defined. These are provided to allow computation of the actual integrations points - as opposed to only their coordinates.

`BEAST.quadrule`

— Function`quadrule(operator,test_refspace,trial_refspace,p,test_element,q_trial_element, qd)`

Returns an object that contains all the dynamic (runtime) information that defines the integration strategy that will be used by `momintegrals!`

to compute the interactions between the local test/trial functions defined on the specified geometric elements. The indices `p`

and `q`

refer to the position of the test and trial elements as encountered during iteration over the output of `geometry`

.

The last argument `qd`

provides access to all precomputed data required for quadrature. For example it might be desirable to precompute all the quadrature points for all possible numerical quadrature schemes that can potentially be required during matrix assembly. This makes sense, since the number of point is order N (where N is the number of faces) but these points will appear in N^2 computations. Precomputation requires some extra memory but can save a lot on computation time.

`BEAST.momintegrals!`

— Functionmomintegrals!(biop, tshs, bshs, tcell, bcell, interactions, strat)

Function for the computation of moment integrals using simple double quadrature.