MOAB: Mesh Oriented datABase  (version 5.3.1)
MBMesquite::DomainClassifier Class Reference

Assign subsets of a mesh do different domains. More...

#include <DomainClassifier.hpp>

+ Inheritance diagram for MBMesquite::DomainClassifier:
+ Collaboration diagram for MBMesquite::DomainClassifier:

Classes

struct  DomainBlock
struct  DomainSet

Public Member Functions

MESQUITE_EXPORT void test_valid_classification (Mesh *mesh, MsqError &err)
 Check that classification maps to B-Rep topology.
MESQUITE_EXPORT DomainClassifier ()
virtual MESQUITE_EXPORT ~DomainClassifier ()
virtual MESQUITE_EXPORT void snap_to (Mesh::VertexHandle entity_handle, Vector3D &coordinate) const
virtual MESQUITE_EXPORT void vertex_normal_at (Mesh::VertexHandle entity_handle, Vector3D &coordinate) const
virtual MESQUITE_EXPORT void element_normal_at (Mesh::ElementHandle entity_handle, Vector3D &coordinate) const
virtual MESQUITE_EXPORT void vertex_normal_at (const Mesh::VertexHandle *handles, Vector3D coordinates[], unsigned count, MsqError &err) const
 evaluate surface normals
virtual MESQUITE_EXPORT void closest_point (Mesh::VertexHandle handle, const Vector3D &position, Vector3D &closest, Vector3D &normal, MsqError &err) const
 evaluate closest point and normal
virtual MESQUITE_EXPORT void domain_DoF (const Mesh::VertexHandle *handle_array, unsigned short *dof_array, size_t num_handles, MsqError &err) const
MESQUITE_EXPORT void clear ()
 Clear all data, including MeshDomain list.
MESQUITE_EXPORT MeshDomainfind_vertex_domain (Mesh::VertexHandle vertex) const
MESQUITE_EXPORT MeshDomainfind_element_domain (Mesh::ElementHandle element) const
MESQUITE_EXPORT void delete_sub_domains (bool yesno)
MESQUITE_EXPORT void delete_all_sub_domains ()

Static Public Member Functions

static MESQUITE_EXPORT void classify_by_tag (DomainClassifier &result, Mesh *mesh, const char *tag_name, MeshDomain **domain_array, const int *id_array, unsigned array_length, MsqError &err)
 Classify mesh entities using tag values.
static MESQUITE_EXPORT void classify_skin_geometrically (DomainClassifier &result, Mesh *mesh, double tolerance, MeshDomain **domain_array, const int *dimension_array, unsigned array_length, MsqError &err)
 Skin mesh and classify skin entities geometrically.
static MESQUITE_EXPORT void classify_geometrically (DomainClassifier &result, Mesh *mesh, double tolerance, MeshDomain **domain_array, const int *dimension_array, unsigned array_length, MsqError &err)
 Classify all mesh entities geometrically.
static MESQUITE_EXPORT void classify_by_handle (DomainClassifier &result, Mesh *mesh, DomainSet *domain_set_array, unsigned array_length, MsqError &err)
 Specify classification explicitly for each entity.

Static Private Member Functions

static MeshDomainfind_domain (Mesh::EntityHandle handle, const std::vector< DomainBlock > &list)

Private Attributes

bool deleteSubDomains
std::vector< DomainBlockvertexList
std::vector< DomainBlockelementList

Detailed Description

Assign subsets of a mesh do different domains.

Provide classification of mesh entities to domain(s). For example, given a domain defined by multiple surfaces and curves of intersection, associate the appropriate elements and vertices to the appropriate geometric entities.

Definition at line 48 of file DomainClassifier.hpp.


Constructor & Destructor Documentation


Member Function Documentation

void MBMesquite::DomainClassifier::classify_by_handle ( DomainClassifier result,
Mesh mesh,
DomainSet domain_set_array,
unsigned  array_length,
MsqError err 
) [static]

Specify classification explicitly for each entity.

For each mesh element and vertex that is classified to a domain, specify that classification explicitly.

Parameters:
resultDomainClassifier to popupate
meshThe MBMesquite::Mesh instance
domain_set_arrayArray of DomainClassifier::DomainSet structs specifying for each MeshDomain: the domain instance and the elements and vertices associated with the domain.
array_lengthLength of 'domain_set_array'

Definition at line 581 of file DomainClassifier.cpp.

References clear(), MBMesquite::DomainClassifier::DomainSet::domain, MBMesquite::DomainClassifier::DomainBlock::domain, elementList, MBMesquite::DomainClassifier::DomainSet::elements, MBMesquite::DomainClassifier::DomainBlock::firstHandle, MBMesquite::Mesh::get_all_elements(), MBMesquite::Mesh::get_all_vertices(), MBMesquite::DomainClassifier::DomainBlock::lastHandle, MSQ_ERRRTN, vertexList, and MBMesquite::DomainClassifier::DomainSet::vertices.

