Classes | Typedefs | Enumerations | Functions | Variables

apf Namespace Reference

All APF symbols are contained in this namespace. More...

Classes

struct  Up
 statically sized container for upward adjacency queries. More...
struct  Copy
 a reference to an object representing the same entity More...
class  Mesh
 Interface to a mesh part. More...
class  Migration
 Migration plan object: local elements to destinations. More...
struct  Sharing
 abstract description of entity copy sharing More...
class  Mesh2
 Extended mesh interface for modification. More...
class  BuildCallback
 User-defined entity creation callback. More...
class  Integrator
 A virtual base for user-defined integrators. More...
struct  Function
 User-defined Analytic Function. More...
struct  Node
 Node identifier. More...
class  EntityShape
 Shape functions over this element. More...
class  FieldShape
 Describes field distribution and shape functions. More...
class  Vector
 template-generic vector of N doubles More...
class  Vector3
 convenience wrapper over apf::Vector<3> More...
class  Matrix
 template-generic matrix of M by N doubles More...
class  Matrix3x3
 convenience wrapper over apf::Matrix<3,3> More...
class  DynamicVector
 A runtime-sized linear algebra vector of doubles. More...
class  DynamicMatrix
 A runtime-sized dense matrix. More...
class  CavityOp
 user-defined mesh cavity operator More...
class  Splitter
 Splits a mesh part into many. More...
class  Balancer
 Load balance over all mesh parts. More...
struct  Remap
 a map from old part ids to new part ids More...
struct  Divide
 divide the part id More...
struct  Multiply
 multiply the part id More...
struct  Modulo
 return part id modulo n More...
struct  Unmodulo
 inverse of apf::Modulo More...
struct  Round
 map to nearest multiple of n More...

Typedefs

typedef NumberingOf< int > Numbering
 Numbering is meant to be a 32-bit local numbering.
typedef NumberingOf< long > GlobalNumbering
 Global numberings use 64-bit integers.
typedef std::map< int,
MeshEntity * > 
Copies
 Remote copy container.
typedef std::set< int > Parts
 Set of unique part ids.
typedef DynamicArray
< MeshEntity * > 
Adjacent
 Set of adjacent mesh entities.
typedef MeshEntity * Downward [12]
 a static array type downward adjacency queries.
typedef Copy Match
 matches are just a special case of copies
typedef DynamicArray< CopyCopyArray
 a set of copies, possibly multiple copies per part
typedef CopyArray Matches
 a set of matched copies
typedef VectorElement MeshElement
 Mesh Elements represent the mesh coordinate vector field.
typedef std::map< int,
MeshEntity * > 
GlobalToVert
 a map from global ids to vertex objects

Enumerations

enum  ValueType {
  SCALAR, VECTOR, MATRIX, PACKED,
  VALUE_TYPES
}
 

The type of value the field stores.

More...
enum  ZoltanMethod {
  RCB, RIB, HYPERGRAPH, PARMETIS,
  GRAPH
}
 

Zoltan partitioning method.

More...
enum  ZoltanApproach {
  PARTITION, REPARTITION, REFINE, PART_KWAY,
  PART_GEOM, PART_GEOM_KWAY, ADAPT_REPART, REFINE_KWAY
}
 

Zoltan partitioning approach.

More...

Functions

void verify (Mesh *m)
 run consistency checks on an apf::Mesh structure
int getDimension (Mesh *m, MeshEntity *e)
 get the dimension of a mesh entity
void unite (Parts &into, Parts const &from)
 unite two sets of unique part ids
void removeTagFromDimension (Mesh *m, MeshTag *tag, int d)
 removes a tag from all entities of dimension (d)
MeshEntity * findUpward (Mesh *m, int type, MeshEntity **down)
 find an entity from one-level downward adjacencies
MeshEntity * findElement (Mesh *m, int type, MeshEntity **verts)
 finds an entity from a set of vertices
MeshEntity * getEdgeVertOppositeVert (Mesh *m, MeshEntity *edge, MeshEntity *v)
 get the other vertex of an edge
void getBridgeAdjacent (Mesh *m, MeshEntity *origin, int bridgeDimension, int targetDimension, Adjacent &result)
 get 2nd-order adjacent entities
int countEntitiesOfType (Mesh *m, int type)
 count all on-part entities of one topological type
bool isSimplex (int type)
 return true if the topological type is a simplex
Vector3 getLinearCentroid (Mesh *m, MeshEntity *e)
 get the average of the entity's vertex coordinates
SharinggetSharing (Mesh *m)
 create a default sharing object for this mesh
void getPeers (Mesh *m, int d, Parts &peers)
 scan the part for [vtx|edge|face]-adjacent part ids
int findIn (MeshEntity **a, int n, MeshEntity *e)
 find pointer (e) in array (a) of length (n)
void findTriDown (Mesh *m, MeshEntity **verts, MeshEntity **down)
 given the vertices of a triangle, find its edges
void changeMeshShape (Mesh *m, FieldShape *newShape, bool project=true)
 deprecated wrapper for apf::Mesh::changeShape
void unfreezeFields (Mesh *m)
 unfreeze all associated fields
int countEntitiesOn (Mesh *m, ModelEntity *me, int dim)
 count the number of mesh entities classified on a model entity
int countOwned (Mesh *m, int dim)
 count the number of owned entities of dimension (dim)
void printStats (Mesh *m)
 print global mesh entity counts per dimension
void warnAboutEmptyParts (Mesh *m)
 print to stderr the number of empty parts, if any
Copy getOtherCopy (Mesh *m, MeshEntity *s)
 given a mesh face, return its remote copy
int getFirstType (Mesh *m, int dim)
 get the type of the first entity in this dimension
void getAlignment (Mesh *m, MeshEntity *elem, MeshEntity *boundary, int &which, bool &flip, int &rotate)
 boundary entity alignment to an element
void migrate (Mesh2 *m, Migration *plan)
 APF's migration function, works on apf::Mesh2.
void setMigrationLimit (size_t maxElements)
 set the maximum elements that apf::migrate moves at once
void displaceMesh (Mesh2 *m, Field *d, double factor=1.0)
 add a field (times a factor) to the mesh coordinates
MeshEntity * makeOrFind (Mesh2 *m, ModelEntity *c, int type, MeshEntity **down, BuildCallback *cb=0)
 like apf::Mesh2::createEntity, but returns already existing entities
MeshEntity * buildElement (Mesh2 *m, ModelEntity *c, int type, MeshEntity **verts, BuildCallback *cb=0)
 build an entity from its vertices
MeshEntity * buildOneElement (Mesh2 *m, ModelEntity *c, int type, Vector3 const *points)
 build a one-element mesh
void initResidence (Mesh2 *m, int dim)
 Set entity residence based on remote copies.
void stitchMesh (Mesh2 *m)
 infer all remote copies from those of vertices
void destroyMesh (Mesh *m)
 Destroys an apf::Mesh.
MeshElementcreateMeshElement (Mesh *m, MeshEntity *e)
 Creates a Mesh Element over an entity.
MeshEntity * getMeshEntity (MeshElement *me)
 Retrieve the mesh entity associated with an apf::MeshElement.
void destroyMeshElement (MeshElement *e)
 Destroys a Mesh Element.
Field * createLagrangeField (Mesh *m, const char *name, int valueType, int order)
 Create an apf::Field using a Lagrange distribution.
Field * createStepField (Mesh *m, const char *name, int valueType)
 Create an apf::Field using a step distribution.
