Experiment

src/CMakeLists.txt

@@ -256,6 +256,7 @@ osh_add_util(osh_reorder)
 osh_add_util(osh_fix)
 osh_add_util(osh_eval_implied)
 osh_add_util(osh_calc)
+osh_add_util(stv)
 if(Omega_h_USE_libMeshb)
   osh_add_util(meshb2osh)
   osh_add_util(osh2meshb)
@@ -276,6 +277,8 @@ osh_add_util(osh_adapt)
 osh_add_util(osh_filesystem)
 osh_add_util(ascii_vtk2osh)
 
+file(COPY stv.osh DESTINATION ${CMAKE_CURRENT_BINARY_DIR})
+
 if(BUILD_TESTING)
   if(Omega_h_USE_MPI)
     get_filename_component(COMPILER_DIR ${CMAKE_CXX_COMPILER} PATH)
@@ -531,6 +534,10 @@ if(Omega_h_USE_STK)
   set(Omega_h_HEADERS ${Omega_h_HEADERS} Omega_h_stk.hpp)
 endif()
 
+if(Omega_h_USE_EGADS)
+  set(Omega_h_HEADERS ${Omega_h_HEADERS} Omega_h_egads.hpp)
+endif()
+
 install(FILES ${Omega_h_HEADERS} DESTINATION include)
 
 if (Omega_h_USE_pybind11)

src/Omega_h_egads.cpp

@@ -182,7 +182,7 @@ void egads_reclassify(Mesh* mesh, Egads* eg) {
     auto ents2faces = mesh->ask_up(dim, FACE);
     auto adj_class_dims = read(unmap(ents2faces.ab2b, face_class_dims, 1));
     auto keep_edges = each_eq_to(adj_class_dims, I8(2));
-    auto ents2eq_faces = filter_graph(ents2faces, keep_edges);
+    auto ents2eq_faces = filter_graph_edges(ents2faces, keep_edges);
     auto adj_eq_face_ids = unmap(ents2eq_faces.ab2b, face_class_ids, 1);
     auto host_a2ab = HostRead<LO>(ents2eq_faces.a2ab);
     auto host_face_ids = HostRead<LO>(adj_eq_face_ids);

src/Omega_h_exodus.cpp

@@ -57,7 +57,7 @@ static void get_elem_type_info(
   } else if (type == "tetra4") {
     *p_dim = 3;
     *p_family = OMEGA_H_SIMPLEX;
-  } else if (type == "TETRA") {
+  } else if (type == "TETRA" || type == "TET" || type == "tet") {
     *p_dim = 3;
     *p_family = OMEGA_H_SIMPLEX;
   } else if (type == "TET4") {

src/Omega_h_metric.cpp

@@ -207,6 +207,12 @@ Reals limit_metric_gradation(
   do {
     values = values2;
     values2 = limit_gradation_once(mesh, values, max_rate);
+    Real max_diff = 0.0;
+    for (LO i=0; i < values.size(); ++i)
+      max_diff = std::max(max_diff, std::abs(values[i] - values2[i]));
+
+    std::cout << "max_diff = " << max_diff << std::endl;
+    std::cout << "tol = " << tol << std::endl;
     ++i;
     if (verbose && can_print(mesh) && i > 0 && i % 50 == 0) {
       std::cout << "warning: gradation limiting is up to step " << i << '\n';

src/stv.cpp

@@ -0,0 +1,60 @@
+#include "Omega_h_mesh.hpp"
+#include "Omega_h_metric.hpp"
+#include "Omega_h_adapt.hpp"
+#include "Omega_h_file.hpp"
+#include "Omega_h_array_ops.hpp"
+#include "Omega_h_quality.hpp"
+
+using Real = Omega_h::Real;
+
+void refine_metric(Omega_h::Write<Real>& metric, Real factor)
+{
+  assert(metric.size() % 6 == 0);
+  int nverts = metric.size()/6;
+  for (int i=0; i < nverts; ++i)
+    for (int j=0; j < 6; ++j)
+      metric[6*i + j] *= factor*factor;
+}
+
+void run_refinement(Omega_h::Mesh* mesh, const std::string& outname_prefix, double refine_factor)
+{
+  std::cout << "number of verts in initial mesh = " << mesh->nverts() << std::endl;
+  std::cout << "number of elements in initial mesh = " << mesh->nelems() << std::endl;
+
+  std::string metric_field_name = "new_metric";
+  Omega_h::Reals orig_metric = Omega_h::get_implied_metrics(mesh);
+  Omega_h::Write<Real> target_metric = Omega_h::deep_copy(orig_metric, "target_metric");
+  refine_metric(target_metric, refine_factor);
+  mesh->add_tag(Omega_h::VERT, metric_field_name, 6, Omega_h::Read<Real>(target_metric));
+
+  Omega_h::vtk::write_vtu(outname_prefix + "_initial_mesh.vtu", mesh);
+
+
+  int num_initial_verts = mesh->nverts();
+  Omega_h::AdaptOpts opts(mesh);
+  opts.min_quality_allowed = 0.15;  // this is needed for the 3X refine case
+  opts.xfer_opts.type_map[metric_field_name] = Omega_h_Transfer::OMEGA_H_METRIC;
+
+  Omega_h::grade_fix_adapt(mesh, opts, target_metric, true);
+
+  std::cout << "number of verts in new mesh = " << mesh->nverts() << std::endl;
+  std::cout << "number of elements in new mesh = " << mesh->nelems() << std::endl;
+
+  Omega_h::Read<Real> target_metric2 = mesh->get_array<Real>(0, metric_field_name);
+  std::cout << "min quality = " << Omega_h::get_min(mesh->comm(), Omega_h::measure_qualities(mesh, target_metric2)) << std::endl;
+
+  Omega_h::vtk::write_vtu(outname_prefix + "_adapted_mesh.vtu", mesh);
+}
+
+int main(int argc, char** argv)
+{
+  std::string fname = "stv.osh";
+
+  auto lib = Omega_h::Library(&argc, &argv);
+  Omega_h::Mesh mesh(&lib);
+  Omega_h::binary::read(fname, lib.world(), &mesh);  
+
+  run_refinement(&mesh, "stv_coarsen2", 0.5);
+
+  return 0;
+}

src/stv.osh/nparts

@@ -0,0 +1 @@
+1

src/stv.osh/version

@@ -0,0 +1 @@
+9