Referenced by classify_by_tag(), classify_geometrically(), classify_skin_geometrically(), get_planar_example(), DomainClassifierTest::test_classify_by_handle(), HigherOrderTest::test_tri_open_domain(), and DomainClassifierTest::test_valid_classification().

{
    result.clear();

    // get all vertices and elements
    std::vector< Mesh::VertexHandle > vertices;
    std::vector< Mesh::ElementHandle > elements;
    mesh->get_all_vertices( vertices, err );MSQ_ERRRTN( err );
    mesh->get_all_elements( elements, err );MSQ_ERRRTN( err );

    // sort all arrays so we can merge
    std::sort( vertices.begin(), vertices.end() );
    std::sort( elements.begin(), elements.end() );
    for( unsigned i = 0; i < array_length; ++i )
    {
        std::sort( sets[i].vertices.begin(), sets[i].vertices.end() );
        std::sort( sets[i].elements.begin(), sets[i].elements.end() );
    }

    // build vertex block list
    std::vector< size_t > indices( array_length, 0 );
    size_t idx = 0;
    while( idx < vertices.size() )
    {
        DomainBlock block;
        block.firstHandle = vertices[idx];
        block.lastHandle  = vertices[idx];

        // find domain
        unsigned dom = 0;
        for( dom = 0; dom < array_length; ++dom )
        {
            const size_t i = indices[dom];
            if( i < sets[dom].vertices.size() && sets[dom].vertices[i] == vertices[idx] )
            {
                ++indices[dom];
                break;
            }
        }
        ++idx;
        if( dom == indices.size() )  // no domain
            continue;

        block.domain = sets[dom].domain;
        while( indices[dom] < sets[dom].vertices.size() && sets[dom].vertices[indices[dom]] == vertices[idx] )
        {
            block.lastHandle = vertices[idx];
            ++idx;
            ++indices[dom];
        }
        result.vertexList.push_back( block );
    }

    // build vertex block list
    indices.clear();
    indices.resize( array_length, 0 );
    idx = 0;
    while( idx < elements.size() )
    {
        DomainBlock block;
        block.firstHandle = elements[idx];
        block.lastHandle  = elements[idx];

        // find domain
        unsigned dom = 0;
        for( dom = 0; dom < array_length; ++dom )
        {
            const size_t i = indices[dom];
            if( i < sets[dom].elements.size() && sets[dom].elements[i] == elements[idx] )
            {
                ++indices[dom];
                break;
            }
        }
        ++idx;
        if( dom == indices.size() )  // no domain
            continue;

        block.domain = sets[dom].domain;
        while( indices[dom] < sets[dom].elements.size() && sets[dom].elements[indices[dom]] == elements[idx] )
        {
            block.lastHandle = elements[idx];
            ++idx;
            ++indices[dom];
        }
        result.elementList.push_back( block );
    }
    /*
      printf("\n");
      for (unsigned i = 0; i < result.vertexList.size(); ++i)
        printf("v %2lu-%2lu @ %p\n", (unsigned long)result.vertexList[i].firstHandle,
                                   (unsigned long)result.vertexList[i].lastHandle,
                                   result.vertexList[i].domain );
      for (unsigned i = 0; i < result.elementList.size(); ++i)
        printf("e %2lu-%2lu @ %p\n", (unsigned long)result.elementList[i].firstHandle,
                                   (unsigned long)result.elementList[i].lastHandle,
                                   result.elementList[i].domain );
      printf("\n");
    */
}
void MBMesquite::DomainClassifier::classify_by_tag ( DomainClassifier result,
Mesh mesh,
const char *  tag_name,
MeshDomain **  domain_array,
const int *  id_array,
unsigned  array_length,
MsqError err 
) [static]

Classify mesh entities using tag values.

Given a list of MeshDomain instances, a unique integer ID for each MeshDomain, and a tag name: classify mesh entities by matching the tag value on each mesh entity to the list of MeshDomain IDs.

Parameters:
resultDomainClassifier to popupate
meshThe MBMesquite::Mesh instance
tag_nameTag containing integer domain ID for each mesh entity.
domain_arrayArray of MeshDomain instances.
id_arrayArray of integer MeshDomain IDs.
array_lengthLength of 'domain_array' and 'id_array'

Definition at line 522 of file DomainClassifier.cpp.

References MBMesquite::arrptr(), classify_by_handle(), MBMesquite::MsqError::clear(), MBMesquite::Mesh::get_all_elements(), MBMesquite::Mesh::get_all_vertices(), MBMesquite::Mesh::INT, MSQ_ERRRTN, MSQ_SETERR, size, MBMesquite::Mesh::tag_get(), MBMesquite::Mesh::tag_get_element_data(), MBMesquite::Mesh::tag_get_vertex_data(), MBMesquite::MsqError::TAG_NOT_FOUND, and MBMesquite::Mesh::tag_properties().

