OOFEM
2.4
OOFEM.org - Object Oriented Finite Element Solver
|
Mesh generator for the PFEM problem, using Bowyer-Watson algorithm of the Delaunay triangulation of a set of nodes (PFEMParticle) creating TR1_2D_PFEM elements. More...
#include <delaunaytriangulator.h>
Public Member Functions | |
DelaunayTriangulator (Domain *d, double setAlpha) | |
Constructor. More... | |
~DelaunayTriangulator () | |
Destructor. More... | |
void | generateMesh () |
Main call. More... | |
Protected Attributes | |
Domain * | domain |
Domain of the PFEM problem containing nodes to be triangulated. More... | |
double | alphaValue |
int | nnode |
Timer | meshingTimer |
Measures overall time of triangulation procedure. More... | |
Timer | searchingTimer |
Measures time needed by searching for non-delaunay triangles. More... | |
Timer | polygonTimer |
Measures time needed for identifying polygon to be retriangulated. More... | |
Timer | creativeTimer |
Measures time needed for creating new Delaunay triangles. More... | |
std::list< DelaunayTriangle * > | generalTriangleList |
Contains all triangles (even not valid) More... | |
std::list< DelaunayTriangle * >::iterator | genIT |
std::list< AlphaEdge2D * > | alphaShapeEdgeList |
Contains resulting alpha-shape in form of a list. More... | |
std::list< AlphaEdge2D * > | edgeList |
contains all edges of the triangulation More... | |
std::list< AlphaEdge2D * >::iterator | elIT |
OctreeSpatialLocalizerT< DelaunayTriangle * > | triangleOctree |
Octree with Delaunay triangles allowing fast search. More... | |
Private Member Functions | |
void | addUniqueEdgeToPolygon (Edge2D *edge, std::list< Edge2D > &polygon) |
Edge is added to the polygon only if it's not contained. Otherwise both are removed (edge shared by two non-Delaunay triangles). More... | |
void | buildInitialBBXMesh (InsertTriangleBasedOnCircumcircle &tInsert) |
Identifies the bounding box of pfemparticles and creates initial triangulation consisting of 2 triangles conecting bounding box nodes. More... | |
void | writeMesh () |
Writes the mesh into the domain by creating new tr1_2d_pfem elements and prescribes zero-pressure boundary condition on alpha-shape nodes. More... | |
void | computeAlphaComplex () |
Reads the triangulation and fills tha edgeList container with alpha-shape edges and set their bounds. More... | |
AlphaEdge2D * | giveBackEdgeIfAlreadyContainedInList (AlphaEdge2D *alphaEdge) |
Fills the edgeList with unique alphaEdges. More... | |
void | giveAlphaShape () |
Iterates through the edgeList container and compares alpha-value with alphaEdge bounds. Alpha shape is stored in the alphaShapeEdgeList. More... | |
void | initializeTimers () |
Initializes Timers and. More... | |
void | findNonDelaunayTriangles (int insertedNode, InsertTriangleBasedOnCircumcircle &tInsert, std::list< Edge2D > &polygon) |
Looks for non-Delaunay triangles in octree and creates a polygon. More... | |
void | meshPolygon (int insertedNode, InsertTriangleBasedOnCircumcircle &tInsert, std::list< Edge2D > &polygon) |
Retriangulates the polygon. More... | |
void | giveTimeReport () |
Prints the time report. More... | |
void | cleanUpTriangleList () |
Iterates through generalTringleList und removes non-valid ones or those containing bounding box nodes. More... | |
void | computeBBXBasedOnNodeData (BoundingBox &BBX) |
Calculates the bounding box base on the domain's nodes. More... | |
Mesh generator for the PFEM problem, using Bowyer-Watson algorithm of the Delaunay triangulation of a set of nodes (PFEMParticle) creating TR1_2D_PFEM elements.
Definition at line 65 of file delaunaytriangulator.h.
oofem::DelaunayTriangulator::DelaunayTriangulator | ( | Domain * | d, |
double | setAlpha | ||
) |
Constructor.
Definition at line 50 of file delaunaytriangulator.C.
References alphaValue, domain, oofem::Domain::giveNumberOfDofManagers(), and nnode.
oofem::DelaunayTriangulator::~DelaunayTriangulator | ( | ) |
Destructor.
Definition at line 61 of file delaunaytriangulator.C.
References edgeList, elIT, generalTriangleList, and genIT.
|
private |
Edge is added to the polygon only if it's not contained. Otherwise both are removed (edge shared by two non-Delaunay triangles).
Definition at line 76 of file delaunaytriangulator.C.
Referenced by findNonDelaunayTriangles().
|
private |
Identifies the bounding box of pfemparticles and creates initial triangulation consisting of 2 triangles conecting bounding box nodes.
Definition at line 500 of file delaunaytriangulator.C.
References oofem::FloatArray::at(), computeBBXBasedOnNodeData(), domain, generalTriangleList, oofem::BoundingBox::giveOrigin(), oofem::BoundingBox::giveSize(), nnode, oofem::Domain::resizeDofManagers(), oofem::Node::setCoordinates(), oofem::Domain::setDofManager(), and triangleOctree.
Referenced by generateMesh().
|
private |
Iterates through generalTringleList und removes non-valid ones or those containing bounding box nodes.
Definition at line 481 of file delaunaytriangulator.C.
References generalTriangleList, genIT, and nnode.
Referenced by generateMesh().
|
private |
Reads the triangulation and fills tha edgeList container with alpha-shape edges and set their bounds.
Definition at line 226 of file delaunaytriangulator.C.
References edgeList, generalTriangleList, genIT, giveBackEdgeIfAlreadyContainedInList(), oofem::AlphaEdge2D::giveInnerAlphaBound(), oofem::AlphaEdge2D::giveOuterAlphaBound(), oofem::min(), oofem::AlphaEdge2D::setHullFlag(), oofem::AlphaEdge2D::setInnerAlphaBound(), oofem::AlphaEdge2D::setOuterAlphaBound(), and oofem::AlphaEdge2D::setSharing().
Referenced by generateMesh().
|
private |
Calculates the bounding box base on the domain's nodes.
Definition at line 549 of file delaunaytriangulator.C.
References oofem::FloatArray::at(), domain, oofem::Domain::giveDofManager(), oofem::Domain::giveNumberOfDofManagers(), oofem::FloatArray::giveSize(), oofem::max(), oofem::min(), nnode, oofem::BoundingBox::setMask(), oofem::BoundingBox::setOrigin(), and oofem::BoundingBox::setSize().
Referenced by buildInitialBBXMesh().
|
private |
Looks for non-Delaunay triangles in octree and creates a polygon.
Definition at line 410 of file delaunaytriangulator.C.
References addUniqueEdgeToPolygon(), domain, oofem::Domain::giveDofManager(), oofem::DelaunayTriangle::giveNode(), oofem::Timer::pauseTimer(), polygonTimer, oofem::Timer::resumeTimer(), searchingTimer, oofem::DelaunayTriangle::setValidFlag(), and triangleOctree.
Referenced by generateMesh().
void oofem::DelaunayTriangulator::generateMesh | ( | ) |
Main call.
Definition at line 100 of file delaunaytriangulator.C.
References alphaValue, buildInitialBBXMesh(), cleanUpTriangleList(), computeAlphaComplex(), domain, findNonDelaunayTriangles(), generalTriangleList, giveAlphaShape(), oofem::Domain::giveDofManager(), initializeTimers(), oofem::PFEMParticle::isActive(), meshPolygon(), nnode, oofem::Domain::resizeDofManagers(), triangleOctree, VERBOSE_PRINT0, and writeMesh().
Referenced by oofem::PFEM::preInitializeNextStep().
|
private |
Iterates through the edgeList container and compares alpha-value with alphaEdge bounds. Alpha shape is stored in the alphaShapeEdgeList.
Definition at line 359 of file delaunaytriangulator.C.
References alphaShapeEdgeList, alphaValue, edgeList, and elIT.
Referenced by generateMesh().
|
private |
Fills the edgeList with unique alphaEdges.
If an edge is already contained a pointer to it is returned and inserted edge is removed.
Definition at line 341 of file delaunaytriangulator.C.
References edgeList, and elIT.
Referenced by computeAlphaComplex().
|
private |
Prints the time report.
Definition at line 460 of file delaunaytriangulator.C.
References creativeTimer, oofem::Timer::getUtime(), meshingTimer, polygonTimer, searchingTimer, and oofem::Timer::stopTimer().
|
private |
Initializes Timers and.
Definition at line 398 of file delaunaytriangulator.C.
References creativeTimer, meshingTimer, oofem::Timer::pauseTimer(), polygonTimer, searchingTimer, and oofem::Timer::startTimer().
Referenced by generateMesh().
|
private |
Retriangulates the polygon.
Definition at line 443 of file delaunaytriangulator.C.
References creativeTimer, domain, generalTriangleList, oofem::Timer::pauseTimer(), oofem::Timer::resumeTimer(), and triangleOctree.
Referenced by generateMesh().
|
private |
Writes the mesh into the domain by creating new tr1_2d_pfem elements and prescribes zero-pressure boundary condition on alpha-shape nodes.
Definition at line 141 of file delaunaytriangulator.C.
References alphaShapeEdgeList, alphaValue, domain, elIT, generalTriangleList, genIT, oofem::PFEM::giveAssociatedCrossSectionNumber(), oofem::PFEM::giveAssociatedMaterialNumber(), oofem::PFEM::giveAssociatedPressureBC(), oofem::Domain::giveDofManager(), oofem::DofManager::giveDofWithID(), oofem::Domain::giveEngngModel(), oofem::Domain::giveNode(), oofem::Domain::giveNumberOfDofManagers(), oofem::Domain::resizeElements(), oofem::Dof::setBcId(), and oofem::Domain::setElement().
Referenced by generateMesh().
|
protected |
Contains resulting alpha-shape in form of a list.
Definition at line 91 of file delaunaytriangulator.h.
Referenced by giveAlphaShape(), and writeMesh().
|
protected |
Definition at line 71 of file delaunaytriangulator.h.
Referenced by DelaunayTriangulator(), generateMesh(), giveAlphaShape(), and writeMesh().
|
protected |
Measures time needed for creating new Delaunay triangles.
Definition at line 85 of file delaunaytriangulator.h.
Referenced by giveTimeReport(), initializeTimers(), and meshPolygon().
|
protected |
Domain of the PFEM problem containing nodes to be triangulated.
Definition at line 69 of file delaunaytriangulator.h.
Referenced by buildInitialBBXMesh(), computeBBXBasedOnNodeData(), DelaunayTriangulator(), findNonDelaunayTriangles(), generateMesh(), meshPolygon(), and writeMesh().
|
protected |
contains all edges of the triangulation
Definition at line 94 of file delaunaytriangulator.h.
Referenced by computeAlphaComplex(), giveAlphaShape(), giveBackEdgeIfAlreadyContainedInList(), and ~DelaunayTriangulator().
|
protected |
Definition at line 95 of file delaunaytriangulator.h.
Referenced by giveAlphaShape(), giveBackEdgeIfAlreadyContainedInList(), writeMesh(), and ~DelaunayTriangulator().
|
protected |
Contains all triangles (even not valid)
Definition at line 87 of file delaunaytriangulator.h.
Referenced by buildInitialBBXMesh(), cleanUpTriangleList(), computeAlphaComplex(), generateMesh(), meshPolygon(), writeMesh(), and ~DelaunayTriangulator().
|
protected |
Definition at line 88 of file delaunaytriangulator.h.
Referenced by cleanUpTriangleList(), computeAlphaComplex(), writeMesh(), and ~DelaunayTriangulator().
|
protected |
Measures overall time of triangulation procedure.
Definition at line 79 of file delaunaytriangulator.h.
Referenced by giveTimeReport(), and initializeTimers().
|
protected |
Definition at line 73 of file delaunaytriangulator.h.
Referenced by buildInitialBBXMesh(), cleanUpTriangleList(), computeBBXBasedOnNodeData(), DelaunayTriangulator(), and generateMesh().
|
protected |
Measures time needed for identifying polygon to be retriangulated.
Definition at line 83 of file delaunaytriangulator.h.
Referenced by findNonDelaunayTriangles(), giveTimeReport(), and initializeTimers().
|
protected |
Measures time needed by searching for non-delaunay triangles.
Definition at line 81 of file delaunaytriangulator.h.
Referenced by findNonDelaunayTriangles(), giveTimeReport(), and initializeTimers().
|
protected |
Octree with Delaunay triangles allowing fast search.
Definition at line 98 of file delaunaytriangulator.h.
Referenced by buildInitialBBXMesh(), findNonDelaunayTriangles(), generateMesh(), and meshPolygon().