Field * createIPField (Mesh *m, const char *name, int valueType, int order)
 Create an apf::Field of integration point data.
Field * createField (Mesh *m, const char *name, int valueType, FieldShape *shape)
 Create a Field from any builtin or user defined FieldShape.
Field * createFieldOn (Mesh *m, const char *name, int valueType)
 Create a field using the mesh's coordinate nodal distribution.
Field * createPackedField (Mesh *m, const char *name, int components, apf::FieldShape *shape=0)
 Create a field of N components without a tensor type.
Field * createGeneralField (Mesh *m, const char *name, int valueType, int components, FieldShape *shape)
 Encompasses both packed and typed fields.
Field * cloneField (Field *f, Mesh *onto)
 Declare a copy of a field on another apf::Mesh.
MeshgetMesh (Field *f)
 Retrieve the Mesh over which a Field is defined.
bool hasEntity (Field *f, MeshEntity *e)
 Returns true iff an entity has data associate with a field.
const char * getName (Field *f)
 Get the name of a Field.
int getValueType (Field *f)
 Retrieve the type of value a field distributes.
void destroyField (Field *f)
 Destroy an apf::Field.
void setScalar (Field *f, MeshEntity *e, int node, double value)
 Set a nodal value of a scalar field.
double getScalar (Field *f, MeshEntity *e, int node)
 Get the node value of a scalar field.
void setVector (Field *f, MeshEntity *e, int node, Vector3 const &value)
 Set the nodal value of a vector field.
void getVector (Field *f, MeshEntity *e, int node, Vector3 &value)
 Get the nodal value of a vector field.
void setMatrix (Field *f, MeshEntity *e, int node, Matrix3x3 const &value)
 Set the nodal value of a matrix field.
void getMatrix (Field *f, MeshEntity *e, int node, Matrix3x3 &value)
 Get the nodal value of a matrix field.
void setComponents (Field *f, MeshEntity *e, int node, double const *components)
 Set the nodal value from an array of component values.
void getComponents (Field *f, MeshEntity *e, int node, double *components)
 Copy the nodal value into an array of component values.
Element * createElement (Field *f, MeshElement *e)
 Create a Field Element from a Mesh Element.
Element * createElement (Field *f, MeshEntity *e)
 Create a Field Element without a parent Mesh Element.
void destroyElement (Element *e)
 Destroy a Field Element.
MeshElementgetMeshElement (Element *e)
 Get the Mesh Element of a Field Element.
MeshEntity * getMeshEntity (Element *e)
 Retrieve the mesh entity of an apf::Element.
double getScalar (Element *e, Vector3 const &param)
 Evaluate a scalar field at a point.
void getGrad (Element *e, Vector3 const &param, Vector3 &grad)
 Get the gradient of a scalar field w.r.t. global coordinates.
void getVector (Element *e, Vector3 const &param, Vector3 &value)
 Evaluate a vector field at a point.
double getDiv (Element *e, Vector3 const &param)
 Evaluate the divergence of a vector field at a point.
void getCurl (Element *e, Vector3 const &param, Vector3 &curl)
 Evaluate the curl of a vector field at a point.
void getVectorGrad (Element *e, Vector3 const &param, Matrix3x3 &deriv)
 Get the gradient of a vector field w.r.t. global coordinates.
void getMatrix (Element *e, Vector3 const &param, Matrix3x3 &value)
 Evaluate the value of a matrix field.
void getComponents (Element *e, Vector3 const &param, double *components)
 Evaluate a field into an array of component values.
int countIntPoints (MeshElement *e, int order)
 Get the number of integration points for an element.
void getIntPoint (MeshElement *e, int order, int point, Vector3 &param)
 Get an integration point in an element.
double getIntWeight (MeshElement *e, int order, int point)
 Get the weight of an integration point in an element.
void mapLocalToGlobal (MeshElement *e, Vector3 const &local, Vector3 &global)
 Map a local coordinate to a global coordinate.
double getDV (MeshElement *e, Vector3 const &param)
 Get the differential volume at a point.
double measure (MeshElement *e)
 Measures the volume, area, or length of a Mesh Element.
double measure (Mesh *m, MeshEntity *e)
 Measures the volume, area, or length of a Mesh Entity.
int getOrder (MeshElement *e)
 Returns the polynomial order of the coordinate field.
void getJacobian (MeshElement *e, Vector3 const &local, Matrix3x3 &j)
 Returns the Jacobian at a local point.
void getJacobianInv (MeshElement *e, Vector3 const &local, Matrix3x3 &jinv)
 Returns the Jacobian inverse at a local point.
int countNodes (Element *e)
 Returns the number of element nodes.
void getScalarNodes (Element *e, NewArray< double > &values)
 Returns the element nodal values for a scalar field.
void getVectorNodes (Element *e, NewArray< Vector3 > &values)
 Returns the element nodal values for a vector field.
void getMatrixNodes (Element *e, NewArray< Matrix3x3 > &values)
 Returns the element nodal values for a matrix field.
void getShapeValues (Element *e, Vector3 const &local, NewArray< double > &values)
 Returns the shape function values at a point.
void getShapeGrads (Element *e, Vector3 const &local, NewArray< Vector3 > &grads)
 Returns the shape function gradients at a point.
FieldShapegetShape (Field *f)
 Retrieve the apf::FieldShape used by a field.
int countComponents (Field *f)
 Count the number of scalar components in the field's value type.
void writeVtkFiles (const char *prefix, Mesh *m)
 Write a set of parallel VTK Unstructured Mesh files from an apf::Mesh with binary (base64) encoding and zlib compression (if LION_COMPRESS=ON).
void writeVtkFiles (const char *prefix, Mesh *m, std::vector< std::string > writeFields)
 Write a set of parallel VTK Unstructured Mesh files from an apf::Mesh with binary (base64) encoding and zlib compression (if LION_COMPRESS=ON).
void writeOneVtkFile (const char *prefix, Mesh *m)
 Output just the .vtu file with ASCII encoding for this part.
void writeASCIIVtkFiles (const char *prefix, Mesh *m)
 Write a set of parallel VTK Unstructured Mesh files from an apf::Mesh with ASCII encoding.
void writeASCIIVtkFiles (const char *prefix, Mesh *m, std::vector< std::string > writeFields)
 Write a set of parallel VTK Unstructured Mesh files from an apf::Mesh with ASCII encoding.
void getGaussPoint (int type, int order, int point, Vector3 &param)
 Return the location of a gaussian integration point.
int countGaussPoints (int type, int order)
 Return the number of Gaussian integration points.
double getJacobianDeterminant (Matrix3x3 const &J, int dimension)
 Return the Jacobian determinant or dimensional equivalent.
int getDimension (MeshElement *me)
 Return the dimension of a MeshElement's MeshEntity.
void synchronize (Field *f, Sharing *shr=0)
 Synchronize field values along partition boundary.
void accumulate (Field *f, Sharing *shr=0)
 Add field values along partition boundary.
void fail (const char *why) __attribute__((noreturn))
 Declare failure of code inside APF.
void freeze (Field *f)
 Convert a Field from Tag to array storage.
void unfreeze (Field *f)
 Convert a Field from array to Tag storage.
bool isFrozen (Field *f)
 Returns true iff the Field uses array storage.
double * getArrayData (Field *f)
 Return the contiguous array storing this field.
void zeroField (Field *f)
 Initialize all nodal values with all-zero components.