Referenced by DomainClassifierTest::test_classify_by_tag().

{
    TagHandle tag = mesh->tag_get( tag_name, err );MSQ_ERRRTN( err );

    std::string name2;
    Mesh::TagType type;
    unsigned size;
    mesh->tag_properties( tag, name2, type, size, err );MSQ_ERRRTN( err );
    if( type != Mesh::INT || size != 1 )
    {
        MSQ_SETERR( err )
        ( MsqError::TAG_NOT_FOUND, "Tag does not contain single integer value: %s\n", tag_name );
        return;
    }

    std::vector< DomainSet > sets( array_length );
    std::map< int, DomainSet* > idmap;
    for( unsigned i = 0; i < array_length; ++i )
    {
        sets[i].domain     = domain_array[i];
        idmap[id_array[i]] = &sets[i];
    }

    std::vector< Mesh::VertexHandle > vertices;
    std::vector< Mesh::ElementHandle > elements;
    mesh->get_all_vertices( vertices, err );MSQ_ERRRTN( err );
    mesh->get_all_elements( elements, err );MSQ_ERRRTN( err );
    if( vertices.empty() ) return;

    std::map< int, DomainSet* >::const_iterator iter;
    std::vector< int > ids( vertices.size() );
    mesh->tag_get_vertex_data( tag, vertices.size(), arrptr( vertices ), arrptr( ids ), err );MSQ_ERRRTN( err );
    for( size_t j = 0; j < vertices.size(); ++j )
    {
        iter = idmap.find( ids[j] );
        if( iter != idmap.end() ) iter->second->vertices.push_back( vertices[j] );
    }

    ids.clear();
    ids.resize( elements.size() );
    if( !elements.empty() )
    {
        mesh->tag_get_element_data( tag, elements.size(), arrptr( elements ), arrptr( ids ), err );
        if( err ) { err.clear(); }
        else
        {
            for( size_t k = 0; k < elements.size(); ++k )
            {
                iter = idmap.find( ids[k] );
                if( iter != idmap.end() ) iter->second->elements.push_back( elements[k] );
            }
        }
    }

    if( !sets.empty() ) classify_by_handle( result, mesh, arrptr( sets ), sets.size(), err );
}
void MBMesquite::DomainClassifier::classify_geometrically ( DomainClassifier result,
Mesh mesh,
double  tolerance,
MeshDomain **  domain_array,
const int *  dimension_array,
unsigned  array_length,
MsqError err 
) [static]

Classify all mesh entities geometrically.

Classify entities by distance from vertex coordiantes to mesh domain instances.

NOTE: Mesquite's limited MeshDomain callback interface is not sufficient for robust geometric classification. If the mesh is not sufficiently refined, this method may produce invalid results.

NOTE: This method is likely to fail for many domains for volume meshes. Use classify_skin_geometrically for volume meshes.

Parameters:
resultDomainClassifier to popupate
meshThe MBMesquite::Mesh instance
toleranceMaximum distance a vertex may deviate from its domain.
domain_arrayArray of MeshDomain instances.
dimension_arrayTopological dimensiono of each MeshDomain
array_lengthLength of 'domain_array' and 'dimension_array'

Definition at line 441 of file DomainClassifier.cpp.

References MBMesquite::arrptr(), classify_by_handle(), MBMesquite::dimension_sort_domains(), domains, MBMesquite::Mesh::elements_get_topologies(), MBMesquite::geom_classify_elements(), MBMesquite::geom_classify_vertices(), MBMesquite::Mesh::get_all_elements(), MBMesquite::Mesh::get_all_vertices(), MSQ_ERRRTN, tolerance, and MBMesquite::vert_classify_elements().

Referenced by PatchDataTest::get_quad8_mesh_and_domain(), process_domain_args(), and DomainClassifierTest::test_classify_by_geometry().

