MOAB: Mesh Oriented datABase
(version 5.4.1)
|
#include "Mesquite.hpp"
#include "MsqError.hpp"
#include "MeshImpl.hpp"
#include "XYRectangle.hpp"
#include "ReferenceMesh.hpp"
#include "RefMeshTargetCalculator.hpp"
#include "TShapeB1.hpp"
#include "TQualityMetric.hpp"
#include "TerminationCriterion.hpp"
#include "QualityAssessor.hpp"
#include "PMeanPTemplate.hpp"
#include "FeasibleNewton.hpp"
#include "ConjugateGradient.hpp"
#include "InstructionQueue.hpp"
#include <iostream>
#include <fstream>
#include <cstdlib>
#include <cstdio>
Go to the source code of this file.
Classes | |
struct | MeshParams |
Struct in which to store mesh description. More... | |
Defines | |
#define | CHECKERR |
Enumerations | |
enum | ExitCodes { NO_ERROR = 0, USAGE_ERROR = 1, NOT_AT_TARGET = 2, FAR_FROM_TARGET = 3, WRONG_DIRECTION = 4, DEGENERATE_ELEMENT = 5, INVERTED_ELEMENT = 6, LAST_EXIT_CODE = INVERTED_ELEMENT } |
enum | ParseState { OPEN, EXPECTING_M, EXPECTING_R, EXPECTING_O, EXPECTING_I } |
Functions | |
std::ostream & | operator<< (std::ostream &str, const MeshParams &p) |
void | usage (const char *argv0, bool brief=true) |
void | create_input_mesh (const MeshParams ¶ms, bool all_fixed, MeshImpl &mesh, MsqError &err) |
void | parse_options (char *argv[], int argc, MeshParams &mesh, MeshParams &ref, std::string &output_file, bool &fixed_boundary, TerminationCriterion::TimeStepFileType &write_timesteps, bool &use_feas_newt, int &num_iterations) |
std::string | base_name (std::string filename) |
int | main (int argc, char *argv[]) |
void | parse_mesh_params (const char *argv, const char *arg, MeshParams &result) |
Variables | |
const char | default_out_file [] = "high_aspect.vtk" |
Default values for parameters. | |
const MeshParams | default_mesh = { 0.6, 0.3, 2.0, 1.0 } |
const MeshParams | default_ref = { 0.5, 0.5, 1.0, 1.0 } |
const int | INNER_ITERATES = 1 |
const int | OUTER_ITERATES = 10 |
const char * | temp_file = "high_aspect_input.vtk" |
#define CHECKERR |
if( err ) \ { \ std::cerr << "Internal error at line " << __LINE__ << ":" << std::endl << ( err ) << std::endl; \ std::exit( LAST_EXIT_CODE + ( err ).error_code() ); \ }
Definition at line 115 of file high_aspect_ratio_test.cpp.
Referenced by main().
enum ExitCodes |
NO_ERROR | |
USAGE_ERROR | |
NOT_AT_TARGET | |
FAR_FROM_TARGET | |
WRONG_DIRECTION | |
DEGENERATE_ELEMENT | |
INVERTED_ELEMENT | |
LAST_EXIT_CODE |
Definition at line 74 of file high_aspect_ratio_test.cpp.
{ NO_ERROR = 0, USAGE_ERROR = 1, NOT_AT_TARGET = 2, FAR_FROM_TARGET = 3, WRONG_DIRECTION = 4, DEGENERATE_ELEMENT = 5, INVERTED_ELEMENT = 6, LAST_EXIT_CODE = INVERTED_ELEMENT };
enum ParseState |
Definition at line 290 of file high_aspect_ratio_test.cpp.
{ OPEN, EXPECTING_M, EXPECTING_R, EXPECTING_O, EXPECTING_I };
std::string base_name | ( | std::string | filename | ) |
void create_input_mesh | ( | const MeshParams & | params, |
bool | all_fixed, | ||
MeshImpl & | mesh, | ||
MsqError & | err | ||
) |
Definition at line 391 of file high_aspect_ratio_test.cpp.
References moab::F, MeshParams::h, MSQ_CHKERR, MBMesquite::MeshImpl::read_vtk(), temp_file, MeshParams::w, MeshParams::x, MeshParams::y, and z.
Referenced by main().
{ const double z = 0; const int F = all_fixed; std::ofstream vtkfile( temp_file ); double bx = all_fixed ? 0.5 * p.w : p.x; double by = all_fixed ? 0.5 * p.h : p.y; vtkfile << "# vtk DataFile Version 3.0" << std::endl << "Mesquite High Aspect Ratio test" << std::endl << "ASCII" << std::endl << "DATASET UNSTRUCTURED_GRID" << std::endl << "POINTS 9 float" << std::endl << 0.0 << ' ' << 0.0 << ' ' << z << std::endl << bx << ' ' << 0.0 << ' ' << z << std::endl << p.w << ' ' << 0.0 << ' ' << z << std::endl << 0.0 << ' ' << by << ' ' << z << std::endl << p.x << ' ' << p.y << ' ' << z << std::endl << p.w << ' ' << by << ' ' << z << std::endl << 0.0 << ' ' << p.h << ' ' << z << std::endl << bx << ' ' << p.h << ' ' << z << std::endl << p.w << ' ' << p.h << ' ' << z << std::endl << "CELLS 4 20" << std::endl << "4 0 1 4 3" << std::endl << "4 1 2 5 4" << std::endl << "4 4 5 8 7" << std::endl << "4 3 4 7 6" << std::endl << "CELL_TYPES 4" << std::endl << "9 9 9 9" << std::endl << "POINT_DATA 9" << std::endl << "SCALARS fixed int" << std::endl << "LOOKUP_TABLE default" << std::endl << "1 " << F << " 1" << std::endl << F << " 0 " << F << std::endl << "1 " << F << " 1" << std::endl; vtkfile.close(); mesh.read_vtk( temp_file, err ); std::remove( temp_file );MSQ_CHKERR( err ); }
int main | ( | int | argc, |
char * | argv[] | ||
) |
Definition at line 154 of file high_aspect_ratio_test.cpp.
References MBMesquite::TerminationCriterion::add_iteration_limit(), MBMesquite::InstructionQueue::add_quality_assessor(), MBMesquite::arrptr(), base_name(), CHECKERR, create_input_mesh(), DEGENERATE_ELEMENT, dof(), MBMesquite::XYRectangle::domain_DoF(), EPS, FAR_FROM_TARGET, MBMesquite::MeshImpl::get_all_vertices(), MBMesquite::QualityAssessor::get_inverted_element_count(), MeshParams::h, init(), MBMesquite::inner(), INNER_ITERATES, MBMesquite::inv(), INVERTED_ELEMENT, MBMesquite::length(), mesh, NOT_AT_TARGET, MBMesquite::TerminationCriterion::NOTYPE, MBMesquite::outer(), parse_options(), MBMesquite::IQInterface::run_instructions(), MBMesquite::QualityImprover::set_inner_termination_criterion(), MBMesquite::InstructionQueue::set_master_quality_improver(), MBMesquite::QualityImprover::set_outer_termination_criterion(), MBMesquite::XYRectangle::setup(), USAGE_ERROR, MBMesquite::PatchSetUser::use_element_on_vertex_patch(), MBMesquite::MeshImpl::vertices_get_coordinates(), MeshParams::w, MBMesquite::TerminationCriterion::write_mesh_steps(), MBMesquite::MeshImpl::write_vtk(), WRONG_DIRECTION, MeshParams::x, and MeshParams::y.
{ MeshParams input_params, reference_params; bool fixed_boundary_vertices, feas_newt_solver; TerminationCriterion::TimeStepFileType write_timestep_files; std::string output_file_name; int num_iterations; parse_options( argv, argc, input_params, reference_params, output_file_name, fixed_boundary_vertices, write_timestep_files, feas_newt_solver, num_iterations ); MsqError err; MeshImpl mesh, refmesh; XYRectangle domain( input_params.w, input_params.h ); create_input_mesh( input_params, fixed_boundary_vertices, mesh, err ); CHECKERR create_input_mesh( reference_params, fixed_boundary_vertices, refmesh, err ); CHECKERR domain.setup( &mesh, err ); CHECKERR ReferenceMesh rmesh( &refmesh ); RefMeshTargetCalculator tc( &rmesh ); TShapeB1 tm; TQualityMetric qm( &tc, &tm ); PMeanPTemplate of( 1.0, &qm ); ConjugateGradient cg( &of ); cg.use_element_on_vertex_patch(); FeasibleNewton fn( &of ); fn.use_element_on_vertex_patch(); VertexMover* solver = feas_newt_solver ? (VertexMover*)&fn : (VertexMover*)&cg; TerminationCriterion inner, outer; inner.add_iteration_limit( INNER_ITERATES ); outer.add_iteration_limit( num_iterations ); if( write_timestep_files != TerminationCriterion::NOTYPE ) outer.write_mesh_steps( base_name( output_file_name ).c_str(), write_timestep_files ); solver->set_inner_termination_criterion( &inner ); solver->set_outer_termination_criterion( &outer ); QualityAssessor qa( &qm ); InstructionQueue q; q.add_quality_assessor( &qa, err ); q.set_master_quality_improver( solver, err ); q.add_quality_assessor( &qa, err ); MeshDomainAssoc mesh_and_domain = MeshDomainAssoc( &mesh, &domain ); q.run_instructions( &mesh_and_domain, err ); CHECKERR mesh.write_vtk( output_file_name.c_str(), err ); CHECKERR // check for inverted elements int inv, unk; qa.get_inverted_element_count( inv, unk, err ); if( inv ) { std::cerr << inv << " inverted elements in final mesh" << std::endl; return INVERTED_ELEMENT; } else if( unk ) { std::cerr << unk << " degenerate elements in final mesh" << std::endl; return DEGENERATE_ELEMENT; } // find the free vertex std::vector< Mesh::VertexHandle > vertices; mesh.get_all_vertices( vertices, err ); if( vertices.empty() ) { std::cerr << "Mesh contains no vertices" << std::endl; return USAGE_ERROR; } std::vector< unsigned short > dof( vertices.size(), 0 ); domain.domain_DoF( arrptr( vertices ), arrptr( dof ), vertices.size(), err ); CHECKERR int idx = std::find( dof.begin(), dof.end(), 2 ) - dof.begin(); const Mesh::VertexHandle free_vertex = vertices[idx]; MsqVertex coords; mesh.vertices_get_coordinates( &free_vertex, &coords, 1, err ); CHECKERR // calculate optimal position for vertex const double xf = reference_params.x / reference_params.w; const double yf = reference_params.y / reference_params.h; Vector3D expect( xf * input_params.w, yf * input_params.h, 0 ); // Test that we aren't further from the expected location // than when we started. const Vector3D init( input_params.x, input_params.y, 0 ); if( ( coords - expect ).length() > ( init - expect ).length() ) { std::cerr << "Vertex moved away from expected optimal position: " << "(" << coords[0] << ", " << coords[1] << std::endl; return WRONG_DIRECTION; } // check if vertex moved MIN_FRACT of the way from the original position // to the desired one in the allowed iterations const double MIN_FRACT = 0.2; // 20% of the way in 10 iterations const double fract = ( coords - init ).length() / ( expect - init ).length(); if( fract < MIN_FRACT ) { std::cerr << "Vertex far from optimimal location" << std::endl << " Expected: (" << expect[0] << ", " << expect[1] << ", " << expect[2] << ")" << std::endl << " Actual: (" << coords[0] << ", " << coords[1] << ", " << coords[2] << ")" << std::endl; return FAR_FROM_TARGET; } // check if vertex is at destired location const double EPS = 5e-2; // allow a little leway if( fabs( coords[0] - expect[0] ) > EPS * input_params.w || fabs( coords[1] - expect[1] ) > EPS * input_params.h || fabs( expect[2] ) > EPS ) { std::cerr << "Vertex not at optimimal location" << std::endl << " Expected: (" << expect[0] << ", " << expect[1] << ", " << expect[2] << ")" << std::endl << " Actual: (" << coords[0] << ", " << coords[1] << ", " << coords[2] << ")" << std::endl; return NOT_AT_TARGET; } return 0; }
std::ostream& operator<< | ( | std::ostream & | str, |
const MeshParams & | p | ||
) |
Definition at line 61 of file high_aspect_ratio_test.cpp.
References MeshParams::h, MeshParams::w, MeshParams::x, and MeshParams::y.
void parse_mesh_params | ( | const char * | argv, |
const char * | arg, | ||
MeshParams & | result | ||
) |
Definition at line 280 of file high_aspect_ratio_test.cpp.
References MeshParams::h, usage, MeshParams::w, MeshParams::x, and MeshParams::y.
Referenced by parse_options().
void parse_options | ( | char * | argv[], |
int | argc, | ||
MeshParams & | mesh, | ||
MeshParams & | ref, | ||
std::string & | output_file, | ||
bool & | fixed_boundary, | ||
TerminationCriterion::TimeStepFileType & | write_timesteps, | ||
bool & | use_feas_newt, | ||
int & | num_iterations | ||
) |
Definition at line 298 of file high_aspect_ratio_test.cpp.
References default_mesh, default_out_file, default_ref, EXPECTING_I, EXPECTING_M, EXPECTING_O, EXPECTING_R, MBMesquite::TerminationCriterion::GNUPLOT, MBMesquite::TerminationCriterion::NOTYPE, OPEN, OUTER_ITERATES, parse_mesh_params(), usage, and MBMesquite::TerminationCriterion::VTK.
Referenced by main().
{ // begin with defaults mesh = default_mesh; ref = default_ref; output_file = default_out_file; fixed_boundary = false; write_timesteps = TerminationCriterion::NOTYPE; feas_newt_solver = false; num_iterations = OUTER_ITERATES; // parse CLI args ParseState state = OPEN; for( int i = 1; i < argc; ++i ) { switch( state ) { case EXPECTING_M: parse_mesh_params( argv[0], argv[i], mesh ); state = OPEN; break; case EXPECTING_R: parse_mesh_params( argv[0], argv[i], ref ); state = OPEN; break; case EXPECTING_O: output_file = argv[i]; state = OPEN; break; case EXPECTING_I: num_iterations = atoi( argv[i] ); state = OPEN; break; case OPEN: if( argv[i][0] != '-' || argv[i][1] == '\0' || argv[i][2] != '\0' ) { std::cerr << "Unexpected argument: \"" << argv[i] << '"' << std::endl; usage( argv[0] ); } switch( argv[i][1] ) { default: usage( argv[0], true ); break; case 'h': usage( argv[0], false ); break; case 'o': state = EXPECTING_O; break; case 'f': fixed_boundary = true; break; case 'F': fixed_boundary = false; break; case 't': write_timesteps = TerminationCriterion::VTK; break; case 'T': write_timesteps = TerminationCriterion::GNUPLOT; break; case 'm': state = EXPECTING_M; break; case 'r': state = EXPECTING_R; break; case 'n': feas_newt_solver = true; break; case 'c': feas_newt_solver = false; break; case 'i': state = EXPECTING_I; break; } break; } } }
void usage | ( | const char * | argv0, |
bool | brief = true |
||
) |
Definition at line 86 of file high_aspect_ratio_test.cpp.
References default_mesh, default_out_file, default_ref, NO_ERROR, OUTER_ITERATES, and USAGE_ERROR.
{ std::ostream& str = brief ? std::cerr : std::cout; str << "Usage: " << argv0 << " [-o <output_file>]" << " [-f|-F] [-t|-T] [-n|-c] [-i <n>]" << " [-m <x>,<y>[,<w>,<h>]]" << " [-r <x>,<y>[,<w>,<h>]]" << std::endl; if( brief ) { str << " " << argv0 << " -h" << std::endl; std::exit( USAGE_ERROR ); } str << " -o Specify output file (default is \"" << default_out_file << "\")" << std::endl << " -f Fixed boundary vertices" << std::endl << " -F Free boundary vertices (default)" << std::endl << " -t Write VTK timesteps" << std::endl << " -T Write GNUPlot timesteps" << std::endl << " -m Specify input mesh parameters (default " << default_mesh << ")" << std::endl << " -r Specify reference mesh parameters (default " << default_ref << ")" << std::endl << " -n Use FeasibleNewton solver" << std::endl << " -c Use ConjugateGradient solver (default)" << std::endl << " -i Specify number of iterations (default:" << OUTER_ITERATES << ")" << std::endl << std::endl; std::exit( NO_ERROR ); }
const MeshParams default_mesh = { 0.6, 0.3, 2.0, 1.0 } |
Definition at line 68 of file high_aspect_ratio_test.cpp.
Referenced by parse_options(), and usage().
const char default_out_file[] = "high_aspect.vtk" |
Default values for parameters.
Definition at line 67 of file high_aspect_ratio_test.cpp.
Referenced by parse_options(), and usage().
const MeshParams default_ref = { 0.5, 0.5, 1.0, 1.0 } |
Definition at line 69 of file high_aspect_ratio_test.cpp.
Referenced by parse_options(), and usage().
const int INNER_ITERATES = 1 |
Definition at line 71 of file high_aspect_ratio_test.cpp.
Referenced by main().
const int OUTER_ITERATES = 10 |
Definition at line 72 of file high_aspect_ratio_test.cpp.
Referenced by parse_options(), and usage().
const char* temp_file = "high_aspect_input.vtk" |
Definition at line 390 of file high_aspect_ratio_test.cpp.
Referenced by create_input_mesh().