MWTree

This is an introduction to the mwtree class. We write a small overarching summary of the class where we define the algorithm/equation/structure reasoning for having this class or where it fits with the rest of the code.

digraph { "MWTree" -> "FunctionTree" "MWTree" -> "OperatorTree" }
template<int D>
class MWTree

Base class for Multiwavelet tree structures, such as FunctionTree and OperatorTree.

The MWTree class is the base class for all tree structures needed for Multiwavelet calculations. The MWTree is a D-dimensional tree structure of MWNodes. The tree starts from a set of root nodes at a common given scale, defining the world box. The most common settings are either a single root node or \( 2^D \) root nodes. Other configurations are however allowed. For example, in 3D one could have a total of 12 root nodes (a 2x2x3 set of root nodes). Once the tree structure is generated, each node will have a parent node (except for the root nodes) and \( 2^D \) child nodes (except for leaf nodes). Most of the methods deal with traversing the tree structure in different ways to fetch specific nodes. Some of them will return a node present in the tree; some other methods will generate the required node on the fly using the MW transform; some methods will return an empty pointer if the node is not present. See specific methods for details.

Subclassed by mrcpp::FunctionTree< 1 >, mrcpp::FunctionTree< 3 >, mrcpp::FunctionTree< D >

Public Functions

MWTree(const MultiResolutionAnalysis<D> &mra, const std::string &n)

MWTree constructor.

Creates an empty tree object, containing only the set of root nodes. The information for the root node configuration to use is in the mra object which is passed to the constructor.

Parameters:
  • mra[in] the multiresolution analysis object

  • n[in] the name of the tree (only for printing purposes)

virtual ~MWTree()

MWTree destructor.

void setZero()

Set the MW coefficients to zero, keeping the same tree structure.

Keeps the node structure of the tree, even though the zero function is representable at depth zero. One should then use cropTree to remove unnecessary nodes.

void clear()

Remove all nodes in the tree.

Leaves the tree in the same state as after construction, i.e. undefined tree structure containing only root nodes without coefficients. The assigned memory, including branch and leaf nodes, (nodeChunks in NodeAllocator) is NOT released, but is immediately available to the new function.

inline double getSquareNorm() const
Returns:

Squared L2 norm of the function

void calcSquareNorm()

Calculate the squared norm \( ||f||^2_{\ldots} \) of a function represented as a tree.

The norm is calculated using endNodes only. The specific type of norm which is computed will depend on the derived class

inline int getNNodes() const
Returns:

the total number of nodes in the tree

int getNNodesAtDepth(int i) const
Parameters:

depth[in] Tree depth (0 depth is the coarsest scale) to count.

Returns:

Total number of nodes in the tree, at given depth (not in use)

int getSizeNodes() const
Returns:

Size of all MW coefs in the tree, in kB

inline NodeBox<D> &getRootBox()
Returns:

void mwTransform(int type, bool overwrite = true)

Full Multiwavelet transform of the tree in either directions.

It performs a Multiwavlet transform of the whole tree. The input parameters will specify the direction (upwards or downwards) and whether the result is added to the coefficients or it overwrites them. See the documentation for the mwTransformUp and mwTransformDown for details.

\[\begin{split} \pmatrix{ s_{nl}\\ d_{nl} } \rightleftarrows \pmatrix{ s_{n+1,2l}\\ s_{n+1,2l+1} } \end{split}\]

Parameters:
  • type[in] TopDown (from roots to leaves) or BottomUp (from leaves to roots) which specifies the direction of the MW transform

  • overwrite[in] if true, the result will overwrite preexisting coefficients.

MWNode<D> *findNode(NodeIndex<D> nIdx)

Finds and returns the node pointer with the given NodeIndex.

Recursive routine to find and return the node with a given NodeIndex. This routine returns the appropriate Node, or a NULL pointer if the node does not exist, or if it is a GenNode. Recursion starts at the appropriate rootNode.

const MWNode<D> *findNode(NodeIndex<D> nIdx) const

Finds and returns the node pointer with the given NodeIndex, const version.

Recursive routine to find and return the node with a given NodeIndex. This routine returns the appropriate Node, or a NULL pointer if the node does not exist, or if it is a GenNode. Recursion starts at the appropriate rootNode.

MWNode<D> &getNode(NodeIndex<D> nIdx)

Finds and returns the node reference with the given NodeIndex.

This routine ALWAYS returns the node you ask for. If the node does not exist, it will be generated by MW transform. Recursion starts at the appropriate rootNode and descends from this.