{
    if( !array_length ) return;

    std::vector< DomainSet > domains;
    int dim_indices[4];
    dimension_sort_domains( domain_array, dimension_array, array_length, domains, dim_indices, err );MSQ_ERRRTN( err );

    // classify vertices by position
    std::vector< Mesh::VertexHandle > vertices;
    mesh->get_all_vertices( vertices, err );MSQ_ERRRTN( err );
    geom_classify_vertices( mesh, arrptr( domains ), dim_indices[3], vertices, tolerance, err );MSQ_ERRRTN( err );

    // get elements and types
    std::vector< Mesh::ElementHandle > elems;
    mesh->get_all_elements( elems, err );MSQ_ERRRTN( err );
    if( elems.empty() ) return;
    std::vector< EntityTopology > types( elems.size() );
    mesh->elements_get_topologies( arrptr( elems ), arrptr( types ), elems.size(), err );MSQ_ERRRTN( err );

    // get rid of elements w/ dimension other than 2
    size_t w = 0;
    for( size_t r = 0; r < elems.size(); ++r )
    {
        if( TopologyInfo::dimension( types[r] ) == 2 )
        {
            elems[w] = elems[r];
            ++w;
        }
    }
    elems.resize( w );

    // classify elements using vertex classification
    std::vector< Mesh::ElementHandle > unknown;
    vert_classify_elements( mesh, arrptr( domains ), domains.size(), dim_indices, elems, unknown, err );MSQ_ERRRTN( err );

    // classify unknown elements geometrically
    elems.swap( unknown );
    unknown.clear();
    geom_classify_elements( mesh, arrptr( domains ), domains.size(), dim_indices, elems, unknown, tolerance, err );

    classify_by_handle( result, mesh, arrptr( domains ), domains.size(), err );
}
void MBMesquite::DomainClassifier::classify_skin_geometrically ( DomainClassifier result,
Mesh mesh,
double  tolerance,
MeshDomain **  domain_array,
const int *  dimension_array,
unsigned  array_length,
MsqError err 
) [static]

Skin mesh and classify skin entities geometrically.

Calculate the boundary of the mesh, and classify boundary entities by calculating the distance from the mesh domain to each vertex.

NOTE: Mesquite's limited MeshDomain callback interface is not sufficient for robust geometric classification. If the mesh is not sufficiently refined, this method may produce invalid results.

NOTE: This method should not be used for surface meshes. Everything in a surface mesh should be classified to a domain. not just the skin. Use classify_geometrically instead.

Parameters:
resultDomainClassifier to popupate
meshThe MBMesquite::Mesh instance
toleranceMaximum distance a vertex may deviate from its domain.
domain_arrayArray of MeshDomain instances.
dimension_arrayTopological dimensiono of each MeshDomain
array_lengthLength of 'domain_array' and 'dimension_array'

Definition at line 487 of file DomainClassifier.cpp.

References MBMesquite::arrptr(), classify_by_handle(), MBMesquite::dimension_sort_domains(), domains, MBMesquite::find_skin(), MBMesquite::geom_classify_elements(), MBMesquite::geom_classify_vertices(), MBMesquite::MsqError::INVALID_MESH, MSQ_ERRRTN, MSQ_SETERR, tolerance, and MBMesquite::vert_classify_elements().

Referenced by get_cut_cube_example(), get_hex_3d_part_example(), get_sphere_cube_example(), get_sphere_cylinder_example(), SlaveBoundaryVerticesTest::make_mesh(), process_domain_args(), and DomainClassifierTest::test_classify_skin().

{
    if( !array_length ) return;

    std::vector< DomainSet > domains;
    int dim_indices[4];
    dimension_sort_domains( domain_array, dimension_array, array_length, domains, dim_indices, err );MSQ_ERRRTN( err );

    std::vector< Mesh::VertexHandle > vertices;
    std::vector< Mesh::ElementHandle > elements;
    find_skin( mesh, vertices, elements, err );MSQ_ERRRTN( err );

    // classify vertices by position
    geom_classify_vertices( mesh, arrptr( domains ), dim_indices[3], vertices, tolerance, err );MSQ_ERRRTN( err );

    // classify elements using vertex classification
    std::vector< Mesh::ElementHandle > unknown;
    vert_classify_elements( mesh, arrptr( domains ), domains.size(), dim_indices, elements, unknown, err );MSQ_ERRRTN( err );

    // classify unknown elements geometrically
    elements.swap( unknown );
    unknown.clear();
    geom_classify_elements( mesh, arrptr( domains ), domains.size(), dim_indices, elements, unknown, tolerance, err );

    if( !unknown.empty() )
    {
        MSQ_SETERR( err )( "Count not classify all skin elements", MsqError::INVALID_MESH );
        return;
    }

    classify_by_handle( result, mesh, arrptr( domains ), domains.size(), err );
}

Clear all data, including MeshDomain list.

Definition at line 190 of file DomainClassifier.hpp.

References elementList, and vertexList.

Referenced by classify_by_handle(), and delete_all_sub_domains().

    {
        vertexList.clear();
        elementList.clear();
    }
void MBMesquite::DomainClassifier::closest_point ( Mesh::VertexHandle  handle,
const Vector3D position,
Vector3D closest,
Vector3D normal,
MsqError err 
) const [virtual]

evaluate closest point and normal

Given a position in space, return the closest position in the domain and the domain normal at that point.

Parameters:
entity_handleEvaluate the subset of the domain contianing this entity
positionInput position for which to evaluate
closestClosest position in the domain.
normalDomain normal at the location of 'closest'

Implements MBMesquite::MeshDomain.

