Add option to specify EB with STL file

src/embedded_boundaries/CMakeLists.txt

@@ -9,6 +9,7 @@ target_sources(incflo
    eb_cyl_tuscan.cpp
    eb_regular.cpp
    eb_sphere.cpp
+   eb_stl.cpp
    eb_chkptfile.cpp
    eb_if.H)
 

src/embedded_boundaries/Make.package

@@ -6,6 +6,7 @@ CEXE_sources += eb_cyl_tuscan.cpp
 CEXE_sources += eb_regular.cpp
 CEXE_sources += eb_sphere.cpp
 CEXE_sources += eb_chkptfile.cpp
+CEXE_sources += eb_stl.cpp
 ifeq ($(DIM), 3)
   CEXE_sources += eb_spherecube.cpp
   CEXE_sources += eb_tuscan.cpp

src/embedded_boundaries/eb_csg.cpp

@@ -16,6 +16,10 @@ void incflo::make_eb_csg(const std::string& geom_file)
     Vector<Real> translation_vec(AMREX_SPACEDIM, 0.0);
 
     ParmParse pp("csg");
+    pp.query("geometry_filename", geom_file);
+    if (geom_file.empty()) {
+        Abort("Missing or invalid input: csg.geometry_filename = ...\n");
+    }
     pp.query("internal_flow", is_internal_flow);
     if(pp.queryarr("scaling_factor", scaling_factor_vec, 0, AMREX_SPACEDIM)) {
       amrex::Print() << "WARNING: The implicit function magnitudes will not be scaled" << "\n";

src/embedded_boundaries/eb_stl.cpp

@@ -0,0 +1,48 @@
+#include <AMReX_EB2.H>
+#include <AMReX_ParmParse.H>
+#include <AMReX_EB2_IndexSpace_STL.H>
+
+#include <incflo.H>
+
+using namespace amrex;
+
+void incflo::make_eb_stl ()
+{
+    bool is_internal_flow = true;
+    Real scaling_factor = 1.0;
+    Vector<Real> translation_vec(3, 0.0);
+    bool use_bvh = true;
+
+    ParmParse pp("stl");
+
+    std::string stl_file;
+    pp.query("geometry_filename", stl_file);
+
+    if (stl_file.empty()) {
+        Abort("\nMissing or invalid input: stl.geometry_filename = ...");
+    } else {
+        Print() << "STL geometry file: " << stl_file << '\n';
+    }
+
+    pp.query("internal_flow", is_internal_flow); // This will flip the normal
+    pp.query("scaling_factor", scaling_factor);
+    pp.queryarr("translation", translation_vec); // This is the stl center
+    pp.query("use_bvh", use_bvh);
+
+    if (is_internal_flow) { Print() << "\n Building geometry for internal flow\n"; }
+    else { Print() << "\n Building geometry for external flow\n"; }
+
+    /************************************************************************
+     *                                                                      *
+     * Build EB levels from STL                                             *
+     *                                                                      *
+     ***********************************************************************/
+    EB2::IndexSpace::push(
+        std::make_unique<EB2::IndexSpaceSTL>
+        (stl_file, scaling_factor,
+         Array<Real,3>{translation_vec[0], translation_vec[1], translation_vec[2]},
+         is_internal_flow, geom[max_level], max_level, 100, /*ngrow*/4,
+         /* build coarse by coarsening*/ true, EB2::ExtendDomainFace(),
+         EB2::NumCoarsenOpt(), use_bvh, /*support mvmc*/false));
+
+}

src/embedded_boundaries/embedded_boundaries.cpp

@@ -19,7 +19,7 @@ void incflo::MakeEBGeometry()
     std::string csg_file;
     pp.query("geometry", geom_type);
     pp.query("geometry_filename", csg_file);
-    amrex::Print() << "incflo.geometry_filename: " << csg_file;
+    amrex::Print() << "incflo.geometry_filename: " << csg_file << "\n";
 
 #ifndef CSG_EB
     AMREX_ALWAYS_ASSERT_WITH_MESSAGE( csg_file.empty(), "CSG Geometry defined in input deck but solver not built with CSG support!");
@@ -33,7 +33,7 @@ void incflo::MakeEBGeometry()
 
     if(geom_type == "cylinder")
     {
-    amrex::Print() << "\n Building cylinder geometry." << "\n";
+        amrex::Print() << "\n Building cylinder geometry." << "\n";
         make_eb_cylinder();
     }
     else if(geom_type == "box")
@@ -44,56 +44,59 @@ void incflo::MakeEBGeometry()
 #if (AMREX_SPACEDIM == 3)
     else if(geom_type == "twocylinders")
     {
-    amrex::Print() << "\n Building twocylinders geometry." << "\n";
+        amrex::Print() << "\n Building twocylinders geometry." << "\n";
         make_eb_twocylinders();
     }
     else if(geom_type == "spherecube")
     {
-    amrex::Print() << "\n Building spherecube geometry." << "\n";
+        amrex::Print() << "\n Building spherecube geometry." << "\n";
         make_eb_spherecube();
     }
     else if(geom_type == "tuscan")
     {
-    amrex::Print() << "\n Building tuscan geometry." << "\n";
+        amrex::Print() << "\n Building tuscan geometry." << "\n";
         make_eb_tuscan();
     }
 #endif
     else if(geom_type == "annulus")
     {
-    amrex::Print() << "\n Building annulus geometry." << "\n";
+        amrex::Print() << "\n Building annulus geometry." << "\n";
         make_eb_annulus();
     }
     else if(geom_type == "sphere")
     {
-    amrex::Print() << "\n Building sphere geometry." << "\n";
+        amrex::Print() << "\n Building sphere geometry." << "\n";
         make_eb_sphere();
     }
     else if(geom_type == "jcap")
     {
-    amrex::Print() << "\n Building JCAP geometry." << "\n";
+        amrex::Print() << "\n Building JCAP geometry." << "\n";
         make_eb_cyl_tuscan();
     }
     else if(geom_type == "chkptfile")
     {
-       make_eb_chkptfile();
+        make_eb_chkptfile();
+    }
+    else if(geom_type == "stl") {
+        make_eb_stl();
     }
 #ifdef CSG_EB
-    else if(!csg_file.empty()) {
-      amrex::Print() << "\n Building geometry from .csg file:  " << csg_file << "\n";
-      make_eb_csg(csg_file);
+    else if(!csg_file.empty() || geom_type == "cgs") {
+        amrex::Print() << "\n Building geometry from .csg file:  " << csg_file << "\n";
+        make_eb_csg(csg_file);
     }
 #endif
     else
     {
-    amrex::Print() << "\n No EB geometry declared in inputs => "
-                   << " Will build all regular geometry." << "\n";
+        amrex::Print() << "\n No EB geometry declared in inputs => "
+                       << " Will build all regular geometry." << "\n";
         make_eb_regular();
     }
     amrex::Print() << "Done making the EB geometry index space.\n" << "\n";
 
     if (m_write_geom_chk) {
-       const auto& is = amrex::EB2::IndexSpace::top();
-       const auto& eb_level = is.getLevel(geom.back());
-       eb_level.write_to_chkpt_file("geom_chk", amrex::EB2::ExtendDomainFace(), amrex::EB2::max_grid_size);
+        const auto& is = amrex::EB2::IndexSpace::top();
+        const auto& eb_level = is.getLevel(geom.back());
+        eb_level.write_to_chkpt_file("geom_chk", amrex::EB2::ExtendDomainFace(), amrex::EB2::max_grid_size);
     }
 }

src/incflo.H

@@ -1004,6 +1004,7 @@ private:
     void make_eb_cyl_tuscan ();
     void make_eb_tuscan ();
     void make_eb_chkptfile ();
+    void make_eb_stl ();
 #ifdef CSG_EB
     void make_eb_csg (const std::string& csg_file);
 #endif