*vfs_DofTabBegin - create an instance of a DofTab object vfs_DofTabEnd - destroy an instance of a DofTab object vfs_DofTabError - return DofTab object error flag
vfs_DofTabDef - define number of nodes and elements Inq - inquire number of nodes and elements vfs_DofTabElemDof - query element global dof vfs_DofTabElemSup - query element global suppressed dof vfs_DofTabElemUse - specify element degree of freedom use vfs_DofTabGetUseTag - query used degree of freedom tags vfs_DofTabInvDof - query inverse degree of freedom information vfs_DofTabMaxElemDof - query for maximum size of element global dof vfs_DofTabNodeDirCos - set direction cosine matrix for preconditioner vfs_DofTabNodeDof - query node global dof vfs_DofTabNodeDofSpec - flag degree of freedom for preconditioner vfs_DofTabNodeSpec - flag node for preconditioner vfs_DofTabNodeSup - query node global suppressed dof vfs_DofTabNodePerm - set nodal degree of freedom constraints vfs_DofTabPointersIndicesPtr - set pointers and indices vfs_DofTabProcess - process and assign global dof vfs_DofTabSetCoords - set the nodal coordinates vfs_DofTabSetFeaType - set the element type vfs_DofTabSetAssemBlock - set running assembling block id vfs_DofTabSetMode - set operating modes vfs_DofTabSetMap - set current element degree of freedom map vfs_DofTabSetMPC - set multipoint constraints vfs_DofTabSetNumThreads - set number of threads vfs_DofTabSetParami - set integer parameters vfs_DofTabSetTopology - set element topology vfs_DofTabThreadAssemNum - query number of assembly groups vfs_DofTabThreadAssemGroups - query assembly groups
Initially all node and element freedoms are assumed to be constrained. Degrees of freedom are activated by declaring their use by a finite element. By looping through all elements and setting the freedoms used by each, all degrees of freedom required by the assembly of all elements will be activated. This is accomplished by setting the degree of freedom map for each element using vfs_DofTabSetMap and then declaring the degrees of freedom used by calling vfs_DofTabElemUse.
After the degree of freedom use for all elements has been established, permanent degree of freedom suppressions at a node are applied using vfs_DofTabNodePerm. These suppressions reflect restraints at nodes. Note that there is no mechanism in DofTab to suppress internal elemental freedoms.
if multipoint constraints (MPCs) are present and they are to be enforced by direct matrix transformation, then the MPCs must be specified using vfs_DofTabSetMPC for each MPC equation.
At this point the DofTab object will have a representation of all active degrees of freedom required by the finite element model. The global active dof can now be numbered by calling vfs_DofTabProcess. The global dof are numbered starting from 1. The total number of global active dof and the total number of global suppressed dof are returned by this function. The function vfs_DofTabProcess also separately numbers the global suppressed dof starting from 1. The degree of freedom numbering is done to minimize the fill for matrix factoring. The current fill reducing algorithm uses Metis.
Once the global dof have been assigned to the active freedoms, the user can query for the global dof numbers at a single node using vfs_DofTabNodeDof or for an entire element using vfs_DofTabElemDof. Calls to vfs_DofTabElemDof must be preceded by vfs_DofTabSetMap. The complete global dof for an element consist of the global dof assigned to the nodal freedoms at each element node plus the global dof assigned to the internal elemental freedoms for the element. A global dof number of 0 returned by these functions represents a suppressed degree of freedom.
The global suppressed dof at a node are returned using vfs_DofTabNodeSup. The global suppressed dof for an element are returned using vfs_DofTabElemSup. A global dof number of 0 returned by this function represents an active degree of freedom. The suppressed dof numbers are useful for assembling reactions.
The function, vfs_DofTabGetUseTag, returns the number and types, or "tags", of all the degrees of freedom used at nodes. For example in a finite element model using only 3D solid elements which use translational degrees of freedom at nodes, the number of degree of freedom types is 3 and the degree of freedom tags are SYS_DOF_TX, SYS_DOF_TY, SYS_DOF_TZ. Use the function vfs_DofTabInvDof to find the node or element ids and degree of freedom tags associated with a set of global degrees of freedom.
Destroy an instance of a DofTab object using vfs_DofTabEnd.
Vint num, hasserial; Vint *numelem, *blockid, *elems; vfs_DofTabThreadAssemNum (doftab,&num,&hasserial); /* allocate arrays */ numelem = (Vint*)malloc(num*sizeof(Vint); blockid = (Vint*)malloc(num*sizeof(Vint); elems = (Vint*)malloc(numel*sizeof(Vint); vfs_DofTabThreadAssemGroups (doftab,&num,numelem,blockid,elems);The number of groups is num. The variable hasserial is non-zero if the last group is meant to be processed serially.
The array numelems contains the number of elements in each group, and the array elems contains the elements in each group. The element numbers in elems are one based. The first numelem[0] elements in elems are in the first group, the next numelem[1] elements are in the second group, etc. Note that if DOFTAB_THREADASSEM has not be enabled, num is returned as 1 and hasserial is enabled.
Assembly of elements can then proceed as follows. In this example 8 threads are used.
Vint n, m; Vint in, istart; Vint numthreads; maxthreads = 8; /* setup SysMatrix for parallel assembly */ vfs_SysMatrixSetNumThreadsAssem (sysmatrix,maxthreads); /* loop over groups */ istart = 0; for(n = 0; n < num; n++) { /* special case of serial last group */ if(blockid[n] <= 0) { numthreads = 1; } else { numthreads = maxthreads; } #pragma omp parallel for num_threads(numthreads) /* loop over elements in group */ for(m = 0; m < numelem[n]; m++) { Vint ithread; Vint in, nedofs, lm[maximum_element_dofs]; Vint nix, ix[maximum_element_connectivity]; Vdouble s[maximum_element_matrix_size]; in = elems[istart+m]; /* get element connectivity, ix[nix] */ ... /* get number of element freedoms and freedom list */ vfs_DofTabElemDof (doftab,in,nix,ix,&nedofs,lm); /* gather/compute element matrix s element in */ ... /* get thread number, increment by one */ ithread = omp_get_thread_num()+1; /* assemble */ vfs_SysMatrixThreadAssem (sysmatrix,ithread,nedofs,lm,s); } istart += numelem[n]; } /* assembly complete, free memory */ free(numelem); free(elems);
In addition vfs_DofTabSetParami must be used to set DOFTAB_REDUCE to DOFTAB_REDUCE_ITER.
A feature of the iterative solver is to flag particular nodes to be given special treatment - the "stiffness" associated with the node shall be represented exactly in the preconditioner at the expense of making the preconditioner somewhat more heavy weight. Use the function vfs_DofTabNodeSpec to specify such nodes. Examples of candidate nodes would be nodes attached to distorted elements, nodes in areas of very high solution gradients, nodes on discrete contacting surfaces, etc. All of this information must be supplied before vfs_DofTabProcess is called.
The iterative solver preconditioner technology has two "aggressiveness" settings, DOFTAB_ITERMEMLEVEL and DOFTAB_ITEREXPLEVEL, which can be changed from the default values if a model is relatively large and contains highly refined regions. These parameters are set using vfs_DofTabSetParami.
*vfs_DofTabBegin - create an instance of a DofTab object
vfs_DofTab *vfs_DofTabBegin ()
None
Destroy an instance of a DofTab object using
void vfs_DofTabEnd (vfs_DofTab *doftab )
Return the current value of a DofTab object error flag using
Vint vfs_DofTabError (vfs_DofTab *doftab )
vfs_DofTabDef - define number of nodes and elements
void vfs_DofTabDef (vfs_DofTab *doftab, Vint numnp, Vint numel)
doftab Pointer to DofTab object. numnp Number of nodes numel Number of elements
None
Inquire of defined numnp and numel as output arguments using
void vfs_DofTabInq (vfs_DofTab *doftab, Vint *numnp, Vint *numel)
vfs_DofTabElemDof - query element global dof
void vfs_DofTabElemDof (vfs_DofTab *doftab, Vint eid, Vint nix, Vint ix[], Vint *nedofs, Vint lm[])
doftab Pointer to DofTab object. eid Element number nix Number of nodes per element ix Vector of element node connectivity
nedofs Number of elements dofs lm Vector of element global dof numbers
vfs_DofTabElemSup - query element global suppressed dof
void vfs_DofTabElemSup (vfs_DofTab *doftab, Vint eid, Vint nix, Vint ix[], Vint *nedofs, Vint lm[])
doftab Pointer to DofTab object. eid Element number nix Number of nodes per element ix Vector of element node connectivity
nedofs Number of element dofs lm Vector of element global dof numbers
vfs_DofTabElemUse - specify element degree of freedom use
void vfs_DofTabElemUse (vfs_DofTab *doftab, Vint eid Vint nix, Vint ix[])
doftab Pointer to DofTab object. eid Element number nix Number of nodes per element ix Vector of element node connectivity
None
vfs_DofTabGetUseTag - query used degree of freedom tags
void vfs_DofTabGetUseTag (vfs_DofTab *doftab, Vint *ntags, Vint tag[])
doftab Pointer to DofTab object.
ntags Number of degrees of freedom types used at nodes tag Vector of degree of freedom types
vfs_DofTabInvDof - query inverse degree of freedom information
void vfs_DofTabInvDof (vfs_DofTab *doftab, Vint ndofs, Vint lm[], Vint ids[], Vint tag[])
doftab Pointer to DofTab object. ndofs Number of degrees of freedom lm Vector of global dof numbers
ids Node or element ids associated with global dof tag Vector of degree of freedom types
vfs_DofTabMaxElemDof - query for maximum size of element global dof
void vfs_DofTabMaxElemDof (vfs_DofTab *doftab, Vint *maxelemdof)
doftab Pointer to DofTab object.
maxelemdof Maximum number of element degrees of freedom
vfs_DofTabNodeDirCos - set direction cosine matrix at nodes for iterative solver
void vfs_DofTabNodeDirCos (vfs_DofTab *doftab, Vint nid, Vdouble tm[3][3])
doftab Pointer to DofTab object. nid Node id tm Direction cosine matrix at node
None
vfs_DofTabNodeDof - query node global dof
void vfs_DofTabNodeDof (vfs_DofTab *doftab, Vint nid, Vint ntags, Vint tag[], Vint lm[])
doftab Pointer to DofTab object. nid Node number ntags Number of degrees of freedom types used at nodes tag Vector of degree of freedom types at nodes
lm Vector of global dof numbers at node
vfs_DofTabNodeDofSpec - flag degree of freedom for preconditioner
void vfs_DofTabNodeDofSpec (vfs_DofTab *doftab, Vint nix, Vint ix[], Vint ntags, Vint tags[])
doftab Pointer to DofTab object. nix Number of nodes being flagged ix List of nodes being flagged ntags Number of nodal degree of freedom tags being flagged tags List of degrees of freedom being flagged
None
vfs_DofTabNodeSpec - flag all of a node's degrees of freedom for preconditioner
void vfs_DofTabNodeSpec (vfs_DofTab *doftab, Vint level, Vint nix, Vint ix[])
doftab Pointer to DofTab object. level Propagation level for flagging node nix Number of nodes being flagged ix List of nodes being flagged
None
vfs_DofTabNodePerm - set nodal degree of freedom constraints
void vfs_DofTabNodePerm (vfs_DofTab *doftab, Vint nid, Vint ntags, Vint tag[])
doftab Pointer to DofTab object. nid Node number ntags Number of degrees of freedom types at node tag Vector of degree of freedom types
None
vfs_DofTabPointersIndicesPtr - set pointers and indices
void vfs_DofTabNodePointersIndicesPtr (vfs_DofTab *doftab, Vlong pointers[], Vint indices[])
doftab Pointer to DofTab object. pointers Location in indices for first row in each column indices List of row indices for each column
None
This function is used in conjunction with vfs_SysMatrixSetValue. In this function, entry a_ij refers to the matrix value for column i, row j.
vfs_DofTabNodeSup - query node global suppressed dof
void vfs_DofTabNodeSup (vfs_DofTab *doftab, Vint nid, Vint ntags, Vint tag[], Vint lm[])
doftab Pointer to DofTab object. nid Node number ntags Number of degrees of freedom types used at nodes tag Vector of degree of freedom types at nodes
lm Vector of global suppressed dof numbers at node
vfs_DofTabProcess - process and assign global dof
void vfs_DofTabProcess (vfs_DofTab *doftab, Vint *neq, Vint *nre)
doftab Pointer to DofTab object.
neq Number of global (active) dof nre Number of global suppressed dof
vfs_DofTabSetCoords - set the nodal coordinates
void vfs_DofTabSetCoords (vfs_DofTab *doftab, Vint nid, Vdouble x[3])
doftab Pointer to DofTab object. nid Node id x Vector of nodal (x,y,z) coordinates
None
vfs_DofTabSetAsseBlock - set running assembling block id
void vfs_DofTabSetAssemBlock (vfs_DofTab *doftab, Vint blockid)
doftab Pointer to DofTab object. blockid Running block id
None
This function is optional, and is only relevant for when users want to assemble elements in vfs_SysMatrix in parallel. When this is the case, DofTab will group together elements that can be assembled in parallel.
Certain applications group elements in a particular way. As a result, the assembly parallelization must be done only for the group of elements available at the time. By assigning a blockid to each element users are guaranteed that a group of elements that can be assembled in parallel belongs to the same blockid.
In general, blockid values are positive. However, users can also specify negative blockid values. When this is the case, it is expected that all elements associated with this block will be assembled in serial mode and need not be grouped.
DofTab guarantees that if both blockid and -blockid are present, then the positive values will be listed immediately before the negative values. Hence, for each positive/negative pair, the elements that can be assembled in parallel will be listed first, immediately followed by the corresponding negative blockid values.
If this function is called for any element, then it must be called for all elements. By default, all elements have blockid set to zero. If any element has an assigned non-zero blockid, then the presence of a zero-valued blockid will generate an error during vfs_DofTabProcess.
vfs_DofTabSetFeaType - set the element type
void vfs_DofTabFeaType (vfs_DofTab *doftab, Vint eid, Vint featype)
doftab Pointer to DofTab object. eid Element id featype The element type
None
vfs_DofTabSetMap - set current element degree of freedom map
void vfs_DofTabSetMap (vfs_DofTab *doftab, Vint nedofs, Vint loc[], Vint tag[])
doftab Pointer to DofTab object. nedofs Number of element degrees of freedom loc Vector of degree of freedom locations tag Vector of degree of freedom types
None
For example, the degree of freedom map for a 4 node tetrahedron, with translational degrees of freedom at each node, is set as follows.
Vint loc[12] = {1,1,1, 2,2,2, 3,3,3, 4,4,4}; Vint tag[12] = {SYS_DOF_TX, SYS_DOF_TY, SYS_DOF_TZ, SYS_DOF_TX, SYS_DOF_TY, SYS_DOF_TZ, SYS_DOF_TX, SYS_DOF_TY, SYS_DOF_TZ, SYS_DOF_TX, SYS_DOF_TY, SYS_DOF_TZ}; vfs_DofTabSetMap (doftab,12,loc,tag);
As another example, the degree of freedom map for a homogeneous multipoint constraint equation of the general form a*tx(node 20) + b*ty(node 20) + c*ty(node 1) = 0. to be enforced using a Lagrange multiplier is set as follows.
Vint loc[4] = {1,1,2,0}; Vint tag[12] = {SYS_DOF_TX, SYS_DOF_TY, SYS_DOF_TY, SYS_DOF_LAGM}; vfs_DofTabSetMap (doftab,4,loc,tag);
vfs_DofTabSetMode - set operating mode
void vfs_DofTabSetMode (vfs_DofTab *doftab, Vint mode, Vint param)
doftab Pointer to DofTab object. mode Type of mode to set =DOFTAB_SHAPEPOINTMODE shape point preconditioner treatment param Mode parameter =SYS_ON Enable mode =SYS_OFF Disable mode
None
Specify the shape point iterative solver preconditioner treatment for shape point elements using DOFTAB_SHAPEPOINTMODE. If enabled shape point element connected nodes will be treated internally as specified nodes using the function vfs_DofTabNodeSpec. If disabled they will not be given this treatment. If a shape point element nodes represent very discrete topologies then enable this mode. If the shape point nodes represent a distributed connection along an edge or particularly a face disable this mode. By default DOFTAB_SHAPEPOINTMODE mode is to SYS_ON.
vfs_DofTabSetMPC - set multipoint constraints
void vfs_DofTabSetMPC (vfs_DofTab *doftab, Vint ndofs, Vint ix[]) Vint tag[]) Vdouble c[], Vdouble rhs)
doftab Pointer to DofTab object. ndofs Number of multipoint degrees of freedom ix Vector of node connectivity tag Vector of degree of freedom types c Vector of multipoint coefficients rhs Right hand side
None
If all coefficients in c are identically zero, the value rhs will still be associated with the first degree of freedom, and forces will be updated with this value during a solution procedure. However, in this case, and in this case only, the degree of freedom will not be considered dependent.
vfs_DofTabSetNumThreads - set number of threads
void vfs_DofTabSetNumThreads (vfs_DofTab *doftab, Vint num)
doftab Pointer to DofTab object. num Number of threads =0 Default, serial execution. >=1 Number of threads for parallel execution
None
vfs_DofTabSetParami - set integer parameters
void vfs_DofTabSetParami (vfs_DofTab *doftab, Vint type, Vint param)
doftab Pointer to DofTab object. type Integer parameter type to set. =DOFTAB_REDUCE Select reduced solver =DOFTAB_ITERMEMLEVEL Iterative low memory level =DOFTAB_ITEREXPLEVEL Iterative expansion level =DOFTAB_THREADASSEM Threaded assembly flag param Specifies the integer value that type will be set to. =DOFTAB_REDUCE_NONE No reduction =DOFTAB_REDUCE_ITER Iterative solver reduction =DOFTAB_REDUCE_RIGID Rigid mode reduction =SYS_OFF Disable =SYS_ON Enable
None
The parameter DOFTAB_REDUCE enables the processing required by the iterative solver or the rigid mode solver. The iterative solver is selected using DOFTAB_REDUCE_ITER, the rigid mode solver is selected using DOFTAB_REDUCE_RIGID. By default DOFTAB_REDUCE is set to DOFTAB_REDUCE_NONE.
The parameter DOFTAB_ITERMEMLEVEL sets the aggressiveness level for the iterative solver used for higher order element contributions to the system matrix. By default DOFTAB_ITERMEMLEVEL is set to zero. A level of 1 would be used for models containing large, highly refined higher order element regions.
The parameter DOFTAB_ITEREXPLEVEL sets the aggressiveness level for the iterative solver used for linear order element contributions to the system matrix. By default DOFTAB_ITEREXPLEVEL is set to 3. A level of 4 or 5 would be used for models containing large, highly refined linear order element regions.
The parameter DOFTAB_THREADASSEM is used to enable the calculation of disjoint element groups using a coloring algorithm. The number of groups and group contents are queried using vfs_DofTabThreadAssemNum and vfs_DofTabThreadAssemGroups. These disjoint groups are to be used with the SysMatrix object using a threaded assembly procedure. By default DOFTAB_THREADASSEM is disabled.
vfs_DofTabSetTopology - set the element topology
void vfs_DofTabSetTopology (vfs_DofTab *doftab, Vint id, Vint shape, Vint maxi, Vint maxj, Vint maxk)
doftab Pointer to DofTab object. id Element id shape Element shape parameter =SYS_SHAPEPOINT Point =SYS_SHAPELINE Line =SYS_SHAPETRI Triangle =SYS_SHAPEQUAD Quadrilateral =SYS_SHAPETET Tetrahedron =SYS_SHAPEWED Wedge =SYS_SHAPEHEX Hexahedron maxi The number of points along the i direction. If maxi = 0 then the linear element form of the specified shape is assumed. maxj The number of points along the j direction. If maxj = 0 and 2 <= maxi <= 4, then a Serendipity finite element is assumed. If 2 <= maxj and 2 <= maxi, then a Lagrange finite element is assumed. For triangle, tetrahedron and wedge shapes, set either maxj = 0 or maxj = maxi. maxk The number of points along the k direction. For tetrahedron shapes, set either maxk = 0 or maxk = maxi For hexahedron and wedge shapes, set either maxk = 2 or maxk = maxi
None
vfs_DofTabSetTopology (doftab,350,SYS_SHAPETET,3,0,0)
vfs_DofTabThreadAssemNum - query number of assembly groups
void vfs_DofTabThreadAssemNum (vfs_DofTab *doftab, Vint *num)
doftab Pointer to DofTab object.
num Number of assembly groups
vfs_DofTabThreadAssemGroups - query assembly groups
void vfs_DofTabThreadAssemGroups (vfs_DofTab *doftab, Vint *num, Vint numelem[], Vint blockid[], Vint elems[])
doftab Pointer to DofTab object.
num Number of assembly groups numelem Vector of number of elements in each group blockid Vector of block ids in each group elems List of elements in each group