Definition at line 1159 of file DomainClassifier.cpp.

References find_vertex_domain().

{
    if( const MeshDomain* dom = find_vertex_domain( handle ) )
        dom->closest_point( handle, position, closest, normal, err );
}

Definition at line 1186 of file DomainClassifier.cpp.

References clear(), domains, elementList, and vertexList.

Referenced by ~DomainClassifier().

{
    // get unique list of domains
    std::set< MeshDomain* > domains;
    std::vector< DomainBlock >::iterator i;
    for( i = vertexList.begin(); i != vertexList.end(); ++i )
        domains.insert( i->domain );
    for( i = elementList.begin(); i != elementList.end(); ++i )
        domains.insert( i->domain );
    std::set< MeshDomain* >::iterator j;
    for( j = domains.begin(); j != domains.end(); ++j )
        delete *j;
    clear();
}
void MBMesquite::DomainClassifier::domain_DoF ( const Mesh::VertexHandle handle_array,
unsigned short *  dof_array,
size_t  num_handles,
MsqError err 
) const [virtual]

Definition at line 1166 of file DomainClassifier.cpp.

References MBMesquite::MeshDomain::domain_DoF(), find_vertex_domain(), and MSQ_ERRRTN.

Referenced by SlaveBoundaryVerticesTest::test_slaved_common().

{
    for( size_t i = 0; i < num_handles; ++i )
    {
        const MeshDomain* dom = find_vertex_domain( handles[i] );
        if( !dom )
            dof_array[i] = 3;
        else
        {
            dom->domain_DoF( handles + i, dof_array + i, 1, err );MSQ_ERRRTN( err );
        }
    }
}
void MBMesquite::DomainClassifier::element_normal_at ( Mesh::ElementHandle  entity_handle,
Vector3D coordinate 
) const [virtual]

Implements MBMesquite::MeshDomain.

Definition at line 1138 of file DomainClassifier.cpp.

References find_element_domain().

{
    if( const MeshDomain* dom = find_element_domain( entity_handle ) )
        dom->element_normal_at( entity_handle, coordinate );
}
MeshDomain * MBMesquite::DomainClassifier::find_domain ( Mesh::EntityHandle  handle,
const std::vector< DomainBlock > &  list 
) [static, private]

Definition at line 1119 of file DomainClassifier.cpp.

Referenced by find_element_domain(), and find_vertex_domain().

{
    std::vector< DomainClassifier::DomainBlock >::const_iterator i;
    i = std::lower_bound( list.begin(), list.end(), handle );
    return ( i != list.end() && i->firstHandle <= handle ) ? i->domain : NULL;
}
void MBMesquite::DomainClassifier::snap_to ( Mesh::VertexHandle  entity_handle,
Vector3D coordinate 
) const [virtual]

Modifies "coordinate" so that it lies on the domain to which "entity_handle" is constrained. The handle determines the domain. The coordinate is the proposed new position on that domain.

Implements MBMesquite::MeshDomain.

Definition at line 1127 of file DomainClassifier.cpp.

References find_vertex_domain().

{
    if( const MeshDomain* dom = find_vertex_domain( entity_handle ) ) dom->snap_to( entity_handle, coordinate );
}

Check that classification maps to B-Rep topology.

Definition at line 782 of file DomainClassifier.cpp.

References MBMesquite::arrptr(), dim, dof(), domains, elementList, MBMesquite::Mesh::elements_get_attached_vertices(), MBMesquite::Mesh::elements_get_topologies(), MBMesquite::faces, MBMesquite::Mesh::get_all_elements(), MBMesquite::Mesh::get_all_vertices(), MBMesquite::MsqError::INVALID_STATE, MSQ_ERRRTN, MSQ_SETERR, n, MBMesquite::next_vertex(), vertexList, MBMesquite::Mesh::vertices_get_attached_elements(), and vtx().

Referenced by DomainClassifierTest::test_valid_classification().

