Add a (claude generated) logger to decrowd outputs / save output and a few other fixes

Used Claude to help generate a logger that would log any of our material model outputs to their own files per mpi rank, output any non-0 rank MFEM logging info to its own MPI rank file, and also duplicate all of the terminal output to its own logging file for future references.
Outside of that had a few bug fixes captured in this that I noticed when running some tougher material models and hitting some interesting edge cases...

Author: rcarson3

src/CMakeLists.txt

@@ -28,6 +28,7 @@ set(EXACONSTIT_HEADERS
     utilities/assembly_ops.hpp
     utilities/rotations.hpp
     utilities/strain_measures.hpp
+    utilities/unified_logger.hpp
     ./TOML_Reader/toml.hpp
     )
 
@@ -59,6 +60,7 @@ set(EXACONSTIT_SOURCES
     solvers/mechanics_solver.cpp
     utilities/dynamic_umat_loader.cpp
     utilities/mechanics_kernels.cpp
+    utilities/unified_logger.cpp
     ./umat_tests/userumat.cxx
     )
 

src/fem_operators/mechanics_operator.cpp

@@ -3,6 +3,7 @@
 #include "models/mechanics_multi_model.hpp"
 #include "utilities/mechanics_kernels.hpp"
 #include "utilities/mechanics_log.hpp"
+#include "utilities/unified_logger.hpp"
 
 #include "mfem/general/forall.hpp"
 #include "RAJA/RAJA.hpp"
@@ -175,7 +176,7 @@ void NonlinearMechOperator::Setup(const mfem::Vector &k) const
    }
    catch(const std::exception &exc) {
       // catch anything thrown within try block that derives from std::exception
-      MFEM_WARNING(exc.what());
+      MFEM_WARNING_0(exc.what());
       succeed = false;
    }
    catch(...) {

src/mechanics_driver.cpp

@@ -78,8 +78,10 @@
 #include "mfem_expt/partial_qfunc.hpp"
 #include "options/option_parser_v2.hpp"
 #include "postprocessing/postprocessing_driver.hpp"
+#include "postprocessing/postprocessing_file_manager.hpp"
 #include "sim_state/simulation_state.hpp"
 #include "utilities/mechanics_log.hpp"
+#include "utilities/unified_logger.hpp"
 
 #include "mfem.hpp"
 #include "mfem/general/forall.hpp"
@@ -129,10 +131,6 @@ int main(int argc, char *argv[])
     * - Enable detailed timing data for strong/weak scaling studies
     */
    double start = MPI_Wtime();
-   // Print MFEM version information for reproducibility and debugging
-   if (myid == 0) {
-      printf("MFEM Version: %d \n", mfem::GetVersion());
-   }
    /**
     * **PHASE 2: COMMAND LINE PROCESSING AND CONFIGURATION**
     */
@@ -157,6 +155,11 @@ int main(int argc, char *argv[])
       return 1;
    }
 
+   // Print MFEM version information for reproducibility and debugging
+   if (myid == 0) {
+      printf("MFEM Version: %d \n", mfem::GetVersion());
+   }
+
    /*
     * Configuration File Parsing:
     * - Load complete simulation configuration from TOML file
@@ -165,6 +168,10 @@ int main(int argc, char *argv[])
     */
    ExaOptions toml_opt;
    toml_opt.parse_options(toml_file, myid);
+
+   exaconstit::UnifiedLogger& logger = exaconstit::UnifiedLogger::getInstance();
+   logger.initialize(toml_opt);
+
    toml_opt.print_options();
 
    /**
@@ -312,7 +319,6 @@ int main(int argc, char *argv[])
        * - Update time-dependent material properties and boundary conditions
        * - Prepare solver state for current time increment
        */
-      const double sim_time = sim_state->getTime();
 
       /*
        * Boundary Condition Change Detection:
@@ -347,7 +353,7 @@ int main(int argc, char *argv[])
        */
       sim_state->finishCycle();
       oper.UpdateModel();