Field * createUserField (Mesh *m, const char *name, int valueType, FieldShape *s, Function *f)
 Create a Field from a user's analytic function.
Field * recoverGradientByVolume (Field *f)
 Compute a nodal gradient field from a nodal input field.
void projectField (Field *to, Field *from)
 Project a field from an existing field.
void getBF (FieldShape *s, MeshElement *e, Vector3 const &p, NewArray< double > &BF)
 Get the basis functions over a mesh element.
void getGradBF (FieldShape *s, MeshElement *e, Vector3 const &p, NewArray< Vector3 > &gradBF)
 Get global gradients of basis functions over a mesh element.
NumberingcreateNumbering (Field *f)
 Create a Numbering of degrees of freedom of a Field.
NumberingcreateNumbering (Mesh *mesh, const char *name, FieldShape *shape, int components)
 Create a generally-defined Numbering.
void destroyNumbering (Numbering *n)
 Destroy a Numbering.
void fix (Numbering *n, MeshEntity *e, int node, int component, bool fixed)
 Set the fixed/free status of a degree of freedom,.
bool isFixed (Numbering *n, MeshEntity *e, int node, int component)
 Check whether a degree of freedom is fixed.
bool isNumbered (Numbering *n, MeshEntity *e, int node, int component)
 Check whether a degree of freedom is numbered.
void number (Numbering *n, MeshEntity *e, int node, int component, int number)
 number a degree of freedom
int getNumber (Numbering *n, MeshEntity *e, int node, int component)
 get a degree of freedom number
Field * getField (Numbering *n)
 get the field being numbered
FieldShapegetShape (Numbering *n)
 get the FieldShape used by a Numbering
const char * getName (Numbering *n)
 get the name of a Numbering
MeshgetMesh (Numbering *n)
 get the mesh associated with a Numbering
int getElementNumbers (Numbering *n, MeshEntity *e, NewArray< int > &numbers)
 returns the node numbers of an element
int countFixed (Numbering *n)
 return the number of fixed degrees of freedom
void synchronize (Numbering *n, Sharing *shr=0)
 numbers non-owned nodes with the values from their owners
NumberingnumberOwnedDimension (Mesh *mesh, const char *name, int dim)
 number the local owned entities of a given dimension
NumberingnumberOverlapDimension (Mesh *mesh, const char *name, int dim)
 number all local entities of a given dimension
NumberingnumberElements (Mesh *mesh, const char *name)
 number the local elements
NumberingnumberOverlapNodes (Mesh *mesh, const char *name, FieldShape *s=0)
 number all local nodes
NumberingnumberOwnedNodes (Mesh *mesh, const char *name, FieldShape *s=0, Sharing *shr=0)
 number the local owned nodes
int countNodes (Numbering *n)
 count the number of nodes that have been numbered
void getNodes (Numbering *n, DynamicArray< Node > &nodes)
 get an array of numbered nodes
void getNodesOnClosure (Mesh *m, ModelEntity *me, DynamicArray< Node > &on, FieldShape *sh=0)
 get nodes on the closure of a model entity
GlobalNumberingcreateGlobalNumbering (Mesh *mesh, const char *name, FieldShape *shape, int components=1)
 create global numbering
MeshgetMesh (GlobalNumbering *n)
 get the mesh associated with a global numbering
int countComponents (GlobalNumbering *n)
 get the components associated with a global numbering
void number (GlobalNumbering *n, Node node, long number, int component=0)
 assign a global number
long getNumber (GlobalNumbering *n, Node node, int component=0)
 get a global number
long getNumber (GlobalNumbering *n, MeshEntity *e, int node, int component=0)
 get a global number
int getElementNumbers (GlobalNumbering *n, MeshEntity *e, NewArray< long > &numbers)
 get an element's global node numbers
GlobalNumberingmakeGlobal (Numbering *n, bool destroy=true)
 converts a local numbering into a global numbering.
void synchronize (GlobalNumbering *n, Sharing *shr=0)
 see the Numbering equivalent and apf::makeGlobal
void destroyGlobalNumbering (GlobalNumbering *n)
 destroy a global numbering
void getNodes (GlobalNumbering *n, DynamicArray< Node > &nodes)
 see the Numbering equivalent
MeshTag * reorder (Mesh *mesh, const char *name)
 Number by adjacency graph traversal.
int NaiveOrder (Numbering *num)
 number all components by simple iteration
int AdjReorder (Numbering *num)
 like apf::reorder, but numbers all free nodal components
void SetNumberingOffset (Numbering *num, int off)
 add an offset to all free nodal component numbers
FieldShapegetLagrange (int order)
 Get the Lagrangian shape function of some polynomial order.
FieldShapegetSerendipity ()
 Get the Serendipity shape functions of second order.
FieldShapegetConstant (int dimension)
 Get the Constant shape function over some dimension.
FieldShapegetIPShape (int dimension, int order)
 Get the Integration Point distribution.
FieldShapegetVoronoiShape (int dimension, int order)
 Get the Voronoi shape function.
FieldShapegetIPFitShape (int dimension, int order)
 Get the IP Fit shape function.
FieldShapegetHierarchic (int order)
 Get the quadratic hierarchic shape function.
void projectHierarchicField (Field *to, Field *from)
 Project a hierarchic field.
int countElementNodes (FieldShape *s, int type)
 count the number of nodes affecting an element type
Vector3 boundaryToElementXi (Mesh *m, MeshEntity *boundary, MeshEntity *element, Vector3 const &xi)
 Reparameterize from boundary entity to element.
void convert (Mesh *in, Mesh2 *out)
 convert one mesh data structure to another
void construct (Mesh2 *m, const int *conn, int nelem, int etype, GlobalToVert &globalToVert)
 construct a mesh from just a connectivity array
void setCoords (Mesh2 *m, const double *coords, int nverts, GlobalToVert &globalToVert)
 Assign coordinates to the mesh.
void destruct (Mesh2 *m, int *&conn, int &nelem, int &etype)
 convert an apf::Mesh2 object into a connectivity array
void extractCoords (Mesh2 *m, double *&coords, int &nverts)
 get a contiguous set of global vertex coordinates
Vector< 3 > cross (Vector< 3 > const &a, Vector< 3 > const &b)
 3D vector cross product
template<std::size_t N>
Vector< N > project (Vector< N > const &a, Vector< N > const &b)
 Returns vector (a) projected onto vector (b).
template<std::size_t N>
Vector< N > reject (Vector< N > const &a, Vector< N > const &b)
 vector rejection
template<std::size_t M, std::size_t N>
Matrix< N, M > transpose (Matrix< M, N > const &m)
 transpose a matrix
template<std::size_t M, std::size_t N>
Matrix< M, N > tensorProduct (Vector< M > const &a, Vector< N > const &b)
 tensor product of two vectors
template<std::size_t M, std::size_t N>
Matrix< M-1, N-1 > getMinor (Matrix< M, N > const &A, std::size_t i, std::size_t j)
 get the minor matrix associated with entry (i,j) of matrix A
template<std::size_t M, std::size_t N>
double getCofactor (Matrix< M, N > const &A, std::size_t i, std::size_t j)
 get the cofactor associated with entry (i,j) of matrix A
template<std::size_t M, std::size_t N>
double getDeterminant (Matrix< M, N > const &A)
 get the determinant of a matrix A
Matrix< 3, 3 > cofactor (Matrix< 3, 3 > const &m)
 get the matrix of cofactors for a given matrix
Matrix< 2, 2 > invert (Matrix< 2, 2 > const &m)
 invert a 2 by 2 matrix
