![]() |
Mesh Oriented datABase
(version 5.4.1)
Array-based unstructured mesh datastructure
|
#include "moab/ParallelComm.hpp"
#include "MBParallelConventions.h"
#include "moab/Core.hpp"
#include "moab/FileOptions.hpp"
#include "ReadParallel.hpp"
#include "Coupler.hpp"
#include "DebugOutput.hpp"
#include "ElemUtil.hpp"
#include <iostream>
#include <iomanip>
#include <sstream>
#include <cstring>
#include <cstdlib>
#include "moab/FindPtFuncs.h"
#include "moab/TupleList.hpp"
#include "moab/gs.hpp"
#include "moab/Types.hpp"
#include "Internals.hpp"
Go to the source code of this file.
Defines | |
#define | IS_BUILDING_MB |
#define | UINT_PER_X(X) ( ( sizeof( X ) + sizeof( uint ) - 1 ) / sizeof( uint ) ) |
#define | UINT_PER_REAL UINT_PER_X( realType ) |
#define | UINT_PER_LONG UINT_PER_X( slong ) |
#define | UINT_PER_ULONG UINT_PER_X( Ulong ) |
#define | UINT_PER_UNSIGNED UINT_PER_X( unsigned ) |
Functions | |
void | get_file_options (int argc, char **argv, std::vector< const char * > &filenames, std::string &norm_tag, std::vector< const char * > &tag_names, std::vector< const char * > &tag_values, std::string &file_opts, int *err) |
void | print_tuples (TupleList *tlp) |
int | print_vertex_fields (Interface *mbi, std::vector< std::vector< EntityHandle > > &groups, Tag &norm_hdl, Coupler::IntegType integ_type) |
double | const_field (double x, double y, double z) |
double | field_1 (double x, double y, double z) |
double | field_2 (double x, double y, double z) |
double | field_3 (double x, double y, double z) |
double | physField (double x, double y, double z) |
ErrorCode | integrate_scalar_field_test () |
int | pack_tuples (TupleList *tl, void **ptr) |
void | unpack_tuples (void *ptr, TupleList **tlp) |
int | main (int argc, char **argv) |
Variables | |
bool | debug = true |
#define IS_BUILDING_MB |
Definition at line 22 of file ssn_test.cpp.
#define UINT_PER_LONG UINT_PER_X( slong ) |
Definition at line 907 of file ssn_test.cpp.
#define UINT_PER_REAL UINT_PER_X( realType ) |
Definition at line 906 of file ssn_test.cpp.
#define UINT_PER_ULONG UINT_PER_X( Ulong ) |
Definition at line 908 of file ssn_test.cpp.
#define UINT_PER_UNSIGNED UINT_PER_X( unsigned ) |
Definition at line 909 of file ssn_test.cpp.
#define UINT_PER_X | ( | X | ) | ( ( sizeof( X ) + sizeof( uint ) - 1 ) / sizeof( uint ) ) |
Definition at line 905 of file ssn_test.cpp.
double const_field | ( | double | x, |
double | y, | ||
double | z | ||
) |
Definition at line 863 of file ssn_test.cpp.
Referenced by integrate_scalar_field_test().
{
// return 5.0/40.0;
return 5.0;
}
double field_1 | ( | double | x, |
double | y, | ||
double | z | ||
) |
Definition at line 870 of file ssn_test.cpp.
Referenced by integrate_scalar_field_test().
{
double f = fabs( x ) + fabs( y ) + fabs( z );
// return f/24.0;
return f;
}
double field_2 | ( | double | x, |
double | y, | ||
double | z | ||
) |
Definition at line 877 of file ssn_test.cpp.
Referenced by integrate_scalar_field_test().
{
double f = x * x + y * y + z * z;
// return f/32.0;
return f;
}
double field_3 | ( | double | x, |
double | y, | ||
double | z | ||
) |
Definition at line 884 of file ssn_test.cpp.
Referenced by integrate_scalar_field_test().
{
double f = 2 * x + 2 * y + 2 * z;
// return f/48.0;
return f;
}
void get_file_options | ( | int | argc, |
char ** | argv, | ||
std::vector< const char * > & | filenames, | ||
std::string & | norm_tag, | ||
std::vector< const char * > & | tag_names, | ||
std::vector< const char * > & | tag_values, | ||
std::string & | file_opts, | ||
int * | err | ||
) |
Definition at line 653 of file ssn_test.cpp.
References moab::debug.
{
int npos = 1;
// get number of files
int nfiles = atoi( argv[npos++] );
// get mesh filenames
filenames.resize( nfiles );
for( int i = 0; i < nfiles; i++ )
filenames[i] = argv[npos++];
// get normTag
if( npos < argc )
normTag = argv[npos++];
else
{
std::cerr << "Insufficient parameters: norm_tag missing" << std::endl;
*err = 1;
return;
}
// get tag selection options
if( npos < argc )
{
char* opts = argv[npos++];
// char sep1[1] = {';'};
// char sep2[1] = {'='};
bool end_vals_seen = false;
std::vector< char* > tmpTagOpts;
// first get the options
for( char* i = strtok( opts, ";" ); i; i = strtok( 0, ";" ) )
{
if( debug ) std::cout << "get_file_options: i=" << i << std::endl;
tmpTagOpts.push_back( i );
}
// parse out the name and val or just name.
for( unsigned int j = 0; j < tmpTagOpts.size(); j++ )
{
char* e = strtok( tmpTagOpts[j], "=" );
if( debug ) std::cout << "get_file_options: name=" << e << std::endl;
tagNames.push_back( e );
e = strtok( 0, "=" );
if( e != NULL )
{
if( debug ) std::cout << "get_file_options: val=" << e << std::endl;
// We have a value
if( end_vals_seen )
{
// ERROR we should not have a value after none are seen
std::cerr << "Incorrect parameters: new value seen after end of values" << std::endl;
*err = 1;
return;
}
// Otherwise get the value string from e and convert it to an int
int* valp = new int;
*valp = atoi( e );
tagValues.push_back( (const char*)valp );
}
else
{
// Otherwise there is no '=' so push a null on the list
end_vals_seen = true;
tagValues.push_back( (const char*)0 );
}
}
}
else
{
std::cerr << "Insufficient parameters: tag_select_opts missing" << std::endl;
*err = 1;
return;
}
// get fileOpts
if( npos < argc )
fileOpts = argv[npos++];
else
{
std::cerr << "Insufficient parameters: file_opts missing" << std::endl;
*err = 1;
return;
}
}
Definition at line 566 of file ssn_test.cpp.
References const_field(), field_1(), field_2(), field_3(), moab::Element::LinearHex::integrate_scalar_field(), MB_CHK_SET_ERR, and MB_SUCCESS.
Referenced by main().
{
// ******************************
std::cout << "********** Test moab::Element::Map::integrate_scalar_field **********" << std::endl;
// Create a simple hex centered at 0,0,0 with sides of length 2.
std::vector< CartVect > biunit_cube( 8 );
biunit_cube[0] = CartVect( -1, -1, -1 );
biunit_cube[1] = CartVect( 1, -1, -1 );
biunit_cube[2] = CartVect( 1, 1, -1 );
biunit_cube[3] = CartVect( -1, 1, -1 );
biunit_cube[4] = CartVect( -1, -1, 1 );
biunit_cube[5] = CartVect( 1, -1, 1 );
biunit_cube[6] = CartVect( 1, 1, 1 );
biunit_cube[7] = CartVect( -1, 1, 1 );
std::vector< CartVect > zerobase_cube( 8 );
zerobase_cube[0] = CartVect( 0, 0, 0 );
zerobase_cube[1] = CartVect( 2, 0, 0 );
zerobase_cube[2] = CartVect( 2, 2, 0 );
zerobase_cube[3] = CartVect( 0, 2, 0 );
zerobase_cube[4] = CartVect( 0, 0, 2 );
zerobase_cube[5] = CartVect( 2, 0, 2 );
zerobase_cube[6] = CartVect( 2, 2, 2 );
zerobase_cube[7] = CartVect( 0, 2, 2 );
// Calculate field values at the corners of both cubes
double bcf[8], bf1[8], bf2[8], bf3[8], zcf[8], zf1[8], zf2[8], zf3[8];
for( int i = 0; i < 8; i++ )
{
bcf[i] = const_field( biunit_cube[i][0], biunit_cube[i][1], biunit_cube[i][2] );
bf1[i] = field_1( biunit_cube[i][0], biunit_cube[i][1], biunit_cube[i][2] );
bf2[i] = field_2( biunit_cube[i][0], biunit_cube[i][1], biunit_cube[i][2] );
bf3[i] = field_3( biunit_cube[i][0], biunit_cube[i][1], biunit_cube[i][2] );
zcf[i] = const_field( zerobase_cube[i][0], zerobase_cube[i][1], zerobase_cube[i][2] );
zf1[i] = field_1( zerobase_cube[i][0], zerobase_cube[i][1], zerobase_cube[i][2] );
zf2[i] = field_2( zerobase_cube[i][0], zerobase_cube[i][1], zerobase_cube[i][2] );
zf3[i] = field_3( zerobase_cube[i][0], zerobase_cube[i][1], zerobase_cube[i][2] );
}
std::cout << "Integrated values:" << std::endl;
try
{
double field_const1, field_const2;
double field_linear1, field_linear2;
double field_quad1, field_quad2;
double field_cubic1, field_cubic2;
int ipoints = 0;
Element::LinearHex biunit_hexMap( biunit_cube );
Element::LinearHex zerobase_hexMap( zerobase_cube );
field_const1 = biunit_hexMap.integrate_scalar_field( bcf );
field_const2 = zerobase_hexMap.integrate_scalar_field( zcf );
std::cout << " binunit_cube, const_field(num_pts=" << ipoints << "): field_val=" << field_const1
<< std::endl;
std::cout << " zerobase_cube, const_field(num_pts=" << ipoints << "): field_val=" << field_const2
<< std::endl;
field_linear1 = biunit_hexMap.integrate_scalar_field( bf1 );
field_linear2 = zerobase_hexMap.integrate_scalar_field( zf1 );
std::cout << " binunit_cube, field_1(num_pts=" << ipoints << "): field_val=" << field_linear1 << std::endl;
std::cout << " zerobase_cube, field_1(num_pts=" << ipoints << "): field_val=" << field_linear2 << std::endl;
field_quad1 = biunit_hexMap.integrate_scalar_field( bf2 );
field_quad2 = zerobase_hexMap.integrate_scalar_field( zf2 );
std::cout << " binunit_cube, field_2(num_pts=" << ipoints << "): field_val=" << field_quad1 << std::endl;
std::cout << " zerobase_cube, field_2(num_pts=" << ipoints << "): field_val=" << field_quad2 << std::endl;
field_cubic1 = biunit_hexMap.integrate_scalar_field( bf3 );
field_cubic2 = zerobase_hexMap.integrate_scalar_field( zf3 );
std::cout << " binunit_cube, field_3(num_pts=" << ipoints << "): field_val=" << field_cubic1 << std::endl;
std::cout << " zerobase_cube, field_3(num_pts=" << ipoints << "): field_val=" << field_cubic2 << std::endl;
}
catch( moab::Element::Map::ArgError )
{
MB_CHK_SET_ERR( MB_FAILURE, "Failed to set vertices on Element::Map." );
}
catch( moab::Element::Map::EvaluationError )
{
MB_CHK_SET_ERR( MB_FAILURE, "Failed to get inverse evaluation of coordinate on Element::Map." );
}
return MB_SUCCESS;
}
int main | ( | int | argc, |
char ** | argv | ||
) |
Definition at line 63 of file ssn_test.cpp.
References moab::Coupler::apply_group_norm_factor(), moab::Range::clear(), moab::Coupler::consolidate_tuples(), moab::Interface::create_meshset(), moab::Coupler::create_tuples(), moab::debug, moab::TupleList::enableWriteAccess(), ErrorCode, moab::Interface::get_entities_by_handle(), moab::Interface::get_entities_by_type_and_tag(), moab::point_locator::io::get_file_options(), moab::Coupler::get_group_integ_vals(), moab::Coupler::get_matching_entities(), moab::TupleList::initialize(), moab::Range::insert(), integrate_scalar_field_test(), moab::Interface::INTERSECT, MB_CHK_SET_ERR, MBENTITYSET, MESHSET_SET, moab::Coupler::normalize_subset(), moab::pack_tuples(), moab::DebugOutput::print(), print_tuples(), print_vertex_fields(), moab::TupleList::set_n(), moab::DebugOutput::set_rank(), moab::DebugOutput::set_verbosity(), moab::Interface::tag_get_data(), moab::Interface::tag_get_handle(), moab::unpack_tuples(), moab::TupleList::vi_wr, moab::TupleList::vl_wr, moab::Coupler::VOLUME, moab::TupleList::vr_wr, and moab::TupleList::vul_wr.
{
// need to init MPI first, to tell how many procs and rank
// Used since Coupler is a parallel code. The Coupler may be run
// in parallel or serial mode but will need MPI either way.
int err = MPI_Init( &argc, &argv );
// Print usage if not enough arguments
if( argc < 3 )
{
std::cerr << "Usage: ";
std::cerr << argv[0] << " ... " << std::endl;
std::cerr << "nfiles : number of mesh files" << std::endl;
std::cerr << "fname1...fnamen : mesh files" << std::endl;
std::cerr << "norm_tag : name of tag to normalize across meshes" << std::endl;
std::cerr << "tag_select_opts : quoted string of tags and values for subset selection, "
"e.g. \"TAG1=VAL1;TAG2=VAL2;TAG3;TAG4\""
<< std::endl;
std::cerr << "file_opts : quoted string of parallel file read options, e.g. "
"\"OPTION1=VALUE1;OPTION2;OPTION3=VALUE3\""
<< std::endl;
err = integrate_scalar_field_test();MB_CHK_SET_ERR( (ErrorCode)err, "Integrate scalar field test failed" );
err = MPI_Finalize();
return err;
}
int nprocs, rank;
err = MPI_Comm_size( MPI_COMM_WORLD, &nprocs );
assert( MPI_SUCCESS == err );
err = MPI_Comm_rank( MPI_COMM_WORLD, &rank );
assert( MPI_SUCCESS == err );
// Create an ofstream to write output. One file each for each proc.
std::stringstream fname;
fname << argv[0] << rank << ".out";
if( !std::freopen( fname.str().c_str(), "a", stdout ) ) return -1;
if( !std::freopen( fname.str().c_str(), "a", stderr ) ) return -1;
// Create the moab instance
Interface* mbi = new Core();
if( NULL == mbi )
{
std::cerr << "MOAB constructor failed" << std::endl;
return 1;
}
// Get the input options
std::cout << "Getting options..." << std::endl;
std::vector< const char* > filenames;
std::vector< const char* > tagNames;
std::vector< const char* > tagValues;
std::string normTag, fileOpts;
get_file_options( argc, argv, filenames, normTag, tagNames, tagValues, fileOpts, &err );MB_CHK_SET_ERR( (ErrorCode)err, "get_file_options failed" );
// Print out the input parameters
std::cout << " Input Parameters - " << std::endl;
std::cout << " Filenames: ";
for( std::vector< const char* >::iterator it = filenames.begin(); it != filenames.end(); ++it )
std::cout << *it << " ";
std::cout << std::endl;
std::cout << " Norm Tag: " << normTag << std::endl;
std::cout << " Selection Data: NumNames=" << tagNames.size() << " NumValues=" << tagValues.size() << std::endl;
std::cout << " TagNames TagValues " << std::endl;
std::cout << " -------------------- --------------------" << std::endl;
std::vector< const char* >::iterator nameIt = tagNames.begin();
std::vector< const char* >::iterator valIt = tagValues.begin();
std::cout << std::setiosflags( std::ios::left );
for( ; nameIt != tagNames.end(); ++nameIt )
{
std::cout << " " << std::setw( 20 ) << *nameIt;
if( *valIt != 0 )
{
std::cout << " " << std::setw( 20 ) << *( (int*)( *valIt ) ) << std::endl;
++valIt;
}
else
std::cout << " NULL " << std::endl;
}
std::cout << std::resetiosflags( std::ios::left );
std::cout << " File Options: " << fileOpts << std::endl;
// Read in mesh(es)
std::cout << "Reading mesh file(s)..." << std::endl;
std::vector< ParallelComm* > pcs( filenames.size() );
std::vector< ReadParallel* > rps( filenames.size() );
// allocate root sets for each mesh for moab
std::vector< EntityHandle > roots( filenames.size() );
ErrorCode result;
for( unsigned int i = 0; i < filenames.size(); i++ )
{
pcs[i] = new ParallelComm( mbi, MPI_COMM_WORLD );
rps[i] = new ReadParallel( mbi, pcs[i] );
result = mbi->create_meshset( MESHSET_SET, roots[i] );
MB_CHK_SET_ERR( result, "Creating root set failed" );
result = rps[i]->load_file( filenames[i], &roots[i], FileOptions( fileOpts.c_str() ) );MB_CHK_SET_ERR( result, "load_file failed" );
}
// Initialize the debug object for Range printing
DebugOutput debugOut( "ssn_test-", std::cerr );
debugOut.set_rank( rank );
debugOut.set_verbosity( 10 );
// Output what is in root sets
for( unsigned int k = 0; k < filenames.size(); k++ )
{
Range rootRg;
result = mbi->get_entities_by_handle( roots[k], rootRg );MB_CHK_SET_ERR( result, "can't get entities" );
debugOut.print( 2, "Root set entities: ", rootRg );
rootRg.clear();
Range partRg;
pcs[k]->get_part_entities( partRg );
debugOut.print( 2, "Partition entities: ", partRg );
partRg.clear();
}
// source is 1st mesh, target is 2nd mesh
Range src_elems, targ_elems;
// ******************************
std::cout << "********** Create Coupler **********" << std::endl;
// Create a coupler
std::cout << "Creating Coupler..." << std::endl;
Coupler mbc( mbi, pcs[0], src_elems, 0 );
// Get tag handles for passed in tags
std::cout << "Getting tag handles..." << std::endl;
int numTagNames = tagNames.size();
std::vector< Tag > tagHandles( numTagNames );
int iTags = 0;
while( iTags < numTagNames )
{
std::cout << "Getting handle for " << tagNames[iTags] << std::endl;
result = mbi->tag_get_handle( tagNames[iTags], tagHandles[iTags] );MB_CHK_SET_ERR( result, "Retrieving tag handles failed" );
iTags++;
}
// ******************************
std::cout << "********** Test create_tuples **********" << std::endl;
// First get some EntitySets for Mesh 1 and Mesh 2
{
Range entsets1, entsets2;
result = mbi->get_entities_by_type_and_tag( roots[0], MBENTITYSET, &tagHandles[0],
(const void* const*)&tagValues[0], tagHandles.size(), entsets1,
Interface::INTERSECT ); // recursive is false
MB_CHK_SET_ERR( result, "sets: get_entities_by_type_and_tag failed on Mesh 1." );
// Create tuple_list for each mesh's
std::cout << "Creating tuples for mesh 1..." << std::endl;
TupleList* m1TagTuples = NULL;
err = mbc.create_tuples( entsets1, &tagHandles[0], tagHandles.size(), &m1TagTuples );MB_CHK_SET_ERR( (ErrorCode)err, "create_tuples failed" );
std::cout << " create_tuples returned" << std::endl;
print_tuples( m1TagTuples );
result = mbi->get_entities_by_type_and_tag( roots[1], MBENTITYSET, &tagHandles[0],
(const void* const*)&tagValues[0], tagHandles.size(), entsets2,
Interface::INTERSECT ); // recursive is false
MB_CHK_SET_ERR( result, "sets: get_entities_by_type_and_tag failed on Mesh 2." );
std::cout << "Creating tuples for mesh 2..." << std::endl;
TupleList* m2TagTuples = NULL;
err = mbc.create_tuples( entsets2, (Tag*)( &tagHandles[0] ), tagHandles.size(), &m2TagTuples );MB_CHK_SET_ERR( (ErrorCode)err, "create_tuples failed" );
std::cout << " create_tuples returned" << std::endl;
print_tuples( m2TagTuples );
// ******************************
std::cout << "********** Test consolidate_tuples **********" << std::endl;
// In this serial version we only have the tuples from Mesh 1 and Mesh 2.
// Just consolidate those for the test.
std::cout << "Consolidating tuple_lists for Mesh 1 and Mesh 2..." << std::endl;
TupleList** tplp_arr = (TupleList**)malloc( 2 * sizeof( TupleList* ) );
TupleList* unique_tpl = NULL;
tplp_arr[0] = m1TagTuples;
tplp_arr[1] = m2TagTuples;
err = mbc.consolidate_tuples( tplp_arr, 2, &unique_tpl );MB_CHK_SET_ERR( (ErrorCode)err, "consolidate_tuples failed" );
std::cout << " consolidate_tuples returned" << std::endl;
print_tuples( unique_tpl );
}
// ******************************
std::cout << "********** Test get_matching_entities **********" << std::endl;
std::vector< std::vector< EntityHandle > > m1EntitySets;
std::vector< std::vector< EntityHandle > > m1EntityGroups;
std::vector< std::vector< EntityHandle > > m2EntitySets;
std::vector< std::vector< EntityHandle > > m2EntityGroups;
// Get matching entities for Mesh 1
std::cout << "Get matching entities for mesh 1..." << std::endl;
err = mbc.get_matching_entities( roots[0], &tagHandles[0], &tagValues[0], tagHandles.size(), &m1EntitySets,
&m1EntityGroups );MB_CHK_SET_ERR( (ErrorCode)err, "get_matching_entities failed" );
std::cout << " get_matching_entities returned " << m1EntityGroups.size() << " entity groups" << std::endl;
// Print out the data in the vector of vectors
std::vector< std::vector< EntityHandle > >::iterator iter_esi;
std::vector< std::vector< EntityHandle > >::iterator iter_egi;
std::vector< EntityHandle >::iterator iter_esj;
std::vector< EntityHandle >::iterator iter_egj;
Range entSetRg;
int icnt;
for( iter_egi = m1EntityGroups.begin(), iter_esi = m1EntitySets.begin(), icnt = 1;
( iter_egi != m1EntityGroups.end() ) && ( iter_esi != m1EntitySets.end() ); ++iter_egi, ++iter_esi, icnt++ )
{
std::cout << " EntityGroup(" << icnt << ") = ";
std::cout.flush();
entSetRg.clear();
for( iter_egj = ( *iter_egi ).begin(); iter_egj != ( *iter_egi ).end(); ++iter_egj )
entSetRg.insert( (EntityHandle)*iter_egj );
debugOut.print( 2, "Mesh1 matching Entities: ", entSetRg );
std::cout.flush();
std::cout << " EntitySet(" << icnt << ") = ";
std::cout.flush();
entSetRg.clear();
for( iter_esj = ( *iter_esi ).begin(); iter_esj != ( *iter_esi ).end(); ++iter_esj )
entSetRg.insert( (EntityHandle)*iter_esj );
debugOut.print( 2, "Mesh1 matching EntitySets: ", entSetRg );
std::cout.flush();
}
// Get matching entities for Mesh 2
std::cout << "Get matching entities for mesh 2..." << std::endl;
err = mbc.get_matching_entities( roots[1], &tagHandles[0], &tagValues[0], tagHandles.size(), &m2EntitySets,
&m2EntityGroups );MB_CHK_SET_ERR( (ErrorCode)err, "get_matching_entities failed" );
std::cout << " get_matching_entities returned " << m2EntityGroups.size() << " entity groups" << std::endl;
for( iter_egi = m2EntityGroups.begin(), iter_esi = m2EntitySets.begin(), icnt = 1;
( iter_egi != m2EntityGroups.end() ) && ( iter_esi != m2EntitySets.end() ); ++iter_egi, ++iter_esi, icnt++ )
{
std::cout << " EntityGroup(" << icnt << ") = ";
std::cout.flush();
entSetRg.clear();
for( iter_egj = ( *iter_egi ).begin(); iter_egj != ( *iter_egi ).end(); ++iter_egj )
entSetRg.insert( (EntityHandle)*iter_egj );
debugOut.print( 2, "Mesh2 matching Entities: ", entSetRg );
std::cout.flush();
std::cout << " EntitySet(" << icnt << ") = ";
std::cout.flush();
entSetRg.clear();
for( iter_esj = ( *iter_esi ).begin(); iter_esj != ( *iter_esi ).end(); ++iter_esj )
entSetRg.insert( (EntityHandle)*iter_esj );
debugOut.print( 2, "Mesh2 matching EntitySets: ", entSetRg );
std::cout.flush();
}
if( debug )
{
// ******************************
std::cout << "********** Test print_tuples **********" << std::endl;
// temporary test function
std::cout << "Testing print_tuples..." << std::endl;
TupleList test_tuple;
int num_ints = 3, num_longs = 2, num_ulongs = 4, num_reals = 6, num_rows = 10;
std::cout << " print of test_tuples zero init..." << std::endl;
test_tuple.initialize( 0, 0, 0, 0, 0 );
test_tuple.enableWriteAccess();
print_tuples( &test_tuple );
std::cout << " print of test_tuples after setting n to 10..." << std::endl;
test_tuple.set_n( 10 );
print_tuples( &test_tuple );
test_tuple.initialize( num_ints, num_longs, num_ulongs, num_reals, num_rows );
std::cout << " print of test_tuples after init..." << std::endl;
print_tuples( &test_tuple );
std::cout << " print of test_tuples after setting n to 10..." << std::endl;
test_tuple.set_n( 10 );
print_tuples( &test_tuple );
for( int i = 0; i < num_rows; i++ )
{
int j;
for( j = 0; j < num_ints; j++ )
test_tuple.vi_wr[i * num_ints + j] = (int)( ( j + 1 ) * ( i + 1 ) );
for( j = 0; j < num_longs; j++ )
test_tuple.vl_wr[i * num_longs + j] = (int)( ( j + 1 ) * ( i + 1 ) );
for( j = 0; j < num_ulongs; j++ )
test_tuple.vul_wr[i * num_ulongs + j] = (int)( ( j + 1 ) * ( i + 1 ) );
for( j = 0; j < num_reals; j++ )
test_tuple.vr_wr[i * num_reals + j] = (int)( ( j + 1 ) * ( i + 1 ) + ( j * 0.01 ) );
}
std::cout << " print of test_tuples after filling with data..." << std::endl;
print_tuples( &test_tuple );
// ******************************
std::cout << "********** Test pack_tuples and unpack_tuples **********" << std::endl;
void* mp_buf;
int buf_sz;
if( rank == 0 )
{
buf_sz = pack_tuples( &test_tuple, &mp_buf );
}
// Send buffer size
err = MPI_Bcast( &buf_sz, 1, MPI_INT, 0, MPI_COMM_WORLD );
if( err != MPI_SUCCESS )
{
std::cerr << "MPI_Bcast of buffer size failed" << std::endl;
return -1;
}
// Allocate a buffer in the other procs
if( rank != 0 )
{
mp_buf = malloc( buf_sz * sizeof( uint ) );
}
err = MPI_Bcast( mp_buf, buf_sz * sizeof( uint ), MPI_UNSIGNED_CHAR, 0, MPI_COMM_WORLD );
if( err != MPI_SUCCESS )
{
std::cerr << "MPI_Bcast of buffer failed" << std::endl;
return -1;
}
TupleList* rcv_tuples;
unpack_tuples( mp_buf, &rcv_tuples );
std::cout << " print of rcv_tuples after unpacking from MPI_Bcast..." << std::endl;
print_tuples( rcv_tuples );
}
err = integrate_scalar_field_test();MB_CHK_SET_ERR( (ErrorCode)err, "Failure in integrating a scalar_field" );
// ******************************
std::cout << "********** Test get_group_integ_vals **********" << std::endl;
std::cout << "Get group integrated field values..." << std::endl;
// print the field values at the vertices before change.
std::cout << " print vertex field values first:" << std::endl;
Tag norm_hdl;
result = mbi->tag_get_handle( normTag.c_str(), norm_hdl );MB_CHK_SET_ERR( (ErrorCode)err, "Failed to get tag handle." );
Coupler::IntegType integ_type = Coupler::VOLUME;
// Mesh 1 field values
std::cout << " Original entity vertex field values (mesh 1): " << std::endl;
print_vertex_fields( mbi, m1EntityGroups, norm_hdl, integ_type );
// Mesh 2 field values
std::cout << " Original entity vertex field values (mesh 2): " << std::endl;
print_vertex_fields( mbi, m2EntityGroups, norm_hdl, integ_type );
// Get the field values
std::vector< double >::iterator iter_ivals;
std::cout << "Get group integrated field values for mesh 1..." << std::endl;
std::vector< double > m1IntegVals( m1EntityGroups.size() );
err = mbc.get_group_integ_vals( m1EntityGroups, m1IntegVals, normTag.c_str(), 4, integ_type );MB_CHK_SET_ERR( (ErrorCode)err, "Failed to get the Mesh 1 group integration values." );
std::cout << "Mesh 1 integrated field values(" << m1IntegVals.size() << "): ";
for( iter_ivals = m1IntegVals.begin(); iter_ivals != m1IntegVals.end(); ++iter_ivals )
{
std::cout << ( *iter_ivals ) << " ";
}
std::cout << std::endl;
std::cout << "Get group integrated field values for mesh 2..." << std::endl;
std::vector< double > m2IntegVals( m2EntityGroups.size() );
err = mbc.get_group_integ_vals( m2EntityGroups, m2IntegVals, normTag.c_str(), 4, integ_type );MB_CHK_SET_ERR( (ErrorCode)err, "Failed to get the Mesh 2 group integration values." );
std::cout << "Mesh 2 integrated field values(" << m2IntegVals.size() << "): ";
for( iter_ivals = m2IntegVals.begin(); iter_ivals != m2IntegVals.end(); ++iter_ivals )
{
std::cout << ( *iter_ivals ) << " ";
}
std::cout << std::endl;
// ******************************
std::cout << "********** Test apply_group_norm_factors **********" << std::endl;
// Make the norm factors by inverting the integration values.
double val;
for( unsigned int i = 0; i < m1IntegVals.size(); i++ )
{
val = m1IntegVals[i];
m1IntegVals[i] = 1 / val;
}
for( unsigned int i = 0; i < m2IntegVals.size(); i++ )
{
val = m2IntegVals[i];
m2IntegVals[i] = 1 / val;
}
std::cout << "Mesh 1 norm factors(" << m1IntegVals.size() << "): ";
for( iter_ivals = m1IntegVals.begin(); iter_ivals != m1IntegVals.end(); ++iter_ivals )
{
std::cout << ( *iter_ivals ) << " ";
}
std::cout << std::endl;
std::cout << "Mesh 2 norm factors(" << m2IntegVals.size() << "): ";
for( iter_ivals = m2IntegVals.begin(); iter_ivals != m2IntegVals.end(); ++iter_ivals )
{
std::cout << ( *iter_ivals ) << " ";
}
std::cout << std::endl;
// Apply the factors and reprint the vertices
err = mbc.apply_group_norm_factor( m1EntitySets, m1IntegVals, normTag.c_str(), integ_type );MB_CHK_SET_ERR( (ErrorCode)err, "Failed to apply norm factors to Mesh 1." );
err = mbc.apply_group_norm_factor( m2EntitySets, m2IntegVals, normTag.c_str(), integ_type );MB_CHK_SET_ERR( (ErrorCode)err, "Failed to apply norm factors to Mesh 2." );
// Get the norm_tag_factor on the EntitySets
// Get the handle for the norm factor tag
Tag norm_factor_hdl;
std::string normFactor = normTag + "_normf";
result = mbi->tag_get_handle( normFactor.c_str(), norm_factor_hdl );MB_CHK_SET_ERR( result, "Failed to get norm factor tag handle." );
// Mesh 1 values
std::cout << "Mesh 1 norm factors per EntitySet...";
for( iter_esi = m1EntitySets.begin(); iter_esi != m1EntitySets.end(); ++iter_esi )
{
for( iter_esj = ( *iter_esi ).begin(); iter_esj != ( *iter_esi ).end(); ++iter_esj )
{
double data = 0;
EntityHandle eh = *iter_esj;
result = mbi->tag_get_data( norm_factor_hdl, &eh, 1, &data );MB_CHK_SET_ERR( result, "Failed to get tag data." );
std::cout << data << ", ";
}
}
std::cout << std::endl;
// Mesh 2 values
std::cout << "Mesh 2 norm factors per EntitySet...";
for( iter_esi = m2EntitySets.begin(); iter_esi != m2EntitySets.end(); ++iter_esi )
{
for( iter_esj = ( *iter_esi ).begin(); iter_esj != ( *iter_esi ).end(); ++iter_esj )
{
double data = 0;
EntityHandle eh = *iter_esj;
result = mbi->tag_get_data( norm_factor_hdl, &eh, 1, &data );MB_CHK_SET_ERR( result, "Failed to get tag data." );
std::cout << data << ", ";
}
}
std::cout << std::endl;
// ******************************
std::cout << "********** Test normalize_subset **********" << std::endl;
// Now call the Coupler::normalize_subset routine and see if we get an error.
std::cout << "Running Coupler::normalize_subset() on mesh 1" << std::endl;
err = mbc.normalize_subset( (EntityHandle)roots[0], normTag.c_str(), &tagNames[0], numTagNames, &tagValues[0],
Coupler::VOLUME, 4 );MB_CHK_SET_ERR( (ErrorCode)err, "Failure in call to Coupler::normalize_subset() on mesh 1" );
// Print the normFactor on each EntitySet after the above call.
// Mesh 1 values
std::cout << "Mesh 1 norm factors per EntitySet...";
for( iter_esi = m1EntitySets.begin(); iter_esi != m1EntitySets.end(); ++iter_esi )
{
for( iter_esj = ( *iter_esi ).begin(); iter_esj != ( *iter_esi ).end(); ++iter_esj )
{
double data = 0;
EntityHandle eh = *iter_esj;
result = mbi->tag_get_data( norm_factor_hdl, &eh, 1, &data );MB_CHK_SET_ERR( result, "Failed to get tag data." );
std::cout << data << ", ";
}
}
std::cout << std::endl;
std::cout << "Running Coupler::normalize_subset() on mesh 2" << std::endl;
err = mbc.normalize_subset( (EntityHandle)roots[1], normTag.c_str(), &tagNames[0], numTagNames, &tagValues[0],
Coupler::VOLUME, 4 );MB_CHK_SET_ERR( (ErrorCode)err, "Failure in call to Coupler::normalize_subset() on mesh 2" );
// Mesh 2 values
std::cout << "Mesh 2 norm factors per EntitySet...";
for( iter_esi = m2EntitySets.begin(); iter_esi != m2EntitySets.end(); ++iter_esi )
{
for( iter_esj = ( *iter_esi ).begin(); iter_esj != ( *iter_esi ).end(); ++iter_esj )
{
double data = 0;
EntityHandle eh = *iter_esj;
result = mbi->tag_get_data( norm_factor_hdl, &eh, 1, &data );MB_CHK_SET_ERR( result, "Failed to get tag data." );
std::cout << data << ", ";
}
}
std::cout << std::endl;
// Done, cleanup
std::cout << "********** ssn_test DONE! **********" << std::endl;
MPI_Finalize();
return 0;
}
int pack_tuples | ( | TupleList * | tl, |
void ** | ptr | ||
) |
double physField | ( | double | x, |
double | y, | ||
double | z | ||
) |
Definition at line 892 of file ssn_test.cpp.
{
double out;
// 1/r^2 decay from {0,0,0}
out = x * x + y * y + z * z;
out += 1e-1; // clamp
out = 1 / out;
return out;
}
void print_tuples | ( | TupleList * | tlp | ) |
Definition at line 748 of file ssn_test.cpp.
References moab::TupleList::get_n(), moab::TupleList::getTupleSize(), moab::TupleList::vi_rd, moab::TupleList::vl_rd, moab::TupleList::vr_rd, and moab::TupleList::vul_rd.
Referenced by main().
{
uint mi, ml, mul, mr;
tlp->getTupleSize( mi, ml, mul, mr );
std::cout << " tuple data: (n=" << tlp->get_n() << ")" << std::endl;
std::cout << " mi:" << mi << " ml:" << ml << " mul:" << mul << " mr:" << mr << std::endl;
std::cout << " [" << std::setw( 11 * mi ) << " int data"
<< " |" << std::setw( 11 * ml ) << " long data"
<< " |" << std::setw( 11 * mul ) << " Ulong data"
<< " |" << std::setw( 11 * mr ) << " real data"
<< " " << std::endl
<< " ";
for( unsigned int i = 0; i < tlp->get_n(); i++ )
{
if( mi > 0 )
{
for( unsigned int j = 0; j < mi; j++ )
{
std::cout << std::setw( 10 ) << tlp->vi_rd[i * mi + j] << " ";
}
}
else
{
std::cout << " ";
}
std::cout << "| ";
if( ml > 0 )
{
for( unsigned int j = 0; j < ml; j++ )
{
std::cout << std::setw( 10 ) << tlp->vl_rd[i * ml + j] << " ";
}
}
else
{
std::cout << " ";
}
std::cout << "| ";
if( mul > 0 )
{
for( unsigned int j = 0; j < mul; j++ )
{
std::cout << std::setw( 10 ) << tlp->vul_rd[i * mul + j] << " ";
}
}
else
{
std::cout << " ";
}
std::cout << "| ";
if( mr > 0 )
{
for( unsigned int j = 0; j < mr; j++ )
{
std::cout << std::setw( 10 ) << tlp->vr_rd[i * mr + j] << " ";
}
}
else
{
std::cout << " ";
}
if( i + 1 < tlp->get_n() ) std::cout << std::endl << " ";
}
std::cout << "]" << std::endl;
}
int print_vertex_fields | ( | Interface * | mbi, |
std::vector< std::vector< EntityHandle > > & | groups, | ||
Tag & | norm_hdl, | ||
Coupler::IntegType | integ_type | ||
) |
Definition at line 819 of file ssn_test.cpp.
References moab::Interface::dimension_from_handle(), ErrorCode, moab::Interface::get_connectivity(), MB_SUCCESS, moab::Interface::tag_get_data(), and moab::Coupler::VOLUME.
Referenced by main().
{
int err = 0;
ErrorCode result;
std::vector< EntityHandle >::iterator iter_j;
for( unsigned int i = 0; i < groups.size(); i++ )
{
std::cout << " Group - " << std::endl << " ";
for( iter_j = groups[i].begin(); iter_j != groups[i].end(); ++iter_j )
{
EntityHandle ehandle = ( *iter_j );
// Check that the entity in iter_j is of the same dimension as the
// integ_type we are performing
int j_type = mbi->dimension_from_handle( ehandle );
if( ( integ_type == Coupler::VOLUME ) && ( j_type != 3 ) ) continue;
// Retrieve the vertices from the element
const EntityHandle* conn = NULL;
int num_verts = 0;
result = mbi->get_connectivity( ehandle, conn, num_verts );
if( MB_SUCCESS != result ) return 1;
std::cout << std::fixed;
for( int iv = 0; iv < num_verts; iv++ )
{
double data = 0;
result = mbi->tag_get_data( norm_hdl, &conn[iv], 1, &data );
if( MB_SUCCESS != result ) return 1;
std::cout << std::setprecision( 8 ) << data << ", ";
}
std::cout << std::endl << " ";
}
std::cout << std::endl;
std::cout.unsetf( std::ios_base::floatfield ); // turn off fixed notation
}
return err;
}
void unpack_tuples | ( | void * | ptr, |
TupleList ** | tlp | ||
) |
bool debug = true |
Definition at line 31 of file ssn_test.cpp.