MOAB: Mesh Oriented datABase
(version 5.4.1)
|
#include <iostream>
#include <cstdlib>
#include "TestUtil.hpp"
#include "Mesquite.hpp"
#include "MeshImpl.hpp"
#include "MsqError.hpp"
#include "Vector3D.hpp"
#include "InstructionQueue.hpp"
#include "PatchData.hpp"
#include "TerminationCriterion.hpp"
#include "QualityAssessor.hpp"
#include "IdealShapeTarget.hpp"
#include "TQualityMetric.hpp"
#include "TInverseMeanRatio.hpp"
#include "ConditionNumberQualityMetric.hpp"
#include "LPtoPTemplate.hpp"
#include "LInfTemplate.hpp"
#include "SteepestDescent.hpp"
#include "ConjugateGradient.hpp"
#include "PlanarDomain.hpp"
Go to the source code of this file.
Functions | |
void | compare_nodes (size_t start_index, size_t end_index, Mesh *mesh1, Mesh *mesh2, MsqError &err) |
InstructionQueue * | create_instruction_queue (MsqError &err) |
int | do_test (bool slave) |
int | do_smooth_ho () |
int | main () |
Variables | |
std::string | LINEAR_INPUT_FILE_NAME = "/linear_input.vtk" |
std::string | QUADRATIC_INPUT_FILE_NAME = "/quadratic_input.vtk" |
std::string | EXPECTED_LINAR_FILE_NAME = "/expected_linear_output.vtk" |
std::string | EXPECTED_QUADRATIC_FILE_NAME = "/expected_quadratic_output.vtk" |
std::string | HOUR_INPUT_FILE_NAME = "/hour-quad8.vtk" |
std::string | OUTPUT_FILE_NAME = "smoothed_qudratic_mesh.vtk" |
const unsigned | NUM_CORNER_VERTICES = 16 |
const unsigned | NUM_MID_NODES = 24 |
const double | SPATIAL_COMPARE_TOLERANCE = 4e-6 |
void compare_nodes | ( | size_t | start_index, |
size_t | end_index, | ||
Mesh * | mesh1, | ||
Mesh * | mesh2, | ||
MsqError & | err | ||
) |
Definition at line 120 of file test/mesquite/higher_order/main.cpp.
References MBMesquite::arrptr(), MBMesquite::MsqError::INTERNAL_ERROR, MBMesquite::length(), MSQ_ERRRTN, MSQ_SETERR, SPATIAL_COMPARE_TOLERANCE, and MBMesquite::Mesh::vertices_get_coordinates().
Referenced by do_test().
{ size_t i, num_verts = end_index - start_index; std::vector< MsqVertex > verts1( num_verts ), verts2( num_verts ); std::vector< Mesh::VertexHandle > handles1( num_verts ), handles2( num_verts ); /* VertexIterator skips higher-order nodes. For now, just assume index == handle std::vector<Mesh::VertexHandle>::iterator handle_iter1, handle_iter2; // Skip start_index vertices VertexIterator* iter1 = mesh1->vertex_iterator( err ); MSQ_ERRRTN(err); VertexIterator* iter2 = mesh2->vertex_iterator( err ); MSQ_ERRRTN(err); for (i = 0; i < start_index; ++i) { if (iter1->is_at_end()) { MSQ_SETERR(err)("start index out of range for first mesh set", MsqError::INVALID_ARG); return; } if (iter2->is_at_end()) { MSQ_SETERR(err)("start index out of range for second mesh set", MsqError::INVALID_ARG); return; } iter1->operator++(); iter2->operator++(); } // Get handles for vertices handle_iter1 = handles1.begin(); handle_iter2 = handles2.begin(); for (i = start_index; i < end_index; ++i) { if (iter1->is_at_end()) { MSQ_SETERR(err)("end index out of range for first mesh set", MsqError::INVALID_ARG); return; } *handle_iter1 = iter1->operator*(); iter1->operator++(); ++handle_iter1; if (iter2->is_at_end()) { MSQ_SETERR(err)("end index out of range for second mesh set", MsqError::INVALID_ARG); return; } *handle_iter2 = iter2->operator*(); iter2->operator++(); ++handle_iter2; } */ for( i = start_index; i < end_index; ++i ) handles1[i - start_index] = handles2[i - start_index] = (void*)i; // Get coordinates from handles mesh1->vertices_get_coordinates( arrptr( handles1 ), arrptr( verts1 ), num_verts, err );MSQ_ERRRTN( err ); mesh2->vertices_get_coordinates( arrptr( handles2 ), arrptr( verts2 ), num_verts, err );MSQ_ERRRTN( err ); // Compare coordinates for( i = 0; i < num_verts; ++i ) { const double diff = ( verts1[i] - verts2[i] ).length(); if( diff > SPATIAL_COMPARE_TOLERANCE ) { MSQ_SETERR( err ) ( MsqError::INTERNAL_ERROR, "%u%s vertices differ. (%f,%f,%f) vs (%f,%f,%f)", (unsigned)( 1 + i ), i % 10 == 0 ? "st" : i % 10 == 1 ? "nd" : i % 10 == 2 ? "rd" : "th", verts1[i][0], verts1[i][1], verts1[i][2], verts2[i][0], verts2[i][1], verts2[i][2] ); return; } } }
InstructionQueue* create_instruction_queue | ( | MsqError & | err | ) |
Definition at line 200 of file test/mesquite/higher_order/main.cpp.
References MBMesquite::TerminationCriterion::add_absolute_vertex_movement(), MBMesquite::TerminationCriterion::add_iteration_limit(), MSQ_ERRZERO, MBMesquite::QualityImprover::set_inner_termination_criterion(), MBMesquite::InstructionQueue::set_master_quality_improver(), MBMesquite::QualityImprover::set_outer_termination_criterion(), and MBMesquite::PatchSetUser::use_global_patch().
Referenced by do_smooth_ho(), and do_test().
{ // creates an intruction queue InstructionQueue* queue1 = new InstructionQueue; // creates a mean ratio quality metric ... // IdealWeightInverseMeanRatio* mean = new IdealWeightInverseMeanRatio(err); MSQ_ERRZERO(err); TargetCalculator* tc = new IdealShapeTarget; TQualityMetric* mean = new TQualityMetric( tc, 0, new TInverseMeanRatio ); LPtoPTemplate* obj_func = new LPtoPTemplate( mean, 1, err ); MSQ_ERRZERO( err ); // creates the optimization procedures // ConjugateGradient* pass1 = new ConjugateGradient( obj_func, err ); SteepestDescent* pass1 = new SteepestDescent( obj_func ); // perform optimization globally pass1->use_global_patch(); // QualityAssessor* mean_qa = new QualityAssessor(mean); //**************Set termination criterion**************** // perform 1 pass of the outer loop (this line isn't essential as it is // the default behavior). TerminationCriterion* tc_outer = new TerminationCriterion; tc_outer->add_iteration_limit( 1 ); pass1->set_outer_termination_criterion( tc_outer ); // perform the inner loop until a certain objective function value is // reached. The exact value needs to be determined (about 18095). // As a safety, also stop if the time exceeds 10 minutes (600 seconds). TerminationCriterion* tc_inner = new TerminationCriterion; tc_inner->add_absolute_vertex_movement( 1e-6 ); pass1->set_inner_termination_criterion( tc_inner ); // adds 1 pass of pass1 to mesh_set1 // queue1->add_quality_assessor(mean_qa,err); MSQ_ERRZERO(err); queue1->set_master_quality_improver( pass1, err ); MSQ_ERRZERO( err ); // queue1->add_quality_assessor(mean_qa,err); MSQ_ERRZERO(err); return queue1; }
int do_smooth_ho | ( | ) |
Definition at line 330 of file test/mesquite/higher_order/main.cpp.
References create_instruction_queue(), geom, HOUR_INPUT_FILE_NAME, MSQ_CHKERR, MBMesquite::MeshImpl::read_vtk(), MBMesquite::IQInterface::run_instructions(), MBMesquite::Settings::set_slaved_ho_node_mode(), MBMesquite::Settings::SLAVE_NONE, MBMesquite::MeshImpl::write_vtk(), and MBMesquite::PlanarDomain::XY.
Referenced by main().
{ MsqPrintError err( cout ); // QuadLagrangeShape quad9; // Create geometry PlanarDomain geom( PlanarDomain::XY ); // Read in one copy of quadratic input mesh cout << "Reading " << HOUR_INPUT_FILE_NAME << endl; MeshImpl* quadratic_in = new MeshImpl; quadratic_in->read_vtk( HOUR_INPUT_FILE_NAME.c_str(), err ); if( MSQ_CHKERR( err ) ) return 1; // Read in expected results // cout << "Reading " << HOUR_EXPECTED_FILE_NAME << endl; // MeshImpl* quadratic_ex = new MeshImpl; // quadratic_ex->read_vtk( results, err ); // if (MSQ_CHKERR(err)) return 1; // Smooth linear mesh and check results cout << "Smoothing higher-order nodes" << endl; InstructionQueue* q1 = create_instruction_queue( err ); if( MSQ_CHKERR( err ) ) return 1; q1->set_slaved_ho_node_mode( Settings::SLAVE_NONE ); // q1->set_mapping_function( &quad9 ); MeshDomainAssoc mesh_and_domain = MeshDomainAssoc( quadratic_in, &geom ); q1->run_instructions( &mesh_and_domain, err ); if( MSQ_CHKERR( err ) ) return 1; cout << "Checking results" << endl; // compare_nodes( 0, NUM_CORNER_VERTICES + NUM_MID_NODES, // quadratic_in, quadratic_ex, err ); if( MSQ_CHKERR( err ) ) { MsqError tmperr; quadratic_in->write_vtk( "bad_mesh.vtk", tmperr ); return 1; } delete q1; quadratic_in->write_vtk( "smooth_ho.vtk", err ); delete quadratic_in; if( MSQ_CHKERR( err ) ) return 1; return 0; }
int do_test | ( | bool | slave | ) |
Definition at line 246 of file test/mesquite/higher_order/main.cpp.
References compare_nodes(), create_instruction_queue(), EXPECTED_LINAR_FILE_NAME, EXPECTED_QUADRATIC_FILE_NAME, geom, LINEAR_INPUT_FILE_NAME, MSQ_CHKERR, NUM_CORNER_VERTICES, NUM_MID_NODES, QUADRATIC_INPUT_FILE_NAME, MBMesquite::MeshImpl::read_vtk(), MBMesquite::IQInterface::run_instructions(), MBMesquite::Settings::set_slaved_ho_node_mode(), MBMesquite::Settings::SLAVE_NONE, MBMesquite::MeshImpl::write_vtk(), and z.
Referenced by main(), IdealTargetTest::test_hex_center(), IdealTargetTest::test_hex_corner(), IdealTargetTest::test_hex_edge(), IdealTargetTest::test_hex_face(), IdealTargetTest::test_tri_center(), IdealTargetTest::test_tri_corner(), and IdealTargetTest::test_tri_edge().
{ MsqPrintError err( cout ); // QuadLagrangeShape quad9; // Create geometry Vector3D z( 0, 0, 1 ), o( 0, 0, 0 ); PlanarDomain geom( z, o ); // Read in linear input mesh cout << "Reading " << LINEAR_INPUT_FILE_NAME << endl; MeshImpl* linear_in = new MeshImpl; linear_in->read_vtk( LINEAR_INPUT_FILE_NAME.c_str(), err ); if( MSQ_CHKERR( err ) ) return 1; // Read in expected linear results cout << "Reading " << EXPECTED_LINAR_FILE_NAME << endl; MeshImpl* linear_ex = new MeshImpl; linear_ex->read_vtk( EXPECTED_LINAR_FILE_NAME.c_str(), err ); if( MSQ_CHKERR( err ) ) return 1; // Read in second copy of quadratic input mesh cout << "Reading " << QUADRATIC_INPUT_FILE_NAME << " again" << endl; MeshImpl* quadratic_in_2 = new MeshImpl; quadratic_in_2->read_vtk( QUADRATIC_INPUT_FILE_NAME.c_str(), err ); if( MSQ_CHKERR( err ) ) return 1; // Read in expected quadratic results cout << "Reading " << EXPECTED_QUADRATIC_FILE_NAME << endl; MeshImpl* quadratic_ex = new MeshImpl; quadratic_ex->read_vtk( EXPECTED_QUADRATIC_FILE_NAME.c_str(), err ); if( MSQ_CHKERR( err ) ) return 1; // Smooth linear mesh and check results cout << "Smoothing linear elements" << endl; InstructionQueue* q1 = create_instruction_queue( err ); if( MSQ_CHKERR( err ) ) return 1; MeshDomainAssoc* mesh_and_domain = new MeshDomainAssoc( linear_in, &geom ); q1->run_instructions( mesh_and_domain, err ); if( MSQ_CHKERR( err ) ) return 1; delete mesh_and_domain; cout << "Checking results" << endl; compare_nodes( 0, NUM_CORNER_VERTICES, linear_in, linear_ex, err ); if( MSQ_CHKERR( err ) ) { MsqError tmperr; linear_in->write_vtk( "bad_mesh.vtk", tmperr ); return 1; } delete q1; delete linear_in; // Smooth corner vertices and adjust mid-side nodes cout << "Smoothing quadratic elements" << endl; InstructionQueue* q3 = create_instruction_queue( err ); if( MSQ_CHKERR( err ) ) return 1; if( !slave ) q3->set_slaved_ho_node_mode( Settings::SLAVE_NONE ); // q3->set_mapping_function( &quad9 ); MeshDomainAssoc* mesh_and_domain2 = new MeshDomainAssoc( quadratic_in_2, &geom ); q3->run_instructions( mesh_and_domain2, err ); if( MSQ_CHKERR( err ) ) return 1; delete mesh_and_domain2; // Make sure corner vertices are the same as in the linear case cout << "Checking results" << endl; compare_nodes( 0, NUM_CORNER_VERTICES, quadratic_in_2, linear_ex, err ); if( MSQ_CHKERR( err ) ) return 1; delete linear_ex; // Make sure mid-side vertices are updated correctly compare_nodes( NUM_CORNER_VERTICES, NUM_CORNER_VERTICES + NUM_MID_NODES, quadratic_in_2, quadratic_ex, err ); if( MSQ_CHKERR( err ) ) { MsqError tmperr; quadratic_in_2->write_vtk( "bad_mesh.vtk", tmperr ); return 1; } delete q3; delete quadratic_in_2; delete quadratic_ex; if( MSQ_CHKERR( err ) ) return 1; return 0; }
int main | ( | ) |
Definition at line 376 of file test/mesquite/higher_order/main.cpp.
References do_smooth_ho(), and do_test().
{ cout << "Running test with all higher-order nodes slaved." << endl; int result1 = do_test( true ); cout << "Running test with no higher-order nodes slaved." << endl; int result2 = do_test( false ); cout << "Running test with only ho-nodes free." << endl; int result3 = do_smooth_ho(); return result1 + result2 + result3; }
std::string EXPECTED_LINAR_FILE_NAME = "/expected_linear_output.vtk" |
Definition at line 112 of file test/mesquite/higher_order/main.cpp.
Referenced by do_test().
std::string EXPECTED_QUADRATIC_FILE_NAME = "/expected_quadratic_output.vtk" |
Definition at line 113 of file test/mesquite/higher_order/main.cpp.
Referenced by do_test().
std::string HOUR_INPUT_FILE_NAME = "/hour-quad8.vtk" |
Definition at line 114 of file test/mesquite/higher_order/main.cpp.
Referenced by do_smooth_ho().
std::string LINEAR_INPUT_FILE_NAME = "/linear_input.vtk" |
Definition at line 110 of file test/mesquite/higher_order/main.cpp.
Referenced by do_test().
const unsigned NUM_CORNER_VERTICES = 16 |
Definition at line 116 of file test/mesquite/higher_order/main.cpp.
Referenced by do_test().
const unsigned NUM_MID_NODES = 24 |
Definition at line 117 of file test/mesquite/higher_order/main.cpp.
Referenced by do_test().
std::string OUTPUT_FILE_NAME = "smoothed_qudratic_mesh.vtk" |
Definition at line 115 of file test/mesquite/higher_order/main.cpp.
std::string QUADRATIC_INPUT_FILE_NAME = "/quadratic_input.vtk" |
Definition at line 111 of file test/mesquite/higher_order/main.cpp.
Referenced by do_test().
const double SPATIAL_COMPARE_TOLERANCE = 4e-6 |
Definition at line 118 of file test/mesquite/higher_order/main.cpp.
Referenced by compare_nodes().