Matrix< 3, 3 > invert (Matrix< 3, 3 > const &m)
 invert a 3 by 3 matrix
template<std::size_t M, std::size_t N>
double getInnerProduct (Matrix< M, N > const &a, Matrix< M, N > const &b)
 get the component-wise inner product of two matrices
Matrix3x3 cross (Vector3 const &u)
 get the skew-symmetric cross product matrix of a vector
Matrix3x3 rotate (Vector3 const &u, double a)
 get the rotation matrix around an axis
Matrix3x3 getFrame (Vector3 const &v)
 derive an orthogonal frame whose x axis is the given vector
int eigen (Matrix3x3 const &A, Vector< 3 > *eigenVectors, double *eigenValues)
 get the eigenvectors and eigenvalues of a 3 by 3 matrix
template<std::size_t N>
DynamicVector fromVector (Vector< N > other)
 convert an apf::Vector into an apf::DynamicVector
void multiply (DynamicMatrix const &a, DynamicVector const &b, DynamicVector &r)
 multiply a DynamicMatrix by a DynamicVector
void multiply (DynamicVector const &b, DynamicMatrix const &a, DynamicVector &r)
 multiply a DynamicVector by a DynamicMatrix
void multiply (DynamicMatrix const &a, DynamicMatrix const &b, DynamicMatrix &r)
 multiply two DynamicMatrix objects
void transpose (DynamicMatrix const &a, DynamicMatrix &r)
 get the transpose of a DynamicMatrix
template<std::size_t N, std::size_t M>
DynamicMatrix fromMatrix (Matrix< N, M > other)
 convert an apf::Matrix into an apf::DynamicMatrix
void remapPartition (apf::Mesh2 *m, Remap &remap)
 remap all part ids in the mesh structure
Mesh2makeEmptyMdsMesh (gmi_model *model, int dim, bool isMatched)
 create an empty MDS part
Mesh2loadMdsMesh (gmi_model *model, const char *meshfile)
 load an MDS mesh from files
Mesh2loadMdsMesh (const char *modelfile, const char *meshfile)
 load an MDS mesh and model from file
Mesh2createMdsMesh (gmi_model *model, Mesh *from)
 create an MDS mesh from an existing mesh
void reorderMdsMesh (Mesh2 *mesh, MeshTag *t=0)
 apply adjacency-based reordering
bool alignMdsMatches (Mesh2 *in)
 align the downward adjacencies of matched entities
bool alignMdsRemotes (Mesh2 *in)
 align the downward adjacencies of remote copies
void deriveMdsModel (Mesh2 *in)
 build a null model such that apf::verify accepts the mesh.
void changeMdsDimension (Mesh2 *in, int d)
 change the dimension of an MDS mesh
int getMdsIndex (Mesh2 *in, MeshEntity *e)
 returns the dimension-unique index for this entity
MeshEntity * getMdsEntity (Mesh2 *in, int dimension, int index)
 retrieve an entity by dimension and index
Mesh2loadMdsFromANSYS (const char *nodefile, const char *elemfile)
 load an MDS mesh from ANSYS .node and .elem files
SplittermakeZoltanSplitter (Mesh *mesh, int method, int approach, bool debug=true, bool sync=true)
 Make a Zoltan Splitter object.
SplittermakeZoltanGlobalSplitter (Mesh *mesh, int method, int approach, bool debug=true)
 Make a Zoltan Splitter object.
BalancermakeZoltanBalancer (Mesh *mesh, int method, int approach, bool debug=true)
 Make a Zoltan Balancer object.
MeshTag * tagOpposites (GlobalNumbering *gn, const char *name)
 Tag global ids of opposite elements to boundary faces.
int * getElementToElement (apf::Mesh *m)
 Get an element-to-element connectivity array.

Variables

int const tri_edge_verts [3][2]
 map from triangle edge order to triangle vertex order
int const quad_edge_verts [4][2]
 map from quad edge order to quad vertex order
int const tet_edge_verts [6][2]
 map from tet edge order to tet vertex order
int const prism_edge_verts [9][2]
 map from prism edge order to prism vertex order
int const pyramid_edge_verts [8][2]
 map from pyramid edge order to pyramid vertex order
int const tet_tri_verts [4][3]
 map from tet triangle order to tet vertex order
int const hex_quad_verts [6][4]
 map from hex quad order to hex vertex order
int const prism_tri_verts [2][3]
 map from prism triangle order to prism vertex order
int const prism_quad_verts [3][4]
 map from prism quad order to prism vertex order
int const pyramid_tri_verts [4][3]
 map from pyramid triangle order to pyramid vertex order
double const pi
 The mathematical constant pi.

Detailed Description

All APF symbols are contained in this namespace.

Wrapping a namespace over everything gives reasonable insurance against future symbol conflicts with other packages while very common names are being used, like Mesh or VECTOR. Users will be able to use the simple names directly with a using statement, likewise all APF code is written without apf::


Typedef Documentation

typedef DynamicArray<MeshEntity*> apf::Adjacent

Set of adjacent mesh entities.

see also apf::Downward and apf::Up

typedef std::map<int,MeshEntity*> apf::Copies

Remote copy container.

the key is the part id, the value is the on-part pointer to the remote copy

typedef MeshEntity* apf::Downward[12]

a static array type downward adjacency queries.

using statically sized arrays saves time by avoiding dynamic allocation, and downward adjacencies have a guaranteed bound.


Enumeration Type Documentation

The type of value the field stores.

The near future may bring more complex tensors.

Enumerator:
SCALAR 

a single scalar value.

VECTOR 

a 3D vector value

MATRIX 

a 3x3 matrix

PACKED 

a user-defined set of components

VALUE_TYPES 

placeholder used to set array sizes

Zoltan partitioning approach.

Enumerator:
PARTITION 

(Hyper)Graph - does not consider the initial distribution

REPARTITION 

(Hyper)Graph - considers the initial distribution

REFINE 

(HYPER)Graph - targets partitions needing only small changes

PART_KWAY 

Graph - multilevel.

PART_GEOM 

Graph - space filling curves.

PART_GEOM_KWAY 

Graph - hybrid method combining PART_KWAY and PART_GEOM.

ADAPT_REPART 

Graph - targets graphs generated from adaptively refined meshes.

REFINE_KWAY 

Graph - targets partitions needing only small changes.

Zoltan partitioning method.

Enumerator:
RCB 

Recursive Coordinate Bisection.

RIB 

Recursive Inertial Bisection.

HYPERGRAPH 

Hyper-graph partitioning.

PARMETIS 

Use ParMetis.

GRAPH 

General graph partitionig.


Function Documentation

void apf::accumulate ( Field *  f,
Sharing *  shr = 0 
)

Add field values along partition boundary.

Using the copies described by an apf::Sharing object, add up the field values of all copies of an entity and assign the sum as the value for all copies.

int apf::AdjReorder ( Numbering *  num  ) 

like apf::reorder, but numbers all free nodal components

Todo:
name should be lower-case
Vector3 apf::boundaryToElementXi ( Mesh *  m,
MeshEntity *  boundary,
MeshEntity *  element,
Vector3 const &  xi 
)

Reparameterize from boundary entity to element.

This function converts a point in the local parametric space of a boundary mesh entity into the equivalent point in the local parametric space of an adjacent element.

MeshEntity* apf::buildElement ( Mesh2 *  m,
ModelEntity *  c,
int  type,
MeshEntity **  verts,
BuildCallback *  cb = 0 
)

