OpenFAST
Wind turbine multiphysics simulator
|
The modules ModMesh and ModMesh_Types provide data structures and subroutines for representing and manipulating meshes and meshed data in the FAST modular framework. More...
Data Types | |
interface | meshconstructelement |
Functions/Subroutines | |
subroutine | meshwrbin (UnIn, M, ErrStat, ErrMsg, FileName) |
This routine writes mesh information in binary form. More... | |
subroutine | meshwrvtkreference (RefPoint, M, FileRootName, ErrStat, ErrMsg) |
This routine writes the reference position and orientations of a mesh in VTK format. More... | |
subroutine | meshwrvtk (RefPoint, M, FileRootName, VTKcount, OutputFieldData, ErrStat, ErrMsg, Twidth, Sib) |
This routine writes mesh information in VTK format. More... | |
subroutine | meshwrvtkfields (Un, M, n) |
This routine writes mesh field information in VTK format. More... | |
subroutine | meshwrvtk_ln2surface (RefPoint, M, FileRootName, VTKcount, OutputFieldData, ErrStat, ErrMsg, Twidth, NumSegments, Radius, verts, Sib) |
This routine writes line2 mesh surface information in VTK format. More... | |
subroutine | meshwrvtk_pointsurface (RefPoint, M, FileRootName, VTKcount, OutputFieldData, ErrStat, ErrMsg, Twidth, NumSegments, Radius, verts, Sib) |
This routine writes point mesh surfaces information in VTK format. More... | |
subroutine | meshprintinfo (U, M, N) |
This routine writes mesh information in text form. More... | |
subroutine | meshcreate (BlankMesh, IOS, Nnodes, ErrStat, ErrMess, Force, Moment, Orientation, TranslationDisp, TranslationVel, RotationVel, TranslationAcc, RotationAcc, nScalars, IsNewSibling) |
Takes a blank, uninitialized instance of Type(MeshType) and defines the number of nodes in the mesh. More... | |
recursive subroutine | meshdestroy (Mesh, ErrStat, ErrMess, IgnoreSibling) |
Destroy the given mesh and deallocate all of its data. More... | |
subroutine | meshpack (Mesh, ReKiBuf, DbKiBuf, IntKiBuf, ErrStat, ErrMess, SizeOnly) |
Given a mesh and allocatable buffers of type INTEGER(IntKi), REAL(ReKi), and REAL(DbKi), return the mesh information compacted into consecutive elements of the corresponding buffers. More... | |
subroutine | meshunpack (Mesh, ReKiBuf, DbKiBuf, IntKiBuf, ErrStat, ErrMess) |
Given a blank, uncreated mesh and buffers of type INTEGER(IntKi), REAL(ReKi), and REAL(DbKi), unpack the mesh information from the buffers. More... | |
subroutine | meshcopy (SrcMesh, DestMesh, CtrlCode, ErrStat, ErrMess, IOS, Force, Moment, Orientation, TranslationDisp, TranslationVel, RotationVel, TranslationAcc, RotationAcc, nScalars) |
Given an existing mesh and a destination mesh, create a completely new copy, a sibling, or update the fields of a second existing mesh from the first mesh. More... | |
subroutine | meshpositionnode (Mesh, Inode, Pos, ErrStat, ErrMess, Orient, Ref) |
For a given node in a mesh, assign the coordinates of the node in the global coordinate space. More... | |
subroutine | meshcommit (Mesh, ErrStat, ErrMess) |
Given a mesh that has been created, spatio-located, and constructed, commit the definition of the mesh, making it ready for initialization and use. More... | |
subroutine | meshconstructelement_1pt (Mesh, Xelement, ErrStat, ErrMess, P1) |
Given a mesh and an element name, construct a point element whose vertex is the node index listed as the remaining argument of the call to MeshConstructElement. More... | |
subroutine | bumpupelementtable_new (Mesh, Xelement, ErrStat, ErrMess) |
subroutine | bumpupelementtable (Mesh, Xelement, ErrStat, ErrMess) |
This subroutine increases the allocated space for MeshElemTable(Xelement)Elements if adding a new element will exceed the pre-allocated space. More... | |
subroutine | meshconstructelement_2pt (Mesh, Xelement, ErrStat, ErrMess, P1, P2) |
Given a mesh and an element name, construct 2-point line (line2) element whose vertices are the node indices listed as the remaining arguments of the call to MeshConstructElement. More... | |
subroutine | meshconstructelement_3pt (Mesh, Xelement, ErrStat, ErrMess, P1, P2, P3) |
added 20130102 as stub for AeroDyn work | |
subroutine | meshconstructelement_4pt (Mesh, Xelement, ErrStat, ErrMess, P1, P2, P3, P4) |
added 20130102 as stub for AeroDyn work | |
subroutine | meshconstructelement_6pt (Mesh, Xelement, ErrStat, ErrMess, P1, P2, P3, P4, P5, P6) |
added 20130102 as stub for AeroDyn work | |
subroutine | meshconstructelement_8pt (Mesh, Xelement, ErrStat, ErrMess, P1, P2, P3, P4, P5, P6, P7, P8) |
added 20130102 as stub for AeroDyn work | |
subroutine | meshconstructelement_10pt (Mesh, Xelement, ErrStat, ErrMess, P1, P2, P3, P4, P5, P6, P7, P8, P9, P10) |
added 20130102 as stub for AeroDyn work | |
subroutine | meshconstructelement_15pt (Mesh, Xelement, ErrStat, ErrMess, P1, P2, P3, P4, P5, P6, P7, P8, P9, P10, P11, P12, P13, P14, P15) |
added 20130102 as stub for AeroDyn work | |
subroutine | meshconstructelement_20pt (Mesh, Xelement, ErrStat, ErrMess, P1, P2, P3, P4, P5, P6, P7, P8, P9, P10, P11, P12, P13, P14, P15, P16, P17, P18, P19, P20) |
added 20130102 as stub for AeroDyn work | |
subroutine | meshsplitelement_2pt (Mesh, Xelement, ErrStat, ErrMess, E1, P1) |
This routine splits a line2 element into two separate elements, using p1 as the new node connecting the two new elements formed from E1. More... | |
subroutine | meshnextelement (Mesh, CtrlCode, ErrStat, ErrMess, Ielement, Xelement, ElemRec) |
Given a control code and a mesh that has been committed, retrieve the next element in the mesh. More... | |
subroutine | packloadmesh_names (M, MeshName, Names, indx_first) |
This subroutine returns the names of the output rows/columns in the Jacobian matrices. More... | |
subroutine | packloadmesh (M, Ary, indx_first) |
This subroutine returns the operating point values of the mesh fields. More... | |
subroutine | packloadmesh_dy (M_p, M_m, dY, indx_first) |
This subroutine computes the differences of two meshes and packs that value into appropriate locations in the dY array. More... | |
subroutine | packmotionmesh_names (M, MeshName, Names, indx_first, FieldMask) |
This subroutine returns the names of rows/columns of motion meshes in the Jacobian matrices. More... | |
subroutine | packmotionmesh (M, Ary, indx_first, FieldMask, UseSmlAngle) |
This subroutine returns the operating point values of the mesh fields. More... | |
subroutine | packmotionmesh_dy (M_p, M_m, dY, indx_first, FieldMask, UseSmlAngle) |
This subroutine computes the differences of two meshes and packs that value into appropriate locations in the dY array. More... | |
subroutine | meshextrapinterp1 (u1, u2, tin, u_out, tin_out, ErrStat, ErrMsg) |
This subroutine calculates a extrapolated (or interpolated) input u_out at time t_out, from previous/future time values of u (which has values associated with times in t). More... | |
subroutine | meshextrapinterp2 (u1, u2, u3, tin, u_out, tin_out, ErrStat, ErrMsg) |
This subroutine calculates a extrapolated (or interpolated) input u_out at time t_out, from previous/future time values of u (which has values associated with times in t). More... | |
Variables | |
integer, parameter, private | bumpup = 64 |
size element list will be increased when adding an element that does not fit in the currently allocated space; do not set to less than 2 | |
character(*), parameter | vtk_aryfmt = '(3(F30.5))' |
text format for triplets written to VTK text files | |
The modules ModMesh and ModMesh_Types provide data structures and subroutines for representing and manipulating meshes and meshed data in the FAST modular framework.
A mesh is comprised of a set of "nodes" (simple points in space) together with information specifying how they are connected to form "elements" representing spatial boundaries between components. ModMesh and ModMesh_Types define point, line, surface, and volume elements in a standard isoparametric mapping from finite element analysis. Currently only points and straight line (line2) elements are implemented.
Associated with a mesh are one or more "fields" that represent the values of variables or "degrees of freedom" at each node. A mesh always has a named "Position" that specifies the location in three-dimensional space as an Xi,Yi,Zi triplet of each node and a field named "RefOrientation" that specifies the orientation (as a direction cosine matrix) of the node. The ModMesh_Types module predefines a number of other fields of triples representing velocities, forces, and moments as well as a field of nine values representing a direction cosine matrix.
The operations on meshes defined in the ModMesh module are creation, spatio-location of nodes, construction, committing the mesh definition, initialization of fields, accessing field data, updating field data, copying, deallocating, and destroying meshes. See https://nwtc.nrel.gov/FAST-Developers and https://nwtc.nrel.gov/system/files/ProgrammingHandbook_Mod20130717.pdf
subroutine modmesh::bumpupelementtable | ( | type(meshtype), intent(inout) | Mesh, |
integer(intki), intent(in) | Xelement, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMess | ||
) |
This subroutine increases the allocated space for MeshElemTable(Xelement)Elements if adding a new element will exceed the pre-allocated space.
[in,out] | mesh | Mesh being constructed |
[in] | xelement | type of element (See Element Names) |
[out] | errstat | Error code |
[out] | errmess | Error message |
subroutine modmesh::meshcommit | ( | type(meshtype), intent(inout) | Mesh, |
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMess | ||
) |
Given a mesh that has been created, spatio-located, and constructed, commit the definition of the mesh, making it ready for initialization and use.
Explicitly committing a mesh provides the opportunity to precompute traversal information, neighbor lists and other information about the mesh. Returns non-zero in value of ErrStat on error.
[in,out] | mesh | Mesh being committed |
[out] | errstat | Error code |
[out] | errmess | Error message |
Check for spatial constraints – can't mix 1D with 2D with 3D
Check that every node is part of an element.
subroutine modmesh::meshconstructelement_1pt | ( | type(meshtype), intent(inout) | Mesh, |
integer(intki), intent(in) | Xelement, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMess, | ||
integer, intent(in) | P1 | ||
) |
Given a mesh and an element name, construct a point element whose vertex is the node index listed as the remaining argument of the call to MeshConstructElement.
Returns a non-zero ErrStat value on error.
[in,out] | mesh | Mesh being constructed |
[in] | xelement | See Element Names |
[out] | errstat | Error code |
[out] | errmess | Error message |
[in] | p1 | node index for this point element |
subroutine modmesh::meshconstructelement_2pt | ( | type(meshtype), intent(inout) | Mesh, |
integer(intki), intent(in) | Xelement, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMess, | ||
integer, intent(in) | P1, | ||
integer, intent(in) | P2 | ||
) |
Given a mesh and an element name, construct 2-point line (line2) element whose vertices are the node indices listed as the remaining arguments of the call to MeshConstructElement.
The adjacency of elements is implied when elements are created that share some of the same nodes. Returns a non-zero value on error.
[in,out] | mesh | Mesh being constructed |
[in] | xelement | See Element Names |
[out] | errstat | Error code |
[out] | errmess | Error message |
[in] | p1 | 1 of 2 points that make a 2-point line element |
[in] | p2 | 1 of 2 points that make a 2-point line element |
subroutine modmesh::meshcopy | ( | type(meshtype), intent(inout), target | SrcMesh, |
type(meshtype), intent(inout), target | DestMesh, | ||
integer(intki), intent(in) | CtrlCode, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMess, | ||
integer(intki), intent(in), optional | IOS, | ||
logical, intent(in), optional | Force, | ||
logical, intent(in), optional | Moment, | ||
logical, intent(in), optional | Orientation, | ||
logical, intent(in), optional | TranslationDisp, | ||
logical, intent(in), optional | TranslationVel, | ||
logical, intent(in), optional | RotationVel, | ||
logical, intent(in), optional | TranslationAcc, | ||
logical, intent(in), optional | RotationAcc, | ||
integer(intki), intent(in), optional | nScalars | ||
) |
Given an existing mesh and a destination mesh, create a completely new copy, a sibling, or update the fields of a second existing mesh from the first mesh.
When CtrlCode is MESH_NEWCOPY or MESH_SIBLING, the destination mesh must be a blank, uncreated mesh.
If CtrlCode is MESH_NEWCOPY, an entirely new copy of the mesh is created, including all fields, with the same data values as the original, but as an entirely separate copy in memory. The new copy is in the same state as the original–if the original has not been committed, neither is the copy; in this case, an all-new copy of the mesh must be committed separately.
If CtrlCode is MESH_SIBLING, the destination mesh is created with the same mesh and position/reference orientation information of the source mesh, and this new sibling is added to the end of the list for the set of siblings. Siblings may have different fields (other than Position and RefOrientation). Therefore, for a sibling, it is necessary, as with MeshCreate, to indicate the fields the sibling will have using optional arguments. Sibling meshes should not be created unless the original mesh has been committed first.
If CtrlCode is MESH_UPDATECOPY, all of the allocatable fields of the destination mesh are updated with the values of the fields in the source. (The underlying mesh is untouched.) The mesh and field definitions of the source and destination meshes must match and both must have been already committed. The destination mesh may be an entirely different copy or it may be a sibling of the source mesh.
[in,out] | srcmesh | Mesh being copied |
[in,out] | destmesh | Copy of mesh |
[in] | ctrlcode | MESH_NEWCOPY, MESH_SIBLING, or MESH_UPDATECOPY |
[out] | errstat | Error code |
[out] | errmess | Error message |
[in] | ios | If present, IOS of new sibling: input (COMPONENT_INPUT), output(COMPONENT_OUTPUT), or state(COMPONENT_STATE) |
[in] | force | If present and true, allocate Force field |
[in] | moment | If present and true, allocate Moment field |
[in] | orientation | If present and true, allocate Orientation field |
[in] | translationdisp | If present and true, allocate TranslationDisp field |
[in] | translationvel | If present and true, allocate TranslationVel field |
[in] | rotationvel | If present and true, allocate RotationVel field |
[in] | translationacc | If present and true, allocate TranslationAcc field |
[in] | rotationacc | If present and true, allocate RotationAcc field |
[in] | nscalars | If present and > 0 , alloc n Scalars |
subroutine modmesh::meshcreate | ( | type(meshtype), intent(inout) | BlankMesh, |
integer, intent(in) | IOS, | ||
integer, intent(in) | Nnodes, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMess, | ||
logical, intent(in), optional | Force, | ||
logical, intent(in), optional | Moment, | ||
logical, intent(in), optional | Orientation, | ||
logical, intent(in), optional | TranslationDisp, | ||
logical, intent(in), optional | TranslationVel, | ||
logical, intent(in), optional | RotationVel, | ||
logical, intent(in), optional | TranslationAcc, | ||
logical, intent(in), optional | RotationAcc, | ||
integer, intent(in), optional | nScalars, | ||
logical, intent(in), optional | IsNewSibling | ||
) |
Takes a blank, uninitialized instance of Type(MeshType) and defines the number of nodes in the mesh.
Optional arguments indicate the fields that will be allocated and associated with the nodes of the mesh. The fields that may be associated with the mesh nodes are Force, Moment, Orientation, Rotation, TranslationDisp, RotationVel, TranslationVel, RotationAcc, TranslationAcc, and an arbitrary number of Scalars. See the definition of ModMeshType for descriptions of these fields.
[in,out] | blankmesh | Mesh to be created |
[in] | ios | input (COMPONENT_INPUT), output(COMPONENT_OUTPUT), or state(COMPONENT_STATE) |
[in] | nnodes | Number of nodes in mesh |
[out] | errstat | error status/level |
[out] | errmess | error message |
[in] | force | If present and true, allocate Force field |
[in] | moment | If present and true, allocate Moment field |
[in] | orientation | If present and true, allocate Orientation field |
[in] | translationdisp | If present and true, allocate TranslationDisp field |
[in] | translationvel | If present and true, allocate TranslationVel field |
[in] | rotationvel | If present and true, allocate RotationVel field |
[in] | translationacc | If present and true, allocate TranslationAcc field |
[in] | rotationacc | If present and true, allocate RotationAcc field |
[in] | nscalars | If present and > 0, allocate nScalars Scalars |
[in] | isnewsibling | If present and true, this is an new sibling so don't allocate new shared fields (RemapFlag, position, RefOrientation, and ElemTable) |
This routine will add any required fields that were not explicitly requested. If the mesh has motion fields and it is an input mesh, it must always have the following fields: TranslationDisp, TranslationVel, TranslationAcc
If the mesh has load fields and it is an input mesh, it must always have the following fields: Moment
recursive subroutine modmesh::meshdestroy | ( | type(meshtype), intent(inout) | Mesh, |
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMess, | ||
logical, intent(in), optional | IgnoreSibling | ||
) |
Destroy the given mesh and deallocate all of its data.
If the optional IgnoreSibling argument is set to TRUE, destroying a sibling in a set has no effect on the other siblings other than to remove the victim from the list of siblings. If IgnoreSibling is omitted or is set to FALSE, all of the other siblings in the set will be destroyed as well.
[in,out] | mesh | Mesh to be vaporized |
[out] | errstat | Error status/code |
[out] | errmess | Error message |
[in] | ignoresibling | if IgnoreSibling is present and true, don't follow the sibling pointers. Instead just unconditionally nullify these. Use this carefully, since it can leave dangling memory if used for a mesh that already exists and has existing siblings. |
subroutine modmesh::meshextrapinterp1 | ( | type(meshtype), intent(in) | u1, |
type(meshtype), intent(in) | u2, | ||
real(dbki), dimension(:), intent(in) | tin, | ||
type(meshtype), intent(inout) | u_out, | ||
real(dbki), intent(in) | tin_out, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMsg | ||
) |
This subroutine calculates a extrapolated (or interpolated) input u_out at time t_out, from previous/future time values of u (which has values associated with times in t).
Order of the interpolation is 1.
[in] | u1 | Inputs at t1 > t2 |
[in] | u2 | Inputs at t2 |
[in] | tin | Times associated with the inputs |
[in,out] | u_out | Inputs at tin_out |
[in] | tin_out | time to be extrap/interp'd to |
[out] | errstat | Error status of the operation |
[out] | errmsg | Error message if ErrStat /= ErrID_None |
subroutine modmesh::meshextrapinterp2 | ( | type(meshtype), intent(in) | u1, |
type(meshtype), intent(in) | u2, | ||
type(meshtype), intent(in) | u3, | ||
real(dbki), dimension(:), intent(in) | tin, | ||
type(meshtype), intent(inout) | u_out, | ||
real(dbki), intent(in) | tin_out, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMsg | ||
) |
This subroutine calculates a extrapolated (or interpolated) input u_out at time t_out, from previous/future time values of u (which has values associated with times in t).
Order of the interpolation is 2.
[in] | u1 | Inputs at t1 > t2 > t3 |
[in] | u2 | Inputs at t2 > t3 |
[in] | u3 | Inputs at t3 |
[in] | tin | Times associated with the inputs |
[in,out] | u_out | Inputs at tin_out |
[in] | tin_out | time to be extrap/interp'd to |
[out] | errstat | Error status of the operation |
[out] | errmsg | Error message if ErrStat /= ErrID_None |
subroutine modmesh::meshnextelement | ( | type(meshtype), intent(inout) | Mesh, |
integer(intki), intent(inout) | CtrlCode, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMess, | ||
integer(intki), intent(out), optional | Ielement, | ||
integer(intki), intent(out), optional | Xelement, | ||
type(elemrectype), intent(inout), optional, pointer | ElemRec | ||
) |
Given a control code and a mesh that has been committed, retrieve the next element in the mesh.
Used to traverse mesh element by element. On entry, the CtrlCode argument contains a control code: zero indicates start from the beginning, an integer between 1 and MeshNelemlist returns that element, and MESH_NEXT means return the next element in traversal. On exit, CtrlCode contains the status of the traversal in (zero or MESH_NOMOREELEMS). The routine optionally outputs the index of the element in the mesh's element list, the name of the element (see "Element Names"), and a pointer to the element.
subroutine modmesh::meshpack | ( | type(meshtype), intent(in) | Mesh, |
real(reki), dimension(:), intent(out), allocatable | ReKiBuf, | ||
real(dbki), dimension(:), intent(out), allocatable | DbKiBuf, | ||
integer(intki), dimension(:), intent(out), allocatable | IntKiBuf, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMess, | ||
logical, intent(in), optional | SizeOnly | ||
) |
Given a mesh and allocatable buffers of type INTEGER(IntKi), REAL(ReKi), and REAL(DbKi), return the mesh information compacted into consecutive elements of the corresponding buffers.
This would be done to allow subsequent writing of the buffers to a file for restarting later. The sense of the name is "pack the data from the mesh into buffers". IMPORTANT: MeshPack allocates the three buffers. It is incumbent upon the calling program to deallocate the buffers when they are no longer needed. For sibling meshes, MeshPack should be called separately for each sibling, because the fields allocated with the siblings are separate and unique to each sibling.
subroutine modmesh::meshpositionnode | ( | type(meshtype), intent(inout) | Mesh, |
integer(intki), intent(in) | Inode, | ||
real(reki), dimension(3), intent(in) | Pos, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMess, | ||
real(r8ki), dimension(3,3), intent(in), optional | Orient, | ||
logical, intent(in), optional | Ref | ||
) |
For a given node in a mesh, assign the coordinates of the node in the global coordinate space.
If an Orient argument is included, the node will also be assigned the specified orientation (orientation is assumed to be the identity matrix if omitted). Returns a non-zero value in ErrStat if Inode is outside the range 1..Nnodes.
[in,out] | mesh | Mesh being spatio-located |
[in] | inode | Number of node being located |
[in] | pos | Xi,Yi,Zi, coordinates of node |
[out] | errstat | Error code |
[out] | errmess | Error message |
[in] | orient | Orientation (direction cosine matrix) of node; identity by default |
subroutine modmesh::meshprintinfo | ( | integer, intent(in) | U, |
type(meshtype), intent(in) | M, | ||
integer, intent(in), optional | N | ||
) |
This routine writes mesh information in text form.
It is used for debugging.
[in] | u | fortran output unit |
[in] | m | mesh to be reported on |
[in] | n | Number to print, default is all nodes |
subroutine modmesh::meshsplitelement_2pt | ( | type(meshtype), intent(inout) | Mesh, |
integer(intki), intent(in) | Xelement, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMess, | ||
integer, intent(in) | E1, | ||
integer, intent(in) | P1 | ||
) |
This routine splits a line2 element into two separate elements, using p1 as the new node connecting the two new elements formed from E1.
[in,out] | mesh | Mesh being constructed |
[in] | xelement | See Element Names |
[out] | errstat | Error code |
[out] | errmess | Error message |
[in] | e1 | number of element in Element Table |
[in] | p1 | node number |
subroutine modmesh::meshunpack | ( | type(meshtype), intent(inout) | Mesh, |
real(reki), dimension(:), intent(in), allocatable | ReKiBuf, | ||
real(dbki), dimension(:), intent(in), allocatable | DbKiBuf, | ||
integer(intki), dimension(:), intent(in), allocatable | IntKiBuf, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMess | ||
) |
Given a blank, uncreated mesh and buffers of type INTEGER(IntKi), REAL(ReKi), and REAL(DbKi), unpack the mesh information from the buffers.
This would be done to recreate a mesh after reading in the buffers on a restart of the program. The sense of the name is "unpack the mesh from buffers." The resulting mesh will be returned in the exact state as when the data in the buffers was packed using MeshPack.
subroutine modmesh::meshwrbin | ( | integer, intent(inout) | UnIn, |
type(meshtype), intent(in) | M, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMsg, | ||
character(*), intent(in), optional | FileName | ||
) |
This routine writes mesh information in binary form.
If UnIn is < 0, it gets a new unit number and opens the file, otherwise the file is appended. It is up to the caller of this routine to close the file when it's finished.
[in,out] | unin | fortran output unit |
[in] | m | mesh to be reported on |
[out] | errstat | Indicates whether an error occurred (see NWTC_Library) |
[out] | errmsg | Error message associated with the ErrStat |
[in] | filename | Name of the file to write the output in |
subroutine modmesh::meshwrvtk | ( | real(siki), dimension(3), intent(in) | RefPoint, |
type(meshtype), intent(in) | M, | ||
character(*), intent(in) | FileRootName, | ||
integer(intki), intent(in) | VTKcount, | ||
logical, intent(in) | OutputFieldData, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMsg, | ||
integer(intki), intent(in) | Twidth, | ||
type(meshtype), intent(in), optional | Sib | ||
) |
This routine writes mesh information in VTK format.
see VTK file information format for XML, here: http://www.vtk.org/wp-content/uploads/2015/04/file-formats.pdf
[in] | refpoint | reference location, normally (0,0,0) |
[in] | m | mesh to be written |
[in] | filerootname | Name of the file to write the output in (excluding extension) |
[in] | vtkcount | Indicates number for VTK output file (when 0, the routine will also write reference information) |
[in] | outputfielddata | flag to determine if we want to output field data or just the absolute position of this mesh |
[out] | errstat | Indicates whether an error occurred (see NWTC_Library) |
[out] | errmsg | Error message associated with the ErrStat |
[in] | twidth | Number of digits in the maximum write-out step (used to pad the VTK write-out in the filename with zeros) |
[in] | sib | "functional" Sibling of M that contains translational displacement information (used to place forces at displaced positions) |
We'll write the mesh reference fields on the first timestep only:
subroutine modmesh::meshwrvtk_ln2surface | ( | real(siki), dimension(3), intent(in) | RefPoint, |
type(meshtype), intent(in) | M, | ||
character(*), intent(in) | FileRootName, | ||
integer(intki), intent(in) | VTKcount, | ||
logical, intent(in) | OutputFieldData, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMsg, | ||
integer(intki), intent(in) | Twidth, | ||
integer(intki), intent(in), optional | NumSegments, | ||
real(siki), dimension(:), intent(in), optional | Radius, | ||
real(siki), dimension(:,:,:), intent(in), optional | verts, | ||
type(meshtype), intent(in), optional | Sib | ||
) |
This routine writes line2 mesh surface information in VTK format.
see VTK file information format for XML, here: http://www.vtk.org/wp-content/uploads/2015/04/file-formats.pdf
[in] | refpoint | reference location, normally (0,0,0) |
[in] | m | mesh to be written |
[in] | filerootname | Name of the file to write the output in (excluding extension) |
[in] | vtkcount | Indicates number for VTK output file (when 0, the routine will also write reference information) |
[in] | outputfielddata | flag to determine if we want to output field data or just the absolute position of this mesh |
[in] | twidth | Number of digits in the maximum write-out step (used to pad the VTK write-out in the filename with zeros) |
[in] | numsegments | Number of segments to split the circle into |
[in] | radius | Radius of each node |
[in] | verts | X-Y verticies (2x{NumSegs}xNNodes) of points that define a shape around each node |
[in] | sib | Sibling of M that contains more field information (used only if OutputFieldData is true, to minimize number of files being written) |
[out] | errstat | Indicates whether an error occurred (see NWTC_Library) |
[out] | errmsg | Error message associated with the ErrStat |
subroutine modmesh::meshwrvtk_pointsurface | ( | real(siki), dimension(3), intent(in) | RefPoint, |
type(meshtype), intent(in) | M, | ||
character(*), intent(in) | FileRootName, | ||
integer(intki), intent(in) | VTKcount, | ||
logical, intent(in) | OutputFieldData, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMsg, | ||
integer(intki), intent(in) | Twidth, | ||
integer(intki), intent(in), optional | NumSegments, | ||
real(siki), intent(in), optional | Radius, | ||
real(siki), dimension(:,:), intent(in), optional | verts, | ||
type(meshtype), intent(in), optional | Sib | ||
) |
This routine writes point mesh surfaces information in VTK format.
see VTK file information format for XML, here: http://www.vtk.org/wp-content/uploads/2015/04/file-formats.pdf
[in] | refpoint | reference location, normally (0,0,0) |
[in] | m | mesh to be written |
[in] | filerootname | Name of the file to write the output in (excluding extension) |
[in] | vtkcount | Indicates number for VTK output file (when 0, the routine will also write reference information) |
[in] | outputfielddata | flag to determine if we want to output field data or just the absolute position of this mesh |
[in] | twidth | Number of digits in the maximum write-out timestep (used to pad the VTK write-out in the filename with zeros) |
[in] | numsegments | Number of segments to split the circle into |
[in] | radius | Radius of each node |
[in] | verts | X-Y-Z verticies (3xn) of points that define a volume around each node |
[in] | sib | Sibling of M that contains more field information (used only if OutputFieldData is true, to minimize number of files being written) |
[out] | errstat | Indicates whether an error occurred (see NWTC_Library) |
[out] | errmsg | Error message associated with the ErrStat |
subroutine modmesh::meshwrvtkfields | ( | integer(intki), intent(in) | Un, |
type(meshtype), intent(in) | M, | ||
integer(intki), intent(in) | n | ||
) |
This routine writes mesh field information in VTK format.
see VTK file information format for XML, here: http://www.vtk.org/wp-content/uploads/2015/04/file-formats.pdf
[in] | un | unit number of already-open vtk file in which to write the field information |
[in] | m | mesh to be written |
[in] | n | number of times to write field value for each mesh node (> 1 when added to surface) |
subroutine modmesh::meshwrvtkreference | ( | real(siki), dimension(3), intent(in) | RefPoint, |
type(meshtype), intent(in) | M, | ||
character(*), intent(in) | FileRootName, | ||
integer(intki), intent(out) | ErrStat, | ||
character(*), intent(out) | ErrMsg | ||
) |
This routine writes the reference position and orientations of a mesh in VTK format.
see VTK file information format for XML, here: http://www.vtk.org/wp-content/uploads/2015/04/file-formats.pdf
[in] | refpoint | reference location, normally (0,0,0) |
[in] | m | mesh to be written |
[in] | filerootname | Name of the file to write the output in (excluding extension) |
[out] | errstat | Indicates whether an error occurred (see NWTC_Library) |
[out] | errmsg | Error message associated with the ErrStat |
subroutine modmesh::packloadmesh | ( | type(meshtype), intent(in) | M, |
real(reki), dimension(:), intent(inout) | Ary, | ||
integer(intki), intent(inout) | indx_first | ||
) |
This subroutine returns the operating point values of the mesh fields.
It assumes both force and moment fields are allocated.
[in] | m | Load mesh |
[in,out] | ary | array to pack this mesh into |
[in,out] | indx_first | index into Ary; gives location of next array position to fill |
subroutine modmesh::packloadmesh_dy | ( | type(meshtype), intent(in) | M_p, |
type(meshtype), intent(in) | M_m, | ||
real(r8ki), dimension(:), intent(inout) | dY, | ||
integer(intki), intent(inout) | indx_first | ||
) |
This subroutine computes the differences of two meshes and packs that value into appropriate locations in the dY array.
Do not change this packing without making sure subroutine aerodyn::init_jacobian is consistant with this routine!
[in] | m_p | AD outputs on given mesh at \( u + \Delta u \) (p=plus) |
[in] | m_m | AD outputs on given mesh at \( u - \Delta u \) (m=minus) |
[in,out] | dy | column of dYdu or dYdz \( \frac{\partial Y}{\partial u_i} = \frac{y_p - y_m}{2 \, \Delta u}\) |
[in,out] | indx_first | index into dY array; gives location of next array position to fill |
subroutine modmesh::packloadmesh_names | ( | type(meshtype), intent(in) | M, |
character(*), intent(in) | MeshName, | ||
character(linchanlen), dimension(:), intent(inout) | Names, | ||
integer(intki), intent(inout) | indx_first | ||
) |
This subroutine returns the names of the output rows/columns in the Jacobian matrices.
It assumes both force and moment fields are allocated.
[in] | m | Load mesh |
[in] | meshname | name of mesh |
[in,out] | names | name of row/column of jacobian |
[in,out] | indx_first | index into Names array; gives location of next array position to fill |
subroutine modmesh::packmotionmesh | ( | type(meshtype), intent(in) | M, |
real(reki), dimension(:), intent(inout) | Ary, | ||
integer(intki), intent(inout) | indx_first, | ||
logical, dimension(fieldmask_size), intent(in), optional | FieldMask, | ||
logical, intent(in), optional | UseSmlAngle | ||
) |
This subroutine returns the operating point values of the mesh fields.
It assumes all fields marked by FieldMask are allocated; Some fields may be allocated by the ModMesh module and not used in the linearization procedure, thus I am not using the check if they are allocated to determine if they should be included.
[in] | m | Motion mesh |
[in,out] | ary | array to pack this mesh into |
[in,out] | indx_first | index into Ary; gives location of next array position to fill |
[in] | fieldmask | flags to determine if this field is part of the packing |
[in] | usesmlangle | flag to determine if the orientation should be packed as a DCM or a log map |
subroutine modmesh::packmotionmesh_dy | ( | type(meshtype), intent(in) | M_p, |
type(meshtype), intent(in) | M_m, | ||
real(r8ki), dimension(:), intent(inout) | dY, | ||
integer(intki), intent(inout) | indx_first, | ||
logical, dimension(fieldmask_size), intent(in), optional | FieldMask, | ||
logical, intent(in), optional | UseSmlAngle | ||
) |
This subroutine computes the differences of two meshes and packs that value into appropriate locations in the dY array.
[in] | m_p | ED outputs on given mesh at \( u + \Delta u \) (p=plus) |
[in] | m_m | ED outputs on given mesh at \( u - \Delta u \) (m=minus) |
[in,out] | dy | column of dYdu \( \frac{\partial Y}{\partial u_i} = \frac{y_p - y_m}{2 \, \Delta u}\) |
[in,out] | indx_first | index into dY array; gives location of next array position to fill |
[in] | fieldmask | flags to determine if this field is part of the packing |
[in] | usesmlangle | flag to determine if the orientation should be packed as a DCM or a log map |
subroutine modmesh::packmotionmesh_names | ( | type(meshtype), intent(in) | M, |
character(*), intent(in) | MeshName, | ||
character(linchanlen), dimension(:), intent(inout) | Names, | ||
integer(intki), intent(inout) | indx_first, | ||
logical, dimension(fieldmask_size), intent(in), optional | FieldMask | ||
) |
This subroutine returns the names of rows/columns of motion meshes in the Jacobian matrices.
It assumes all fields marked by FieldMask are allocated; Some fields may be allocated by the ModMesh module and not used in the linearization procedure, thus I am not using the check if they are allocated to determine if they should be included.
[in] | m | Motion mesh |
[in] | meshname | name of mesh |
[in,out] | names | name of row/column of jacobian |
[in,out] | indx_first | index into Names array; gives location of next array position to fill |
[in] | fieldmask | flags to determine if this field is part of the packing |