The InterpolationMatrix
Class¶
-
class
stat_fem.InterpolationMatrix.
InterpolationMatrix
(function_space, coords)¶ Class representing an interpolation matrix
This class holds the sparse matrix needed to interpolate from the mesh to the data. It also handles much of the parallel communication between processes to enable serial linear algebra on the data while maintaining parallel FEM solves that are distributed across processes.
In most cases, the user should not need to create
InterpolationMatrix
objects directly – these can be constructed from the information available in other objects used in thestat-fem
solves.Variables: - function_space – The FEM function space for the solution
- coords – Spatial locations where data observations are available. Must be a 2D numpy array, where the first index is the index representing the different sensor locations and the second index represents the spatial dimensions.
- n_data – Number of sensor locations
- dataspace_distrib – Distributed PETSc Vector holding data items over all MPI processes involved
- n_data_local – Number of data items held on the local process
- dataspace_gathered – PETSc Vector with all sensor data collected on the root process.
- petsc_scatter – PETSc Scatter object used to transfer data between
dataspace_gathered
anddataspace_distrib
- meshspace_vector – Firedrake Vector holding FEM mesh data
- n_mesh – Number of mesh DOFs in the Firedrake Function
- n_mesh_local – Number of local mesh DOFs on the local process
- interp – PETSc Sparse Matrix for interpolating between the data space and mesh space solutions
- is_assembled – Boolean indicating if
interp
has been assembled
-
__init__
(function_space, coords)¶ create and assemble interpolation matrix
-
assemble
()¶ Compute values and assemble interpolation matrix
Assembly function to compute the nonzero sparse matrix entries and assemble the sparse matrix. Should only need to be called once for each analysis and thus will return if
is_assembled
isTrue
without re-assembling the matrix.
-
destroy
()¶ Deallocate memory for PETSc vectors and matrix
This function deallocates the PETSc objects with a PETSc-like interface. It will call the
destroy
method on all allocated PETSc objects within theInterpolationMatrix
class.
-
get_meshspace_column_vector
(idx)¶ Returns distributed meshspace column vector for a given data point
The FEM solve requires extracting the column vectors (distributed over the full mesh) for each sensor location. This function checks that the index is valid and returns the appropriate column vector as a Firedrake Vector.
Parameters: idx (int) – Index of desired sensor location. Must be a non-negative integer less than the total number of sensors. Returns: Firedrake Vector holding the appropriate column of the interpolation matrix. Return type: Firedrake Vector
-
interp_data_to_mesh
(data_array)¶ Interpolate a vector of sensor values to the FEM mesh
This function takes a gathered numpy array in the data space that is held on the root process and interpolates it to a distributed mesh vector. The provided sensor values are first scattered to a distributed PETSc Vector and then interpolation is done by multiplying the vector by the sparse interpolation matrix.
The numpy array must be defined only on the root process (all other processes must have an array length of zero). This is checked against the expected sizes of the
dataspace_gathered
Vector on each process, an an assertion error will be triggered if these do not match.Returns a Firedrake Function holding the interpolated values on the FEM mesh.
Parameters: data_array (ndarray) – Numpy array holding the sensor values. Must be a 1D array with the full data on the root process and empty (length zero) arrays on all other processes. Returns: Data interpolated to the FEM mesh DOFs. Return type: Firedrake Vector
-
interp_mesh_to_data
(input_mesh_vector)¶ Function to interpolate from the FEM mesh to the sensor locations
This function takes a distributed mesh vector and interpolates it to a gathered numpy array on the root process (all other processes will return an empty array). Input must be a Firedrake Vector defined on the FEM mesh.
Parameters: input_mesh_vector (Firedrake Vector) – Firedrake Vector to be interpolated to the mesh locations. Returns: Numpy array holding interpolated solution at the sensor location. Root process will hold all data and return a 1D array of length (n_data,)
, while all other processes will return an array of zero length.Return type: ndarray