build an entity from its vertices

any missing intermediate entities will also be built, and all will be classified to (c). If a non-zero BuildCallback pointer is given, it will be called for each entity created, including intermediate ones.

MeshEntity* apf::buildOneElement ( Mesh2 *  m,
ModelEntity *  c,
int  type,
Vector3 const *  points 
)

build a one-element mesh

this is mostly useful for debugging

Todo:
this doesn't get used much, maybe remove it
void apf::changeMdsDimension ( Mesh2 *  in,
int  d 
)

change the dimension of an MDS mesh

this should be called before adding entities of dimension higher than the previous mesh dimension (when building a higher dimensional mesh from a lower one), or after removing all entities of higher dimension (when reducing a high dimensional mesh to a lower one)

Field* apf::cloneField ( Field *  f,
Mesh *  onto 
)

Declare a copy of a field on another apf::Mesh.

This will just make a Field object with the same properties, but not fill in any data.

void apf::construct ( Mesh2 *  m,
const int *  conn,
int  nelem,
int  etype,
GlobalToVert &  globalToVert 
)

construct a mesh from just a connectivity array

this function is here to interface with very simple mesh formats. Given a set of elements described only in terms of the ordered global ids of their vertices, this function builds a reasonable apf::Mesh2 structure and as a side effect returns a map from global ids to local vertices.

This is a fully scalable parallel mesh construction algorithm, no processor incurs memory or runtime costs proportional to the global mesh size.

Note that all vertices will have zero coordinates, so it is often good to use apf::setCoords after this.

void apf::convert ( Mesh *  in,
Mesh2 *  out 
)

convert one mesh data structure to another

this function will fill in a structure that fully implements apf::Mesh2 by using information from an implementation of apf::Mesh. This is a fully scalable parallel mesh conversion tool.

int apf::countElementNodes ( FieldShape *  s,
int  type 
)

count the number of nodes affecting an element type

Parameters:
type select from apf::Mesh::Type
int apf::countIntPoints ( MeshElement *  e,
int  order 
)

Get the number of integration points for an element.

Parameters:
order the polynomial order of accuracy desired for the integration (not to be confused with the polynomial order of shape functions).
int apf::countNodes ( Element *  e  ) 

Returns the number of element nodes.

This is the number of nodes affecting an element, as opposed to the nodes tagged to an entity.

Element* apf::createElement ( Field *  f,
MeshElement *  e 
)

Create a Field Element from a Mesh Element.

A Field Element object caches elemental data for use in evaluation, mapping, and integration. use destroyElement to free this data.

Parameters:
f The field which the Element will represent
e An existing MeshElement for the desired entity
Returns:
The new field Element
Element* apf::createElement ( Field *  f,
MeshEntity *  e 
)

Create a Field Element without a parent Mesh Element.

Warning: most users should call the version which takes a MeshElement as input. Only call this function if you know the other one isn't right for you.

GlobalNumbering* apf::createGlobalNumbering ( Mesh *  mesh,
const char *  name,
FieldShape *  shape,
int  components = 1 
)

create global numbering

see apf::createNumbering. so far global numberings have one component

Field* apf::createIPField ( Mesh *  m,
const char *  name,
int  valueType,
int  order 
)

Create an apf::Field of integration point data.

Parameters:
m the mesh over which the field is defined
name a unique name for this field
valueType the type of field data
order polynomial order of accuracy
Field* apf::createLagrangeField ( Mesh *  m,
const char *  name,
int  valueType,
int  order 
)

Create an apf::Field using a Lagrange distribution.

Parameters:
m the mesh over which the field is defined
name a unique name for this field
valueType the type of field data
order the polynomial order of the shape functions (so far 1 or 2)
Mesh2* apf::createMdsMesh ( gmi_model model,
Mesh *  from 
)

create an MDS mesh from an existing mesh

Parameters:
from the mesh to copy

this function uses apf::convert to copy any apf::Mesh

MeshElement* apf::createMeshElement ( Mesh *  m,
MeshEntity *  e 
)

Creates a Mesh Element over an entity.

A Mesh Element allows queries to the coordinate field, including mapping, differential and total volume, as well as gauss integration point data. A Mesh Element is also required to build a Field Element.

Numbering* apf::createNumbering ( Mesh *  mesh,
const char *  name,
FieldShape *  shape,
int  components 
)

Create a generally-defined Numbering.

This numbering will be available via mesh->findNumbering, etc. The shape determines where the nodes are, and the component count determines how many integers there are per node.

Field* apf::createPackedField ( Mesh *  m,
const char *  name,
int  components,
apf::FieldShape shape = 0 
)

Create a field of N components without a tensor type.

Packed fields are used to interface with applications whose fields are not easily categorized as tensors of order 0,1,2. They contain enough information to interpolate values in an element and such, but some higher-level functionality is left out.

Field* apf::createStepField ( Mesh *  m,
const char *  name,
int  valueType 
)

Create an apf::Field using a step distribution.

A step-wise distribution is a C-1 continuous field defined by one node at each element, with the field value being constant over the element and discontinuous at element boundaries

Parameters:
m the mesh over which the field is defined
name a unique name for this field
valueType the type of field data
Field* apf::createUserField ( Mesh *  m,
const char *  name,
int  valueType,
FieldShape *  s,
Function *  f 
)

Create a Field from a user's analytic function.

This field will use no memory, has values on all nodes, and calls the user Function for all value queries. Writing to this field does nothing.

void apf::deriveMdsModel ( Mesh2 *  in  ) 

build a null model such that apf::verify accepts the mesh.

given an MDS mesh that is (wrongly) classified on a null model, this algorithm will classify all interior entities onto a model region and all boundary entities onto a boundary model entity, as defined by mesh upward adjacencies.

void apf::destroyField ( Field *  f  ) 

Destroy an apf::Field.

This function will also remove any field data that this field attached to its Mesh domain.

void apf::destroyMesh ( Mesh *  m  ) 

Destroys an apf::Mesh.

This only destroys the apf::Mesh object, the underlying mesh database is unaffected. Mesh objects are wrappers over mesh databases, and are created by functions provided outside the APF core.

void apf::destroyMeshElement ( MeshElement *  e  ) 

Destroys a Mesh Element.

This only destroys the apf::MeshElement object, the underlying mesh entity and field data are unaffected.

void apf::destruct ( Mesh2 *  m,
int *&  conn,
int &  nelem,
int &  etype 
)

convert an apf::Mesh2 object into a connectivity array

this is useful for debugging the apf::convert function

void apf::displaceMesh ( Mesh2 *  m,
Field *  d,
double  factor = 1.0 
)

add a field (times a factor) to the mesh coordinates

this is useful in mechanical deformation problems to transform the mesh from reference space to deformed space. Setting the factor to -1 will undo the deformation

void apf::extractCoords ( Mesh2 *  m,
double *&  coords,
int &  nverts 
)

get a contiguous set of global vertex coordinates

this is used for debugging apf::setCoords

void apf::fail ( const char *  why  ) 

Declare failure of code inside APF.

This function prints the string as an APF failure to stderr and then calls abort. It can be called from code that is part of the apf namespace, but not outside of that.

int apf::findIn ( MeshEntity **  a,
int  n,
MeshEntity *  e 
)

find pointer (e) in array (a) of length (n)

Returns:
-1 if not found, otherwise i such that a[i] = e
void apf::findTriDown ( Mesh *  m,
MeshEntity **  verts,
MeshEntity **  down 
)

