[MOAB-dev] commit/MOAB: 3 new changesets
commits-noreply at bitbucket.org
commits-noreply at bitbucket.org
Mon Feb 3 16:18:59 CST 2014
3 new commits in MOAB:
https://bitbucket.org/fathomteam/moab/commits/3569b60c9bf1/
Changeset: 3569b60c9bf1
Branch: None
User: vijaysm
Date: 2014-01-11 07:46:55
Summary: Merged fathomteam/moab into master
Affected #: 7 files
diff --git a/config/compiler.m4 b/config/compiler.m4
index ce83a5f..4787100 100644
--- a/config/compiler.m4
+++ b/config/compiler.m4
@@ -548,6 +548,8 @@ case "$cc_compiler:$host_cpu" in
bgq)
FATHOM_CC_32BIT=-q32
FATHOM_CC_64BIT=-q64
+ FATHOM_CC_SPECIAL=-qarch=qp
+ FATHOM_CXX_SPECIAL="-qarch=qp -qpic=large -qmaxmem=-1"
AR="ar"
NM="nm -B"
;;
diff --git a/src/SpatialLocator.cpp b/src/SpatialLocator.cpp
index f765ab9..36016c1 100644
--- a/src/SpatialLocator.cpp
+++ b/src/SpatialLocator.cpp
@@ -72,7 +72,7 @@ namespace moab
locTable.enableWriteAccess();
// pass storage directly into locate_points, since we know those arrays are contiguous
- ErrorCode rval = locate_points(pos, num_points, locTable.vul_wr, locTable.vr_wr, NULL, rel_iter_tol, abs_iter_tol,
+ ErrorCode rval = locate_points(pos, num_points, (EntityHandle*)locTable.vul_wr, locTable.vr_wr, NULL, rel_iter_tol, abs_iter_tol,
inside_tol);
std::fill(locTable.vi_wr, locTable.vi_wr+num_points, 0);
locTable.set_n(num_points);
diff --git a/tools/mbcslam/CslamUtils.cpp b/tools/mbcslam/CslamUtils.cpp
index e806d61..43668a6 100644
--- a/tools/mbcslam/CslamUtils.cpp
+++ b/tools/mbcslam/CslamUtils.cpp
@@ -906,16 +906,34 @@ void velocity_case1(CartVect & arrival_point, double t, CartVect & velo)
//
ErrorCode enforce_convexity(Interface * mb, EntityHandle lset, int my_rank)
{
- // look at each quad; compute all 4 angles; if one is reflex, break along that diagonal
+ // look at each polygon; compute all angles; if one is reflex, break that angle with
+ // the next triangle; put the 2 new polys in the set;
+ // still look at the next poly
// replace it with 2 triangles, and remove from set;
// it should work for all polygons / tested first for case 1, with dt 0.5 (too much deformation)
// get all entities of dimension 2
// then get the connectivity, etc
+
Range inputRange;
ErrorCode rval = mb->get_entities_by_dimension(lset, 2, inputRange);
if (MB_SUCCESS != rval)
return rval;
+ Tag corrTag=0;
+ EntityHandle dumH=0;
+ rval = mb->tag_get_handle(CORRTAGNAME,
+ 1, MB_TYPE_HANDLE, corrTag,
+ MB_TAG_DENSE, &dumH);
+ if(rval==MB_TAG_NOT_FOUND)
+ corrTag = 0;
+
+ Tag gidTag;
+ rval = mb->tag_get_handle("GLOBAL_ID", 1, MB_TYPE_INTEGER,
+ gidTag, MB_TAG_DENSE);
+
+ if(rval!=MB_SUCCESS)
+ return rval;
+
std::vector<double> coords;
coords.resize(3*MAXEDGES); // at most 10 vertices per polygon
// we should create a queue with new polygons that need processing for reflex angles
@@ -942,6 +960,17 @@ ErrorCode enforce_convexity(Interface * mb, EntityHandle lset, int my_rank)
rval = mb->get_connectivity(eh, verts, num_nodes);
if (MB_SUCCESS != rval)
return rval;
+ EntityHandle corrHandle=0;
+ if (corrTag)
+ {
+ rval = mb->tag_get_data(corrTag, &eh, 1, &corrHandle);
+ if (MB_SUCCESS != rval)
+ return rval;
+ }
+ int gid=0;
+ rval = mb->tag_get_data(gidTag, &eh, 1, &gid);
+ if (MB_SUCCESS != rval)
+ return rval;
coords.resize(3 * num_nodes);
if (num_nodes < 4)
continue; // if already triangles, don't bother
@@ -998,6 +1027,15 @@ ErrorCode enforce_convexity(Interface * mb, EntityHandle lset, int my_rank)
rval = mb->add_entities(lset, &newElement, 1);
if (MB_SUCCESS != rval)
return rval;
+ if (corrTag)
+ {
+ rval = mb->tag_set_data(corrTag, &newElement, 1, &corrHandle);
+ if (MB_SUCCESS != rval)
+ return rval;
+ }
+ rval = mb->tag_set_data(gidTag, &newElement, 1, &gid);
+ if (MB_SUCCESS != rval)
+ return rval;
if (num_nodes == 4)
{
// create another triangle
@@ -1017,6 +1055,15 @@ ErrorCode enforce_convexity(Interface * mb, EntityHandle lset, int my_rank)
rval = mb->add_entities(lset, &newElement, 1);
if (MB_SUCCESS != rval)
return rval;
+ if (corrTag)
+ {
+ rval = mb->tag_set_data(corrTag, &newElement, 1, &corrHandle);
+ if (MB_SUCCESS != rval)
+ return rval;
+ }
+ rval = mb->tag_set_data(gidTag, &newElement, 1, &gid);
+ if (MB_SUCCESS != rval)
+ return rval;
mb->remove_entities(lset, &eh, 1);
brokenPolys++;
/*std::cout<<"remove: " ;
@@ -1029,7 +1076,7 @@ ErrorCode enforce_convexity(Interface * mb, EntityHandle lset, int my_rank)
}
}
}
- std::cout << "on rank " << my_rank << " " << brokenPolys << " concave polygons were decomposed in convex ones \n";
+ std::cout << "on local process " << my_rank << " " << brokenPolys << " concave polygons were decomposed in convex ones \n";
return MB_SUCCESS;
}
ErrorCode create_span_quads(Interface * mb, EntityHandle euler_set, int rank)
diff --git a/tools/mbcslam/Intx2Mesh.cpp b/tools/mbcslam/Intx2Mesh.cpp
index 625a442..99d7ca3 100644
--- a/tools/mbcslam/Intx2Mesh.cpp
+++ b/tools/mbcslam/Intx2Mesh.cpp
@@ -800,6 +800,12 @@ ErrorCode Intx2Mesh::create_departure_mesh_2nd_alg(EntityHandle & euler_set, Ent
ERRORR(rval, "can't set global id tag on new el");
}
}
+ // order the remote cells tuple list, with the global id, because we will search in it
+ //remote_cells->print("remote_cells before sorting");
+ moab::TupleList::buffer sort_buffer;
+ sort_buffer.buffer_init(n);
+ remote_cells->sort(1, &sort_buffer);
+ sort_buffer.reset();
return MB_SUCCESS;
}
diff --git a/tools/mbcslam/create_dp.cpp b/tools/mbcslam/create_dp.cpp
index bd2459b..43e1735 100644
--- a/tools/mbcslam/create_dp.cpp
+++ b/tools/mbcslam/create_dp.cpp
@@ -185,10 +185,11 @@ int main(int argc, char **argv)
if (argc < 3)
{
- std::cout << " usage: create_dp <input><output> -t <time> -dt <delta_t> -h \n";
+ std::cout << " usage: create_dp <input><output> -t <time> -dt <delta_t> [-skipdp] -f <field> -h \n";
return 1;
}
+ bool skip= false;
double dt=0.1;
double t=0.1; // corresponding to diffusion first step
@@ -208,9 +209,18 @@ int main(int argc, char **argv)
if (!strcmp(argv[index], "-h"))
{
- std::cout << " usage: create_dp <input><output> -t <time> -dt <delta_t> -h \n";
+ std::cout << " usage: create_dp <input><output> -t <time> -dt <delta_t> [-skipdp] -f <field> -h \n";
return 1;
}
+ if (!strcmp(argv[index], "-f")) // delete partition sets
+ {
+ field_type = atoi(argv[++index]);
+ }
+
+ if (!strcmp(argv[index], "-skipdp") )
+ {
+ skip = true;
+ }
index++;
}
@@ -224,43 +234,50 @@ int main(int argc, char **argv)
std::cout << " -t " << t << " -dt " << dt << " input: " << input_mesh1 <<
" output: " << output << "\n";
- Range verts;
- rval = mb.get_entities_by_dimension(0, 0, verts);
- if (MB_SUCCESS != rval)
- return 1;
-
- double *x_ptr, *y_ptr, *z_ptr;
- int count;
- rval = mb.coords_iterate(verts.begin(), verts.end(), x_ptr, y_ptr, z_ptr, count);
- if (MB_SUCCESS != rval)
+ // skip if we need for DP (already existing)
+ if (skip)
+ {
+ std::cout<<" do not add DP tag \n";
+ }
+ else
+ {
+ Range verts;
+ rval = mb.get_entities_by_dimension(0, 0, verts);
+ if (MB_SUCCESS != rval)
return 1;
- assert(count == (int) verts.size()); // should end up with just one contiguous chunk of vertices
- Tag tagh = 0;
- std::string tag_name("DP");
- rval = mb.tag_get_handle(tag_name.c_str(), 3, MB_TYPE_DOUBLE, tagh, MB_TAG_DENSE | MB_TAG_CREAT);
- CHECK_ERR(rval);
- void *data; // pointer to the LOC in memory, for each vertex
- int count_tag;
+ double *x_ptr, *y_ptr, *z_ptr;
+ int count;
+ rval = mb.coords_iterate(verts.begin(), verts.end(), x_ptr, y_ptr, z_ptr, count);
+ if (MB_SUCCESS != rval)
+ return 1;
+ assert(count == (int) verts.size()); // should end up with just one contiguous chunk of vertices
- rval = mb.tag_iterate(tagh, verts.begin(), verts.end(), count_tag, data);
- CHECK_ERR(rval);
- // here we are checking contiguity
- assert(count_tag == (int) verts.size());
- double * ptr_DP=(double*)data;
+ Tag tagh = 0;
+ std::string tag_name("DP");
+ rval = mb.tag_get_handle(tag_name.c_str(), 3, MB_TYPE_DOUBLE, tagh, MB_TAG_DENSE | MB_TAG_CREAT);
+ CHECK_ERR(rval);
+ void *data; // pointer to the LOC in memory, for each vertex
+ int count_tag;
- for (int v = 0; v < count; v++) {
- //EntityHandle v = verts[v];
- CartVect pos( x_ptr[v], y_ptr[v] , z_ptr[v]);
- CartVect newPos;
- departure_point_case1(pos, t, dt, newPos);
- ptr_DP[0]=newPos[0];
- ptr_DP[1]=newPos[1];
- ptr_DP[2]=newPos[2];
- ptr_DP+=3; // increment to the next vertex
+ rval = mb.tag_iterate(tagh, verts.begin(), verts.end(), count_tag, data);
+ CHECK_ERR(rval);
+ // here we are checking contiguity
+ assert(count_tag == (int) verts.size());
+ double * ptr_DP=(double*)data;
+
+ for (int v = 0; v < count; v++) {
+ //EntityHandle v = verts[v];
+ CartVect pos( x_ptr[v], y_ptr[v] , z_ptr[v]);
+ CartVect newPos;
+ departure_point_case1(pos, t, dt, newPos);
+ ptr_DP[0]=newPos[0];
+ ptr_DP[1]=newPos[1];
+ ptr_DP[2]=newPos[2];
+ ptr_DP+=3; // increment to the next vertex
+ }
}
-
rval = add_field_value(mb);
mb.write_file(output);
diff --git a/tools/mbcslam/intx_mpas.cpp b/tools/mbcslam/intx_mpas.cpp
index 570c97f..d737156 100644
--- a/tools/mbcslam/intx_mpas.cpp
+++ b/tools/mbcslam/intx_mpas.cpp
@@ -216,7 +216,7 @@ int main(int argc, char **argv)
if (Verbose)
{
std::stringstream lagrIni;
- lagrIni<<"def0" << rank<<".h5m";
+ lagrIni<<"lagr0" << rank<<".h5m";
rval = mb.write_file(lagrIni.str().c_str(), 0, 0, &covering_lagr_set, 1);
}
@@ -226,7 +226,7 @@ int main(int argc, char **argv)
if (Verbose)
{
std::stringstream ste;
- ste<<"lagr0" << rank<<".h5m";
+ ste<<"euler0" << rank<<".h5m";
rval = mb.write_file(ste.str().c_str(), 0, 0, &euler_set, 1);
}
diff --git a/tools/mbcslam/wrap_intx.cpp b/tools/mbcslam/wrap_intx.cpp
index b973acf..878a832 100644
--- a/tools/mbcslam/wrap_intx.cpp
+++ b/tools/mbcslam/wrap_intx.cpp
@@ -17,7 +17,7 @@
using namespace moab;
double radius = 1.;
double gtol = 1.e-9;
-bool debug = false;
+bool debug = true;
#ifdef __cplusplus
extern "C" {
@@ -52,6 +52,10 @@ void update_tracer( iMesh_Instance instance, iBase_EntitySetHandle imesh_euler_s
ERRORV(rval , "can't write covering set ");
}
+ //
+ rval = enforce_convexity(mb, covering_lagr_set);
+ ERRORV(rval , "can't write covering set ");
+
EntityHandle outputSet;
rval = mb->create_meshset(MESHSET_SET, outputSet);
ERRORV(rval , "can't create output set ");
https://bitbucket.org/fathomteam/moab/commits/963a275cea8b/
Changeset: 963a275cea8b
Branch: None
User: vijaysm
Date: 2014-02-03 23:11:59
Summary: Use the fortran linker for fortran examples.
Affected #: 1 file
diff --git a/examples/makefile b/examples/makefile
index a6e3d1f..6c1ed9a 100644
--- a/examples/makefile
+++ b/examples/makefile
@@ -36,7 +36,7 @@ DirectAccessNoHoles: DirectAccessNoHoles.o ${MOAB_LIBDIR}/libMOAB.la
${MOAB_CXX} -o $@ $< ${MOAB_LIBS_LINK}
DirectAccessNoHolesF90: DirectAccessNoHolesF90.o ${MOAB_LIBDIR}/libMOAB.la
- ${MOAB_CXX} -o $@ $< ${IMESH_LIBS}
+ ${MOAB_FC} -o $@ $< ${IMESH_LIBS}
ReduceExchangeTags : ReduceExchangeTags.o ${MOAB_LIBDIR}/libMOAB.la
${MOAB_CXX} -o $@ $< ${MOAB_LIBS_LINK}
@@ -51,7 +51,7 @@ point_in_elem_search: point_in_elem_search.o ${MOAB_LIBDIR}/libMOAB.la
${MOAB_CXX} -o $@ $< ${MOAB_LIBS_LINK}
PushParMeshIntoMoabF90: PushParMeshIntoMoabF90.o
- ${MOAB_CXX} -o $@ $< ${IMESH_LIBS} -lgfortran -L/usr/lib/openmpi/lib -lmpi_f90 -lmpi_f77 -lmpi -lopen-rte -lopen-pal -ldl
+ ${MOAB_FC} -o $@ $< ${IMESH_LIBS}
DeformMeshRemap: DeformMeshRemap.o ${MOAB_LIBDIR}/libMOAB.la
${MOAB_CXX} -o $@ $< ${MOAB_LIBS_LINK} -lmbcoupler ${MOAB_LIBS_LINK}
https://bitbucket.org/fathomteam/moab/commits/be8a392b1845/
Changeset: be8a392b1845
Branch: master
User: vijaysm
Date: 2014-02-03 23:13:58
Summary: Merge remote-tracking branch 'upstream/master'
Affected #: 2 files
diff --git a/tools/mbcslam/diffusion.cpp b/tools/mbcslam/diffusion.cpp
index 378473e..be70cd6 100644
--- a/tools/mbcslam/diffusion.cpp
+++ b/tools/mbcslam/diffusion.cpp
@@ -397,17 +397,28 @@ ErrorCode compute_tracer_case1(Interface * mb, Intx2MeshOnSphere & worker, Entit
// if in parallel, we have to move some elements to another proc, and receive other cells
// from other procs
// lagr and euler are preserved
+ if (writeFiles) // so if need to write lagr files too
+ {
+
+ }
EntityHandle covering_set;
rval = worker.create_departure_mesh_3rd_alg(lagr_set, covering_set);
if (writeFiles) // so if write
{
+ std::stringstream departureMesh;
+ departureMesh << "Departure" << rank << "_" << tStep << ".vtk";
+ rval = mb->write_file(departureMesh.str().c_str(), 0, 0, &lagr_set, 1);
+ CHECK_ERR(rval);
+
std::stringstream newTracer;
newTracer << "Tracer" << rank << "_" << tStep << ".vtk";
rval = mb->write_file(newTracer.str().c_str(), 0, 0, &euler_set, 1);
+ CHECK_ERR(rval);
std::stringstream lagr_cover;
lagr_cover << "Cover" << rank << "_" << tStep << ".vtk";
rval = mb->write_file(lagr_cover.str().c_str(), 0, 0, &covering_set, 1);
+ CHECK_ERR(rval);
}
// so we have now the departure at the previous time
diff --git a/tools/mbcslam/intx_mpas.cpp b/tools/mbcslam/intx_mpas.cpp
index d737156..b3d49d9 100644
--- a/tools/mbcslam/intx_mpas.cpp
+++ b/tools/mbcslam/intx_mpas.cpp
@@ -22,6 +22,7 @@
#include "moab/ParallelComm.hpp"
#include "CslamUtils.hpp"
+#include <time.h>
#ifdef MESHDIR
std::string TestDir( STRINGIFY(MESHDIR) );
@@ -43,6 +44,7 @@ double radius = 1.;// in m: 6371220.
double t = 0.1, delta_t = 0.05; // check the script
bool Verbose = false;
+double rot= M_PI/10;
ErrorCode manufacture_lagrange_mesh_on_sphere(Interface * mb,
EntityHandle euler_set)
@@ -84,7 +86,7 @@ ErrorCode manufacture_lagrange_mesh_on_sphere(Interface * mb,
//int vix=0; // vertex index in new array
double T=5;// check the script
- double rot= M_PI/10;
+
for (Range::iterator vit=connecVerts.begin();vit!=connecVerts.end(); vit++ )
{
EntityHandle oldV=*vit;
@@ -152,12 +154,23 @@ int main(int argc, char **argv)
if (!strcmp(argv[index], "-FF"))
{
flux_form= true;
- index++;
}
if (!strcmp(argv[index], "-v"))
{
Verbose = true;
- index++;
+ }
+ if (!strcmp(argv[index], "-t"))
+ {
+ t = atof(argv[++index]);
+ }
+ if (!strcmp(argv[index], "-t"))
+ {
+ t = atof(argv[++index]);
+ }
+ if (!strcmp(argv[index], "-rot"))
+ {
+ rot = atof(argv[++index]);
+ rot = M_PI/rot; // so rot 50 means rotate with M_PI/50 radians
}
index++;
@@ -165,7 +178,7 @@ int main(int argc, char **argv)
}
// start copy
std::string opts = std::string("PARALLEL=READ_PART;PARTITION_METHOD=RCBZOLTAN")+
- std::string(";PARALLEL_RESOLVE_SHARED_ENTS;VARIABLE=;")+extra_read_opts;
+ std::string(";PARALLEL_RESOLVE_SHARED_ENTS;VARIABLE=;NO_EDGES;")+extra_read_opts;
Core moab;
Interface & mb = moab;
EntityHandle euler_set;
@@ -173,20 +186,30 @@ int main(int argc, char **argv)
rval = mb.create_meshset(MESHSET_SET, euler_set);
CHECK_ERR(rval);
+ clock_t tt = clock();
rval = mb.load_file(filename_mesh1, &euler_set, opts.c_str());
+ CHECK_ERR(rval);
ParallelComm* pcomm = ParallelComm::get_pcomm(&mb, 0);
CHECK_ERR(rval);
- rval = pcomm->check_all_shared_handles();
- CHECK_ERR(rval);
+ /*rval = pcomm->check_all_shared_handles();
+ CHECK_ERR(rval);*/
// end copy
int rank = pcomm->proc_config().proc_rank();
+ int procs = pcomm->proc_config().proc_size();
if (0==rank)
std::cout << " case 1: use -gtol " << gtol << " -dt " << delta_t <<
- " -R " << radius << " -input " << filename_mesh1 << "\n";
+ " -R " << radius << " -input " << filename_mesh1 << " -t " << t << " -rot " << rot << "\n";
+
+ if (0==rank)
+ {
+ std::cout << "load mesh from " << filename_mesh1 << "\n on " << procs << " processors in "
+ << (clock() - tt) / (double) CLOCKS_PER_SEC << " seconds" << std::endl;
+ tt = clock();
+ }
rval = manufacture_lagrange_mesh_on_sphere(&mb, euler_set);
if (MB_SUCCESS != rval)
@@ -208,11 +231,23 @@ int main(int argc, char **argv)
//double radius = 1.; // input
worker.SetRadius(radius);
-
+ if (0==rank)
+ {
+ std::cout << "manufacture departure mesh " << filename_mesh1 << "\n on " << procs << " processors in "
+ << (clock() - tt) / (double) CLOCKS_PER_SEC << " seconds" << std::endl;
+ tt = clock();
+ }
worker.SetErrorTolerance(gtol);
rval = worker.create_departure_mesh_2nd_alg(euler_set, covering_lagr_set);
CHECK_ERR(rval);
+ if (0==rank)
+ {
+ std::cout << "communicate covering mesh on " << procs << " processors in "
+ << (clock() - tt) / (double) CLOCKS_PER_SEC << " seconds" << std::endl;
+ tt = clock();
+ }
+
if (Verbose)
{
std::stringstream lagrIni;
@@ -240,8 +275,13 @@ int main(int argc, char **argv)
rval = worker.intersect_meshes(covering_lagr_set, euler_set, outputSet);
if (MB_SUCCESS != rval)
return 1;
-
- if (rank<=4)
+ if (0==rank)
+ {
+ std::cout << "intersect meshes in " << procs << " processors in "
+ << (clock() - tt) / (double) CLOCKS_PER_SEC << " seconds" << std::endl;
+ tt = clock();
+ }
+ if (Verbose && rank<=4)
{
std::string opts_write("");
std::stringstream outf;
@@ -250,11 +290,15 @@ int main(int argc, char **argv)
if (MB_SUCCESS != rval)
std::cout << "can't write output\n";
}
- double intx_area = area_on_sphere_lHuiller(&mb, outputSet, radius);
- double arrival_area = area_on_sphere_lHuiller(&mb, euler_set, radius);
- std::cout << " Arrival area: " << arrival_area
- << " intersection area:" << intx_area << " rel error: "
- << fabs((intx_area - arrival_area) / arrival_area) << "\n";
+
+ if (rank <= 4)
+ {
+ double intx_area = area_on_sphere_lHuiller(&mb, outputSet, radius);
+ double arrival_area = area_on_sphere_lHuiller(&mb, euler_set, radius);
+ std::cout << "On proc " << rank << " arrival area: " << arrival_area
+ << " intersection area:" << intx_area << " rel error: "
+ << fabs((intx_area - arrival_area) / arrival_area) << "\n";
+ }
MPI_Finalize();
if (MB_SUCCESS != rval)
Repository URL: https://bitbucket.org/fathomteam/moab/
--
This is a commit notification from bitbucket.org. You are receiving
this because you have the service enabled, addressing the recipient of
this email.
More information about the moab-dev
mailing list