MWNode<D> &getNodeOrEndNode(NodeIndex<D> nIdx)

Finds and returns the node with the given NodeIndex.

This routine returns the Node you ask for, or the EndNode on the path to the requested node, if the requested one is deeper than the leaf node ancestor. It will never create or return GenNodes. Recursion starts at the appropriate rootNode and decends from this.

const MWNode<D> &getNodeOrEndNode(NodeIndex<D> nIdx) const

Finds and returns the node reference with the given NodeIndex. Const version.

This routine ALWAYS returns the node you ask for. If the node does not exist, it will be generated by MW transform. Recursion starts at the appropriate rootNode and decends from this.

MWNode<D> &getNode(Coord<D> r, int depth = -1)

Finds and returns the node at a given depth that contains a given coordinate.

This routine ALWAYS returns the node you ask for, and will generate nodes that do not exist. Recursion starts at the appropriate rootNode and decends from this.

Parameters:
  • depth[in] requested node depth from root scale.

  • r[in] coordinates of an arbitrary point in space

MWNode<D> &getNodeOrEndNode(Coord<D> r, int depth = -1)

Finds and returns the node at a given depth that contains a given coordinate.

This routine returns the Node you ask for, or the EndNode on the path to the requested node, and will never create or return GenNodes. Recursion starts at the appropriate rootNode and decends from this.

Parameters:
  • depth[in] requested node depth from root scale.

  • r[in] coordinates of an arbitrary point in space

const MWNode<D> &getNodeOrEndNode(Coord<D> r, int depth = -1) const

Finds and returns the node at a given depth that contains a given coordinate. Const version.

This routine returns the Node you ask for, or the EndNode on the path to the requested node, and will never create or return GenNodes. Recursion starts at the appropriate rootNode and decends from this.

Parameters:
  • depth[in] requested node depth from root scale.

  • r[in] coordinates of an arbitrary point in space

MWNodeVector<D> *copyEndNodeTable()

Returns the list of all EndNodes.

copies the list of all EndNode pointers into a new vector and retunrs it.

void deleteRootNodes()

Deletes all the nodes in the tree.

This method will recursively delete all the nodes, including the root nodes. Derived classes will call this method when the object is deleted.

void resetEndNodeTable()

Recreate the endNodeTable.

the endNodeTable is first deleted and then rebuilt from scratch. It makes use of the TreeIterator to traverse the tree.

int getIx(NodeIndex<D> nIdx)

gives serialIx of a node from its NodeIndex

Peter will document this!

void makeMaxSquareNorms()

sets values for maxSquareNorm in all nodes

it defines the upper bound of the squared norm \( ||f||^2_{\ldots} \) in this node or its descendents

Public Members

MWNodeVector<D> endNodeTable

Final projected nodes.

Protected Functions

virtual void mwTransformDown(bool overwrite)

Regenerates all scaling coeffs by MW transformation of existing s/w-coeffs on coarser scales.

The transformation starts at the rootNodes and proceeds recursively all the way to the leaf nodes. The existing scaling coefficeints will either be overwritten or added to. The latter operation is generally used after the operator application.

Parameters:

overwrite[in] if true the preexisting coefficients are overwritten

virtual void mwTransformUp()

Regenerates all s/d-coeffs by backtransformation.

It starts at the bottom of the tree (scaling coefficients of the leaf nodes) and it generates the scaling and wavelet coefficients if the parent node. It then proceeds recursively all the way up to the root nodes. This is generally used after a function projection to purify the coefficients obtained by quadrature at coarser scales which are therefore not precise enough.

void incrementNodeCount(int scale)

Increments node counter by one for non-GenNodes.

TO BE DOCUMENTED

Warning

: This routine is not thread safe, and must NEVER be called outside a critical region in parallel. It’s way. way too expensive to lock the tree, so don’t even think about it.

void decrementNodeCount(int scale)

Decrements node counter by one for non-GenNodes.

TO BE DOCUMENTED

Warning

: This routine is not thread safe, and must NEVER be called outside a critical region in parallel. It’s way. way too expensive to lock the tree, so don’t even think about it.

virtual std::ostream &print(std::ostream &o) const

Prints a summary of the tree structure on the output file.

Protected Attributes

NodeBox<D> rootBox

The actual container of nodes.

std::vector<int> nodesAtDepth

Node counter.

std::vector<int> nodesAtNegativeDepth

Node counter.