given the vertices of a triangle, find its edges

Parameters:
down the resulting array of edges
MeshEntity* apf::findUpward ( Mesh *  m,
int  type,
MeshEntity **  down 
)

find an entity from one-level downward adjacencies

this function ignores the ordering of adjacent entities

void apf::fix ( Numbering *  n,
MeshEntity *  e,
int  node,
int  component,
bool  fixed 
)

Set the fixed/free status of a degree of freedom,.

must be called prior to making any isFixed() calls on the same node/component.

Parameters:
n the numbering object
e the mesh entity with which the node is associated
node the node number withing the mesh entity
component the component number within the nodal tensor
void apf::getAlignment ( Mesh *  m,
MeshEntity *  elem,
MeshEntity *  boundary,
int &  which,
bool &  flip,
int &  rotate 
)

boundary entity alignment to an element

Parameters:
m the mesh
elem the element
boundary an entity on the boundary of elem
which index of (boundary) in getDownward((elem)...)
flip true iff orientation of (boundary) is opposite canonical
rotate position of canonical vertex 0 in boundary vertices, or in boundary vertices reversed if (flip)==true
double* apf::getArrayData ( Field *  f  ) 

Return the contiguous array storing this field.

This function is only defined for fields which are using array storage, for which apf::isFrozen returns true.

void apf::getBF ( FieldShape *  s,
MeshElement *  e,
Vector3 const &  p,
NewArray< double > &  BF 
)

Get the basis functions over a mesh element.

Parameters:
s the field shape that defines the basis functions
e the mesh element over which the basis functions are defined
p the local coordinates at which the basis functions are evaluated
BF nodal array of basis functions evaluated at p
template<std::size_t M, std::size_t N>
double apf::getCofactor ( Matrix< M, N > const &  A,
std::size_t  i,
std::size_t  j 
)

get the cofactor associated with entry (i,j) of matrix A

this is only instantiated for square matrices up to 4 by 4

void apf::getComponents ( Field *  f,
MeshEntity *  e,
int  node,
double *  components 
)

Copy the nodal value into an array of component values.

the output array must already be allocated, apf::countComponents can help with this.

FieldShape* apf::getConstant ( int  dimension  ) 

Get the Constant shape function over some dimension.

this pseudo-shape function places a node on every element of the given dimension. Dimensions up to 3 are available

void apf::getCurl ( Element *  e,
Vector3 const &  param,
Vector3 &  curl 
)

Evaluate the curl of a vector field at a point.

Parameters:
param The local coordinates in the element.
curl The curl vector at that point.
template<std::size_t M, std::size_t N>
double apf::getDeterminant ( Matrix< M, N > const &  A  ) 

get the determinant of a matrix A

this is only instantiated for square matrices up to 4 by 4

double apf::getDiv ( Element *  e,
Vector3 const &  param 
)

Evaluate the divergence of a vector field at a point.

Parameters:
param The local coordinates in the element.
Returns:
The divergence at that point.
double apf::getDV ( MeshElement *  e,
Vector3 const &  param 
)

Get the differential volume at a point.

This function is meant to provide the differential measure of an element at a point, based on the Jacobian determinant in the case of regions, and equivalent terms for lower dimensions.

Returns:
The differential volume
int apf::getElementNumbers ( Numbering *  n,
MeshEntity *  e,
NewArray< int > &  numbers 
)

returns the node numbers of an element

numbers are returned in the standard element node ordering for its shape functions

Returns:
the number of element nodes
int* apf::getElementToElement ( apf::Mesh m  ) 

Get an element-to-element connectivity array.

this function assumes the mesh has one element type. the resulting array is created with new int[nelements * nsides]. nsides is the number of faces of an element. entry [i * nsides + j] is the global id of the j'th adjacent element to local element i, which can be -1 for a geometric boundary.

Matrix3x3 apf::getFrame ( Vector3 const &  v  ) 

derive an orthogonal frame whose x axis is the given vector

this is a robust algorithm for choosing an arbitrary frame that is aligned with a given vector while avoiding the numerical issues that arise when the vector is near global axes.

Note, however, that the resulting frame is not normalized. some uses of this function require that the x basis vector be exactly the input, so that is the default behavior. It is the user's responsibility to normalize the result if desired

void apf::getGaussPoint ( int  type,
int  order,
int  point,
Vector3 &  param 
)

Return the location of a gaussian integration point.

Parameters:
type the element type, from apf::Mesh::getType
order the order of the integration rule
point which point of the integration rule
param the resulting parent element coordinates
void apf::getGrad ( Element *  e,
Vector3 const &  param,
Vector3 &  grad 
)

Get the gradient of a scalar field w.r.t. global coordinates.

Parameters:
param The local coordinates in the element.
grad The gradient vector at that point.
void apf::getGradBF ( FieldShape *  s,
MeshElement *  e,
Vector3 const &  p,
NewArray< Vector3 > &  gradBF 
)

Get global gradients of basis functions over a mesh element.

Parameters:
gradBF nodal array of gradients of basis functions evaluated at p
FieldShape* apf::getHierarchic ( int  order  ) 

Get the quadratic hierarchic shape function.

only first and second order so far

void apf::getIntPoint ( MeshElement *  e,
int  order,
int  point,
Vector3 &  param 
)

Get an integration point in an element.

Parameters:
order The polynomial order of accuracy.
point The integration point number.
param The resulting local coordinate of the integration point.
double apf::getIntWeight ( MeshElement *  e,
int  order,
int  point 
)

Get the weight of an integration point in an element.

All integration point tables in APF are scaled such that the sum of the weights equals the area of of the parent element.

FieldShape* apf::getIPFitShape ( int  dimension,
int  order 
)

Get the IP Fit shape function.

the IP Fit FieldShape is equivalent to the IPShape except that it is capable of evaluating as a shape function whose value at any point in the element is a polynomial fit to the integration point data in that element.

FieldShape* apf::getIPShape ( int  dimension,
int  order 
)

Get the Integration Point distribution.

Parameters:
dimension The dimensionality of the elements
order The order of accuracy, determines the integration points

This allows users to create a field which has values at the integration points of elements. Orders 1 to 3 for dimension 2 or 3 are available

double apf::getJacobianDeterminant ( Matrix3x3 const &  J,
int  dimension 
)

Return the Jacobian determinant or dimensional equivalent.

this is a special function taking into account the various formulae for differential volume, area, lenght, etc.

Parameters:
J Jacobian matrix, vector gradient of coordinate field with respect to parent element coordinates
dimension spacial dimension of the entity being measured
void apf::getJacobianInv ( MeshElement *  e,
Vector3 const &  local,
Matrix3x3 &  jinv 
)

Returns the Jacobian inverse at a local point.

computes the Moore-Penrose psuedo-inverse of J. This is needed to handle non-square Jacobians that arise from edges embedded in 2D or higher and faces embeded in 3D. If J is invertible, the Moore-Penrose inverse equals the inverse.

FieldShape* apf::getLagrange ( int  order  ) 

Get the Lagrangian shape function of some polynomial order.

we have only first and second order so far

Vector3 apf::getLinearCentroid ( Mesh *  m,
MeshEntity *  e 
)

get the average of the entity's vertex coordinates

this also works if given just a vertex, so its a convenient way to get the centroid of any entity

void apf::getMatrix ( Element *  e,
Vector3 const &  param,
Matrix3x3 &  value 
)

Evaluate the value of a matrix field.