-      post_process.Update(ti, sim_time);
+      post_process.Update(ti, sim_state->getTrueCyleTime());
    } // end loop over time steps
 
    /**
@@ -371,7 +377,7 @@ int main(int argc, char *argv[])
    if (myid == 0) {
       printf("The process took %lf seconds to run\n", (avg_sim_time / world_size));
    }
-
+   logger.shutdown();
 } // End of main simulation scope for proper resource cleanup
 
    /*

src/models/mechanics_ecmech.cpp

@@ -2,6 +2,7 @@
 #include "models/mechanics_model.hpp"
 #include "utilities/mechanics_log.hpp"
 #include "utilities/mechanics_kernels.hpp"
+#include "utilities/unified_logger.hpp"
 
 #include "mfem.hpp"
 #include "mfem/general/forall.hpp"
@@ -424,6 +425,11 @@ void ExaCMechModel::ModelSetup(const int nqpts, const int nelems, const int /*sp
                                const int nnodes, const mfem::Vector &jacobian,
                                const mfem::Vector &loc_grad, const mfem::Vector &vel)
 {
+
+   auto& logger = exaconstit::UnifiedLogger::getInstance();
+   std::string material_log = logger.getMaterialLogFilename("exacmech", m_region);
+   exaconstit::UnifiedLogger::ScopedCapture capture(material_log);
+
    const int nstatev = numStateVars;
 
    const double *jacobian_array = jacobian.Read();

src/models/mechanics_multi_model.cpp

@@ -5,6 +5,7 @@
 #include "mfem_expt/partial_qfunc.hpp"
 #include "utilities/mechanics_log.hpp"
 #include "utilities/dynamic_umat_loader.hpp"
+#include "utilities/unified_logger.hpp"
 
 #include <stdexcept>
 #include <filesystem>
@@ -235,11 +236,11 @@ bool MultiExaModel::SetupChildModel(int region_idx, const int nqpts, const int n
         return true;
     }
     catch (const std::exception& e) {
-        MFEM_WARNING("Region " + std::to_string(actual_region_id) + " failed: " + e.what());
+        MFEM_WARNING_0("[Cycle " << std::to_string(m_sim_state->getSimulationCycle() + 1) <<  " ]Region " + std::to_string(actual_region_id) + " failed: " + e.what());
         return false;
     }
     catch (...) {
-        MFEM_WARNING("Region " + std::to_string(actual_region_id) + " failed with unknown error");
+        MFEM_WARNING_0("Region " + std::to_string(actual_region_id) + " failed with unknown error");
         return false;
     }
 }

src/models/mechanics_umat.cpp

@@ -2,6 +2,7 @@
 #include "boundary_conditions/BCManager.hpp"
 #include "utilities/assembly_ops.hpp"
 #include "utilities/strain_measures.hpp"
+#include "utilities/unified_logger.hpp"
 
 #include "RAJA/RAJA.hpp"
 #include "mfem/fem/qfunction.hpp"
@@ -350,6 +351,9 @@ void AbaqusUmatModel::ModelSetup(const int nqpts, const int nelems, const int sp
                                  const int /*nnodes*/, const mfem::Vector &jacobian,
                                  const mfem::Vector & /*loc_grad*/, const mfem::Vector &/*vel*/)
 {
+   auto& logger = exaconstit::UnifiedLogger::getInstance();
+   std::string material_log = logger.getMaterialLogFilename("umat", m_region);
+   exaconstit::UnifiedLogger::ScopedCapture capture(material_log);
     // Load UMAT library if using on-demand loading
    if (use_dynamic_loading_ && load_strategy_ == DynamicUmatLoader::LoadStrategy::LOAD_ON_SETUP) {
       if (!LoadUmatLibrary()) {

src/postprocessing/projection_class.cpp

@@ -1,6 +1,7 @@
 #include "projection_class.hpp"
 
 #include "utilities/rotations.hpp"
+#include "utilities/unified_logger.hpp"
 
 #include "ECMech_const.h"
 #include "SNLS_linalg.h"
@@ -187,11 +188,11 @@ StateVariableProjection::Execute(std::shared_ptr<SimulationState> sim_state,
     m_component_length = (m_component_length == -1) ? vdim : m_component_length;
 
     if ((m_component_length + m_component_index) > vdim) {
-        MFEM_ABORT("StateVariableProjection provided a length and index that pushes us past the state variable length");
+        MFEM_ABORT_0("StateVariableProjection provided a length and index that pushes us past the state variable length");
     };
 
     if (m_component_length > state_gf->VectorDim()) {
-        MFEM_ABORT("StateVariableProjection provided length is greater than the gridfunction vector length");
+        MFEM_ABORT_0("StateVariableProjection provided length is greater than the gridfunction vector length");
     };
 
     const auto l2g = qpts2mesh.Read();
@@ -265,11 +266,11 @@ ElasticStrainProjection::Execute(std::shared_ptr<SimulationState> sim_state,
     m_component_length = (m_component_length == -1) ? vdim : m_component_length;
 
     if ((m_component_length + m_component_index) > vdim) {
-        MFEM_ABORT("ElasticStrainProjection provided a length and index that pushes us past the state variable length");
+        MFEM_ABORT_0("ElasticStrainProjection provided a length and index that pushes us past the state variable length");
     };
 
     if (m_component_length > elastic_strain_gf->VectorDim()) {
-        MFEM_ABORT("ElasticStrainProjection provided length is greater than the gridfunction vector length");
+        MFEM_ABORT_0("ElasticStrainProjection provided length is greater than the gridfunction vector length");
     };
 
     const int estrain_ind = sim_state->GetQuadratureFunctionStatePair("elastic_strain", region).first;

src/postprocessing/projection_class.hpp

@@ -710,7 +710,9 @@ class ElasticStrainProjection final : public StateVariableProjection {
                             int component_index,
                             int component_length,
                             const std::string& display_name)
-                            : StateVariableProjection(state_var_name, component_index, component_length, display_name, ptmc::EXACMECH_ONLY) {}
+                            : StateVariableProjection(state_var_name, component_index, component_length, display_name, ptmc::EXACMECH_ONLY) {
+                                m_component_length = 6;
+                            }
     /**
      * @brief Execute elastic strain projection with coordinate transformation
      * 

src/sim_state/simulation_state.cpp

@@ -288,6 +288,12 @@ TimeManagement::updateDeltaTime(const int nr_steps, const bool success) {
         updateTime();
         num_failures++;
         num_sub_steps = 1;
+        if (internal_tracker == TimeStep::FINAL) {
+            const double tf_dt = std::abs(time - time_final);
+            if (tf_dt > std::abs(1e-3 * dt)) {
+                internal_tracker = TimeStep::RETRIAL;
+            } 
+        }
         // If we've failed too many times just give up at this point
         if (num_failures > max_failures) {
             return TimeStep::FAILED;
@@ -298,6 +304,8 @@ TimeManagement::updateDeltaTime(const int nr_steps, const bool success) {
         }
     }
 
+    old_time = time;
+
     if (internal_tracker == TimeStep::FINAL) {
         internal_tracker = TimeStep::FINISHED;
         return TimeStep::FINISHED;

src/sim_state/simulation_state.hpp

@@ -55,6 +55,9 @@ class TimeManagement {
 private:
     /** @brief Current simulation time */
     double time = 0.0;
+
+    /** @brief Old simulation time */
+    double old_time = 0.0;
 
     /** @brief Final simulation time (target end time) */
     double time_final = 0.0;
@@ -140,6 +143,13 @@ class TimeManagement {
      * @return Current time value
      */
     double getTime() const { return time; }
+
+    /**
+     * @brief Get actual simulation time if auto-time stepping used
+     * 
+     * @return Actual time step value for a step
+     */
+    double getTrueCyleTime() const { return old_time; }
 
     /**
      * @brief Get current time step size
@@ -248,6 +258,22 @@ class TimeManagement {
         dt = dt_restart;
     }
 
+    /**
+     * @brief Print retrial diagnostic information
+     * 
+     * @details Outputs detailed information about cycle time step info including:
+     * - Original time step size before we reduced things down
+     * - Current time
+     * - Current cycle
+     * - Current time step size
+     * 
+     * Used for debugging convergence issues and understanding when/why
+     * retrying a time step is required.
+     */
+    void printRetrialStats() const {
+        std::cout << "[Cycle: "<< (simulation_cycle + 1)  << " , time: " << time << "] Previous attempts to converge failed step: dt old was " << dt_orig << " new dt is " << dt  << std::endl;
+    }
+
     /**
      * @brief Print sub-stepping diagnostic information
      * 
@@ -260,7 +286,7 @@ class TimeManagement {
      * sub-stepping is being triggered.
      */
     void printSubStepStats() const {
-        std::cout << "Previous attempts to converge failed but now starting sub-stepping of our desired time step: desired dt old was " << dt_orig << " sub-stepping dt is " << dt << " and number of sub-steps required is " << required_num_sub_steps << std::endl;
+        std::cout << "[Cycle: "<< (simulation_cycle + 1)  << " , time: " << time << "] Previous attempts to converge failed but now starting sub-stepping of our desired time step: desired dt old was " << dt_orig << " sub-stepping dt is " << dt << " and number of sub-steps required is " << required_num_sub_steps << std::endl;
     }
 
     /**
@@ -935,6 +961,8 @@ class SimulationState
         return GetRegionRootRank(region_id) == my_id;
     }
 
+    size_t GetMPIID() const { return my_id; }
+
     // =========================================================================
     // SOLUTION FIELD ACCESS
     // =========================================================================
@@ -978,6 +1006,13 @@ class SimulationState
      */
     double getTime() const { return m_time_manager.getTime(); }
 
+    /**
+     * @brief Get actual simulation time for a given cycle as auto-time step might have changed things
+     * 
+     * @return Current time value from TimeManagement
+     */
+    double getTrueCyleTime() const { return m_time_manager.getTrueCyleTime(); }
+
     /**
      * @brief Get current time step size
      * 
@@ -1027,6 +1062,15 @@ class SimulationState
      */
     void printTimeStats() const { m_time_manager.printTimeStats(); }
 
+    /**
+     * @brief Print retrial time step statistics
+     * 
+     * @details Outputs current time and time step information for monitoring
+     * adaptive time step behavior. Delegates to TimeManagement.
+     */
+    void printRetrialTimeStats() const { m_time_manager.printRetrialStats(); }
+
+
 private:
     /** @brief Create MPI communicators for each region containing only ranks with that region
      *  @details This prevents deadlocks in collective operations when some ranks have no