Introduction
Intrepid is a library of interoperable tools for compatible discretizations of Partial Differential Equations (PDEs). Included with the Trilinos 10.0 release is the "expert version" of Intrepid. This version is intended primarily for application developers who want to reuse large parts of their existing code frameworks such as I/O, data structures, assembly routines, etc. while gaining access to advanced discretization capabilities provided by Intrepid. In such cases the bulk of the data is owned and managed by the user rather than by Intrepid. To avoid unnecessary and possibly expensive copying of data to and from Intrepid, the expert version of the package comprises of mostly stateless classes operating on user-owned data. Virtually all numerical data required by PDE codes can be represented as a multi-dimensional array of scalar values. For this reason, and to enhance interoprability, Intrepid classes are templated on generic multi-dimensional arrays. The Shards package provides an implementation of a multi-dimensional array that can be used for that purpose, or users can write their own multi-dimensional arrays as long as a minimal interface is supported.
Overview
Current release of Intrepid includes the following features:
- Default finite element basis functions for H(grad), H(curl), H(div) and L2 spaces of orders up to 2 on standard cell topologies in 1D, 2D and 3D
- High-order (up to 10) basis functions for H(grad), H(curl), H(div) and L2 spaces on select cell topologies
- Pullbacks (transformations) from reference coordinate frame of H(grad), H(curl), H(div) and L2 fields
- Pullbacks of gradient, curl and divergence of H(grad), H(curl), H(div) fields
- Cubature rules of orders up to 20 on most standard 1D, 2D and 3D cell topologies
- Implementation of multi-diumensional arrays and algebraic operations on them
- Examples showing solution of basic 2nd order elliptic boundary value problems (Poisson, div-curl, and curl-curl systems) using Intrepid
Quick Start
Familiarity with with the following concepts, objects, and tools is required:
The following example demonstrates, in 7 steps, the computation of finite element stiffness matrices on a set of tetrahedral cells using a piecewise linear basis and an appropriate integration rule.
Step 1: Select a cell topology
shards::CellTopology cellType = shards::getCellTopologyData< shards::Tetrahedron<> >();
int spaceDim = cellType->getDimension();
int numNodes = cellType->getNodeCount();
We additionally set the number of computational cells numCells
.
Step 2: Select integration (cubature) rule
DefaultCubatureFactory<double> cubFactory;
int cubDegree = 2;
Teuchos::RCP<Cubature<double> > myCub = cubFactory.create(cellType, cubDegree);
int numCubPoints = myCub->getNumPoints();
Step 3: Select discrete basis
Basis_HGRAD_TET_C1_FEM<double, FieldContainer<double> > tetBasis;
int numFields = tetBasis.getCardinality();
Step 4: Format multi-dimensional arrays
FieldContainer<double> cub_points(numCubPoints, spaceDim);
FieldContainer<double> cub_weights(numCubPoints);
FieldContainer<double> cell_nodes(numCells, numNodes, spaceDim);
FieldContainer<double> jacobian(numCells, numCubPoints, spaceDim, spaceDim);
FieldContainer<double> jacobian_inv(numCells, numCubPoints, spaceDim, spaceDim);
FieldContainer<double> jacobian_det(numCells, numCubPoints);
FieldContainer<double> weighted_measure(numCells, numCubPoints);
FieldContainer<double> grad_at_cub_points(numFields, numCubPoints, spaceDim);
FieldContainer<double> transformed_grad_at_cub_points(numCells, numFields, numCubPoints, spaceDim);
FieldContainer<double> weighted_transformed_grad_at_cub_points(numCells, numFields, numCubPoints, spaceDim);
FieldContainer<double> stiffness_matrices(numCells, numFields, numFields);
We assume that the array cell_nodes
is filled with nodes defining a set of computational (physical) cells.
Step 5: Evaluate differential operator applied to basis at cubature points
myCub->getCubature(cub_points, cub_weights);
tetBasis.getValues(grad_at_cub_points, cub_points, OPERATOR_GRAD);
Step 6: Apply cell tools
CellTools<double>::setJacobian(jacobian, cub_points, cell_nodes, cellType);
CellTools<double>::setJacobianInv(jacobian_inv, jacobian);
CellTools<double>::setJacobianDet(jacobian_det, jacobian);
Step 7: Apply function space tools
FunctionSpaceTools::computeCellMeasure<double>(weighted_measure,
jacobian_det,
cub_weights);
FunctionSpaceTools::HGRADtransformGRAD<double>(transformed_grad_at_cub_points,
jacobian_inv,
grad_at_cub_points);
FunctionSpaceTools::multiplyMeasure<double>(weighted_transformed_grad_at_cub_points,
weighted_measure,
transformed_grad_at_cub_points);
FunctionSpaceTools::integrate<double>(stiffness_matrices,
transformed_grad_at_cub_points,
weighted_transformed_grad_at_cub_points,
COMP_CPP);
The computed (local) stiffness matrices can now be used in the assembly of a (global) discrete differential operator, e.g. a discrete Laplacian.
Done!
Development Plans
The next release of Intrepid is expected to support Finite Difference and Finite Volume discretizations on standard and non-standard (polygon and polyhedron) cell topologies. A "user-friendly" version for rapid development of PDE codes is also under development.