Parameters:
param The local coordinates in the element.
value The field value at that point.
MeshEntity* apf::getMdsEntity ( Mesh2 *  in,
int  dimension,
int  index 
)

retrieve an entity by dimension and index

indices follow iteration order, so this function is equivalent to iterating (index) times, but is actually much faster than that. this function only works when the arrays have no gaps, so call apf::reorderMdsMesh after any mesh modification.

int apf::getMdsIndex ( Mesh2 *  in,
MeshEntity *  e 
)

returns the dimension-unique index for this entity

this function only works when the arrays have no gaps, so call apf::reorderMdsMesh after any mesh modification.

MeshElement* apf::getMeshElement ( Element *  e  ) 

Get the Mesh Element of a Field Element.

Each apf::Element operates over an apf::MeshElement, which must be maintained as long as the apf::Element exists. Multiple apf::Elements may share an apf::MeshElement.

template<std::size_t M, std::size_t N>
Matrix<M-1,N-1> apf::getMinor ( Matrix< M, N > const &  A,
std::size_t  i,
std::size_t  j 
)

get the minor matrix associated with entry (i,j) of matrix A

this is only instantiated for square matrices up to 4 by 4

const char* apf::getName ( Field *  f  ) 

Get the name of a Field.

Both for use convenience and for technical reasons related to tagging, each Field should have a unique name.

void apf::getNodes ( Numbering *  n,
DynamicArray< Node > &  nodes 
)

get an array of numbered nodes

the array is sorted by dimension first and then by mesh iterator order

void apf::getNodesOnClosure ( Mesh *  m,
ModelEntity *  me,
DynamicArray< Node > &  on,
FieldShape *  sh = 0 
)

get nodes on the closure of a model entity

all local nodes associated with mesh entities classified on the closure of the model entity are returned. this is useful for applying boundary conditions, and correctly handles cases when a part has, for example, vertices classified on the edge of a model face, but none of the triangles classified on the model face itself.

double apf::getScalar ( Field *  f,
MeshEntity *  e,
int  node 
)

Get the node value of a scalar field.

Returns:
The field value at that node.
double apf::getScalar ( Element *  e,
Vector3 const &  param 
)

Evaluate a scalar field at a point.

Parameters:
param The local coordinates in the element.
Returns:
The field value at that point.
void apf::getShapeGrads ( Element *  e,
Vector3 const &  local,
NewArray< Vector3 > &  grads 
)

Returns the shape function gradients at a point.

these are gradients with respect to global coordinates.

Sharing* apf::getSharing ( Mesh *  m  ) 

create a default sharing object for this mesh

for normal meshes, the sharing object just describes remote copies. For matched meshes, the sharing object describes matches for matched entities and remote copies for other entities

void apf::getVector ( Element *  e,
Vector3 const &  param,
Vector3 &  value 
)

Evaluate a vector field at a point.

Parameters:
param The local coordinates in the element.
value The field value at that point.
void apf::getVectorGrad ( Element *  e,
Vector3 const &  param,
Matrix3x3 &  deriv 
)

Get the gradient of a vector field w.r.t. global coordinates.

Parameters:
param The local coordinates in the element.
deriv The gradient matrix at that point.
FieldShape* apf::getVoronoiShape ( int  dimension,
int  order 
)

Get the Voronoi shape function.

the Voronoi FieldShape is equivalent to the IPShape except that it is capable of evaluating as a shape function whose value at any point in the element is the value of the closest integration point in that element.

void apf::initResidence ( Mesh2 *  m,
int  dim 
)

Set entity residence based on remote copies.

this function acts on all entities of one dimension

Mesh2* apf::loadMdsFromANSYS ( const char *  nodefile,
const char *  elemfile 
)

load an MDS mesh from ANSYS .node and .elem files

this call takes two filenames, one for a .node and another for a .elem file.

the resulting MDS mesh will be constructed with a null geometric model via gmi_load(".null"), so be sure to call gmi_register_null before this function.

currently, ANSYS element types SOLID72 and SOLID92 are supported, which become linear and quadratic tetrahedra, respectively.

Mesh2* apf::loadMdsMesh ( gmi_model model,
const char *  meshfile 
)

load an MDS mesh from files

Parameters:
model the geometric model interface
meshfile The path to an SMB format mesh. If the path is "something"".smb", then the file "somethingN.smb" will be loaded where N is the part number. If the path is "something/", then the file "something/N.smb" will be loaded. For both of these cases, if the path is prepended with "bz2:", then it will be uncompressed using PCU file IO functions. Calling apf::Mesh::writeNative on the resulting object will do the same in reverse.
Mesh2* apf::loadMdsMesh ( const char *  modelfile,
const char *  meshfile 
)

load an MDS mesh and model from file

Parameters:
modelfile will be passed to gmi_load to get the model
Mesh2* apf::makeEmptyMdsMesh ( gmi_model model,
int  dim,
bool  isMatched 
)

create an empty MDS part

Parameters:
model the geometric model interface
dim the eventual mesh dimension. MDS needs to allocate arrays based on this before users add entities.
isMatched whether or not there will be matched entities
GlobalNumbering* apf::makeGlobal ( Numbering *  n,
bool  destroy = true 
)

converts a local numbering into a global numbering.

Parameters:
destroy Should the input Numbering* be destroyed?

the original local numbering is destroyed. All local numbers are increased by an offset; the offset on part P is the sum of the numbered nodes on parts [0,P-1].

this is done in O(log N) time in parallel, where N is the part count.

this function has no intrinsic knowledge of ownership, it operates simply on nodes which have been explicitly numbered. the input to this function is usually a numbering produced by a numberOwned* function, and the result is a global numbering of all the owned nodes. subsequently calling synchronize on the global numbering completes the typical process which leaves all nodes (owned and not) with a global number attached.

Balancer* apf::makeZoltanBalancer ( Mesh *  mesh,
int  method,
int  approach,
bool  debug = true 
)

Make a Zoltan Balancer object.

this Balancer will apply Zoltan to the global mesh to improve load balance.

Also note that this Balancer will create a Zoltan edge between two elements that share matched faces.

Parameters:
method select from apf::ZoltanMethod
approach select from apf::ZoltanApproach
debug print the full Zoltan configuration
Splitter* apf::makeZoltanGlobalSplitter ( Mesh *  mesh,
int  method,
int  approach,
bool  debug = true 
)

Make a Zoltan Splitter object.

the resulting splitter will apply Zoltan to the global mesh part to break it into several new parts.

Parameters:
method select from apf::ZoltanMethod
approach select from apf::ZoltanApproach
debug print the full Zoltan configuration when splitting
Splitter* apf::makeZoltanSplitter ( Mesh *  mesh,
int  method,
int  approach,
bool  debug = true,
bool  sync = true 
)

Make a Zoltan Splitter object.

the resulting splitter will apply Zoltan to the local mesh part to break it into several new parts.

Parameters:
method select from apf::ZoltanMethod
approach select from apf::ZoltanApproach
debug print the full Zoltan configuration when splitting
sync all parts are splitting by the same factor, multiply the part ids in the resulting apf::Migration accordingly
double apf::measure ( MeshElement *  e  ) 

Measures the volume, area, or length of a Mesh Element.

By integrating the differential volume over the element, a general measure is obtained. This correctly measures curved meshes.

Returns:
The measure of the element
double apf::measure ( Mesh *  m,
MeshEntity *  e 
)

Measures the volume, area, or length of a Mesh Entity.

