From 8ee92d0ee0d6126464f4abb64a98c230f9fe4b7c Mon Sep 17 00:00:00 2001
From: psakiev <psakiev@sandia.gov>
Date: Mon, 23 Oct 2023 16:45:27 -0600
Subject: [PATCH] Add option to dump mesh on failed jacobian check

---
 include/Realm.h                          |  3 ++-
 src/Realm.C                              |  9 +++++----
 src/ngp_algorithms/GeometryInteriorAlg.C | 12 ++++++++++--
 3 files changed, 17 insertions(+), 7 deletions(-)

diff --git a/include/Realm.h b/include/Realm.h
index 343e5d5757..df9325f975 100644
--- a/include/Realm.h
+++ b/include/Realm.h
@@ -183,7 +183,7 @@ class Realm
   void compute_vrtm(const std::string& = "velocity");
   void compute_l2_scaling();
   void output_converged_results();
-  void provide_output();
+  void provide_output(bool forcedOutput=false);
   void provide_restart_output();
 
   void register_interior_algorithm(stk::mesh::Part* part);
@@ -481,6 +481,7 @@ class Realm
 
   // check if there are negative Jacobians
   bool checkJacobians_;
+  bool outputFailedJacobians_;
 
   // types of physics
   bool isothermalFlow_;
diff --git a/src/Realm.C b/src/Realm.C
index 69743afcb3..8b1cfec7e8 100644
--- a/src/Realm.C
+++ b/src/Realm.C
@@ -240,6 +240,7 @@ Realm::Realm(Realms& realms, const YAML::Node& node)
     edgesPart_(0),
     checkForMissingBcs_(false),
     checkJacobians_(false),
+    outputFailedJacobians_(false),
     isothermalFlow_(true),
     uniformFlow_(true),
     provideEntityCount_(false),
@@ -670,6 +671,7 @@ Realm::load(const YAML::Node& node)
 
   // check for bad Jacobians in the mesh
   get_if_present(node, "check_jacobians", checkJacobians_, checkJacobians_);
+  get_if_present(node, "output_on_failed_jacobian_check", outputFailedJacobians_, outputFailedJacobians_);
 
   // entity count
   get_if_present(
@@ -3152,7 +3154,7 @@ Realm::overset_field_update(
 //-------- provide_output --------------------------------------------------
 //--------------------------------------------------------------------------
 void
-Realm::provide_output()
+Realm::provide_output(bool forcedOutput)
 {
   stk::diag::TimeBlock mesh_output_timeblock(Simulation::outputTimer());
   const double start_time = NaluEnv::self().nalu_time();
@@ -3169,7 +3171,6 @@ Realm::provide_output()
     const int modStep = timeStepCount - outputInfo_->outputStart_;
 
     // check for elapsed WALL time threshold
-    bool forcedOutput = false;
     if (outputInfo_->userWallTimeResults_.first) {
       const double elapsedWallTime = stk::wall_time() - wallTimeStart_;
       // find the max over all core
@@ -3179,8 +3180,7 @@ Realm::provide_output()
         1);
       // convert to hours
       g_elapsedWallTime /= 3600.0;
-      // only force output the first time the timer is exceeded
-      if (g_elapsedWallTime > outputInfo_->userWallTimeResults_.second) {
+      if (g_elapsedWallTime > outputInfo_->userWallTimeResults_.second || forcedOutput) {
         forcedOutput = true;
         outputInfo_->userWallTimeResults_.first = false;
         NaluEnv::self().naluOutputP0() << "Realm::provide_output()::Forced "
@@ -3234,6 +3234,7 @@ Realm::provide_output()
   }
 }
 
+
 //--------------------------------------------------------------------------
 //-------- provide_restart_output ------------------------------------------
 //--------------------------------------------------------------------------
diff --git a/src/ngp_algorithms/GeometryInteriorAlg.C b/src/ngp_algorithms/GeometryInteriorAlg.C
index 33b210648b..3df57ec29a 100644
--- a/src/ngp_algorithms/GeometryInteriorAlg.C
+++ b/src/ngp_algorithms/GeometryInteriorAlg.C
@@ -57,8 +57,16 @@ template <typename AlgTraits>
 void
 GeometryInteriorAlg<AlgTraits>::execute()
 {
-  if (realm_.checkJacobians_)
-    impl_negative_jacobian_check();
+  if (realm_.checkJacobians_){
+    try{
+      impl_negative_jacobian_check();
+    }
+    catch (const std::exception& e){
+      // dump exodus file if the user enabled this feature then rethrow
+      realm_.provide_output(realm_.outputFailedJacobians_);
+      throw e;
+    }
+  }
 
   impl_compute_dual_nodal_volume();