prisms-center · landinjm · Jan 6, 2025 · Jan 6, 2025
diff --git a/applications/CHiMaD_benchmarks/CHiMaD_benchmark1c/customPDE.h b/applications/CHiMaD_benchmarks/CHiMaD_benchmark1c/customPDE.h
@@ -1,5 +1,3 @@
-#include <deal.II/base/exceptions.h>
-
 #include <core/matrixFreePDE.h>
 
 using namespace dealii;
@@ -102,6 +100,8 @@ class customPDE : public MatrixFreePDE<dim, degree>
 
 #include <deal.II/grid/grid_generator.h>
 
+#include <core/exceptions.h>
+
 template <int dim, int degree>
 void
 customPDE<dim, degree>::create_triangulation(

diff --git a/applications/CHiMaD_benchmarks/CHiMaD_benchmark6b/customPDE.h b/applications/CHiMaD_benchmarks/CHiMaD_benchmark6b/customPDE.h
@@ -106,6 +106,8 @@ class customPDE : public MatrixFreePDE<dim, degree>
 
 #include <deal.II/grid/grid_generator.h>
 
+#include <core/exceptions.h>
+
 template <int dim, int degree>
 void
 customPDE<dim, degree>::create_triangulation(

diff --git a/applications/allenCahn_conserved/customPDE.h b/applications/allenCahn_conserved/customPDE.h
@@ -106,6 +106,8 @@ class customPDE : public MatrixFreePDE<dim, degree>
 // solve each time increment
 #include <deal.II/lac/solver_cg.h>
 
+#include <core/exceptions.h>
+
 template <int dim, int degree>
 void
 customPDE<dim, degree>::solveIncrement(bool skip_time_dependent)
@@ -298,10 +300,10 @@ customPDE<dim, degree>::solveIncrement(bool skip_time_dependent)
                             }
                           else
                             {
-                              std::cerr << "PRISMS-PF Error: Nonlinear solver "
-                                           "tolerance types other than ABSOLUTE_CHANGE "
-                                           "have yet to be implemented."
-                                        << std::endl;
+                              AssertThrow(
+                                false,
+                                FeatureNotImplemented(
+                                  "Nonlinear solver tolerances besides ABSOLUTE_CHANGE"));
                             }
                         }
                     }

diff --git a/include/core/boundary_conditions/nonUniformDirichletBC.h b/include/core/boundary_conditions/nonUniformDirichletBC.h
@@ -24,7 +24,8 @@ class NonUniformDirichletBC : public dealii::Function<dim>
 
   // IC for scalar values
   [[nodiscard]] double