Returns:
The measure of the element
void apf::migrate ( Mesh2 *  m,
Migration *  plan 
)

APF's migration function, works on apf::Mesh2.

if your database implements apf::Mesh2 (and residence is separate from remote copies) then you may use this to implement most of apf::Mesh::migrate. Users of APF are encouraged to call apf::Mesh::migrate instead of calling this directly

int apf::NaiveOrder ( Numbering *  num  ) 

number all components by simple iteration

Todo:
name should be lower-case
Numbering* apf::numberOverlapNodes ( Mesh *  mesh,
const char *  name,
FieldShape *  s = 0 
)

number all local nodes

Parameters:
s if non-zero, use nodes from this FieldShape, otherwise use the mesh's coordinate nodes
Numbering* apf::numberOwnedNodes ( Mesh *  mesh,
const char *  name,
FieldShape *  s = 0,
Sharing *  shr = 0 
)

number the local owned nodes

Parameters:
s if non-zero, use nodes from this FieldShape, otherwise use the mesh's coordinate nodes
shr if non-zero, use this Sharing to determine ownership, otherwise call apf::getSharing
Field* apf::recoverGradientByVolume ( Field *  f  ) 

Compute a nodal gradient field from a nodal input field.

given a nodal field, compute approximate nodal gradient values by giving each node a volume-weighted average of the gradients computed at each element around it.

void apf::remapPartition ( apf::Mesh2 m,
Remap &  remap 
)

remap all part ids in the mesh structure

when using sub-group partitioning schemes or splitting meshes (see Parma_SplitPartition or apf::Splitter), it is useful to be able to update all partition model structures in a mesh to reflect a transition from one partitioning scheme to the next.

this function applies the given map to all part ids in the remote copies, resident sets, and matching using the apf::Mesh2 interface

MeshTag* apf::reorder ( Mesh *  mesh,
const char *  name 
)

Number by adjacency graph traversal.

a plain single-integer tag is used to number the vertices and elements of a mesh

void apf::reorderMdsMesh ( Mesh2 *  mesh,
MeshTag *  t = 0 
)

apply adjacency-based reordering

Parameters:
t Optional user-defined ordering of the vertices. Set this to NULL to use the internal ordering system. Otherwise, attach a unique integer to each vertex in the range [0, vertices). this will indicate the order in which they appear after reordering.

similar to the algorithm for apf::reorder, this function will traverse adjacencies to reorder each topological type. Then all MDS arrays are re-formed in this new order. An important side effect of this function is that there are no gaps in the MDS arrays after this

Matrix3x3 apf::rotate ( Vector3 const &  u,
double  a 
)

get the rotation matrix around an axis

Parameters:
u the axis to rotate around
a the amount of rotation in radians
void apf::setCoords ( Mesh2 *  m,
const double *  coords,
int  nverts,
GlobalToVert &  globalToVert 
)

Assign coordinates to the mesh.

Each peer provides a set of the coordinates. The coords most be ordered according to the global ids of the vertices. Peer 0 provides the coords for vertices 0 to m-1, peer to for m to n-1, ... After this call, all vertices in the apf::Mesh2 object have correct coordinates assigned.

void apf::setMatrix ( Field *  f,
MeshEntity *  e,
int  node,
Matrix3x3 const &  value 
)

Set the nodal value of a matrix field.

Parameters:
value The vector value.
void apf::setMigrationLimit ( size_t  maxElements  ) 

set the maximum elements that apf::migrate moves at once

apf::migrate implements gradual limited migration in an effort to help applications keep memory use to a minimum. This function globally sets the limit on migration, which causes any migration requests greater than the limit to be performed as several consecutive migrations.

void apf::SetNumberingOffset ( Numbering *  num,
int  off 
)

add an offset to all free nodal component numbers

Todo:
name should be lower-case
void apf::setScalar ( Field *  f,
MeshEntity *  e,
int  node,
double  value 
)

Set a nodal value of a scalar field.

Parameters:
node The node number within the entity. So far, it is just 0 for vertices and for edges in 2nd order. higher order will bring more nodes per edge and onto faces and such.
void apf::setVector ( Field *  f,
MeshEntity *  e,
int  node,
Vector3 const &  value 
)

Set the nodal value of a vector field.

Parameters:
value The vector value.
void apf::stitchMesh ( Mesh2 *  m  ) 

infer all remote copies from those of vertices

given that the remote copies of the vertices are set up correctly, this function will synchronize the remote copies and resident part sets for all other entities correctly.

void apf::synchronize ( Field *  f,
Sharing *  shr = 0 
)

Synchronize field values along partition boundary.

Using the ownership and copies described by an apf::Sharing object, copy values from the owned nodes to their copies, possibly assigning them values for the first time.

void apf::synchronize ( Numbering *  n,
Sharing *  shr = 0 
)

numbers non-owned nodes with the values from their owners

Works even if the non-owned nodes have no number currently assigned, after this they are all numbered

Parameters:
shr if non-zero, use this Sharing model to determine ownership and copies, otherwise call apf::getSharing
MeshTag* apf::tagOpposites ( GlobalNumbering *  gn,
const char *  name 
)

Tag global ids of opposite elements to boundary faces.

this function creates a LONG tag of one value and attaches to all partition boundary faces the global id of the element on the other side.

Parameters:
gn global element numbering
name the name of the resulting tag
void apf::unfreezeFields ( Mesh *  m  ) 

unfreeze all associated fields

see apf::unfreezeField

void apf::unite ( Parts &  into,
Parts const &  from 
)

unite two sets of unique part ids

Parameters:
into becomes the union
void apf::verify ( Mesh *  m  ) 

run consistency checks on an apf::Mesh structure

this can be used to implement apf::Mesh::verify. Other implementations may define their own.

void apf::writeASCIIVtkFiles ( const char *  prefix,
Mesh *  m,
std::vector< std::string >  writeFields 
)

Write a set of parallel VTK Unstructured Mesh files from an apf::Mesh with ASCII encoding.

Only fields whose name appears in the vector writeFields will be output. Nodal fields whose shape differs from the mesh shape will not be output. Fields with incomplete data will not be output.

void apf::writeASCIIVtkFiles ( const char *  prefix,
Mesh *  m 
)

Write a set of parallel VTK Unstructured Mesh files from an apf::Mesh with ASCII encoding.

Nodal fields whose shape differs from the mesh shape will not be output. Fields with incomplete data will not be output.

void apf::writeOneVtkFile ( const char *  prefix,
Mesh *  m 
)

Output just the .vtu file with ASCII encoding for this part.

this function is useful for debugging large parallel meshes.

void apf::writeVtkFiles ( const char *  prefix,
Mesh *  m 
)

Write a set of parallel VTK Unstructured Mesh files from an apf::Mesh with binary (base64) encoding and zlib compression (if LION_COMPRESS=ON).

Nodal fields whose shape differs from the mesh shape will not be output. Fields with incomplete data will not be output.

void apf::writeVtkFiles ( const char *  prefix,
Mesh *  m,
std::vector< std::string >  writeFields 
)

Write a set of parallel VTK Unstructured Mesh files from an apf::Mesh with binary (base64) encoding and zlib compression (if LION_COMPRESS=ON).

Only fields whose name appears in the vector writeFields will be output. Nodal fields whose shape differs from the mesh shape will not be output. Fields with incomplete data will not be output.


Variable Documentation

double const apf::pi

The mathematical constant pi.

although it doesn't fit perfectly in this header, there is no more appropriate place for this in APF.

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Defines