SimplexTemplateRefiner.cpp

Go to the documentation of this file.

#include "SimplexTemplateRefiner.hpp"

#include "EdgeSizeEvaluator.hpp"
#include "RefinerTagManager.hpp"

#include "moab/Interface.hpp"

#include <iostream>
#include <stack>

// Static arrays holding parametric coordinates of element vertices
static double MBVertexParametric[] = { 0., 0., 0. };
static double MBEdgeParametric[]   = { 0., 0., 0., 1., 0., 0. };
static double MBTriParametric[]    = { 0., 0., 0., 1., 0., 0., 0., 1., 0. };
static double MBTetParametric[]    = { 0., 0., 0., 1., 0., 0., 0., 1., 0., 0., 0., 1. };

#ifdef MB_DEBUG_TESSELLATOR
#define MB_TESSELLATOR_INCR_CASE_COUNT( cs )        this->case_counts[cs]++
#define MB_TESSELLATOR_INCR_SUBCASE_COUNT( cs, sc ) this->subcase_counts[cs][sc]++
#else  // MB_DEBUG_TESSELLATOR
#define MB_TESSELLATOR_INCR_CASE_COUNT( cs )
#define MB_TESSELLATOR_INCR_SUBCASE_COUNT( cs, sc )
#endif  // MB_DEBUG_TESSELLATOR