-  value(const dealii::Point<dim> &p, const unsigned int component = 0) const override
+  value(const dealii::Point<dim>           &p,
+        [[maybe_unused]] const unsigned int component = 0) const override
   {
     double                 scalar_BC = 0.0;
     dealii::Vector<double> vector_BC(dim);

diff --git a/include/core/exceptions.h b/include/core/exceptions.h
@@ -0,0 +1,21 @@
+#ifndef EXCEPTIONS_H
+#define EXCEPTIONS_H
+
+#include <deal.II/base/exceptions.h>
+
+using namespace dealii;
+
+/**
+ * Exception for parts of the library that have yet to be implemented yet. The argument is
+ * used to provide context for the feature that has yet to be implemented.
+ */
+DeclException1(
+  FeatureNotImplemented,
+  std::string,
+  << "The following feature has yet to be implemented in PRISMS-PF:\n"
+  << arg1
+  << "\nCheck the issues section of PRISMS-PF's github to see if this feature is under "
+     "development. Additionally, please considering provided a patch to PRISMS-PF if you "
+     "feel that feature is worthwhile for yourself and others.");
+
+#endif
diff --git a/include/core/initial_conditions/initialConditions.h b/include/core/initial_conditions/initialConditions.h
@@ -1,10 +1,3 @@
-/*
- * initialConditions.h
- *
- *  Created on: Feb 27, 2017
- *      Author: stephendewitt
- */
-
 #ifndef INCLUDE_INITIALCONDITIONS_H_
 #define INCLUDE_INITIALCONDITIONS_H_
 
@@ -32,7 +25,8 @@ class InitialCondition : public dealii::Function<dim>
 
   // IC for scalar values
   [[nodiscard]] double
-  value(const dealii::Point<dim> &p, const unsigned int component = 0) const override
+  value(const dealii::Point<dim>           &p,
+        [[maybe_unused]] const unsigned int component = 0) const override
   {
     double                 scalar_IC = 0.0;
     dealii::Vector<double> vector_IC(dim);
@@ -78,4 +72,4 @@ class InitialConditionVector : public dealii::Function<dim>
   MatrixFreePDE<dim, degree> *matrix_free_pde;
 };
 
-#endif /* INCLUDE_INITIALCONDITIONS_H_ */
+#endif
diff --git a/include/core/matrixFreePDE.h b/include/core/matrixFreePDE.h
@@ -3,7 +3,6 @@
 
 // dealii headers
 #include <deal.II/base/config.h>
-#include <deal.II/base/exceptions.h>
 #include <deal.II/base/quadrature.h>
 #include <deal.II/base/timer.h>
 #include <deal.II/distributed/solution_transfer.h>

diff --git a/include/core/userInputParameters.h b/include/core/userInputParameters.h
@@ -5,7 +5,6 @@
 #define INCLUDE_USERINPUTPARAMETERS_H_
 
 #include <deal.II/base/conditional_ostream.h>
-#include <deal.II/base/exceptions.h>
 #include <deal.II/base/point.h>
 #include <deal.II/lac/la_parallel_vector.h>
 #include <deal.II/lac/vector.h>
@@ -16,6 +15,7 @@
 #include <boost/variant.hpp>
 
 #include <core/boundary_conditions/varBCs.h>
+#include <core/exceptions.h>
 #include <core/inputFileReader.h>
 #include <core/model_variables.h>
 #include <core/refinement/RefinementCriterion.h>

diff --git a/include/core/variableContainer.h b/include/core/variableContainer.h
@@ -3,7 +3,6 @@
 #ifndef VARIBLECONTAINER_H
 #define VARIBLECONTAINER_H
 
-#include <deal.II/base/exceptions.h>
 #include <deal.II/lac/vector.h>
 #include <deal.II/matrix_free/evaluation_flags.h>
 #include <deal.II/matrix_free/fe_evaluation.h>

diff --git a/include/core/variableValueContainer.h b/include/core/variableValueContainer.h
@@ -43,26 +43,4 @@ class variableValueContainer
   unsigned int              num_entries;
 };
 
-// void variableValueContainer::set(unsigned int global_variable_index, double
-// variable_value){
-//     var_index.push_back(global_variable_index);
-//     value.push_back(variable_value);
-//     num_entries++;
-// }
-//
-// double variableValueContainer::operator()(unsigned int
-// global_variable_index){
-//     for (unsigned int i=0; i < num_entries; i++){
-//         if (global_variable_index == var_index[i]){
-//             return value[i];
-//         }
-//     }
-//
-//     // If this point is reached, the index isn't allowed so an error should
-//     be thrown std::cerr << "PRISMS-PF Error: Attempted access of a variable
-//     value that was not marked as needed in 'parameters.in'. Double-check the
-//     indices in user functions where a variable value is requested." <<
-//     std::endl; abort();
-// }
-
 #endif
diff --git a/include/field_input/IntegrationTools/pfield/Mesh.hh b/include/field_input/IntegrationTools/pfield/Mesh.hh
@@ -1,17 +1,17 @@
 #ifndef Mesh_HH
 #define Mesh_HH
 
-#include <deal.II/base/exceptions.h>
-
 #include <algorithm>
 #include <cmath>
+#include <core/exceptions.h>
 #include <cstdlib>
 #include <field_input/IntegrationTools/datastruc/Bin.hh>
 #include <field_input/IntegrationTools/pfield/interpolation/Hexahedron.hh>
 #include <field_input/IntegrationTools/pfield/interpolation/Interpolator.hh>
 #include <field_input/IntegrationTools/pfield/interpolation/Quad.hh>
 #include <field_input/IntegrationTools/pfunction/PFuncBase.hh>
 #include <fstream>
