example_adaptivemesht.cc

This is an example of how to use the Adaptive Mesh class together with items, stellar and adaptive policy classes.

See also:

adaptivemesht.hh

#include <iostream>

#include <cstring>

#include <mesh.hh>

typedef unsigned int uint;

template< class Traits >
class MyVertex : public a48::AdaptiveVertexT< Traits > {
public: uint id; float pos[2]; 
};

class MyTraits {
public:
    typedef MyVertex< MyTraits > vertex_type;
    typedef a48::AdaptiveHalfedgeT< MyTraits > halfedge_type;
    typedef a48::AdaptiveFaceT< MyTraits > face_type;
};

typedef a48::AdaptiveMeshT< MyTraits > MyMesh; 

template< class Traits >
class MyItemsPolicy : public virtual a48::BaseItemsPolicy< Traits > {

    typedef typename Traits::vertex_type vertex_type;
    typedef typename Traits::halfedge_type halfedge_type;
    typedef typename Traits::face_type face_type;

public:

    uint num_verts, num_faces;
    float *vertices;
    uint *faces;

    MyItemsPolicy() {
        num_verts = 4; num_faces = 2;
        const float _vertices[] = {
            0.f, 0.f,
            10.f, 0.f,
            10.f, 10.f,
            0.f, 10.f
        };
        const uint _faces[] = {
            2, 0, 1,
            0, 2, 3
        };
        vertices = new float[ num_verts * 2 ];
        faces = new uint[ num_faces * 3 ];
        memcpy( vertices, _vertices, num_verts * 2 * sizeof(float) );
        memcpy( faces, _faces, num_faces * 3 * sizeof(uint) );
    }

    ~MyItemsPolicy() {
        delete [] vertices;
        delete [] faces;
    }

    void set_vertex_attributes(const uint& i, vertex_type *v) {
        v->id = i;
        v->pos[0] = vertices[ i*2 + 0 ];
        v->pos[1] = vertices[ i*2 + 1 ];
    }

    void set_vertex_attributes(vertex_type *v1, const vertex_type *v0) { }
    void set_halfedge_attributes(halfedge_type *h1, const halfedge_type *h0) { }
    void set_halfedge_attributes(halfedge_type *h, halfedge_type *o) { }
    void set_halfedge_attributes(const unsigned int& i, const unsigned int& j, const unsigned int& hi, halfedge_type *h) { }
    void set_face_attributes(face_type *f1, const face_type *f0) { }
    void set_face_attributes(const unsigned int& i, face_type *f) { }

};

template< class Traits >
class MyStellarPolicy : public virtual a48::BaseStellarPolicy< Traits > {

    typedef typename Traits::vertex_type vertex_type;
    typedef typename Traits::halfedge_type halfedge_type;
    typedef typename Traits::face_type face_type;

public:

    void sample_vertex(const halfedge_type *h, vertex_type *v) {
        v->pos[0] = ( h->vertex()->pos[0] + h->from_vertex()->pos[0] ) / 2.f;
        v->pos[1] = ( h->vertex()->pos[1] + h->from_vertex()->pos[1] ) / 2.f;
        v->id = ip.num_verts;
        ++ip.num_verts;
    }

    void sample_vertex(const face_type *f, vertex_type *v) {
        v->pos[0] = ( f->halfedge()->vertex()->pos[0] +
                  f->halfedge()->next()->vertex()->pos[0] +
                  f->halfedge()->previous()->vertex()->pos[0] ) / 3.f;
        v->pos[1] = ( f->halfedge()->vertex()->pos[1] +
                  f->halfedge()->next()->vertex()->pos[1] +
                  f->halfedge()->previous()->vertex()->pos[1] ) / 3.f;
        v->id = ip.num_verts;
        ++ip.num_verts;
    }

    void remove_vertex(const halfedge_type *h, const vertex_type *v) {
        --ip.num_verts;
    }

    void remove_vertex(const face_type *f, const vertex_type *v) {
        --ip.num_verts;
    }

    void flip_halfedge(const halfedge_type *h) { }

    MyItemsPolicy< Traits > ip;

};

template< class Traits >
class MyAdaptivePolicy : public virtual a48::BaseAdaptivePolicy< Traits > {

    typedef typename Traits::vertex_type vertex_type;
    typedef typename Traits::halfedge_type halfedge_type;
    typedef typename Traits::face_type face_type;

    inline int rank(vertex_type *v) { return (int)(maximum_level-v->level()); }
    inline int rank(halfedge_type *h) { return (int)(maximum_level-h->level()); }

    class element_s {
    public: vertex_type *v; int r;
        element_s(vertex_type *_v = 0, const int& _r = 0) : v(_v), r(_r) { }
        bool operator == (const element_s& e) { return this->v == e.v; }
        bool operator < (const element_s& e) { return this->r < e.r; }
    };