namespace moab
{

/// Construct a template refiner.
SimplexTemplateRefiner::SimplexTemplateRefiner()
{
    this->tag_manager  = 0;
    this->tag_assigner = new SimplexTemplateTagAssigner( this );
    this->corner_coords.resize( 6 * 8 );  // Hex has 8 verts w/ 6 coordinates each
    this->corner_tags.resize( 8 );        // Hex has 8 verts (this is a pointer, not the actual tag data)
    this->corner_handles.resize( 8 );     // Hex has 8 verts (this is a pointer, not actual hash data)
    this->input_is_output = false;        // Until we know better
}

/// Empty destructor for good form.
SimplexTemplateRefiner::~SimplexTemplateRefiner()
{
    delete this->tag_assigner;
}

/**\brief Stream a single mesh entity through the refiner.
 */
bool SimplexTemplateRefiner::refine_entity( EntityType etyp, EntityHandle entity )
{
    this->reset_heap_pointers();
    bool rval = true;
    const EntityHandle* conn;
    int num_nodes;
    Interface* imesh = this->tag_manager->get_input_mesh();
    if( imesh->get_connectivity( entity, conn, num_nodes ) != MB_SUCCESS )
    {
        return false;
    }
    this->corner_coords.resize( 6 * num_nodes );
    this->corner_tags.resize( num_nodes );
    this->corner_handles.resize( num_nodes );

    // Have to make num_nodes calls to get_coords() because we need xyz interleaved with rst coords.
    Tag tag_handle;
    int tag_offset;
    void* tag_data;
    for( int n = 0; n < num_nodes; ++n )
    {
        if( imesh->get_coords( &conn[n], 1, &corner_coords[6 * n + 3] ) != MB_SUCCESS )
        {
            return false;
        }
        tag_data = this->heap_tag_storage();
        for( int i = 0; i < this->tag_manager->get_number_of_vertex_tags(); ++i )
        {
            this->tag_manager->get_input_vertex_tag( i, tag_handle, tag_offset );
            if( imesh->tag_get_data( tag_handle, &conn[n], 1, (char*)tag_data + tag_offset ) != MB_SUCCESS )
            {
                return false;
            }
        }
        if( this->input_is_output )
        {
            this->corner_handles[n] = conn[n];
        }
        else
        {
            this->corner_handles[n] = this->output_functor->map_vertex( conn[n], &corner_coords[6 * n], tag_data );
#if 0
      std::cout << "#+# " << this->corner_handles[n] << " < "
        << corner_coords[ 6 * n + 3 ] << ", "
        << corner_coords[ 6 * n + 4 ] << ", "
        << corner_coords[ 6 * n + 5 ] << ">\n";
#endif  // 0
        }
        this->corner_tags[n] = tag_data;
    }

    switch( etyp )
    {
        case MBVERTEX:
            this->assign_parametric_coordinates( 1, MBVertexParametric, &this->corner_coords[0] );
            this->refine_0_simplex( &this->corner_coords[0], this->corner_tags[0], this->corner_handles[0] );
            rval = false;
            break;
        case MBEDGE:
            this->assign_parametric_coordinates( 2, MBEdgeParametric, &this->corner_coords[0] );
            rval = this->refine_1_simplex( this->maximum_number_of_subdivisions, &this->corner_coords[0],
                                           this->corner_tags[0], this->corner_handles[0], &this->corner_coords[6],
                                           this->corner_tags[1], this->corner_handles[1] );
            break;
        case MBTRI:
            this->assign_parametric_coordinates( 3, MBTriParametric, &this->corner_coords[0] );
            rval = this->refine_2_simplex( this->maximum_number_of_subdivisions, 7, &this->corner_coords[0],
                                           this->corner_tags[0], this->corner_handles[0], &this->corner_coords[6],
                                           this->corner_tags[1], this->corner_handles[1], &this->corner_coords[12],
                                           this->corner_tags[2], this->corner_handles[2] );
            break;
        case MBQUAD:
            std::cerr << "Quadrilaterals not handled yet\n";
            rval = false;
            break;
        case MBPOLYGON:
            std::cerr << "Polygons not handled yet\n";
            rval = false;
            break;
        case MBTET:
            this->assign_parametric_coordinates( 4, MBTetParametric, &this->corner_coords[0] );
            rval = this->refine_3_simplex( this->maximum_number_of_subdivisions, &this->corner_coords[0],
                                           this->corner_tags[0], this->corner_handles[0], &this->corner_coords[6],
                                           this->corner_tags[1], this->corner_handles[1], &this->corner_coords[12],
                                           this->corner_tags[2], this->corner_handles[2], &this->corner_coords[18],
                                           this->corner_tags[3], this->corner_handles[3] );
            break;
        case MBPYRAMID:
            std::cerr << "Pyramid schemes not handled yet\n";
            rval = false;
            break;
        case MBPRISM:
            std::cerr << "Prisms not handled yet\n";
            rval = false;
            break;
        case MBKNIFE:
            std::cerr << "Knifahedra not handled yet\n";
            rval = false;
            break;
        case MBHEX:
            std::cerr << "Hexahedra not handled yet\n";
            rval = false;
            break;
        case MBPOLYHEDRON:
            std::cerr << "Polyhedra not handled yet\n";
            rval = false;
            break;
        case MBENTITYSET:
            std::cerr << "How should entity sets be handled? We might iterate over its entries or "
                         "skip it."
                         "Must coordinate with MeshRefiner's loop over entities...\n";
            rval = false;
            break;
        case MBMAXTYPE:
            rval = false;
            break;
    }
    return rval;
}

/**\brief Set the function object used to decide which tag values an edge or face midpoint is
 * assigned.
 *
 * This will change the tag assigner's edge size evaluator to match the refiner's.
 * @param[in] ta The new tag assigner. This must be non-NULL.
 * @retval True if the tag assigner was changed and false otherwise.
 */
bool SimplexTemplateRefiner::set_tag_assigner( SimplexTemplateTagAssigner* ta )
{
    if( !ta || ta == this->tag_assigner ) return false;

    this->tag_assigner = ta;
    return true;
}

/**\brief Set the function object used to decide whether an edge is subdivided or not.
 */
bool SimplexTemplateRefiner::prepare( RefinerTagManager* tmgr, EntityRefinerOutputFunctor* ofunc )
{
    this->tag_manager = tmgr;
    this->tag_assigner->set_tag_manager( tmgr );
    this->input_is_output = ( this->tag_manager->get_input_mesh() == this->tag_manager->get_output_mesh() );
    // this->tag_assigner->set_edge_size_evaluator( this->edge_size_evaluator );
    return this->EntityRefiner::prepare( tmgr, ofunc );
}

/**\fn unsigned long SimplexTemplateRefiner::get_heap_size_bound( int max_recursions ) const
 *\brief Bound on the number of new vertices used to allocate the heap.
 *
 * This bound is based on a hexahedron that is divided into 48 tetrahedra (a point is added to
 * the center of each face so that compatible boundaries are guaranteed on neighboring hexahedra),
 * each of which has 4 edges.
 */

/**\brief "Refine" a vertex by passing it through to the output.
 *
 * FIXME: There is some question as to whether this should pass vertices through
 * since there does not appear to be a distinction between vertices as points
 * in space and vertices as degrees-of-freedom in a mesh (i.e. a vertex that is
 * treated as a lumped-parameter model).
 */
void SimplexTemplateRefiner::refine_0_simplex( const double* v0, const void* t0, EntityHandle h0 )
{
    // Ignore these arguments... the caller is responsible for calling the output functor to ensure
    // the vertex exists on the output.
    (void)v0;
    (void)t0;

    ( *this->output_functor )( h0 );
    ( *this->output_functor )( MBVERTEX );
}

/**\brief Refine an edge.
 */
bool SimplexTemplateRefiner::refine_1_simplex( int max_depth,
                                               const double* v0,
                                               const void* t0,
                                               EntityHandle h0,
                                               const double* v1,
                                               const void* t1,
                                               EntityHandle h1 )
{
    int edge_code = 0;

    double* midptc;
    void* midptt;
    // NB: If support for multiple recursive subdivisions (without a merge in between) is required,
    //     this will have to change as it assumes that h0 and h1 have only a single entry.
    EntityHandle midpth = 0;

    if( max_depth-- > 0 )
    {
        midptc = this->heap_coord_storage();
        midptt = this->heap_tag_storage();
        int i;
        // make valgrind happy
        // vtkstd::fill( midpt0, midpt0 + 6, 0. );
        for( i = 0; i < 6; i++ )
            midptc[i] = ( v0[i] + v1[i] ) / 2.;

        ( *this->tag_assigner )( v0, t0, h0, midptc, midptt, v1, t1, h1 );
        if( this->edge_size_evaluator->evaluate_edge( v0, t0, midptc, midptt, v1, t1 ) )
        {
            edge_code = 1;
            midpth    = ( *this->output_functor )( h0, h1, midptc, midptt );
        }
    }

    switch( edge_code )
    {
            // No edges to subdivide
        case 0:
            ( *this->output_functor )( h0 );
            ( *this->output_functor )( h1 );
            ( *this->output_functor )( MBEDGE );
            break;

            // One edge to subdivide
        case 1:
            this->refine_1_simplex( max_depth, v0, t0, h0, midptc, midptt, midpth );
            this->refine_1_simplex( max_depth, midptc, midptt, midpth, v1, t1, h1 );
            break;
    }

    return ( 0 == edge_code ) ? false : true;
}

/**\brief Refine a triangle.
 */
bool SimplexTemplateRefiner::refine_2_simplex( int max_depth,
                                               int move,
                                               const double* v0,
                                               const void* t0,
                                               EntityHandle h0,
                                               const double* v1,
                                               const void* t1,
                                               EntityHandle h1,
                                               const double* v2,
                                               const void* t2,
                                               EntityHandle h2 )
{
    int edge_code = 0;

    double* midpt0c;
    double* midpt1c;
    double* midpt2c;
    void* midpt0t;
    void* midpt1t;
    void* midpt2t;
    EntityHandle midpt0h = 0;
    EntityHandle midpt1h = 0;
    EntityHandle midpt2h = 0;

    if( max_depth-- > 0 )
    {
        int i;
        midpt0c = this->heap_coord_storage();
        midpt1c = this->heap_coord_storage();
        midpt2c = this->heap_coord_storage();
        midpt0t = this->heap_tag_storage();
        midpt1t = this->heap_tag_storage();
        midpt2t = this->heap_tag_storage();
        midpt0h = 0;
        midpt1h = 0;
        midpt2h = 0;
        for( i = 0; i < 6; ++i )
        {
            midpt0c[i] = ( v0[i] + v1[i] ) / 2.;
            midpt1c[i] = ( v1[i] + v2[i] ) / 2.;
            midpt2c[i] = ( v2[i] + v0[i] ) / 2.;
        }
        ( *this->tag_assigner )( v0, t0, h0, midpt0c, midpt0t, v1, t1, h1 );
        ( *this->tag_assigner )( v1, t1, h1, midpt1c, midpt1t, v2, t2, h2 );
        ( *this->tag_assigner )( v2, t2, h2, midpt2c, midpt2t, v0, t0, h0 );
        if( ( move & 1 ) && this->edge_size_evaluator->evaluate_edge( v0, t0, midpt0c, midpt0t, v1, t1 ) )
        {
            edge_code += 1;
            midpt0h = ( *this->output_functor )( h0, h1, midpt0c, midpt0t );
        }
        if( ( move & 2 ) && this->edge_size_evaluator->evaluate_edge( v1, t1, midpt1c, midpt1t, v2, t2 ) )
        {
            edge_code += 2;
            midpt1h = ( *this->output_functor )( h1, h2, midpt1c, midpt1t );
        }
        if( ( move & 4 ) && this->edge_size_evaluator->evaluate_edge( v2, t2, midpt2c, midpt2t, v0, t0 ) )
        {
            edge_code += 4;
            midpt2h = ( *this->output_functor )( h2, h0, midpt2c, midpt2t );
        }
    }

    switch( edge_code )
    {
            // No edges to subdivide
        case 0:
            ( *this->output_functor )( h0 );
            ( *this->output_functor )( h1 );
            ( *this->output_functor )( h2 );
            ( *this->output_functor )( MBTRI );
            break;

            // One edge to subdivide
        case 1:
            this->refine_2_simplex( max_depth, move | 2, v0, t0, h0, midpt0c, midpt0t, midpt0h, v2, t2, h2 );
            this->refine_2_simplex( max_depth, move | 4, midpt0c, midpt0t, midpt0h, v1, t1, h1, v2, t2, h2 );
            break;
        case 2:
            this->refine_2_simplex( max_depth, move | 4, v0, t0, h0, v1, t1, h1, midpt1c, midpt1t, midpt1h );
            this->refine_2_simplex( max_depth, move | 1, v0, t0, h0, midpt1c, midpt1t, midpt1h, v2, t2, h2 );
            break;
        case 4:
            this->refine_2_simplex( max_depth, move | 2, v0, t0, h0, v1, t1, h1, midpt2c, midpt2t, midpt2h );
            this->refine_2_simplex( max_depth, move | 1, midpt2c, midpt2t, midpt2h, v1, t1, h1, v2, t2, h2 );
            break;

            // Two edges to subdivide
        case 3:
            this->refine_2_simplex( max_depth, move | 4, midpt0c, midpt0t, midpt0h, v1, t1, h1, midpt1c, midpt1t,
                                    midpt1h );
            if( this->compare_Hopf_cross_string_dist( v2, midpt0c, v0, midpt1c ) )
            {
                this->refine_2_simplex( max_depth, move | 5, midpt0c, midpt0t, midpt0h, midpt1c, midpt1t, midpt1h, v2,
                                        t2, h2 );
                this->refine_2_simplex( max_depth, move | 2, v0, t0, h0, midpt0c, midpt0t, midpt0h, v2, t2, h2 );
            }
            else
            {
                this->refine_2_simplex( max_depth, move | 6, v0, t0, h0, midpt0c, midpt0t, midpt0h, midpt1c, midpt1t,
                                        midpt1h );
                this->refine_2_simplex( max_depth, move | 1, v0, t0, h0, midpt1c, midpt1t, midpt1h, v2, t2, h2 );
            }
            break;
        case 5:
            this->refine_2_simplex( max_depth, move | 2, v0, t0, h0, midpt0c, midpt0t, midpt0h, midpt2c, midpt2t,
                                    midpt2h );
            if( this->compare_Hopf_cross_string_dist( v2, midpt0c, v1, midpt2c ) )
            {
                this->refine_2_simplex( max_depth, move | 4, midpt0c, midpt0t, midpt0h, v1, t1, h1, v2, t2, h2 );
                this->refine_2_simplex( max_depth, move | 3, midpt2c, midpt2t, midpt2h, midpt0c, midpt0t, midpt0h, v2,
                                        t2, h2 );
            }
            else
            {
                this->refine_2_simplex( max_depth, move | 6, midpt0c, midpt0t, midpt0h, v1, t1, h1, midpt2c, midpt2t,
                                        midpt2h );
                this->refine_2_simplex( max_depth, move | 1, midpt2c, midpt2t, midpt2h, v1, t1, h1, v2, t2, h2 );
            }
            break;
        case 6:
            this->refine_2_simplex( max_depth, move | 1, midpt2c, midpt2t, midpt2h, midpt1c, midpt1t, midpt1h, v2, t2,
                                    h2 );
            if( this->compare_Hopf_cross_string_dist( v0, midpt1c, v1, midpt2c ) )
            {
                this->refine_2_simplex( max_depth, move | 3, v0, t0, h0, midpt1c, midpt1t, midpt1h, midpt2c, midpt2t,
                                        midpt2h );
                this->refine_2_simplex( max_depth, move | 4, v0, t0, h0, v1, t1, h1, midpt1c, midpt1t, midpt1h );
            }
            else
            {
                this->refine_2_simplex( max_depth, move | 2, v0, t0, h0, v1, t1, h1, midpt2c, midpt2t, midpt2h );
                this->refine_2_simplex( max_depth, move | 5, midpt2c, midpt2t, midpt2h, v1, t1, h1, midpt1c, midpt1t,
                                        midpt1h );
            }
            break;

            // Three edges to subdivide
        case 7:
            this->refine_2_simplex( max_depth, 7, midpt0c, midpt0t, midpt0h, midpt1c, midpt1t, midpt1h, midpt2c,
                                    midpt2t, midpt2h );
            this->refine_2_simplex( max_depth, move | 2, v0, t0, h0, midpt0c, midpt0t, midpt0h, midpt2c, midpt2t,
                                    midpt2h );
            this->refine_2_simplex( max_depth, move | 4, midpt0c, midpt0t, midpt0h, v1, t1, h1, midpt1c, midpt1t,
                                    midpt1h );
            this->refine_2_simplex( max_depth, move | 1, midpt2c, midpt2t, midpt2h, midpt1c, midpt1t, midpt1h, v2, t2,
                                    h2 );
            break;
    }

    return true;
}

/**\brief Refine a tetrahedron.
 */
bool SimplexTemplateRefiner::refine_3_simplex( int max_depth,
                                               double* v0,
                                               void* t0,
                                               EntityHandle h0,
                                               double* v1,
                                               void* t1,
                                               EntityHandle h1,
                                               double* v2,
                                               void* t2,
                                               EntityHandle h2,
                                               double* v3,
                                               void* t3,
                                               EntityHandle h3 )
{
    int edge_code = 0;

    double* midpt0c;
    double* midpt1c;
    double* midpt2c;
    double* midpt3c;
    double* midpt4c;
    double* midpt5c;

    void* midpt0t;
    void* midpt1t;
    void* midpt2t;
    void* midpt3t;
    void* midpt4t;
    void* midpt5t;

    EntityHandle midpt0h = 0;
    EntityHandle midpt1h = 0;
    EntityHandle midpt2h = 0;
    EntityHandle midpt3h = 0;
    EntityHandle midpt4h = 0;
    EntityHandle midpt5h = 0;

    if( max_depth-- > 0 )
    {
        midpt0c = this->heap_coord_storage();
        midpt1c = this->heap_coord_storage();
        midpt2c = this->heap_coord_storage();
        midpt3c = this->heap_coord_storage();
        midpt4c = this->heap_coord_storage();
        midpt5c = this->heap_coord_storage();

        midpt0t = this->heap_tag_storage();
        midpt1t = this->heap_tag_storage();
        midpt2t = this->heap_tag_storage();
        midpt3t = this->heap_tag_storage();
        midpt4t = this->heap_tag_storage();
        midpt5t = this->heap_tag_storage();

        for( int i = 0; i < 6; ++i )
        {
            midpt0c[i] = ( v0[i] + v1[i] ) * .5;
            midpt1c[i] = ( v1[i] + v2[i] ) * .5;
            midpt2c[i] = ( v2[i] + v0[i] ) * .5;
            midpt3c[i] = ( v0[i] + v3[i] ) * .5;
            midpt4c[i] = ( v1[i] + v3[i] ) * .5;
            midpt5c[i] = ( v2[i] + v3[i] ) * .5;
        }

        ( *this->tag_assigner )( v0, t0, h0, midpt0c, midpt0t, v1, t1, h1 );
        ( *this->tag_assigner )( v1, t1, h1, midpt1c, midpt1t, v2, t2, h2 );
        ( *this->tag_assigner )( v2, t2, h2, midpt2c, midpt2t, v0, t0, h0 );
        ( *this->tag_assigner )( v0, t0, h0, midpt3c, midpt3t, v3, t3, h3 );
        ( *this->tag_assigner )( v1, t1, h1, midpt4c, midpt4t, v3, t3, h3 );
        ( *this->tag_assigner )( v2, t2, h2, midpt5c, midpt5t, v3, t3, h3 );

        if( this->edge_size_evaluator->evaluate_edge( v0, t0, midpt0c, midpt0t, v1, t1 ) )
        {
            edge_code |= 1;
            midpt0h = ( *this->output_functor )( h0, h1, midpt0c, midpt0t );
        }
        if( this->edge_size_evaluator->evaluate_edge( v1, t1, midpt1c, midpt1t, v2, t2 ) )
        {
            edge_code |= 2;
            midpt1h = ( *this->output_functor )( h1, h2, midpt1c, midpt1t );
        }
        if( this->edge_size_evaluator->evaluate_edge( v2, t2, midpt2c, midpt2t, v0, t0 ) )
        {
            edge_code |= 4;
            midpt2h = ( *this->output_functor )( h2, h0, midpt2c, midpt2t );
        }
        if( this->edge_size_evaluator->evaluate_edge( v0, t0, midpt3c, midpt3t, v3, t3 ) )
        {
            edge_code |= 8;
            midpt3h = ( *this->output_functor )( h0, h3, midpt3c, midpt3t );
        }
        if( this->edge_size_evaluator->evaluate_edge( v1, t1, midpt4c, midpt4t, v3, t3 ) )
        {
            edge_code |= 16;
            midpt4h = ( *this->output_functor )( h1, h3, midpt4c, midpt4t );
        }
        if( this->edge_size_evaluator->evaluate_edge( v2, t2, midpt5c, midpt5t, v3, t3 ) )
        {
            edge_code |= 32;
            midpt5h = ( *this->output_functor )( h2, h3, midpt5c, midpt5t );
        }
    }

    double edge_length2[6];
    for( int ei = 0; ei < 6; ++ei )
        edge_length2[ei] = 0.;

    for( int c = 0; c < 3; ++c )
    {
        double tmp;
        tmp = v1[c] - v0[c];
        edge_length2[0] += tmp * tmp;
        tmp = v2[c] - v1[c];
        edge_length2[1] += tmp * tmp;
        tmp = v2[c] - v0[c];
        edge_length2[2] += tmp * tmp;
        tmp = v3[c] - v0[c];
        edge_length2[3] += tmp * tmp;
        tmp = v3[c] - v1[c];
        edge_length2[4] += tmp * tmp;
        tmp = v3[c] - v2[c];
        edge_length2[5] += tmp * tmp;
    }

    if( 0 == edge_code )
    {
        // No edges to subdivide
        ( *this->output_functor )( h0 );
        ( *this->output_functor )( h1 );
        ( *this->output_functor )( h2 );
        ( *this->output_functor )( h3 );
        ( *this->output_functor )( MBTET );

        return false;
    }

    double* facept0c;
    double* facept1c;
    double* facept2c;
    double* facept3c;
    facept0c                  = this->heap_coord_storage();
    facept1c                  = this->heap_coord_storage();
    facept2c                  = this->heap_coord_storage();
    facept3c                  = this->heap_coord_storage();
    double* vertex_coords[14] = { v0,      v1,      v2,      v3,       midpt0c,  midpt1c,  midpt2c,
                                  midpt3c, midpt4c, midpt5c, facept0c, facept1c, facept2c, facept3c };

    void* facept0t        = this->heap_tag_storage();
    void* facept1t        = this->heap_tag_storage();
    void* facept2t        = this->heap_tag_storage();
    void* facept3t        = this->heap_tag_storage();
    void* vertex_tags[14] = { t0,      t1,      t2,      t3,       midpt0t,  midpt1t,  midpt2t,
                              midpt3t, midpt4t, midpt5t, facept0t, facept1t, facept2t, facept3t };

    EntityHandle vertex_hash[14] = { h0, h1, h2, h3, midpt0h, midpt1h, midpt2h, midpt3h, midpt4h, midpt5h, 0, 0, 0, 0 };

    // Generate tetrahedra that are compatible except when edge
    // lengths are equal on indeterminately subdivided faces.
    double* permuted_coords[14];
    void* permuted_tags[14];
    EntityHandle permuted_hash[14];
    double permlen[6];  // permuted edge lengths
    int C = SimplexTemplateRefiner::template_index[edge_code][0];
    int P = SimplexTemplateRefiner::template_index[edge_code][1];

    // 1. Permute the tetrahedron into our canonical configuration
    for( int i = 0; i < 14; ++i )
    {
        permuted_coords[i] = vertex_coords[SimplexTemplateRefiner::permutations_from_index[P][i]];
        permuted_tags[i]   = vertex_tags[SimplexTemplateRefiner::permutations_from_index[P][i]];
        permuted_hash[i]   = vertex_hash[SimplexTemplateRefiner::permutations_from_index[P][i]];
    }

    for( int i = 4; i < 10; ++i )
    {
        // permute edge lengths too
        permlen[i - 4] = edge_length2[SimplexTemplateRefiner::permutations_from_index[P][i] - 4];
    }

    int comparison_bits;
    std::stack< int* > output_tets;
    std::stack< int* > output_perm;
    std::stack< int > output_sign;

    // cout << "Case " << C << "  Permutation " << P << endl;
    // 2. Generate tetrahedra based on the configuration.
    //    Note that case 0 is handled above (edgeCode == 0).

    switch( C )
    {
        case 1:  // Ruprecht-Müller Case 1
            MB_TESSELLATOR_INCR_CASE_COUNT( 0 );
            output_tets.push( SimplexTemplateRefiner::templates + 0 );
            output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
            output_sign.push( 1 );
            MB_TESSELLATOR_INCR_SUBCASE_COUNT( 0, 0 );
            break;
        case 2:  // Ruprecht-Müller Case 2a
            comparison_bits = ( permlen[0] <= permlen[1] ? 1 : 0 ) | ( permlen[0] >= permlen[1] ? 2 : 0 ) | 0;
            if( ( comparison_bits & 3 ) == 3 )
            {
                // Compute face point and tag
                for( int i = 0; i < 6; ++i )
                {
                    permuted_coords[10][i] =
                        ( permuted_coords[0][i] + permuted_coords[2][i] ) * .375 + permuted_coords[1][i] * .25;
                }
                ( *this->tag_assigner )( permuted_tags[0], permuted_tags[2], permuted_tags[1], permuted_tags[10] );
                permuted_hash[10] = ( *this->output_functor )( permuted_hash[0], permuted_hash[2], permuted_hash[1],
                                                               permuted_coords[10], permuted_tags[10] );
            }
            MB_TESSELLATOR_INCR_CASE_COUNT( 1 );
            output_tets.push( SimplexTemplateRefiner::templates + 9 );
            output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
            output_sign.push( 1 );
            MB_TESSELLATOR_INCR_SUBCASE_COUNT( 1, 0 );
            switch( comparison_bits )
            {
                case 2:  // 0>1
                    output_tets.push( SimplexTemplateRefiner::templates + 14 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 1, 1 );
                    break;
                case 1:  // 1>0
                    output_tets.push( SimplexTemplateRefiner::templates + 14 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[13] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 1, 2 );
                    break;
                case 3:  // 0=1
                    output_tets.push( SimplexTemplateRefiner::templates + 23 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 1, 3 );
                    break;
            }
            break;
        case 3:  // Ruprecht-Müller Case 2b
            MB_TESSELLATOR_INCR_CASE_COUNT( 2 );
            output_tets.push( SimplexTemplateRefiner::templates + 40 );
            output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
            output_sign.push( 1 );
            MB_TESSELLATOR_INCR_SUBCASE_COUNT( 2, 0 );
            break;
        case 4:  // Ruprecht-Müller Case 3a
            comparison_bits = ( permlen[0] <= permlen[3] ? 1 : 0 ) | ( permlen[0] >= permlen[3] ? 2 : 0 ) |
                              ( permlen[2] <= permlen[3] ? 4 : 0 ) | ( permlen[2] >= permlen[3] ? 8 : 0 ) |
                              ( permlen[0] <= permlen[2] ? 16 : 0 ) | ( permlen[0] >= permlen[2] ? 32 : 0 ) | 0;
            if( ( comparison_bits & 3 ) == 3 )
            {
                // Compute face point and tag
                for( int i = 0; i < 6; ++i )
                {
                    permuted_coords[11][i] =
                        ( permuted_coords[1][i] + permuted_coords[3][i] ) * .375 + permuted_coords[0][i] * .25;
                }
                ( *this->tag_assigner )( permuted_tags[1], permuted_tags[3], permuted_tags[0], permuted_tags[11] );
                permuted_hash[11] = ( *this->output_functor )( permuted_hash[1], permuted_hash[3], permuted_hash[0],
                                                               permuted_coords[11], permuted_tags[11] );
            }
            if( ( comparison_bits & 12 ) == 12 )
            {
                // Compute face point and tag
                for( int i = 0; i < 6; ++i )
                {
                    permuted_coords[13][i] =
                        ( permuted_coords[2][i] + permuted_coords[3][i] ) * .375 + permuted_coords[0][i] * .25;
                }
                ( *this->tag_assigner )( permuted_tags[2], permuted_tags[3], permuted_tags[0], permuted_tags[13] );
                permuted_hash[13] = ( *this->output_functor )( permuted_hash[2], permuted_hash[3], permuted_hash[0],
                                                               permuted_coords[13], permuted_tags[13] );
            }
            if( ( comparison_bits & 48 ) == 48 )
            {
                // Compute face point and tag
                for( int i = 0; i < 6; ++i )
                {
                    permuted_coords[10][i] =
                        ( permuted_coords[1][i] + permuted_coords[2][i] ) * .375 + permuted_coords[0][i] * .25;
                }
                ( *this->tag_assigner )( permuted_tags[1], permuted_tags[2], permuted_tags[0], permuted_tags[10] );
                permuted_hash[10] = ( *this->output_functor )( permuted_hash[1], permuted_hash[2], permuted_hash[0],
                                                               permuted_coords[10], permuted_tags[10] );
            }
            MB_TESSELLATOR_INCR_CASE_COUNT( 3 );
            output_tets.push( SimplexTemplateRefiner::templates + 57 );
            output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
            output_sign.push( 1 );
            MB_TESSELLATOR_INCR_SUBCASE_COUNT( 3, 0 );
            switch( comparison_bits )
            {
                case 42:  // 0>2>3<0
                    output_tets.push( SimplexTemplateRefiner::templates + 62 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 3, 1 );
                    break;
                case 25:  // 2>3>0<2
                    output_tets.push( SimplexTemplateRefiner::templates + 62 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[11] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 3, 2 );
                    break;
                case 37:  // 3>0>2<3
                    output_tets.push( SimplexTemplateRefiner::templates + 62 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[3] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 3, 3 );
                    break;
                case 21:  // 3>2>0<3
                    output_tets.push( SimplexTemplateRefiner::templates + 62 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[22] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 3, 4 );
                    break;
                case 26:  // 2>0>3<2
                    output_tets.push( SimplexTemplateRefiner::templates + 62 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[12] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 3, 5 );
                    break;
                case 38:  // 0>3>2<0
                    output_tets.push( SimplexTemplateRefiner::templates + 62 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[15] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 3, 6 );
                    break;
                case 58:  // 0=2>3<0
                    output_tets.push( SimplexTemplateRefiner::templates + 75 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 3, 7 );
                    break;
                case 29:  // 2=3>0<2
                    output_tets.push( SimplexTemplateRefiner::templates + 75 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[11] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 3, 8 );
                    break;
                case 39:  // 0=3>2<0
                    output_tets.push( SimplexTemplateRefiner::templates + 75 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[3] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 3, 9 );
                    break;
                case 53:  // 3>0=2<3
                    output_tets.push( SimplexTemplateRefiner::templates + 96 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 3, 10 );
                    break;
                case 46:  // 0>2=3<0
                    output_tets.push( SimplexTemplateRefiner::templates + 96 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[11] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 3, 11 );
                    break;
                case 27:  // 2>0=3<2
                    output_tets.push( SimplexTemplateRefiner::templates + 96 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[3] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 3, 12 );
                    break;
                case 63:  // 0=2=3=0
                    output_tets.push( SimplexTemplateRefiner::templates + 117 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 3, 13 );
                    break;
            }
            break;
        case 5:  // Ruprecht-Müller Case 3b
            MB_TESSELLATOR_INCR_CASE_COUNT( 4 );
            output_tets.push( SimplexTemplateRefiner::templates + 162 );
            output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
            output_sign.push( 1 );
            MB_TESSELLATOR_INCR_SUBCASE_COUNT( 4, 0 );
            break;
        case 6:  // Ruprecht-Müller Case 3c
            comparison_bits = ( permlen[0] <= permlen[1] ? 1 : 0 ) | ( permlen[0] >= permlen[1] ? 2 : 0 ) |
                              ( permlen[0] <= permlen[3] ? 4 : 0 ) | ( permlen[0] >= permlen[3] ? 8 : 0 ) | 0;
            if( ( comparison_bits & 3 ) == 3 )
            {
                // Compute face point and tag
                for( int i = 0; i < 6; ++i )
                {
                    permuted_coords[10][i] =
                        ( permuted_coords[0][i] + permuted_coords[2][i] ) * .375 + permuted_coords[1][i] * .25;
                }
                ( *this->tag_assigner )( permuted_tags[0], permuted_tags[2], permuted_tags[1], permuted_tags[10] );
                permuted_hash[10] = ( *this->output_functor )( permuted_hash[0], permuted_hash[2], permuted_hash[1],
                                                               permuted_coords[10], permuted_tags[10] );
            }
            if( ( comparison_bits & 12 ) == 12 )
            {
                // Compute face point and tag
                for( int i = 0; i < 6; ++i )
                {
                    permuted_coords[11][i] =
                        ( permuted_coords[1][i] + permuted_coords[3][i] ) * .375 + permuted_coords[0][i] * .25;
                }
                ( *this->tag_assigner )( permuted_tags[1], permuted_tags[3], permuted_tags[0], permuted_tags[11] );
                permuted_hash[11] = ( *this->output_functor )( permuted_hash[1], permuted_hash[3], permuted_hash[0],
                                                               permuted_coords[11], permuted_tags[11] );
            }
            MB_TESSELLATOR_INCR_CASE_COUNT( 5 );
            switch( comparison_bits )
            {
                case 10:  // 0>1,0>3
                    output_tets.push( SimplexTemplateRefiner::templates + 179 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 5, 0 );
                    break;
                case 5:  // 1>0,3>0
                    output_tets.push( SimplexTemplateRefiner::templates + 200 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 5, 1 );
                    break;
                case 6:  // 0>1,3>0
                    output_tets.push( SimplexTemplateRefiner::templates + 221 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 5, 2 );
                    break;
                case 9:  // 1>0,0>3
                    output_tets.push( SimplexTemplateRefiner::templates + 242 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 5, 3 );
                    break;
                case 11:  // 0=1,0>3
                    output_tets.push( SimplexTemplateRefiner::templates + 263 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 5, 4 );
                    break;
                case 14:  // 0=3,0>1
                    output_tets.push( SimplexTemplateRefiner::templates + 263 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[5] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 5, 5 );
                    break;
                case 7:  // 3>0,0=1
                    output_tets.push( SimplexTemplateRefiner::templates + 292 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 5, 6 );
                    break;
                case 13:  // 1>0,0=3
                    output_tets.push( SimplexTemplateRefiner::templates + 292 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[5] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 5, 7 );
                    break;
                case 15:  // 0=1,0=3
                    output_tets.push( SimplexTemplateRefiner::templates + 321 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 5, 8 );
                    break;
            }
            break;
        case 7:  // Ruprecht-Müller Case 3d
            comparison_bits = ( permlen[0] <= permlen[2] ? 1 : 0 ) | ( permlen[0] >= permlen[2] ? 2 : 0 ) |
                              ( permlen[0] <= permlen[4] ? 4 : 0 ) | ( permlen[0] >= permlen[4] ? 8 : 0 ) | 0;
            if( ( comparison_bits & 3 ) == 3 )
            {
                // Compute face point and tag
                for( int i = 0; i < 6; ++i )
                {
                    permuted_coords[10][i] =
                        ( permuted_coords[1][i] + permuted_coords[2][i] ) * .375 + permuted_coords[0][i] * .25;
                }
                ( *this->tag_assigner )( permuted_tags[1], permuted_tags[2], permuted_tags[0], permuted_tags[10] );
                permuted_hash[10] = ( *this->output_functor )( permuted_hash[1], permuted_hash[2], permuted_hash[0],
                                                               permuted_coords[10], permuted_tags[10] );
            }
            if( ( comparison_bits & 12 ) == 12 )
            {
                // Compute face point and tag
                for( int i = 0; i < 6; ++i )
                {
                    permuted_coords[11][i] =
                        ( permuted_coords[0][i] + permuted_coords[3][i] ) * .375 + permuted_coords[1][i] * .25;
                }
                ( *this->tag_assigner )( permuted_tags[0], permuted_tags[3], permuted_tags[1], permuted_tags[11] );
                permuted_hash[11] = ( *this->output_functor )( permuted_hash[0], permuted_hash[3], permuted_hash[1],
                                                               permuted_coords[11], permuted_tags[11] );
            }
            MB_TESSELLATOR_INCR_CASE_COUNT( 6 );
            switch( comparison_bits )
            {
                case 10:  // 0>4,0>2
                    output_tets.push( SimplexTemplateRefiner::templates + 362 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 6, 0 );
                    break;
                case 5:  // 4>0,2>0
                    output_tets.push( SimplexTemplateRefiner::templates + 383 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 6, 1 );
                    break;
                case 9:  // 0>4,2>0
                    output_tets.push( SimplexTemplateRefiner::templates + 404 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 6, 2 );
                    break;
                case 6:  // 4>0,0>2
                    output_tets.push( SimplexTemplateRefiner::templates + 425 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 6, 3 );
                    break;
                case 14:  // 0=4,0>2
                    output_tets.push( SimplexTemplateRefiner::templates + 446 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 6, 4 );
                    break;
                case 11:  // 0=2,0>4
                    output_tets.push( SimplexTemplateRefiner::templates + 446 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[5] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 6, 5 );
                    break;
                case 13:  // 2>0,0=4
                    output_tets.push( SimplexTemplateRefiner::templates + 475 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 6, 6 );
                    break;
                case 7:  // 4>0,0=2
                    output_tets.push( SimplexTemplateRefiner::templates + 475 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[5] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 6, 7 );
                    break;
                case 15:  // 0=4,0=2
                    output_tets.push( SimplexTemplateRefiner::templates + 504 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 6, 8 );
                    break;
            }
            break;
        case 8:  // Ruprecht-Müller Case 4a
            comparison_bits = ( permlen[4] <= permlen[5] ? 1 : 0 ) | ( permlen[4] >= permlen[5] ? 2 : 0 ) |
                              ( permlen[3] <= permlen[4] ? 4 : 0 ) | ( permlen[3] >= permlen[4] ? 8 : 0 ) | 0;
            if( ( comparison_bits & 3 ) == 3 )
            {
                // Compute face point and tag
                for( int i = 0; i < 6; ++i )
                {
                    permuted_coords[12][i] =
                        ( permuted_coords[1][i] + permuted_coords[2][i] ) * .375 + permuted_coords[3][i] * .25;
                }
                ( *this->tag_assigner )( permuted_tags[1], permuted_tags[2], permuted_tags[3], permuted_tags[12] );
                permuted_hash[12] = ( *this->output_functor )( permuted_hash[1], permuted_hash[2], permuted_hash[3],
                                                               permuted_coords[12], permuted_tags[12] );
            }
            if( ( comparison_bits & 12 ) == 12 )
            {
                // Compute face point and tag
                for( int i = 0; i < 6; ++i )
                {
                    permuted_coords[11][i] =
                        ( permuted_coords[0][i] + permuted_coords[1][i] ) * .375 + permuted_coords[3][i] * .25;
                }
                ( *this->tag_assigner )( permuted_tags[0], permuted_tags[1], permuted_tags[3], permuted_tags[11] );
                permuted_hash[11] = ( *this->output_functor )( permuted_hash[0], permuted_hash[1], permuted_hash[3],
                                                               permuted_coords[11], permuted_tags[11] );
            }
            MB_TESSELLATOR_INCR_CASE_COUNT( 7 );
            output_tets.push( SimplexTemplateRefiner::templates + 545 );
            output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
            output_sign.push( 1 );
            MB_TESSELLATOR_INCR_SUBCASE_COUNT( 7, 0 );
            switch( comparison_bits )
            {
                case 5:  // 5>4>3
                    output_tets.push( SimplexTemplateRefiner::templates + 554 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 7, 1 );
                    break;
                case 10:  // 3>4>5
                    output_tets.push( SimplexTemplateRefiner::templates + 554 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[13] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 7, 2 );
                    break;
                case 6:  // 3<4>5
                    output_tets.push( SimplexTemplateRefiner::templates + 571 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 7, 3 );
                    break;
                case 9:  // 3>4<5
                    output_tets.push( SimplexTemplateRefiner::templates + 588 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 7, 4 );
                    break;
                case 14:  // 3=4>5
                    output_tets.push( SimplexTemplateRefiner::templates + 605 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 7, 5 );
                    break;
                case 7:  // 4=5,4>3
                    output_tets.push( SimplexTemplateRefiner::templates + 605 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[13] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 7, 6 );
                    break;
                case 13:  // 5>4,3=4
                    output_tets.push( SimplexTemplateRefiner::templates + 630 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 7, 7 );
                    break;
                case 11:  // 3>4=5
                    output_tets.push( SimplexTemplateRefiner::templates + 630 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[13] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 7, 8 );
                    break;
                case 15:  // 3=4=5
                    output_tets.push( SimplexTemplateRefiner::templates + 655 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 7, 9 );
                    break;
            }
            break;
        case 9:  // Ruprecht-Müller Case 4b
            comparison_bits = ( permlen[1] <= permlen[2] ? 1 : 0 ) | ( permlen[1] >= permlen[2] ? 2 : 0 ) |
                              ( permlen[2] <= permlen[3] ? 4 : 0 ) | ( permlen[2] >= permlen[3] ? 8 : 0 ) |
                              ( permlen[3] <= permlen[4] ? 16 : 0 ) | ( permlen[3] >= permlen[4] ? 32 : 0 ) |
                              ( permlen[1] <= permlen[4] ? 64 : 0 ) | ( permlen[1] >= permlen[4] ? 128 : 0 ) | 0;
            if( ( comparison_bits & 3 ) == 3 )
            {
                // Compute face point and tag
                for( int i = 0; i < 6; ++i )
                {
                    permuted_coords[10][i] =
                        ( permuted_coords[1][i] + permuted_coords[0][i] ) * .375 + permuted_coords[2][i] * .25;
                }
                ( *this->tag_assigner )( permuted_tags[1], permuted_tags[0], permuted_tags[2], permuted_tags[10] );
                permuted_hash[10] = ( *this->output_functor )( permuted_hash[1], permuted_hash[0], permuted_hash[2],
                                                               permuted_coords[10], permuted_tags[10] );
            }
            if( ( comparison_bits & 12 ) == 12 )
            {
                // Compute face point and tag
                for( int i = 0; i < 6; ++i )
                {
                    permuted_coords[13][i] =
                        ( permuted_coords[2][i] + permuted_coords[3][i] ) * .375 + permuted_coords[0][i] * .25;
                }
                ( *this->tag_assigner )( permuted_tags[2], permuted_tags[3], permuted_tags[0], permuted_tags[13] );
                permuted_hash[13] = ( *this->output_functor )( permuted_hash[2], permuted_hash[3], permuted_hash[0],
                                                               permuted_coords[13], permuted_tags[13] );
            }
            if( ( comparison_bits & 48 ) == 48 )
            {
                // Compute face point and tag
                for( int i = 0; i < 6; ++i )
                {
                    permuted_coords[11][i] =
                        ( permuted_coords[0][i] + permuted_coords[1][i] ) * .375 + permuted_coords[3][i] * .25;
                }
                ( *this->tag_assigner )( permuted_tags[0], permuted_tags[1], permuted_tags[3], permuted_tags[11] );
                permuted_hash[11] = ( *this->output_functor )( permuted_hash[0], permuted_hash[1], permuted_hash[3],
                                                               permuted_coords[11], permuted_tags[11] );
            }
            if( ( comparison_bits & 192 ) == 192 )
            {
                // Compute face point and tag
                for( int i = 0; i < 6; ++i )
                {
                    permuted_coords[12][i] =
                        ( permuted_coords[2][i] + permuted_coords[3][i] ) * .375 + permuted_coords[1][i] * .25;
                }
                ( *this->tag_assigner )( permuted_tags[2], permuted_tags[3], permuted_tags[1], permuted_tags[12] );
                permuted_hash[12] = ( *this->output_functor )( permuted_hash[2], permuted_hash[3], permuted_hash[1],
                                                               permuted_coords[12], permuted_tags[12] );
            }
            MB_TESSELLATOR_INCR_CASE_COUNT( 8 );
            switch( comparison_bits )
            {
                case 85:  // 2>1,3>2,4>3,4>1
                    output_tets.push( SimplexTemplateRefiner::templates + 688 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 0 );
                    break;
                case 102:  // 1>2,3>2,3>4,4>1
                    output_tets.push( SimplexTemplateRefiner::templates + 688 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[14] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 1 );
                    break;
                case 170:  // 1>2,2>3,3>4,1>4
                    output_tets.push( SimplexTemplateRefiner::templates + 688 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[15] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 2 );
                    break;
                case 153:  // 2>1,2>3,4>3,1>4
                    output_tets.push( SimplexTemplateRefiner::templates + 688 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[5] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 3 );
                    break;
                case 90:  // 1>2,2>3,4>3,4>1
                    output_tets.push( SimplexTemplateRefiner::templates + 688 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[9] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 4 );
                    break;
                case 105:  // 2>1,2>3,3>4,4>1
                    output_tets.push( SimplexTemplateRefiner::templates + 688 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[7] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 5 );
                    break;
                case 165:  // 2>1,3>2,3>4,1>4
                    output_tets.push( SimplexTemplateRefiner::templates + 688 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[19] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 6 );
                    break;
                case 150:  // 1>2,3>2,4>3,1>4
                    output_tets.push( SimplexTemplateRefiner::templates + 688 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[23] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 7 );
                    break;
                case 101:  // 2>1,3>2,3>4,4>1
                {
                    int alternates[] = { 713, 738, -1 };
                    output_tets.push( SimplexTemplateRefiner::templates +
                                      this->best_tets( alternates, permuted_coords, 0, 1 ) );
                }
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 8 );
                    break;
                case 86:  // 1>2,3>2,4>3,4>1
                {
                    int alternates[] = { 713, 738, -1 };
                    output_tets.push( SimplexTemplateRefiner::templates +
                                      this->best_tets( alternates, permuted_coords, 14, -1 ) );
                }
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[14] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 9 );
                    break;
                case 154:  // 1>2,2>3,4>3,1>4
                {
                    int alternates[] = { 713, 738, -1 };
                    output_tets.push( SimplexTemplateRefiner::templates +
                                      this->best_tets( alternates, permuted_coords, 5, 1 ) );
                }
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[5] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 10 );
                    break;
                case 169:  // 2>1,2>3,3>4,1>4
                {
                    int alternates[] = { 713, 738, -1 };
                    output_tets.push( SimplexTemplateRefiner::templates +
                                      this->best_tets( alternates, permuted_coords, 15, -1 ) );
                }
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[15] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 11 );
                    break;
                case 89:  // 2>1,2>3,4>3,4>1
                    output_tets.push( SimplexTemplateRefiner::templates + 763 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 12 );
                    break;
                case 166:  // 1>2,3>2,3>4,1>4
                    output_tets.push( SimplexTemplateRefiner::templates + 763 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[15] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 13 );
                    break;
                case 103:  // 1=2,3>2,3>4,4>1
                    output_tets.push( SimplexTemplateRefiner::templates + 788 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 14 );
                    break;
                case 87:  // 1=2,3>2,4>3,4>1
                    output_tets.push( SimplexTemplateRefiner::templates + 788 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[14] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 15 );
                    break;
                case 185:  // 2>1,2>3,3=4,1>4
                    output_tets.push( SimplexTemplateRefiner::templates + 788 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[15] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 16 );
                    break;
                case 186:  // 1>2,2>3,3=4,1>4
                    output_tets.push( SimplexTemplateRefiner::templates + 788 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[5] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 17 );
                    break;
                case 158:  // 1>2,2=3,4>3,1>4
                    output_tets.push( SimplexTemplateRefiner::templates + 788 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[9] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 18 );
                    break;
                case 229:  // 2>1,3>2,3>4,1=4
                    output_tets.push( SimplexTemplateRefiner::templates + 788 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[7] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 19 );
                    break;
                case 233:  // 2>1,2>3,3>4,1=4
                    output_tets.push( SimplexTemplateRefiner::templates + 788 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[19] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 20 );
                    break;
                case 94:  // 1>2,2=3,4>3,4>1
                    output_tets.push( SimplexTemplateRefiner::templates + 788 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[23] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 21 );
                    break;
                case 155:  // 1=2,2>3,4>3,1>4
                    output_tets.push( SimplexTemplateRefiner::templates + 825 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 22 );
                    break;
                case 171:  // 1=2,2>3,3>4,1>4
                    output_tets.push( SimplexTemplateRefiner::templates + 825 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[14] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 23 );
                    break;
                case 118:  // 1>2,3>2,3=4,4>1
                    output_tets.push( SimplexTemplateRefiner::templates + 825 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[15] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 24 );
                    break;
                case 117:  // 2>1,3>2,3=4,4>1
                    output_tets.push( SimplexTemplateRefiner::templates + 825 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[5] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 25 );
                    break;
                case 109:  // 2>1,2=3,3>4,4>1
                    output_tets.push( SimplexTemplateRefiner::templates + 825 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[9] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 26 );
                    break;
                case 218:  // 1>2,2>3,4>3,1=4
                    output_tets.push( SimplexTemplateRefiner::templates + 825 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[7] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 27 );
                    break;
                case 214:  // 1>2,3>2,4>3,1=4
                    output_tets.push( SimplexTemplateRefiner::templates + 825 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[19] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 28 );
                    break;
                case 173:  // 2>1,2=3,3>4,1>4
                    output_tets.push( SimplexTemplateRefiner::templates + 825 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[23] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 29 );
                    break;
                case 91:  // 1=2,2>3,4>3,4>1
                    output_tets.push( SimplexTemplateRefiner::templates + 862 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 30 );
                    break;
                case 167:  // 1=2,3>2,3>4,1>4
                    output_tets.push( SimplexTemplateRefiner::templates + 862 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[14] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 31 );
                    break;
                case 182:  // 1>2,3>2,3=4,1>4
                    output_tets.push( SimplexTemplateRefiner::templates + 862 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[15] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 32 );
                    break;
                case 121:  // 2>1,2>3,3=4,4>1
                    output_tets.push( SimplexTemplateRefiner::templates + 862 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[5] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 33 );
                    break;
                case 93:  // 2>1,2=3,4>3,4>1
                    output_tets.push( SimplexTemplateRefiner::templates + 862 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[9] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 34 );
                    break;
                case 217:  // 2>1,2>3,4>3,1=4
                    output_tets.push( SimplexTemplateRefiner::templates + 862 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[7] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 35 );
                    break;
                case 230:  // 1>2,3>2,3>4,1=4
                    output_tets.push( SimplexTemplateRefiner::templates + 862 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[19] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 36 );
                    break;
                case 174:  // 1>2,2=3,3>4,1>4
                    output_tets.push( SimplexTemplateRefiner::templates + 862 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[23] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 37 );
                    break;
                case 119:  // 1=2,3>2,3=4,4>1
                    output_tets.push( SimplexTemplateRefiner::templates + 899 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 38 );
                    break;
                case 187:  // 1=2>3=4,1>4
                    output_tets.push( SimplexTemplateRefiner::templates + 899 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[15] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 39 );
                    break;
                case 222:  // 1>2,2=3,4>3,1=4
                    output_tets.push( SimplexTemplateRefiner::templates + 899 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[9] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 40 );
                    break;
                case 237:  // 2>1,2=3,3>4,1=4
                    output_tets.push( SimplexTemplateRefiner::templates + 899 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[7] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 41 );
                    break;
                case 95:  // 4>1=2=3,4>3
                    output_tets.push( SimplexTemplateRefiner::templates + 944 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 42 );
                    break;
                case 231:  // 1=2,3>2,3>4,1=4
                    output_tets.push( SimplexTemplateRefiner::templates + 944 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[14] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 43 );
                    break;
                case 190:  // 1>2=3=4,1>4
                    output_tets.push( SimplexTemplateRefiner::templates + 944 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[15] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 44 );
                    break;
                case 249:  // 2>1,2>3,3=4,1=4
                    output_tets.push( SimplexTemplateRefiner::templates + 944 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[5] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 45 );
                    break;
                case 175:  // 1=2=3>4,1>4
                    output_tets.push( SimplexTemplateRefiner::templates + 993 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 46 );
                    break;
                case 219:  // 1=2>3,4>3,1=4
                    output_tets.push( SimplexTemplateRefiner::templates + 993 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[14] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 47 );
                    break;
                case 125:  // 2>1,2=3=4>1
                    output_tets.push( SimplexTemplateRefiner::templates + 993 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[15] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 48 );
                    break;
                case 246:  // 1>2,3>2,3=4=1
                    output_tets.push( SimplexTemplateRefiner::templates + 993 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[5] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 49 );
                    break;
                case 255:  // 1=2=3=4=1
                    output_tets.push( SimplexTemplateRefiner::templates + 1042 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 8, 50 );
                    break;
            }
            break;
        case 10:  // Ruprecht-Müller Case 5
            comparison_bits = ( permlen[1] <= permlen[2] ? 1 : 0 ) | ( permlen[1] >= permlen[2] ? 2 : 0 ) |
                              ( permlen[3] <= permlen[4] ? 4 : 0 ) | ( permlen[3] >= permlen[4] ? 8 : 0 ) | 0;
            if( ( comparison_bits & 3 ) == 3 )
            {
                // Compute face point and tag
                for( int i = 0; i < 6; ++i )
                {
                    permuted_coords[10][i] =
                        ( permuted_coords[1][i] + permuted_coords[0][i] ) * .375 + permuted_coords[2][i] * .25;
                }
                ( *this->tag_assigner )( permuted_tags[1], permuted_tags[0], permuted_tags[2], permuted_tags[10] );
                permuted_hash[10] = ( *this->output_functor )( permuted_hash[1], permuted_hash[0], permuted_hash[2],
                                                               permuted_coords[10], permuted_tags[10] );
            }
            if( ( comparison_bits & 12 ) == 12 )
            {
                // Compute face point and tag
                for( int i = 0; i < 6; ++i )
                {
                    permuted_coords[11][i] =
                        ( permuted_coords[0][i] + permuted_coords[1][i] ) * .375 + permuted_coords[3][i] * .25;
                }
                ( *this->tag_assigner )( permuted_tags[0], permuted_tags[1], permuted_tags[3], permuted_tags[11] );
                permuted_hash[11] = ( *this->output_functor )( permuted_hash[0], permuted_hash[1], permuted_hash[3],
                                                               permuted_coords[11], permuted_tags[11] );
            }
            MB_TESSELLATOR_INCR_CASE_COUNT( 9 );
            output_tets.push( SimplexTemplateRefiner::templates + 1107 );
            output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
            output_sign.push( 1 );
            MB_TESSELLATOR_INCR_SUBCASE_COUNT( 9, 0 );
            switch( comparison_bits )
            {
                case 10:  // 1>2,3>4
                    output_tets.push( SimplexTemplateRefiner::templates + 1116 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 9, 1 );
                    break;
                case 5:  // 2>1,4>3
                    output_tets.push( SimplexTemplateRefiner::templates + 1116 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[14] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 9, 2 );
                    break;
                case 6:  // 1>2,4>3
                {
                    int alternates[] = { 1137, 1158, -1 };
                    output_tets.push( SimplexTemplateRefiner::templates +
                                      this->best_tets( alternates, permuted_coords, 0, 1 ) );
                }
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 9, 3 );
                    break;
                case 9:  // 2>1,3>4
                {
                    int alternates[] = { 1137, 1158, -1 };
                    output_tets.push( SimplexTemplateRefiner::templates +
                                      this->best_tets( alternates, permuted_coords, 14, -1 ) );
                }
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[14] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 9, 4 );
                    break;
                case 11:  // 1=2,3>4
                {
                    int alternates[] = { 1179, 1212, 1245, -1 };
                    output_tets.push( SimplexTemplateRefiner::templates +
                                      this->best_tets( alternates, permuted_coords, 0, 1 ) );
                }
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 9, 5 );
                    break;
                case 7:  // 1=2,4>3
                {
                    int alternates[] = { 1179, 1212, 1245, -1 };
                    output_tets.push( SimplexTemplateRefiner::templates +
                                      this->best_tets( alternates, permuted_coords, 14, -1 ) );
                }
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[14] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 9, 6 );
                    break;
                case 14:  // 3=4,1>2
                {
                    int alternates[] = { 1179, 1212, 1245, -1 };
                    output_tets.push( SimplexTemplateRefiner::templates +
                                      this->best_tets( alternates, permuted_coords, 5, 1 ) );
                }
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[5] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 9, 7 );
                    break;
                case 13:  // 3=4,2>1
                {
                    int alternates[] = { 1179, 1212, 1245, -1 };
                    output_tets.push( SimplexTemplateRefiner::templates +
                                      this->best_tets( alternates, permuted_coords, 15, -1 ) );
                }
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[15] );
                    output_sign.push( -1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 9, 8 );
                    break;
                case 15:  // 1=2,3=4
                    output_tets.push( SimplexTemplateRefiner::templates + 1278 );
                    output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
                    output_sign.push( 1 );
                    MB_TESSELLATOR_INCR_SUBCASE_COUNT( 9, 9 );
                    break;
            }
            break;
        case 11:  // Ruprecht-Müller Case 6
            MB_TESSELLATOR_INCR_CASE_COUNT( 10 );
            output_tets.push( SimplexTemplateRefiner::templates + 1319 );
            output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
            output_sign.push( 1 );
            MB_TESSELLATOR_INCR_SUBCASE_COUNT( 10, 0 );
            {
                int alternates[] = { 1336, 1353, 1370, -1 };
                output_tets.push( SimplexTemplateRefiner::templates +
                                  this->best_tets( alternates, permuted_coords, 0, 1 ) );
            }
            output_perm.push( SimplexTemplateRefiner::permutations_from_index[0] );
            output_sign.push( 1 );
            MB_TESSELLATOR_INCR_SUBCASE_COUNT( 10, 1 );
            break;
    }

    int* tets;
    int ntets;
    int* perm;
    int sgn;