{
    size_t i;

    // Get all mesh entities
    std::vector< Mesh::VertexHandle > vertices;
    std::vector< Mesh::ElementHandle > elements;
    mesh->get_all_vertices( vertices, err );MSQ_ERRRTN( err );
    mesh->get_all_elements( elements, err );MSQ_ERRRTN( err );
    std::sort( vertices.begin(), vertices.end() );
    std::sort( elements.begin(), elements.end() );

    // Get contents of each domain.
    std::map< MeshDomain*, int >::iterator iter;
    std::map< MeshDomain*, int > idxmap;
    for( i = 0; i < vertexList.size(); ++i )
        idxmap[vertexList[i].domain] = 0;
    for( i = 0; i < elementList.size(); ++i )
        idxmap[elementList[i].domain] = 0;
    std::vector< DomainSet > domains( idxmap.size() );
    int idx = 0;
    for( iter = idxmap.begin(); iter != idxmap.end(); ++iter )
    {
        iter->second        = idx;
        domains[idx].domain = iter->first;
        ++idx;
    }
    for( i = 0; i < vertexList.size(); ++i )
    {
        std::vector< Mesh::VertexHandle >::const_iterator s, e;
        s              = std::lower_bound( vertices.begin(), vertices.end(), vertexList[i].firstHandle );
        e              = std::upper_bound( vertices.begin(), vertices.end(), vertexList[i].lastHandle );
        DomainSet& set = domains[idxmap[vertexList[i].domain]];
        std::copy( s, e, std::back_inserter( set.vertices ) );
    }
    for( i = 0; i < elementList.size(); ++i )
    {
        std::vector< Mesh::ElementHandle >::const_iterator s, e;
        s              = std::lower_bound( elements.begin(), elements.end(), elementList[i].firstHandle );
        e              = std::upper_bound( elements.begin(), elements.end(), elementList[i].lastHandle );
        DomainSet& set = domains[idxmap[elementList[i].domain]];
        std::copy( s, e, std::back_inserter( set.elements ) );
    }

    // guess geometric dimension for each domain
    std::map< DomainSet*, int > dimmap;
    std::vector< unsigned short > dof;
    for( i = 0; i < domains.size(); ++i )
    {
        if( !domains[i].elements.empty() )
            dimmap[&domains[i]] = 2;
        else
        {
            dof.resize( domains[i].vertices.size() );
            domains[i].domain->domain_DoF( &( domains[i].vertices[0] ), arrptr( dof ), dof.size(), err );MSQ_ERRRTN( err );
            unsigned short dim  = *std::max_element( dof.begin(), dof.end() );
            dimmap[&domains[i]] = dim;
        }
    }

    // group domains by dimension
    std::vector< DomainSet* > points, curves, surfaces;
    for( std::map< DomainSet*, int >::iterator it = dimmap.begin(); it != dimmap.end(); ++it )
    {
        switch( it->second )
        {
            case 0:
                points.push_back( it->first );
                break;
            case 1:
                curves.push_back( it->first );
                break;
            case 2:
                surfaces.push_back( it->first );
                break;
            default:
                MSQ_SETERR( err )
                ( MsqError::INVALID_STATE, "Invalid domain dimension: %d\n", it->second );
                break;
        }
    }

    // check that each point-domain has a single vertex
    for( i = 0; i < points.size(); ++i )
    {
        if( points[i]->vertices.size() != 1 || !points[i]->elements.empty() )
        {
            MSQ_SETERR( err )( "Point domain mesh not a single vertex.", MsqError::INVALID_STATE );
            return;
        }
    }

    // check that each curve domain has a chain of connected edges
    std::set< Mesh::VertexHandle > unseen;
    for( i = 0; i < curves.size(); ++i )
    {
        if( !curves[i]->elements.empty() )
        {
            MSQ_SETERR( err )( "Elements associated with 1D domain.", MsqError::INVALID_STATE );
            return;
        }

        unseen.clear();
        std::copy( curves[i]->vertices.begin(), curves[i]->vertices.end(), std::inserter( unseen, unseen.begin() ) );

        const Mesh::VertexHandle first_vtx = *unseen.begin();
        unseen.erase( unseen.begin() );
        Mesh::VertexHandle vtx = first_vtx;
        // find chain of vertices
        while( next_vertex( mesh, vtx, unseen, err ) )
            MSQ_ERRRTN( err );
        // search again from the starting vertex because
        // it probably wasn't the first one on the curve
        vtx = first_vtx;
        while( next_vertex( mesh, vtx, unseen, err ) )
            MSQ_ERRRTN( err );
        // were all vertices in a chain?
        if( !unseen.empty() )
        {
            MSQ_SETERR( err )
            ( "Curve domain contains vertices not in a simply connected chain.", MsqError::INVALID_STATE );
            return;
        }
    }

    std::set< Mesh::VertexHandle > seen;
    std::set< Mesh::ElementHandle > remaining;
    std::vector< Mesh::VertexHandle > verts, verts2;
    std::vector< Mesh::ElementHandle > stack, vert_elems;
    std::vector< EntityTopology > types;
    std::vector< size_t > junk;
    // if surface contains elements...
    for( i = 0; i < surfaces.size(); ++i )
    {
        if( surfaces[i]->elements.empty() ) continue;

        // Check that any vertices on surface are contained in an
        // element on the surface.
        verts.clear();
        mesh->elements_get_attached_vertices( &( surfaces[i]->elements[0] ), surfaces[i]->elements.size(), verts, junk,
                                              err );MSQ_ERRRTN( err );
        seen.clear();
        std::copy( verts.begin(), verts.end(), std::inserter( seen, seen.begin() ) );

        std::vector< Mesh::VertexHandle >::const_iterator v;
        for( v = surfaces[i]->vertices.begin(); v != surfaces[i]->vertices.end(); ++v )
        {
            std::set< Mesh::VertexHandle >::iterator j = seen.find( *v );
            if( j == seen.end() )
            {
                MSQ_SETERR( err )
                ( "Vertex on surface domain not in any element.", MsqError::INVALID_STATE );
                return;
            }
        }

        // check that elements form 2D patch
        stack.clear();
        remaining.clear();
        std::copy( surfaces[i]->elements.begin(), surfaces[i]->elements.end(),
                   std::inserter( remaining, remaining.begin() ) );
        stack.push_back( *remaining.begin() );
        remaining.erase( remaining.begin() );
        while( !stack.empty() )
        {
            Mesh::ElementHandle elem = stack.back();
            stack.pop_back();
            verts.clear();
            mesh->elements_get_attached_vertices( &elem, 1, verts, junk, err );MSQ_ERRRTN( err );
            // for each edge
            for( size_t j = 0; j < verts.size(); ++j )
            {
                Mesh::VertexHandle v1 = verts[j], v2 = verts[( j + 1 ) % verts.size()];
                vert_elems.clear();
                mesh->vertices_get_attached_elements( &v1, 1, vert_elems, junk, err );MSQ_ERRRTN( err );
                types.resize( vert_elems.size() );
                if( !vert_elems.empty() )
                {
                    mesh->elements_get_topologies( arrptr( vert_elems ), arrptr( types ), vert_elems.size(), err );MSQ_ERRRTN( err );
                }
                while( !vert_elems.empty() )
                {
                    Mesh::ElementHandle e2 = vert_elems.back();
                    EntityTopology type    = types.back();
                    vert_elems.pop_back();
                    types.pop_back();
                    if( TopologyInfo::dimension( type ) != 2 ) continue;
                    verts2.clear();
                    mesh->elements_get_attached_vertices( &e2, 1, verts2, junk, err );MSQ_ERRRTN( err );
                    size_t jdx = std::find( verts2.begin(), verts2.end(), v1 ) - verts2.begin();
                    if( verts2[( jdx + 1 ) % verts2.size()] != v2 &&
                        verts2[( jdx + verts2.size() - 1 ) % verts2.size()] != v2 )
                        continue;
                    std::set< Mesh::ElementHandle >::iterator r = remaining.find( e2 );
                    if( r == remaining.end() ) continue;

                    stack.push_back( *r );
                    remaining.erase( r );
                }
            }
        }

        if( !remaining.empty() )
        {
            MSQ_SETERR( err )
            ( "Surface mesh not a single, simply-connected patch", MsqError::INVALID_STATE );
            return;
        }
    }

    // check that sides of volume elements that are on surface
    // form simply connected patch
    for( i = 0; i < surfaces.size(); ++i )
    {
        // build list of sides known to be on surface
        std::sort( surfaces[i]->vertices.begin(), surfaces[i]->vertices.end() );
        std::set< std::pair< Mesh::ElementHandle, int > > sides;
        for( size_t j = 0; j < surfaces[i]->vertices.size(); ++j )
        {
            Mesh::VertexHandle v = surfaces[i]->vertices[j];
            vert_elems.clear();
            mesh->vertices_get_attached_elements( &v, 1, vert_elems, junk, err );MSQ_ERRRTN( err );
            types.resize( vert_elems.size() );
            if( !vert_elems.empty() )
            {
                mesh->elements_get_topologies( arrptr( vert_elems ), arrptr( types ), vert_elems.size(), err );MSQ_ERRRTN( err );
            }
            while( !vert_elems.empty() )
            {
                Mesh::ElementHandle e = vert_elems.back();
                EntityTopology type   = types.back();
                vert_elems.pop_back();
                types.pop_back();
                if( TopologyInfo::dimension( type ) != 3 ) continue;
                verts.clear();
                mesh->elements_get_attached_vertices( &e, 1, verts, junk, err );MSQ_ERRRTN( err );

                for( unsigned s = 0; s < TopologyInfo::faces( type ); ++s )
                {
                    unsigned n;
                    const unsigned* si = TopologyInfo::face_vertices( type, s, n );
                    unsigned ns        = 0;
                    for( unsigned k = 0; k < n; ++k )
                    {
                        if( std::binary_search( surfaces[i]->vertices.begin(), surfaces[i]->vertices.end(),
                                                verts[si[k]] ) )
                            ++ns;
                    }
                    if( ns >= 3 ) sides.insert( std::pair< Mesh::ElementHandle, int >( e, s ) );
                }
            }
        }

        std::vector< std::pair< Mesh::ElementHandle, int > > sstack;
        sstack.push_back( *sides.begin() );
        sides.erase( sides.begin() );
        while( !sstack.empty() )
        {
            Mesh::ElementHandle e = sstack.back().first;
            int s                 = sstack.back().second;
            sstack.pop_back();

            verts.clear();
            mesh->elements_get_attached_vertices( &e, 1, verts, junk, err );MSQ_ERRRTN( err );
            EntityTopology type;
            mesh->elements_get_topologies( &e, &type, 1, err );MSQ_ERRRTN( err );
            unsigned n;
            const unsigned* si = TopologyInfo::face_vertices( type, s, n );

            // for each edge
            for( unsigned j = 0; j < n; ++j )
            {
                Mesh::VertexHandle v1 = verts[si[j]], v2 = verts[si[( j + 1 ) % n]];
                vert_elems.clear();
                mesh->vertices_get_attached_elements( &v1, 1, vert_elems, junk, err );MSQ_ERRRTN( err );
                types.resize( vert_elems.size() );
                if( !vert_elems.empty() )
                {
                    mesh->elements_get_topologies( arrptr( vert_elems ), arrptr( types ), vert_elems.size(), err );MSQ_ERRRTN( err );
                }
                while( !vert_elems.empty() )
                {
                    Mesh::ElementHandle e2 = vert_elems.back();
                    EntityTopology type2   = types.back();
                    vert_elems.pop_back();
                    types.pop_back();
                    if( TopologyInfo::dimension( type ) != 3 ) continue;
                    verts2.clear();
                    mesh->elements_get_attached_vertices( &e2, 1, verts2, junk, err );MSQ_ERRRTN( err );
                    // for each face
                    for( unsigned s2 = 0; s2 < TopologyInfo::faces( type2 ); ++s2 )
                    {
                        std::set< std::pair< Mesh::ElementHandle, int > >::iterator side;
                        side = sides.find( std::pair< Mesh::ElementHandle, int >( e2, s2 ) );
                        if( side == sides.end() ) continue;

                        unsigned n2;
                        const unsigned* si2 = TopologyInfo::face_vertices( type2, s2, n2 );
                        unsigned jdx;
                        for( jdx = 0; jdx < n2; ++jdx )
                            if( verts2[si2[jdx]] == v1 ) break;
                        assert( jdx < n2 );

                        if( verts2[si2[( jdx + 1 ) % n2]] == v2 || verts2[si2[( jdx + n2 - 1 ) % n2]] == v2 )
                        {
                            sstack.push_back( *side );
                            sides.erase( side );
                        }
                    }
                }
            }
        }

        if( !sides.empty() )
        {
            MSQ_SETERR( err )
            ( "Surface mesh not a single, simply-connected patch", MsqError::INVALID_STATE );
            return;
        }
    }
}
void MBMesquite::DomainClassifier::vertex_normal_at ( Mesh::VertexHandle  entity_handle,
Vector3D coordinate 
) const [virtual]