    class element_r {
    public: halfedge_type *h; int r;
        element_r(halfedge_type *_h = 0, const int& _r = 0) : h(_h), r(_r) { }
        bool operator == (const element_r& e) { return this->h == e.h or (!this->h->is_boundary() and this->h->opposite() == e.h); }
        bool operator < (const element_r& e) { return this->r < e.r; }
    };

    void insert_s(vertex_type *v) {
        element_s es( v, rank(v) );
        if( std::find( s.begin(), s.end(), es ) != s.end() ) return;
        s.push_back( es );
        std::push_heap( s.begin(), s.end() );
    }

    void insert_r(halfedge_type *h) {
        element_r er( h, rank(h) );
        if( std::find( r.begin(), r.end(), er ) != r.end() ) return;
        r.push_back( er );
        std::push_heap( r.begin(), r.end() );
    }

    void remove_s(vertex_type *v) { s.resize( std::remove( s.begin(), s.end(), element_s(v) ) - s.begin() ); }
    void remove_r(halfedge_type *h) { r.resize( std::remove( r.begin(), r.end(), element_r(h) ) - r.begin() ); }

    typedef typename std::vector< element_s > simplification_list;
    typedef typename simplification_list::iterator simplification_iterator;

    typedef typename std::vector< element_r > refinement_list;
    typedef typename refinement_list::iterator refinement_iterator;

    MyMesh *m;
    uint maximum_level;
    simplification_list s;
    refinement_list r;

public:

    MyAdaptivePolicy(MyMesh *_m, const uint& _ml = 4) : m(_m), maximum_level(_ml) {
        for (MyMesh::face_iterator fit =  m->faces_begin();
             fit != m->faces_end(); ++fit) {
            insert_s( (*fit)->weld_vertex() );
            insert_r( (*fit)->split_halfedge() );
        }
    }

    void print_simp(void) {
        std::cout << "[simplification]" << std::flush;
        for (simplification_iterator sit = s.begin();
             sit != s.end(); ++sit) {
            if( (*sit).v == 0 ) break;
            std::cout << " [" << (*sit).v->id << "]" << std::flush;
        }
        std::cout << std::endl;
    }

    void print_ref(void) {
        std::cout << "[refinement]" << std::flush;
        for (refinement_iterator rit = r.begin();
             rit != r.end(); ++rit) {
            if( (*rit).h == 0 ) break;
            std::cout << " [" << (*rit).h->from_vertex()->id << ":"
                  << (*rit).h->vertex()->id << "]" << std::flush;
        }
        std::cout << std::endl;
    }

    void print(void) { print_simp(); print_ref(); }

    void update_simplification(vertex_type *v) {
        remove_s( v );
        halfedge_type *h = v->halfedge();
        do {
            remove_r( h->previous() );
            h = h->next()->opposite();
        } while( h and h != v->halfedge() );
        std::make_heap( s.begin(), s.end() );
        std::make_heap( r.begin(), r.end() );
    }

    void update_simplification(halfedge_type *h) {
        insert_r( h );
        if( h->face()->weld_vertex() == h->next()->vertex() )
            insert_s( h->next()->vertex() );
        if( !h->is_boundary() and h->opposite()->face()->weld_vertex() == h->opposite()->next()->vertex() )
            insert_s( h->opposite()->next()->vertex() );
    }

    void update_simplification(void) {
        for (simplification_iterator sit = s.begin(); sit != s.end(); ++sit) (*sit).r = rank( (*sit).v );
        std::make_heap( s.begin(), s.end() );
    }

    void update_refinement(halfedge_type *h) {
        remove_r( h );
        if( !h->is_boundary() ) {
            remove_s( h->next()->vertex() );
            remove_s( h->opposite()->next()->vertex() );
        }
        std::make_heap( s.begin(), s.end() );
        std::make_heap( r.begin(), r.end() );
    }

    void update_refinement(vertex_type *v) {
        insert_s( v );
        halfedge_type *h = v->halfedge();
        do {
            if( h->previous()->is_split() ) insert_r( h->previous() );
            h = h->next()->opposite();
        } while( h and h != v->halfedge() );        
    }

    void update_refinement(void) {
        for (refinement_iterator rit = r.begin(); rit != r.end(); ++rit) (*rit).r = rank( (*rit).h );
        std::make_heap( r.begin(), r.end() );
    }

    vertex_type* next_simplify_vertex(void) {
        if( s.empty() ) return 0;
        vertex_type *vpop = s.front().v;
        std::pop_heap( s.begin(), s.end() ); s.pop_back();
        return vpop;
    }

    halfedge_type* next_refine_halfedge(void) {
        if( r.empty() ) return 0;
        halfedge_type *hpop = r.front().h;
        if( hpop->level() >= maximum_level-2 ) return 0;
        std::pop_heap( r.begin(), r.end() ); r.pop_back();
        return hpop;
    }

};