#ifdef MB_DEBUG_TESSELLATOR
    if( output_tets.empty() )
    {
        cout << "Argh! Case " << C << " Perm " << P << " has no output!" << endl;
    }
#endif  // MB_DEBUG_TESSELLATOR
    while( !output_tets.empty() )
    {
        tets  = output_tets.top();
        ntets = *tets;
        tets++;
        perm = output_perm.top();
        sgn  = output_sign.top();

        output_tets.pop();
        output_perm.pop();
        output_sign.pop();

        int t;
        if( sgn > 0 )
        {
            for( t = 0; t < ntets; ++t )
            {
                this->refine_3_simplex( max_depth, permuted_coords[perm[tets[0]]], permuted_tags[perm[tets[0]]],
                                        permuted_hash[perm[tets[0]]], permuted_coords[perm[tets[1]]],
                                        permuted_tags[perm[tets[1]]], permuted_hash[perm[tets[1]]],
                                        permuted_coords[perm[tets[2]]], permuted_tags[perm[tets[2]]],
                                        permuted_hash[perm[tets[2]]], permuted_coords[perm[tets[3]]],
                                        permuted_tags[perm[tets[3]]], permuted_hash[perm[tets[3]]] );
                tets += 4;
            }
        }
        else
        {
            // we have an inverted tet... reverse the first 2 vertices
            // so the orientation is positive.
            for( t = 0; t < ntets; ++t )
            {
                this->refine_3_simplex( max_depth, permuted_coords[perm[tets[1]]], permuted_tags[perm[tets[1]]],
                                        permuted_hash[perm[tets[1]]], permuted_coords[perm[tets[0]]],
                                        permuted_tags[perm[tets[0]]], permuted_hash[perm[tets[0]]],
                                        permuted_coords[perm[tets[2]]], permuted_tags[perm[tets[2]]],
                                        permuted_hash[perm[tets[2]]], permuted_coords[perm[tets[3]]],
                                        permuted_tags[perm[tets[3]]], permuted_hash[perm[tets[3]]] );
                tets += 4;
            }
        }
    }

    return true;
}

/**\brief This is used by refine_entity to assign parametric coordinates to corners of each element.
 */
void SimplexTemplateRefiner::assign_parametric_coordinates( int num_nodes, const double* src, double* tgt )
{
    for( int i = 0; i < num_nodes; ++i, src += 3, tgt += 6 )
        for( int j = 0; j < 3; ++j )
            tgt[j] = src[j];
}

/**\brief Returns true if || a0a1 || < || b0b1 ||
 *
 * We use this to test which triangulation has the best
 * aspect ratio when there are 2 to choose from.
 */
bool SimplexTemplateRefiner::compare_Hopf_cross_string_dist( const double* a0,
                                                             const double* a1,
                                                             const double* b0,
                                                             const double* b1 )
{
    double sq_mag_a = 0.;
    double sq_mag_b = 0.;
    for( int i = 0; i < 3; ++i )
    {
        double tmp;
        tmp = a0[i] - a1[i];
        sq_mag_a += tmp * tmp;
        tmp = b0[i] - b1[i];
        sq_mag_b += tmp * tmp;
    }
    return sq_mag_a < sq_mag_b;
}

/*
 * The array below is indexed by the edge code for a tetrahedron.
 * Looking up a row with a tet's edge code will return C and P.
 * C is a configuration number and P is a permutation index.
 *
 * C is based on the case number from Ruprecht and
 * Müller's (1998) paper on adaptive tetrahedra. (The case
 * numbers are shown to the left of the row in the column
 * labeled case. The only difference is that we introduce
 * a case 3d which is part of case 3c in the paper.)
 *
 * P is an index into the permutations_from_index array below,
 * and is used to transform the current tetrahedron into
 * the canonical configuration associated with C.
 *
 * The 6-digit binary number to the left (which is shown in
 * the horribly UNconventional LSB->MSB order) is the edge
 * code for the row. The 6 digits correspond to the 6 edges
 * of the tetrahedron; a '0' implies no subdivision while
 * a '1' implies subdivision should occur. The ordering of
 * the bits is
 *
 * Edge 0-1, Edge 1-2, Edge 2-0, Edge 0-3, Edge 1-3, Edge 2-3,
 *
 * where the numbers are vertices of the tetrahedron 0-1-2-3.
 * Note that Tet 0-1-2-3 must be positive (i.e., the plane
 * specified by Triangle 0-1-2 must have a normal pointing
 * towards vertex 3, and Triangle 0-1-2's normal must be
 * calculated using the cross-product (Edge 0-1) x (Edge 0-2)).
 *
 * ===========
 * References:
 * (Ruprect and Müller, 1998) A Scheme for Edge-based Adaptive
 *   Tetrahedron Subdivision, Mathematical Visualization (eds.
 *   Hege and Polthier), pp. 61--70. Springer-Verlag. 1998.
 */
int SimplexTemplateRefiner::template_index[64][2] = {
    /*      code case      C    P */
    /* 000000  0  0  */ { 0, 0 },
    /* 100000  1  1  */ { 1, 0 },
    /* 010000  2  1  */ { 1, 1 },
    /* 110000  3  2a */ { 2, 0 },
    /* 001000  4  1  */ { 1, 2 },
    /* 101000  5  2a */ { 2, 2 },
    /* 011000  6  2a */ { 2, 1 },
    /* 111000  7  3b */ { 5, 11 },
    /* 000100  8  1  */ { 1, 10 },
    /* 100100  9  2a */ { 2, 5 },
    /* 010100 10  2b */ { 3, 1 },
    /* 110100 11  3c */ { 6, 0 },
    /* 001100 12  2a */ { 2, 10 },
    /* 101100 13  3a */ { 4, 0 },
    /* 011100 14  3d */ { 7, 2 },
    /* 111100 15  4a */ { 8, 6 },
    /* 000010 16  1  */ { 1, 6 },
    /* 100010 17  2a */ { 2, 4 },
    /* 010010 18  2a */ { 2, 8 },
    /* 110010 19  3a */ { 4, 1 },
    /* 001010 20  2b */ { 3, 2 },
    /* 101010 21  3d */ { 7, 0 },
    /* 011010 22  3c */ { 6, 1 },
    /* 111010 23  4a */ { 8, 9 },
    /* 000110 24  2a */ { 2, 3 },
    /* 100110 25  3b */ { 5, 0 },
    /* 010110 26  3d */ { 7, 4 },
    /* 110110 27  4a */ { 8, 11 },
    /* 001110 28  3c */ { 6, 10 },
    /* 101110 29  4a */ { 8, 7 },
    /* 011110 30  4b */ { 9, 0 },
    /* 111110 31  5  */ { 10, 7 },
    /* 000001 32  1  */ { 1, 7 },
    /* 100001 33  2b */ { 3, 0 },
    /* 010001 34  2a */ { 2, 7 },
    /* 110001 35  3d */ { 7, 1 },
    /* 001001 36  2a */ { 2, 11 },
    /* 101001 37  3c */ { 6, 2 },
    /* 011001 38  3a */ { 4, 2 },
    /* 111001 39  4a */ { 8, 3 },
    /* 000101 40  2a */ { 2, 9 },
    /* 100101 41  3d */ { 7, 10 },
    /* 010101 42  3c */ { 6, 7 },
    /* 110101 43  4b */ { 9, 2 },
    /* 001101 44  3b */ { 5, 7 },
    /* 101101 45  4a */ { 8, 8 },
    /* 011101 46  4a */ { 8, 4 },
    /* 111101 47  5  */ { 10, 6 },
    /* 000011 48  2a */ { 2, 6 },
    /* 100011 49  3c */ { 6, 4 },
    /* 010011 50  3b */ { 5, 1 },
    /* 110011 51  4a */ { 8, 10 },
    /* 001011 52  3d */ { 7, 7 },
    /* 101011 53  4b */ { 9, 1 },
    /* 011011 54  4a */ { 8, 5 },
    /* 111011 55  5  */ { 10, 10 },
    /* 000111 56  3a */ { 4, 10 },
    /* 100111 57  4a */ { 8, 1 },
    /* 010111 58  4a */ { 8, 2 },
    /* 110111 59  5  */ { 10, 2 },
    /* 001111 60  4a */ { 8, 0 },
    /* 101111 61  5  */ { 10, 1 },
    /* 011111 62  5  */ { 10, 0 },
    /* 111111 63  6  */ { 11, 0 },
};

/* Does this mean anything? If so, then you are either
 * superstitious or much more clever than I (or both?).
 */
/* permutation index, P:  0  1  2  3  4  5  6  7  8  9 10 11 */
/* number of references: 12  9  9  3  4  2  5  6  2  3  7  2 */

/*
 * The array below is a list of all the _positive_
 * permutations of Tetrahedron 0-1-2-3. Given a
 * permutation index, it returns a row of 14 values:
 * these are the vertex numbers of the permuted
 * tetrahedron. The first 4 values are the permuted
 * corner indices, the next 6 values are the
 * permuted edge midpoint indices, and the final
 * entries reference mid-face points inserted
 * to maintain a compatible tetrahedralization.
 *
 * There are 24 entries, 6 for each of the 4 faces of
 * the tetrahedron.
 */
int SimplexTemplateRefiner::permutations_from_index[24][14] = {
    /* corners      midpoints          face points   */
    /* POSITIVE ARRANGEMENTS                         */
    { 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13 }, /* Face 0-1-2 */
    { 1, 2, 0, 3, 5, 6, 4, 8, 9, 7, 10, 12, 13, 11 },
    { 2, 0, 1, 3, 6, 4, 5, 9, 7, 8, 10, 13, 11, 12 },

    { 0, 3, 1, 2, 7, 8, 4, 6, 9, 5, 11, 13, 12, 10 }, /* Face 0-3-1 */
    { 3, 1, 0, 2, 8, 4, 7, 9, 5, 6, 11, 12, 10, 13 },
    { 1, 0, 3, 2, 4, 7, 8, 5, 6, 9, 11, 10, 13, 12 },

    { 1, 3, 2, 0, 8, 9, 5, 4, 7, 6, 12, 11, 13, 10 }, /* Face 1-3-2 */
    { 3, 2, 1, 0, 9, 5, 8, 7, 6, 4, 12, 13, 10, 11 },
    { 2, 1, 3, 0, 5, 8, 9, 6, 4, 7, 12, 10, 11, 13 },

    { 2, 3, 0, 1, 9, 7, 6, 5, 8, 4, 13, 12, 11, 10 }, /* Face 2-3-0 */
    { 3, 0, 2, 1, 7, 6, 9, 8, 4, 5, 13, 11, 10, 12 },
    { 0, 2, 3, 1, 6, 9, 7, 4, 5, 8, 13, 10, 12, 11 },

    /* NEGATIVE ARRANGEMENTS                         */
    { 0, 2, 1, 3, 6, 5, 4, 7, 9, 8, 10, 13, 12, 11 }, /* Face 0-1-2 */
    { 2, 1, 0, 3, 5, 4, 6, 9, 8, 7, 10, 12, 11, 13 },
    { 1, 0, 2, 3, 4, 6, 5, 8, 7, 9, 10, 11, 13, 12 },

    { 0, 1, 3, 2, 4, 8, 7, 6, 5, 9, 11, 10, 12, 13 }, /* Face 0-3-1 */
    { 1, 3, 0, 2, 8, 7, 4, 5, 9, 6, 11, 12, 13, 10 },
    { 3, 0, 1, 2, 7, 4, 8, 9, 6, 5, 11, 13, 10, 12 },

    { 1, 2, 3, 0, 5, 9, 8, 4, 6, 7, 12, 10, 13, 11 }, /* Face 1-3-2 */
    { 2, 3, 1, 0, 9, 8, 5, 6, 7, 4, 12, 13, 11, 10 },
    { 3, 1, 2, 0, 8, 5, 9, 7, 4, 6, 12, 11, 10, 13 },

    { 2, 0, 3, 1, 6, 7, 9, 5, 4, 8, 13, 10, 11, 12 }, /* Face 2-3-0 */
    { 0, 3, 2, 1, 7, 9, 6, 4, 8, 5, 13, 11, 12, 10 },
    { 3, 2, 0, 1, 9, 6, 7, 8, 5, 4, 13, 12, 10, 11 } };

/*
 * Below is a list of output tetrahedra. The array is
 * generated by TessellatorGenerator.py
 * which also generates the code that references it.
 * Each set of tetrahedra begins with a single integer
 * that is the number of tetrahedra for that particular
 * case. It is followed by 5 integers for each output
 * tetrahedron; the first four numbers on each row are
 * indices of the output tetrahedron. The final number
 * is a bit vector specifying which edges of the
 * tetrahedron are internal to the parent tetrahedron
 * being decomposed.
 *
 * Multiple lists of output tetrahedra may be
 * combined to create the tessellation of a single
 * input tetrahedron.
 */

int SimplexTemplateRefiner::templates[] = {
    // case 1_0
    2,
    0,
    4,
    2,
    3,
    4,
    1,
    2,
    3,

    // case 2a_0
    1,
    3,
    4,
    5,
    1,

    // case 2a, 0>1
    2,
    0,
    4,
    2,
    3,
    4,
    5,
    2,
    3,

    // case 2a, 0=1
    4,
    10,
    3,
    0,
    4,
    10,
    3,
    4,
    5,
    10,
    3,
    5,
    2,
    10,
    3,
    2,
    0,

    // case 2b_0
    4,
    0,
    4,
    9,
    3,
    4,
    1,
    9,
    3,
    0,
    4,
    2,
    9,
    4,
    1,
    2,
    9,

    // case 3a_0
    1,
    4,
    7,
    6,
    0,

    // case 3a, 0>2>3<0
    3,
    1,
    3,
    2,
    4,
    4,
    6,
    3,
    2,
    4,
    6,
    7,
    3,

    // case 3a, 0=2>3<0
    5,
    4,
    6,
    7,
    3,
    10,
    1,
    2,
    3,
    10,
    2,
    6,
    3,
    10,
    6,
    4,
    3,
    10,
    4,
    1,
    3,

    // case 3a, 3>0=2<3
    5,
    1,
    3,
    2,
    7,
    10,
    1,
    2,
    7,
    10,
    2,
    6,
    7,
    10,
    6,
    4,
    7,
    10,
    4,
    1,
    7,

    // case 3a, 0=2=3=0
    11,
    2,
    6,
    10,
    13,
    3,
    7,
    13,
    11,
    4,
    1,
    10,
    11,
    11,
    6,
    10,
    4,
    11,
    6,
    13,
    10,
    11,
    6,
    7,
    13,
    11,
    6,
    4,
    7,
    2,
    10,
    11,
    13,
    1,
    10,
    11,
    2,
    2,
    11,
    3,
    13,
    3,
    2,
    1,
    11,

    // case 3b_0
    4,
    0,
    7,
    4,
    2,
    4,
    7,
    8,
    2,
    4,
    8,
    1,
    2,
    7,
    3,
    8,
    2,

    // case 3c, 0>1,0>3
    5,
    4,
    2,
    7,
    5,
    4,
    2,
    0,
    7,
    4,
    3,
    1,
    5,
    4,
    3,
    5,
    7,
    3,
    5,
    7,
    2,

    // case 3c, 1>0,3>0
    5,
    0,
    5,
    2,
    7,
    0,
    5,
    7,
    4,
    7,
    1,
    4,
    5,
    7,
    1,
    5,
    3,
    3,
    5,
    7,
    2,

    // case 3c, 0>1,3>0
    5,
    4,
    2,
    7,
    5,
    4,
    2,
    0,
    7,
    7,
    1,
    4,
    5,
    7,
    1,
    5,
    3,
    3,
    5,
    7,
    2,

    // case 3c, 1>0,0>3
    5,
    0,
    5,
    2,
    7,
    0,
    5,
    7,
    4,
    4,
    3,
    1,
    5,
    4,
    3,
    5,
    7,
    3,
    5,
    7,
    2,

    // case 3c, 0=1,0>3
    7,
    4,
    1,
    5,
    3,
    10,
    0,
    4,
    7,
    10,
    2,
    0,
    7,
    10,
    7,
    4,
    3,
    10,
    2,
    7,
    3,
    10,
    5,
    2,
    3,
    10,
    4,
    5,
    3,

    // case 3c, 3>0,0=1
    7,
    7,
    1,
    5,
    3,
    7,
    5,
    2,
    3,
    10,
    0,
    4,
    7,
    10,
    2,
    0,
    7,
    10,
    5,
    2,
    7,
    10,
    4,
    5,
    7,
    1,
    5,
    4,
    7,

    // case 3c, 0=1,0=3
    10,
    4,
    1,
    5,
    11,
    11,
    1,
    5,
    3,
    10,
    0,
    4,
    7,
    10,
    2,
    0,
    7,
    10,
    5,
    2,
    3,
    10,
    2,
    7,
    3,
    10,
    7,
    4,
    11,
    10,
    7,
    11,
    3,
    10,
    4,
    5,
    11,
    10,
    11,
    5,
    3,

    // case 3d, 0>4,0>2
    5,
    4,
    3,
    6,
    0,
    4,
    3,
    8,
    6,
    4,
    2,
    8,
    1,
    4,
    2,
    6,
    8,
    2,
    3,
    6,
    8,

    // case 3d, 4>0,2>0
    5,
    8,
    0,
    6,
    4,
    8,
    0,
    3,
    6,
    6,
    1,
    8,
    4,
    6,
    1,
    2,
    8,
    2,
    3,
    6,
    8,

    // case 3d, 0>4,2>0
    5,
    4,
    3,
    6,
    0,
    4,
    3,
    8,
    6,
    6,
    1,
    8,
    4,
    6,
    1,
    2,
    8,
    2,
    3,
    6,
    8,

    // case 3d, 4>0,0>2
    5,
    8,
    0,
    6,
    4,
    8,
    0,
    3,
    6,
    4,
    2,
    8,
    1,
    4,
    2,
    6,
    8,
    2,
    3,
    6,
    8,

    // case 3d, 0=4,0>2
    7,
    4,
    1,
    2,
    8,
    11,
    4,
    0,
    6,
    11,
    0,
    3,
    6,
    11,
    2,
    4,
    6,
    11,
    3,
    2,
    6,
    11,
    3,
    8,
    2,
    11,
    8,
    4,
    2,

    // case 3d, 2>0,0=4
    7,
    6,
    2,
    8,
    1,
    6,
    8,
    2,
    3,
    11,
    4,
    0,
    6,
    11,
    0,
    3,
    6,
    8,
    11,
    3,
    6,
    8,
    4,
    11,
    6,
    1,
    6,
    4,
    8,

    // case 3d, 0=4,0=2
    10,
    4,
    1,
    10,
    8,
    10,
    2,
    8,
    1,
    11,
    4,
    0,
    6,
    11,
    0,
    3,
    6,
    11,
    3,
    8,
    2,
    11,
    3,
    2,
    6,
    11,
    10,
    4,
    6,
    11,
    10,
    6,
    2,
    8,
    4,
    11,
    10,
    11,
    10,
    2,
    8,

    // case 4a_0
    2,
    7,
    8,
    9,
    3,
    7,
    9,
    8,
    6,

    // case 4a, 5>4>3
    4,
    8,
    0,
    6,
    1,
    8,
    0,
    7,
    6,
    9,
    1,
    6,
    2,
    9,
    1,
    8,
    6,

    // case 4a, 3<4>5
    4,
    8,
    0,
    6,
    1,
    8,
    0,
    7,
    6,
    8,
    2,
    6,
    9,
    8,
    2,
    1,
    6,

    // case 4a, 3>4<5
    4,
    6,
    9,
    8,
    1,
    6,
    9,
    1,
    2,
    6,
    7,
    0,
    1,
    6,
    7,
    1,
    8,

    // case 4a, 3=4>5
    6,
    6,
    7,
    0,
    11,
    6,
    0,
    1,
    11,
    6,
    7,
    11,
    8,
    6,
    11,
    1,
    8,
    1,
    2,
    6,
    8,
    2,
    6,
    8,
    9,

    // case 4a, 5>4,3=4
    6,
    6,
    7,
    0,
    11,
    6,
    0,
    1,
    11,
    6,
    7,
    11,
    8,
    6,
    11,
    1,
    8,
    1,
    2,
    6,
    9,
    1,
    6,
    8,
    9,

    // case 4a, 3=4=5
    8,
    6,
    7,
    0,
    11,
    6,
    0,
    1,
    11,
    6,
    7,
    11,
    8,
    6,
    11,
    1,
    8,
    6,
    1,
    2,
    12,
    6,
    2,
    9,
    12,
    6,
    9,
    8,
    12,
    6,
    8,
    1,
    12,

    // case 4b, 2>1,3>2,4>3,4>1
    6,
    6,
    8,
    1,
    5,
    6,
    8,
    0,
    1,
    6,
    8,
    7,
    0,
    6,
    8,
    2,
    7,
    7,
    8,
    2,
    3,
    6,
    8,
    5,
    2,

    // case 4b, 2>1,3>2,3>4,4>1
    6,
    6,
    8,
    1,
    5,
    6,
    8,
    7,
    1,
    6,
    7,
    0,
    1,
    8,
    7,
    3,
    2,
    6,
    8,
    5,
    2,
    6,
    8,
    2,
    7,

    // case 4b, 2>1,3>2,3>4,4>1, a
    6,
    7,
    8,
    1,
    5,
    6,
    5,
    7,
    1,
    6,
    7,
    0,
    1,
    8,
    7,
    3,
    2,
    7,
    8,
    5,
    2,
    6,
    5,
    2,
    7,

    // case 4b, 2>1,2>3,4>3,4>1
    6,
    6,
    8,
    5,
    2,
    6,
    8,
    2,
    3,
    6,
    8,
    3,
    7,
    6,
    8,
    7,
    0,
    6,
    8,
    0,
    1,
    6,
    8,
    1,
    5,

    // case 4b, 1=2,3>2,3>4,4>1
    9,
    10,
    6,
    0,
    7,
    10,
    1,
    5,
    8,
    10,
    0,
    1,
    7,
    10,
    7,
    1,
    8,
    6,
    7,
    10,
    8,
    6,
    10,
    5,
    8,
    6,
    2,
    7,
    8,
    6,
    5,
    2,
    8,
    7,
    8,
    2,
    3,

    // case 4b, 1=2,2>3,4>3,1>4
    9,
    10,
    6,
    0,
    7,
    10,
    1,
    5,
    8,
    10,
    0,
    1,
    8,
    10,
    7,
    0,
    8,
    6,
    7,
    10,
    8,
    6,
    10,
    5,
    8,
    6,
    3,
    7,
    8,
    6,
    5,
    3,
    8,
    6,
    5,
    2,
    3,

    // case 4b, 1=2,2>3,4>3,4>1
    9,
    10,
    6,
    0,
    7,
    10,
    1,
    5,
    8,
    10,
    0,
    1,
    8,
    10,
    7,
    0,
    8,
    6,
    7,
    10,
    8,
    6,
    10,
    5,
    8,
    6,
    3,
    7,
    8,
    6,
    5,
    2,
    8,
    6,
    2,
    3,
    8,

    // case 4b, 1=2,3>2,3=4,4>1
    11,
    10,
    6,
    0,
    7,
    10,
    1,
    5,
    8,
    10,
    0,
    1,
    11,
    10,
    11,
    1,
    8,
    10,
    0,
    11,
    7,
    10,
    7,
    11,
    8,
    6,
    7,
    10,
    8,
    6,
    10,
    5,
    8,
    6,
    2,
    7,
    8,
    6,
    5,
    2,
    8,
    7,
    8,
    2,
    3,

    // case 4b, 4>1=2=3,4>3
    12,
    10,
    6,
    0,
    7,
    10,
    1,
    5,
    8,
    10,
    0,
    1,
    8,
    10,
    7,
    0,
    8,
    13,
    6,
    2,
    5,
    13,
    3,
    7,
    8,
    13,
    2,
    3,
    8,
    13,
    2,
    8,
    5,
    6,
    7,
    10,
    8,
    6,
    10,
    5,
    8,
    6,
    13,
    7,
    8,
    6,
    5,
    13,
    8,

    // case 4b, 1=2=3>4,1>4
    12,
    10,
    6,
    0,
    7,
    10,
    1,
    5,
    8,
    10,
    0,
    1,
    7,
    10,
    7,
    1,
    8,
    13,
    6,
    2,
    5,
    13,
    3,
    7,
    8,
    13,
    2,
    3,
    5,
    13,
    3,
    8,
    5,
    6,
    7,
    10,
    8,
    6,
    10,
    5,
    8,
    6,
    13,
    7,
    8,
    6,
    5,
    13,
    8,

    // case 4b, 1=2=3=4=1
    16,
    10,
    6,
    0,
    7,
    10,
    1,
    5,
    8,
    10,
    0,
    1,
    11,
    10,
    11,
    1,
    8,
    10,
    0,
    11,
    7,
    10,
    7,
    11,
    8,
    13,
    6,
    2,
    5,
    13,
    3,
    7,
    8,
    13,
    2,
    3,
    12,
    13,
    2,
    12,
    5,
    13,
    12,
    3,
    8,
    13,
    12,
    5,
    8,
    6,
    7,
    10,
    8,
    6,
    10,
    5,
    8,
    6,
    5,
    13,
    8,
    6,
    13,
    7,
    8,

    // case 5_0
    2,
    7,
    8,
    9,
    3,
    6,
    5,
    2,
    9,

    // case 5, 1>2,3>4
    5,
    5,
    7,
    1,
    8,
    5,
    7,
    0,
    1,
    5,
    7,
    6,
    0,
    5,
    7,
    9,
    6,
    5,
    7,
    8,
    9,

    // case 5, 1>2,4>3
    5,
    0,
    5,
    6,
    7,
    0,
    5,
    7,
    8,
    0,
    5,
    8,
    1,
    5,
    7,
    9,
    6,
    5,
    7,
    8,
    9,

    // case 5, 1>2,4>3, a
    5,
    0,
    5,
    6,
    8,
    0,
    6,
    7,
    8,
    0,
    5,
    8,
    1,
    5,
    8,
    9,
    6,
    6,
    7,
    8,
    9,

    // case 5, 1=2,3>4
    8,
    10,
    6,
    0,
    7,
    10,
    1,
    5,
    8,
    10,
    0,
    1,
    7,
    10,
    7,
    1,
    8,
    10,
    8,
    5,
    9,
    10,
    6,
    7,
    9,
    10,
    7,
    8,
    9,
    10,
    5,
    6,
    9,

    // case 5, 1=2,3>4, a
    8,
    10,
    6,
    0,
    7,
    10,
    1,
    5,
    8,
    10,
    0,
    1,
    7,
    10,
    7,
    1,
    8,
    7,
    8,
    5,
    9,
    10,
    6,
    7,
    5,
    10,
    7,
    8,
    5,
    5,
    9,
    6,
    7,

    // case 5, 1=2,3>4, b
    8,
    10,
    6,
    0,
    7,
    10,
    1,
    5,
    8,
    10,
    0,
    1,
    7,
    10,
    7,
    1,
    8,
    6,
    8,
    5,
    9,
    10,
    6,
    7,
    8,
    10,
    6,
    8,
    5,
    8,
    9,
    6,
    7,

    // case 5, 1=2,3=4
    10,
    10,
    6,
    0,
    7,
    10,
    1,
    5,
    8,
    10,
    0,
    1,
    11,
    10,
    11,
    1,
    8,
    10,
    0,
    11,
    7,
    10,
    7,
    11,
    8,
    10,
    8,
    5,
    9,
    10,
    6,
    7,
    9,
    10,
    7,
    8,
    9,
    10,
    5,
    6,
    9,

    // case 6_0
    4,
    7,
    8,
    9,
    3,
    6,
    5,
    2,
    9,
    4,
    1,
    5,
    8,
    0,
    4,
    6,
    7,

    // case 6_1
    4,
    6,
    4,
    5,
    8,
    6,
    5,
    9,
    8,
    6,
    9,
    7,
    8,
    6,
    7,
    4,
    8,

    // case 6_1, a
    4,
    5,
    8,
    9,
    7,
    5,
    9,
    6,
    7,
    5,
    6,
    4,
    7,
    5,
    4,
    8,
    7,

    // case 6_1, b
    4,
    4,
    5,
    6,
    9,
    4,
    6,
    7,
    9,
    4,
    7,
    8,
    9,
    4,
    8,
    5,
    9,

};

}  // namespace moab