CIBC:Development:virtual:Basis

From SCIRun Documentation Wiki

Contents 1 Improvements to Basis functions 1.1 Expansion of the topological elements 1.2 Access to field data (constant - linear - quadratic - cubic) 1.3 None isoparametric fields 1.4 Cubic interpolation

Improvements to Basis functions

Using the non-linear functionality which has been added to SCIRun, a few improvements are added to this branch to prepare for a better integration of this kind of functionality into SCIRun.

Expansion of the topological elements

The current topological pieces in SCIRun consist of Node, Edge, Face and Cell. To incorporate quadratic meshes this needds to be extended with ENode. The latter denotes the additional node at each edge. To generalize programming the ENode concept should be present in each mesh no matter whether it is linear or quadratic. The reason for this is that the field data can be quadratic while the mesh type is linear. Hence in the virtual branch each mesh type has an ENode. In case of a linear one it is in the center of the edge in case of quadratic mesh one can set the actual location.

To support this concept I added in:

• VMesh::ENode::iterator
• VMesh::ENode::index_type
• VMesh::ENode::size_type
• VMesh::ENode::array_type

Although the concept is often equal to the Edge, their are conceptual differences, e.g. one cannot set a position of an edge or one cannot assign a value to an edge.

Similarly I added the following concepts:

• get_enodes(array,edge)
• get_enodes(array,face)
• get_enodes(array,cell)
• set_enode(point,idx)
• get_enode(point,idx)

In the field concept I added in:

• set_evalue(val,idx)
• get_evalue(val,idx)
• set_evalues(ptr,size)
• get_evalues(ptr,size)

Note: as the edge numbering changes when new nodes are added the final numbering of the edge nodes can only be used when the mesh is finalized. If one adds an element the ordering is changed.

Access to field data (constant - linear - quadratic - cubic)

In order to access values on the enodes and the additional derivatives needed for the cubic approximations, in the virtual interface we grouped the vectors for the fielddata, the quadratic nodes, and the cubic derivatives. All of these have the same datatype for a field and hence are combined into one object called VFData.

This class has access points to every of these arrays and hence only needs a similar structure like the ElemData to do the interpolations.

The VFData currently has support for:

• constant data
• linear data
• quadratic data

Functions:

data access functions:

• get_value(val&, idx)
• set_value(val, idx)
• get_evalue(val&, idx) quadratic interface
• set_evalue(val, idx) quadratic interface

data interpolation functions:

• interpolate(elem_interpolator,value&)
• gradient(elem_interpolator,value&)

data copying:

• copy_value(src_field,idx_src,idx_dst)
• copy_evalue(src_field,idx_src,idx_dst)

get size of arrays:

• num_values()
• num_evalues()

resize everything consistently:

• resize_fdata() -> this one will now resize all arrays

get all values at once:

• get_values(ptr,size)
• set_values(ptr,size)
• get_evalues(ptr,size)
• set_evalues(ptr,size)

I am still looking for an elegant way to describe the derivatives:

Possibly

• get_dvalue(val&, idx, didx) // dixd is the number of the derivative
• set_dvalue(val, idx, didx)

etc..

None isoparametric fields

SCIRun allows the field data to have a different order from the mesh element description. This leads to the problem of having two basis classes one for the data and one for the mesh. The trouble starts by the notion that the interpolation of data in fact depends on both bases, especially with derivatives one needs to consider both.

Hence it is advantegeous for the mesh to know about all the different basis orders even though it may be linear.

Therefor I reorganized the basis classes:

Each basis will have the following eight functions:

• get_constant_weights(coords,elem,weights)
• get_linear_weights(coords,elem,weights)
• get_quadratic_weights(coords,elem,weights)
• get_cubic_weights(coords,elem,weights)

• get_derivate_constant_weights(coords,elem,weights)
• get_derivate_linear_weights(coords,elem,weights)
• get_derivate_quadratic_weights(coords,elem,weights)
• get_derivate_cubic_weights(coords,elem,weights)

In combination in VFData having all the data stored, one can combine both to generate an interpolation function. In the mesh get the weights and indices for the element interpolation of the order of the field data and hen use the data container to get the final value.

This way interpolation may not be as fast as the elem data path but it can be made completely virtual and separates the actual interpolation code from the datatype: allowing to add more datatypes without too much trouble.

Each get_weights() and get_derivate_weights() in each of the basis classes will be based on the ones above. The latter are concentrated in a separate class.

Problem: currently the basis class depends on the order, whereas we should have combined all the different orders into one basis class that only depends on the element type. The Basis class should have been something called Element. The current implementation allows for a cleaner interface and is still compatible with the old Basis classes.

Cubic interpolation

The cubic interpolation has not fully been implemented as one piece of the functionality is missing:

The cubic version makes use of derivatives, however as the interpolation is done in local coordinates there is a problem: the space from a mesh viewpoint is only C0, meaning in a local coordinate space the value on one edge of the element is equal to the value at the other side, but not the derivative. As the derivatives depend on the size and the deformation of the element.

The current implementation assigns scaling factors. These however need to be computed. A more elegant way is computing the scaling factors based on the mesh. Hence the interpolation depends on values (data), the basis function interpolation functions, and the scaling factors (from mesh). Hence to do it properly one should add code to the mesh to compute the factors. This should be done automatically and should not depend on external assignment.

This has currently not yet been implemented, but we should think about this as it would complete the basis support and will render a good solution for the scaling factors (these classes are ugly and only a few people know how that stuff actually works)