{
int err = MPI_Init(&argc, &argv);
if (argc < 3) {
std::cerr << "Usage: ";
std::cerr << argv[0] << " <nfiles> <fname1> ... <fnamen> <norm_tag> <tag_select_opts> <file_opts>" << 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 = MPI_Finalize();
return 1;
}
int nprocs, rank;
err = MPI_Comm_size(MPI_COMM_WORLD, &nprocs);
err = MPI_Comm_rank(MPI_COMM_WORLD, &rank);
std::stringstream fname;
fname << argv[0] << rank << ".out";
if (!std::freopen(fname.str().c_str(), "a", stdout)) return false;
if (!std::freopen(fname.str().c_str(), "a", stderr)) return false;
Interface *mbi = new Core();
if (NULL == mbi) {
std::cerr << "MOAB constructor failed" << std::endl;
return 1;
}
iMesh_Instance iMeshInst = reinterpret_cast<iMesh_Instance>(mbi);
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);
CHKERR(err, "get_file_options failed");
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;
std::cout << "Reading mesh file(s)..." << std::endl;
std::vector<ParallelComm *> pcs(filenames.size());
std::vector<ReadParallel *> rps(filenames.size());
iBase_EntitySetHandle *roots = (iBase_EntitySetHandle *) malloc(sizeof(iBase_EntitySetHandle) * 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]);
iMesh_createEntSet(iMeshInst, 0, &(roots[i]), &err);
CHKERR(err, "Creating root set failed");
result = rps[i]->load_file(filenames[i], (EntityHandle *)&roots[i], FileOptions(fileOpts.c_str()));
CHKERR(result, "iMeshInstance::load_file failed");
}
DebugOutput debugOut("ssn_test-", std::cerr);
debugOut.set_rank(rank);
debugOut.set_verbosity(10);
for (unsigned int k = 0; k < filenames.size(); k++) {
iBase_EntityHandle *rootEnts = NULL;
int rootEntsAlloc = 0;
int rootEntsSize = 0;
err = 0;
iMesh_getEntities(iMeshInst, roots[k], iBase_ALL_TYPES, iMesh_ALL_TOPOLOGIES,
&rootEnts, &rootEntsAlloc, &rootEntsSize, &err);
Range rootRg;
for (int j = 0; j < rootEntsSize; j++)
rootRg.insert((EntityHandle) rootEnts[j]);
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();
}
Range src_elems, targ_elems;
std::cout << "********** Create Coupler **********" << std::endl;
std::cout << "Creating Coupler..." << std::endl;
Coupler mbc(mbi, pcs[0], src_elems, 0);
std::cout << "Getting tag handles..." << std::endl;
int numTagNames = tagNames.size();
err = iBase_SUCCESS;
std::vector<iBase_TagHandle> tagHandles(numTagNames);
int iTags = 0;
while (iTags < numTagNames) {
std::cout << "Getting handle for " << tagNames[iTags] << std::endl;
iMesh_getTagHandle(iMeshInst, tagNames[iTags], &tagHandles[iTags], &err, strlen(tagNames[iTags]));
CHKERR(err, "Retrieving tag handles failed");
iTags++;
}
std::cout << "********** Test create_tuples **********" << std::endl;
iBase_EntitySetHandle *m1EntSets = NULL;
iBase_EntitySetHandle *m2EntSets = NULL;
int m1EntSetsSize=0;
int m2EntSetsSize=0;
int m1EntSetsAlloc=0;
int m2EntSetsAlloc=0;
iMesh_getEntSetsByTagsRec(iMeshInst, roots[0], &tagHandles[0],
&tagValues[0], tagHandles.size(), 0,
&m1EntSets, &m1EntSetsAlloc, &m1EntSetsSize, &err);
CHKERR(err, "iMesh_getEntSetsByTagsRec failed on Mesh 1.");
iMesh_getEntSetsByTagsRec(iMeshInst, roots[1], &tagHandles[0],
&tagValues[0], tagHandles.size(), 0,
&m2EntSets, &m2EntSetsAlloc, &m2EntSetsSize, &err);
CHKERR(err, "iMesh_getEntSetsByTagsRec failed on Mesh 2.");
std::cout << "Creating tuples for mesh 1..." << std::endl;
TupleList *m1TagTuples = NULL;
err = mbc.create_tuples(m1EntSets, m1EntSetsSize,
&tagHandles[0], tagHandles.size(), &m1TagTuples);
CHKERR(err, "create_tuples failed");
std::cout << " create_tuples returned" << std::endl;
print_tuples(m1TagTuples);
std::cout << "Creating tuples for mesh 2..." << std::endl;
TupleList *m2TagTuples = NULL;
err = mbc.create_tuples(m2EntSets, m2EntSetsSize,
&tagHandles[0], tagHandles.size(), &m2TagTuples);
CHKERR(err, "create_tuples failed");
std::cout << " create_tuples returned" << std::endl;
print_tuples(m2TagTuples);
std::cout << "********** Test consolidate_tuples **********" << std::endl;
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);
CHKERR(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<iBase_EntitySetHandle> > m1EntitySets;
std::vector< std::vector<iBase_EntityHandle> > m1EntityGroups;
std::vector< std::vector<iBase_EntitySetHandle> > m2EntitySets;
std::vector< std::vector<iBase_EntityHandle> > m2EntityGroups;
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);
CHKERR(err, "get_matching_entities failed");
std::cout << " get_matching_entities returned " << m1EntityGroups.size() << " entity groups" << std::endl;
std::vector< std::vector<iBase_EntitySetHandle> >::iterator iter_esi;
std::vector< std::vector<iBase_EntityHandle> >::iterator iter_egi;
std::vector<iBase_EntitySetHandle>::iterator iter_esj;
std::vector<iBase_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();
}
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);
CHKERR(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;
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);
}
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;
}
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);
}
std::cout << "********** Test moab::Element::Map::integrate_scalar_field **********" << std::endl;
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 );
double field_val = 0.0;
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;
Element::LinearHex hexMap;
try {
for (int i = 1; i <=5; i++) {
hexMap.set_vertices(biunit_cube);
field_val = hexMap.integrate_scalar_field(bcf);
std::cout << " binunit_cube, const_field(num_pts=" << i << "): field_val=" << field_val << std::endl;
field_val = 0.0;
field_val = hexMap.integrate_scalar_field(bf1);
std::cout << " binunit_cube, field_1(num_pts=" << i << "): field_val=" << field_val << std::endl;
field_val = 0.0;
field_val = hexMap.integrate_scalar_field(bf2);
std::cout << " binunit_cube, field_2(num_pts=" << i << "): field_val=" << field_val << std::endl;
field_val = 0.0;
field_val = hexMap.integrate_scalar_field(bf3);
std::cout << " binunit_cube, field_3(num_pts=" << i << "): field_val=" << field_val << std::endl;
field_val = 0.0;
hexMap.set_vertices(zerobase_cube);
field_val = hexMap.integrate_scalar_field(zcf);
std::cout << " zerobase_cube, const_field(num_pts=" << i << "): field_val=" << field_val << std::endl;
field_val = 0.0;
field_val = hexMap.integrate_scalar_field(zf1);
std::cout << " zerobase_cube, field_1(num_pts=" << i << "): field_val=" << field_val << std::endl;
field_val = 0.0;
field_val = hexMap.integrate_scalar_field(zf2);
std::cout << " zerobase_cube, field_2(num_pts=" << i << "): field_val=" << field_val << std::endl;
field_val = 0.0;
field_val = hexMap.integrate_scalar_field(zf3);
std::cout << " zerobase_cube, field_3(num_pts=" << i << "): field_val=" << field_val << std::endl;
field_val = 0.0;
}
}
catch (moab::Element::Map::ArgError) {
CHKERR(iBase_FAILURE, "Failed to set vertices on Element::Map.");
}
catch (moab::Element::Map::EvaluationError) {
CHKERR(iBase_FAILURE, "Failed to get inverse evaluation of coordinate on Element::Map.");
}
std::cout << "********** Test get_group_integ_vals **********" << std::endl;
std::cout << "Get group integrated field values..." << std::endl;
std::cout << " print vertex field values first:" << std::endl;
iBase_TagHandle norm_hdl;
iMesh_getTagHandle(iMeshInst, normTag.c_str(), &norm_hdl, &err, strlen(normTag.c_str()));
CHKERR(err, "Failed to get tag handle.");
Coupler::IntegType integ_type = Coupler::VOLUME;
std::cout << " Original entity vertex field values (mesh 1): " << std::endl;
print_vertex_fields(mbi, iMeshInst, m1EntityGroups, norm_hdl, integ_type);
std::cout << " Original entity vertex field values (mesh 2): " << std::endl;
print_vertex_fields(mbi, iMeshInst, m2EntityGroups, norm_hdl, integ_type);
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);
CHKERR(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);
CHKERR(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;
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;
err = mbc.apply_group_norm_factor(m1EntitySets, m1IntegVals, normTag.c_str(), integ_type);
CHKERR(err, "Failed to apply norm factors to Mesh 1.");
err = mbc.apply_group_norm_factor(m2EntitySets, m2IntegVals, normTag.c_str(), integ_type);
CHKERR(err, "Failed to apply norm factors to Mesh 2.");
iBase_TagHandle norm_factor_hdl;
std::string normFactor = normTag + "_normf";
iMesh_getTagHandle(iMeshInst, normFactor.c_str(), &norm_factor_hdl, &err, strlen(normFactor.c_str()));
CHKERR(err, "Failed to get norm factor tag handle.");
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;
iMesh_getEntSetDblData(iMeshInst, (iBase_EntitySetHandle) *iter_esj, norm_factor_hdl, &data, &err);
CHKERR(err, "Failed to get tag data.");
std::cout << data << ", ";
}
}
std::cout << std::endl;
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;
iMesh_getEntSetDblData(iMeshInst, (iBase_EntitySetHandle) *iter_esj, norm_factor_hdl, &data, &err);
CHKERR(err, "Failed to get tag data.");
std::cout << data << ", ";
}
}
std::cout << std::endl;
std::cout << "********** Test normalize_subset **********" << std::endl;
std::cout << "Running Coupler::normalize_subset() on mesh 1" << std::endl;
err = mbc.normalize_subset(roots[0],
normTag.c_str(),
&tagNames[0],
numTagNames,
&tagValues[0],
Coupler::VOLUME,
4);
CHKERR(err, "Failure in call to Coupler::normalize_subset() on mesh 1");
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;
iMesh_getEntSetDblData(iMeshInst, (iBase_EntitySetHandle) *iter_esj, norm_factor_hdl, &data, &err);
CHKERR(err, "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(roots[1],
normTag.c_str(),
&tagNames[0],
numTagNames,
&tagValues[0],
Coupler::VOLUME,
4);
CHKERR(err, "Failure in call to Coupler::normalize_subset() on mesh 2");
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;
iMesh_getEntSetDblData(iMeshInst, (iBase_EntitySetHandle) *iter_esj, norm_factor_hdl, &data, &err);
CHKERR(err, "Failed to get tag data.");
std::cout << data << ", ";
}
}
std::cout << std::endl;
std::cout << "********** ssn_test DONE! **********" << std::endl;
MPI_Finalize();
return 0;
}