void show_mesh(const MyMesh& m) {

    std::cout << "_mesh_ Vertices: " << m.size_of_vertices() << std::endl;

    for (MyMesh::vertex_const_iterator vit = m.vertices_begin();
         vit != m.vertices_end(); ++vit) {

        std::cout << "\t(" << (*vit)->id << std::flush;
        std::cout << ") [" << (*vit)->pos[0] << " " << (*vit)->pos[1] << std::flush;
        std::cout << " " << ((*vit)->is_boundary() ? "b" : "i");
        std::cout << " d=" << (*vit)->degree();
        std::cout << " l=" << (*vit)->level();
        std::cout << "]\n";

    }

    std::cout << "_mesh_ Halfedges: " << m.size_of_halfedges() << std::endl;

    for (MyMesh::halfedge_const_iterator hit = m.halfedges_begin();
         hit != m.halfedges_end(); ++hit) {

        std::cout << "\t[ " << (*hit)->from_vertex()->id
              << ":" << (*hit)->vertex()->id
              << " " << ((*hit)->is_boundary() ? "b" : "i")
              << " (" << (*hit)->face()->halfedge()->previous()->vertex()->id
              << ", " << (*hit)->face()->halfedge()->vertex()->id
              << ", " << (*hit)->face()->halfedge()->next()->vertex()->id
              << ") ]";

        if( !(*hit)->is_boundary() ) {
            std::cout << " / [ "
                  << (*hit)->opposite()->from_vertex()->id
                  << ":" << (*hit)->opposite()->vertex()->id
                  << " " << ((*hit)->opposite()->is_boundary() ? "b" : "i")
                  << "]";
        }

        std::cout << "\n";
        
    }

    std::cout << "_mesh_ Faces: " << m.size_of_faces() << std::endl;

    for (MyMesh::face_const_iterator fit = m.faces_begin();
         fit != m.faces_end(); ++fit) {

        std::cout << "\t[ " << (*fit)->halfedge()->previous()->vertex()->id
              << ", " << (*fit)->halfedge()->vertex()->id
              << ", " << (*fit)->halfedge()->next()->vertex()->id
              << " " << ((*fit)->is_boundary() ? "b" : "i")
              << " ] / [" << (*fit)->split_halfedge()->from_vertex()->id
              << ":" << (*fit)->split_halfedge()->vertex()->id << "]\n";

    }

    std::cout << "_mesh_ Vertices stars:" << std::endl;

    for (MyMesh::vertex_const_iterator vit = m.vertices_begin();
         vit != m.vertices_end(); ++vit) {

        std::cout << "\t(" << (*vit)->id << ")" << std::flush;

        for (const MyTraits::halfedge_type *hcurr = (*vit)->halfedge();
             hcurr != 0; hcurr = (*vit)->star_next(hcurr)) {

            std::cout << " [" << hcurr->from_vertex()->id
                  << ":" << hcurr->vertex()->id
                  << "]" << std::flush;

        }

        std::cout << "\n";

    }

}

void refine_mesh(MyMesh& m, const uint& id0, const uint& id1) {

    for (MyMesh::halfedge_const_iterator hit = m.halfedges_begin();
         hit != m.halfedges_end(); ++hit) {

        if( (*hit)->from_vertex()->id == id0 and
            (*hit)->vertex()->id == id1 ) { m.refine( *hit ); break; }

    }

}

void simplify_mesh(MyMesh& m, const uint& id) {

    for (MyMesh::vertex_const_iterator vit = m.vertices_begin();
         vit != m.vertices_end(); ++vit) {

        if( (*vit)->id == id ) { m.simplify( *vit ); break; }

    }

}

int main(int argc, char* argv[]) {

    MyStellarPolicy< MyTraits > sp;

    MyMesh m( &sp.ip, &sp );

    m.set_base( sp.ip.num_verts, sp.ip.num_faces, sp.ip.faces );

    std::cout << "[test] 16-faces Refinement ...\n";

    refine_mesh( m, 0, 2 );
    refine_mesh( m, 0, 1 );
    refine_mesh( m, 1, 2 );
    refine_mesh( m, 2, 3 );
    refine_mesh( m, 3, 0 );
    refine_mesh( m, 4, 0 );
    refine_mesh( m, 4, 1 );
    refine_mesh( m, 4, 2 );
    refine_mesh( m, 4, 3 );

    show_mesh( m );

    std::cout << "[test] Creating adaptive policy ...\n";

    MyAdaptivePolicy< MyTraits > ap( &m );

    ap.print();

    m.set_adaptive_policy( &ap );

    std::cout << "[test] Simplifying Vertex 4\n";

    simplify_mesh( m, 4 );

    show_mesh( m );

    ap.print();

    std::cout << "[test] Adapt Refine\n";

    m.adapt_refine();

    show_mesh( m );

    std::cout << "[test] Adapt Simplify\n";

    m.adapt_simplify();

    show_mesh( m );

    return 0;

}