+#include <math.h>
 #include <sstream>
 
 namespace PRISMS
@@ -329,15 +329,20 @@ namespace PRISMS
       std::cout << "Read rectilinear file and create mesh" << std::endl;
 
       bool               mesh_as_points = true;
-      std::vector<float> x_coord, y_coord, z_coord;
+      std::vector<float> x_coord;
+      std::vector<float> y_coord;
+      std::vector<float> z_coord;
 
       std::istringstream ss;
-      std::string        line, str, type;
+      std::string        line;
+      std::string        str;
+      std::string        type;
 
-      unsigned int uli_dummy;
-      double       d_dummy;
+      unsigned int uli_dummy = 0;
 
-      unsigned long int         Npoints, Ncells, Ncell_numbers, u;
+      unsigned long int         Npoints = 0;
+      unsigned long int         Ncells  = 0;
+      unsigned long int         u       = 0;
       std::vector<unsigned int> cell_node;
 
       PRISMS::Coordinate<dim> _coord;
@@ -364,7 +369,7 @@ namespace PRISMS
                   std::cout << "  reserve OK" << std::endl;
                   for (unsigned int i = 0; i < Npoints; i++)
                     {
-                      float temp_coord;
+                      float temp_coord = NAN;
 
                       infile >> temp_coord;
 
@@ -389,7 +394,7 @@ namespace PRISMS
                   std::cout << "  reserve OK" << std::endl;
                   for (unsigned int i = 0; i < Npoints; i++)
                     {
-                      float temp_coord;
+                      float temp_coord = NAN;
 
                       infile >> temp_coord;
 
@@ -414,7 +419,7 @@ namespace PRISMS
                   std::cout << "  reserve OK" << std::endl;
                   for (unsigned int i = 0; i < Npoints; i++)
                     {
-                      float temp_coord;
+                      float temp_coord = NAN;
 
                       infile >> temp_coord;
 
@@ -440,7 +445,9 @@ namespace PRISMS
             }
 
           // interpolated coordinates for each node of a cell
-          std::vector<float> COORD_X(Npoints), COORD_Y(Npoints), COORD_Z(Npoints);
+          std::vector<float> COORD_X(Npoints);
+          std::vector<float> COORD_Y(Npoints);
+          std::vector<float> COORD_Z(Npoints);
 
           u = 0; //(defined at the beginning of read_vtk)
           if (dim > 2)
@@ -527,7 +534,9 @@ namespace PRISMS
                 }
 
               for (int m = 0; m < dim; m++)
-                add_once(value[m], hist[m], _coord[m]);
+                {
+                  add_once(value[m], hist[m], _coord[m]);
+                }
 
               _node.push_back(_coord);
             }
@@ -560,11 +569,15 @@ namespace PRISMS
             }
           std::cout << "  Min Coordinate: ";
           for (int j = 0; j < dim; j++)
-            std::cout << _min[j] << " ";
+            {
+              std::cout << _min[j] << " ";
+            }
           std::cout << std::endl;
           std::cout << "  Max Coordinate: ";
           for (int j = 0; j < dim; j++)
-            std::cout << _max[j] << " ";
+            {
+              std::cout << _max[j] << " ";
+            }
           std::cout << std::endl;
 
           std::cout << "  done" << std::endl;

diff --git a/include/field_input/IntegrationTools/pfield/interpolation/Interpolator.hh b/include/field_input/IntegrationTools/pfield/interpolation/Interpolator.hh
@@ -58,14 +58,14 @@ namespace PRISMS
     }
 
     virtual bool
-    is_in_range(const Coordinate &coord)
+    is_in_range([[maybe_unused]] const Coordinate &coord)
     {
       undefined("bool is_in_range(Coordinate coord) const");
       return false;
     }
 
     virtual double