Returns the normal of the domain to which "entity_handle" is constrained. For non-planar surfaces, the normal is calculated at the point on the domain that is closest to the passed in value of "coordinate". If the domain does not have a normal, or the normal cannot be determined, "coordinate" is set to (0,0,0). Otherwise, "coordinate" is set to the domain's normal at the appropriate point. In summary, the handle determines the domain. The coordinate determines the point of interest on that domain.

User should see also PatchData::get_domain_normal_at_vertex and PatchData::get_domain_normal_at_element .

Implements MBMesquite::MeshDomain.

Definition at line 1132 of file DomainClassifier.cpp.

References find_vertex_domain().

{
    if( const MeshDomain* dom = find_vertex_domain( entity_handle ) )
        dom->vertex_normal_at( entity_handle, coordinate );
}
void MBMesquite::DomainClassifier::vertex_normal_at ( const Mesh::VertexHandle handles,
Vector3D  coordinates[],
unsigned  count,
MsqError err 
) const [virtual]

evaluate surface normals

Returns normals for a domain.

Parameters:
handlesThe domain evaluated is the one in which this mesh entity is constrained.
coordinatesAs input, a list of positions at which to evaluate the domain. As output, the resulting domain normals.
countThe length of the coordinates array.

Implements MBMesquite::MeshDomain.

Definition at line 1144 of file DomainClassifier.cpp.

References find_vertex_domain(), MBMesquite::MsqError::INVALID_ARG, MSQ_ERRRTN, MSQ_SETERR, and MBMesquite::MeshDomain::vertex_normal_at().

{
    for( unsigned i = 0; i < count; ++i )
    {
        const MeshDomain* dom = find_vertex_domain( handles[i] );
        if( !dom )
        {
            MSQ_SETERR( err )( MsqError::INVALID_ARG );
            return;
        }
        dom->vertex_normal_at( handles + i, coordinates + i, 1, err );MSQ_ERRRTN( err );
    }
}

Member Data Documentation

List of all members.


The documentation for this class was generated from the following files:
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines