MOAB: Mesh Oriented datABase  (version 5.4.1)
ParShapeImprover::ParShapeImprovementWrapper Class Reference
+ Inheritance diagram for ParShapeImprover::ParShapeImprovementWrapper:
+ Collaboration diagram for ParShapeImprover::ParShapeImprovementWrapper:

Public Member Functions

 ParShapeImprovementWrapper (int inner_iterations=100, double cpu_time=0.0, double grad_norm=1.e-8, int parallel_iterations=10)

Public Attributes

bool m_do_untangle_only

Protected Member Functions

void run_wrapper (MeshDomainAssoc *mesh_and_domain, ParallelMesh *pmesh, Settings *settings, QualityAssessor *qa, MsqError &err)

Private Attributes

int innerIter
double maxTime
double gradNorm
const double untBeta
const double successiveEps
int parallelIterations

Detailed Description

Definition at line 122 of file par_hex_untangle_shape.cpp.


Constructor & Destructor Documentation

ParShapeImprover::ParShapeImprovementWrapper::ParShapeImprovementWrapper ( int  inner_iterations = 100,
double  cpu_time = 0.0,
double  grad_norm = 1.e-8,
int  parallel_iterations = 10 
) [inline]

Definition at line 127 of file par_hex_untangle_shape.cpp.

            : innerIter( inner_iterations ), maxTime( cpu_time ), gradNorm( grad_norm ), untBeta( DEF_UNT_BETA ),
              successiveEps( DEF_SUC_EPS ), parallelIterations( parallel_iterations ), m_do_untangle_only( false )
        {
        }

Member Function Documentation

void ParShapeImprover::ParShapeImprovementWrapper::run_wrapper ( MeshDomainAssoc mesh_and_domain,
ParallelMesh pmesh,
Settings settings,
QualityAssessor quality_assessor,
MsqError err 
) [protected, virtual]

Function that each wrapper must implement

term_outer.add_relative_successive_improvement( successiveEps );

Implements MBMesquite::Wrapper.

Definition at line 160 of file par_hex_untangle_shape.cpp.

References MBMesquite::TerminationCriterion::add_absolute_gradient_L2_norm(), MBMesquite::TerminationCriterion::add_absolute_quality_improvement(), MBMesquite::TerminationCriterion::add_absolute_vertex_movement(), MBMesquite::TerminationCriterion::add_iteration_limit(), MBMesquite::QualityAssessor::add_quality_assessment(), MBMesquite::InstructionQueue::add_quality_assessor(), MBMesquite::MeshDomainAssoc::get_domain(), MBMesquite::MeshDomainAssoc::get_mesh(), MBMesquite::QualityMetric::LINEAR, mesh, MPI_COMM_WORLD, MSQ_DBGOUT, MSQ_ERRRTN, rank, MBMesquite::InstructionQueue::run_common(), MBMesquite::IQInterface::run_instructions(), MBMesquite::AveragingQM::set_averaging_method(), MBMesquite::QualityImprover::set_inner_termination_criterion(), MBMesquite::InstructionQueue::set_master_quality_improver(), MBMesquite::QualityImprover::set_outer_termination_criterion(), MBMesquite::Timer::since_birth(), and MBMesquite::PatchSetUser::use_global_patch().

{
    int rank, nprocs;
    MPI_Comm_rank( MPI_COMM_WORLD, &rank );
    MPI_Comm_size( MPI_COMM_WORLD, &nprocs );

    // Define an untangler
    // UntangleBetaQualityMetric untangle_metric( untBeta );
    UntangleBetaQualityMetric untangle_metric( 1.e-6 );

    bool check_untangle = true;
    if( check_untangle )
    {
        std::cout << "\nP[" << rank << "]  ParShapeImprover.... running QA with untangle_metric before... "
                  << std::endl;
        InstructionQueue q1;
        QualityAssessor qa_untangle( &untangle_metric );
        q1.add_quality_assessor( &qa_untangle, err );MSQ_ERRRTN( err );
        q1.run_common( mesh_and_domain, pmesh, settings, err );
        std::cout << "\nP[" << rank << "]  ParShapeImprover.... running QA with untangle_metric... before... done "
                  << std::endl;
    }

    LPtoPTemplate untangle_func( 2, &untangle_metric );
    ConjugateGradient untangle_solver( &untangle_func );
    // untangle_solver.set_debugging_level(3);

    // SteepestDescent untangle_solver( &untangle_func );
    TerminationCriterion untangle_inner( "<type:untangle_inner>" ), untangle_outer( "<type:untangle_outer>" );
    untangle_solver.use_global_patch();

    // untangle_inner.add_absolute_gradient_L2_norm( gradNorm );
    // untangle_inner.add_absolute_successive_improvement( successiveEps );
    // untangle_inner.add_relative_successive_improvement( 1.e-6 );
    // untangle_inner.add_untangled_mesh();

    untangle_inner.write_iterations( "untangle.gpt", err );

    // For parallel runs, we generally need to have the inner and outer TerminationCriterion
    //   have the same criteria else we can get an infinite loop (see VertexMover::loop_over_mesh)
    untangle_inner.add_absolute_quality_improvement( 0.0 );
    untangle_inner.add_iteration_limit( 20 );

    untangle_outer.add_absolute_quality_improvement( 0.0 );
    untangle_outer.add_iteration_limit( pmesh ? parallelIterations : 1 );

    untangle_solver.set_inner_termination_criterion( &untangle_inner );
    untangle_solver.set_outer_termination_criterion( &untangle_outer );

    // define shape improver
    IdealWeightInverseMeanRatio inverse_mean_ratio;
    inverse_mean_ratio.set_averaging_method( QualityMetric::LINEAR );
    LPtoPTemplate obj_func( 2, &inverse_mean_ratio );

    ConjugateGradient shape_solver( &obj_func );
    TerminationCriterion term_inner( "<type:shape_inner>" ), term_outer( "<type:shape_outer>" );
    term_inner.write_iterations( "shape.gpt", err );

    shape_solver.use_global_patch();
    qa->add_quality_assessment( &inverse_mean_ratio );

    // For parallel runs, we generally need to have the inner and outer TerminationCriterion
    //   have the same criteria else we can get an infinite loop (see VertexMover::loop_over_mesh)
    term_inner.add_absolute_gradient_L2_norm( gradNorm );
    term_inner.add_absolute_vertex_movement( 0.0 );
    term_inner.add_iteration_limit( innerIter );

    term_outer.add_absolute_gradient_L2_norm( gradNorm );
    term_outer.add_absolute_vertex_movement( 0.0 );
    term_outer.add_iteration_limit( pmesh ? parallelIterations : 1 );

    // term_outer.add_absolute_quality_improvement( 1.e-6 );
    //! term_outer.add_relative_successive_improvement( successiveEps );

    shape_solver.set_inner_termination_criterion( &term_inner );
    shape_solver.set_outer_termination_criterion( &term_outer );

    // Apply CPU time limit to untangler
    if( maxTime > 0.0 ) untangle_inner.add_cpu_time( maxTime );

    Timer totalTimer;

    // Run untangler
    std::cout << "\nP[" << rank << "] "
              << " ParShapeImprovementWrapper: running untangler...\n " << std::endl;
    bool use_untangle_wrapper = false;
    if( use_untangle_wrapper )
    {
        UntangleWrapper uw;
        // uw.set_untangle_metric(UntangleWrapper::BETA);
        uw.run_instructions( mesh_and_domain, err );
    }
    else
    {
        InstructionQueue q1;
        QualityAssessor qa_untangle( &untangle_metric );
        q1.add_quality_assessor( &qa_untangle, err );MSQ_ERRRTN( err );
        q1.set_master_quality_improver( &untangle_solver, err );MSQ_ERRRTN( err );
        q1.add_quality_assessor( &qa_untangle, err );MSQ_ERRRTN( err );
        q1.run_common( mesh_and_domain, pmesh, settings, err );
    }
    std::cout << "\nP[" << rank << "] "
              << " ParShapeImprovementWrapper: running untangler... done\n " << std::endl;
    std::cout << "\nP[" << rank << "] "
              << " ParShapeImprovementWrapper: MsqError after untangler: " << err << std::endl;

    bool check_quality_after_untangler = true;
    if( check_quality_after_untangler )
    {
        Mesh* mesh         = mesh_and_domain->get_mesh();
        MeshDomain* domain = mesh_and_domain->get_domain();
        int num_invalid    = count_invalid_elements( *mesh, domain );
        std::cout << "\nP[" << rank << "] "
                  << " ParShapeImprover num_invalid after untangler= " << num_invalid << " "
                  << ( num_invalid ? " ERROR still have invalid elements after MBMesquite untangle"
                                   : " SUCCESS: untangled invalid elements " )
                  << std::endl;

        if( check_untangle )
        {
            std::cout << "\nP[" << rank << "]  ParShapeImprover.... running QA with untangle_metric " << std::endl;
            InstructionQueue q1;
            QualityAssessor qa_untangle( &untangle_metric );
            q1.add_quality_assessor( &qa_untangle, err );MSQ_ERRRTN( err );
            q1.run_common( mesh_and_domain, pmesh, settings, err );
            std::cout << "\nP[" << rank << "]  ParShapeImprover.... running QA with untangle_metric... done "
                      << std::endl;
        }

        if( num_invalid ) return;
    }
    if( m_do_untangle_only ) return;MSQ_ERRRTN( err );

    // If limited by CPU time, limit next step to remaning time
    if( maxTime > 0.0 )
    {
        double remaining = maxTime - totalTimer.since_birth();
        if( remaining <= 0.0 )
        {
            MSQ_DBGOUT( 2 ) << "Optimization is terminating without perfoming shape improvement." << std::endl;
            remaining = 0.0;
        }
        term_inner.add_cpu_time( remaining );
    }

    // Run shape improver
    InstructionQueue q2;
    std::cout << "\nP[" << rank << "] "
              << " ParShapeImprovementWrapper: running shape improver... \n"
              << std::endl;
    q2.add_quality_assessor( qa, err );MSQ_ERRRTN( err );
    q2.set_master_quality_improver( &shape_solver, err );MSQ_ERRRTN( err );
    q2.add_quality_assessor( qa, err );MSQ_ERRRTN( err );
    q2.run_common( mesh_and_domain, pmesh, settings, err );
    std::cout << "\nP[" << rank << "] "
              << " ParShapeImprovementWrapper: running shape improver... done \n"
              << std::endl;MSQ_ERRRTN( err );
}

Member Data Documentation

Definition at line 147 of file par_hex_untangle_shape.cpp.

List of all members.


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