-    operator()(const Coordinate &coord)
+    operator()([[maybe_unused]] const Coordinate &coord)
     {
       undefined("double operator()(Coordinate coord)");
       return double();

diff --git a/include/grains/OrderParameterRemapper.h b/include/grains/OrderParameterRemapper.h
@@ -25,8 +25,7 @@ class OrderParameterRemapper
     std::vector<SimplifiedGrainRepresentation<dim>> &grain_representations,
     std::vector<dealii::LinearAlgebra::distributed::Vector<double> *> &solution_fields,
     dealii::DoFHandler<dim>                                           &dof_handler,
-    unsigned int                                                       dofs_per_cell,
-    double                                                             buffer);
+    unsigned int                                                       dofs_per_cell);
 
   /**
    * This method does the core work of the class to reassign grains across
@@ -39,8 +38,7 @@ class OrderParameterRemapper
     const dealii::LinearAlgebra::distributed::Vector<double> *grain_index_field,
     std::vector<dealii::LinearAlgebra::distributed::Vector<double> *> &solution_fields,
     dealii::DoFHandler<dim>                                           &dof_handler,
-    unsigned int                                                       dofs_per_cell,
-    double                                                             buffer);
+    unsigned int                                                       dofs_per_cell);
 
 protected:
 };

diff --git a/include/utilities/utilities.h b/include/utilities/utilities.h
@@ -4,6 +4,7 @@
 #include <deal.II/base/tensor.h>
 #include <deal.II/base/vectorization.h>
 
+#include <core/exceptions.h>
 #include <vector>
 
 /**

diff --git a/src/core/checkpoint.cc b/src/core/checkpoint.cc
@@ -1,8 +1,7 @@
-#include <deal.II/base/exceptions.h>
-
 #include <boost/archive/binary_iarchive.hpp>
 #include <boost/archive/binary_oarchive.hpp>
 
+#include <core/exceptions.h>
 #include <core/matrixFreePDE.h>
 #include <filesystem>
 #include <fstream>

diff --git a/src/core/initial_conditions/initialConditions.cc b/src/core/initial_conditions/initialConditions.cc
@@ -1,6 +1,5 @@
-// methods to apply initial conditions
-
 #include <cmath>
+#include <core/exceptions.h>
 #include <core/initial_conditions/initialConditions.h>
 #include <core/matrixFreePDE.h>
 #include <field_input/IntegrationTools/PField.hh>
@@ -229,8 +228,7 @@ MatrixFreePDE<dim, degree>::applyInitialConditions()
         &grain_index_field,
         solutionSet,
         *dofHandlersSet_nonconst.at(scalar_field_index),
-        FESet.at(scalar_field_index)->dofs_per_cell,
-        userInputs.buffer_between_grains);
+        FESet.at(scalar_field_index)->dofs_per_cell);
 
       // Smooth the order parameters according to Fick's 2nd Law
       // In the time cycle below, we evolve the weak form of Eq.:
@@ -423,29 +421,3 @@ MatrixFreePDE<dim, degree>::applyInitialConditions()
         }
     }
 }
-
-// =================================================================================
-
-// I don't think vector fields are implemented in PFields yet
-// template <int dim>
-// class InitialConditionPFieldVec : public Function<dim>
-//{
-// public:
-//  unsigned int index;
-//  Vector<double> values;
-//  typedef PRISMS::PField<double*, double, 2> ScalarField2D;
-//  ScalarField2D &inputField;
-//
-//  InitialConditionPFieldVec (const unsigned int _index, ScalarField2D
-//  &_inputField) : Function<dim>(1), index(_index), inputField(_inputField) {}
-//
-//  void vector_value (const Point<dim> &p,Vector<double> &vector_IC) const
-//  {
-//	  double coord[dim];
-//	  for (unsigned int i = 0; i < dim; i++){
-//		  coord[i] = p(i);
-//	  }
-//
-//	  vector_IC = inputField(coord);
-//  }
-//};
diff --git a/src/core/inputFileReader.cc b/src/core/inputFileReader.cc
@@ -1,6 +1,7 @@
 #include <deal.II/base/mpi.h>
 #include <deal.II/base/utilities.h>
 
+#include <core/exceptions.h>
 #include <core/inputFileReader.h>
 #include <core/refinement/RefinementCriterion.h>
 #include <fstream>