This page was generated from docs\source\notebooks/approximation.ipynb.

Interpolation and extrapolation#

Every Lagrangian cell comes with the necessary tools to interpolate or extrapolate multidimensional data. This can be done manually using the shape functions and shape function matrices of the cells, but there is an easier and quicker way that uses these in the background for us. Take the 8-noded hexahedron as an example.

[1]:
from sigmaepsilon.mesh.cells import H8

master_coordinates = H8.Geometry.master_coordinates()

Every cell is equipped an interpolator, which is a factory function that returns another callable, that we use for interpolation or extrapolation of known data at known locations. Note that the mechanism uses only class level information, therefore there is no need to create an instance.

[2]:
from sigmaepsilon.mesh import LagrangianCellApproximator

approximator = H8.Geometry.approximator()
isinstance(approximator, LagrangianCellApproximator)
[2]:
True

We have to feed the interpolator with the locations of the known data, the knowd data itself and the locations where we want to know the data:

[3]:
source_coordinates = master_coordinates / 2
source_values = [1, 2, 3, 4, 5, 6, 7, 8]
target_coordinates = master_coordinates

approximator(source=source_coordinates, values=source_values, target=target_coordinates)
[3]:
array([-3.5,  0.5,  0.5,  4.5,  4.5,  8.5,  8.5, 12.5])

It is possible to pass the source coordinates to the factory function. This is useful if we plan to reuse the approximator with the same source points and can save a little time. In this case only the source values and the target points need to be provided.

[4]:
source_coordinates = master_coordinates / 2
source_values = [1, 2, 3, 4, 5, 6, 7, 8]
target_coordinates = master_coordinates

approximator = H8.Geometry.approximator(source_coordinates)
approximator(values=source_values, target=target_coordinates)
[4]:
array([-3.5,  0.5,  0.5,  4.5,  4.5,  8.5,  8.5, 12.5])

As noted in the documentation, if the number of source coorindates does not match the number of nodes (and hence the number of shape functions) of the master element of the class, the approximation is gonna be under or overdetermined and the operation involves the calculation of a generalized inverse. This means, that for instance feeding a 4-noded quadrilateral with 9 source points and data values is more information than what the class is normally able to precisely handle and the resulting approximator will represent a fitting function. In that case, if you want a precise approximation, you would want to use a 9-node quadrilateral, or accept the loss of information.

[5]:
from sigmaepsilon.mesh.cells import Q4, Q9

master_coordinates = Q9.Geometry.master_coordinates()
source_coordinates = master_coordinates / 2
source_values = [i + 1 for i in range(9)]
target_coordinates = master_coordinates

approximator = Q4.Geometry.approximator()
approximator(source=source_coordinates, values=source_values, target=target_coordinates)
[5]:
array([1.66666667, 4.33333333, 4.33333333, 9.66666667, 3.        ,
       4.33333333, 7.        , 5.66666667, 5.        ])

Usiing the 9-noded quadrilateral is a better choice here and you can have an exact interpolation or extrapolation.

[6]:
approximator = Q9.Geometry.approximator()
approximator(source=source_coordinates, values=source_values, target=target_coordinates)
[6]:
array([-45., -29., -29., -37.,  -5.,  -1.,  -1.,   3.,   9.])

All is the same for one dimensional cells:

[7]:
from sigmaepsilon.mesh.cells import L3

master_coordinates = L3.Geometry.master_coordinates()
source_coordinates = master_coordinates / 2
source_values = [i + 1 for i in range(3)]
target_coordinates = master_coordinates

approximator = L3.Geometry.approximator()
approximator(source=source_coordinates, values=source_values, target=target_coordinates)
[7]:
array([0., 2., 4.])

Multidimensional data#

Set up an approximator:

[8]:
import numpy as np

approximator = H8.Geometry.approximator()
master_coordinates = H8.Geometry.master_coordinates()

source_coordinates = master_coordinates / 2
target_coordinates = master_coordinates * 2

By default, multidimensional data is expected such that the last axis goes along source points (it must have the same length than the number of source points), and this layout lieds to the best performance.

[9]:
source_values = np.random.rand(10, 2, 8)
approximator(
    source=source_coordinates,
    values=source_values,
    target=target_coordinates[:3]
).shape
[9]:
(10, 2, 3)

You can use the ‘axis’ parameter to indicate that this axis is not the last, but you would probably have to accept a loss in performance (probably alongside a warning), since the memory layout of your array is not optimal.

[10]:
source_values = np.random.rand(8, 2, 10)
approximator(
    source=source_coordinates,
    values=source_values,
    target=target_coordinates[:3],
    axis=0
).shape
[10]:
(3, 2, 10)

The workaround here is to use numpy.ascontiguousarray and reordering the input data.

Custom cells#

Getting an approximator for a custom cell goes the same way, after the cell have been properly set up.

[11]:
from sigmaepsilon.mesh.geometry import PolyCellGeometry1d

Custom1dCell: PolyCellGeometry1d = PolyCellGeometry1d.generate_class(number_of_nodes=4)

master_coordinates = Custom1dCell.master_coordinates()
source_coordinates = master_coordinates / 2
source_values = [i + 1 for i in range(Custom1dCell.number_of_nodes)]
target_coordinates = master_coordinates

approximator = Custom1dCell.approximator()
approximator(source=source_coordinates, values=source_values, target=target_coordinates)
[11]:
array([-0.5,  1.5,  3.5,  5.5])