constant relaxation for density. Refs enrico-dev#9

include/enrico/coupled_driver.h

@@ -105,6 +105,10 @@ class CoupledDriver {
   //! relaxation aplied to the heat source if not set
   double alpha_T_{alpha_};
 
+  //! Constant relaxation factor for the density, defaults to the
+  //! relaxation applied to the heat source if not set
+  double alpha_d_{alpha_};
+
   //! Enumeration of available temperature initial condition specifications.
   //! 'neutronics' sets temperature condition from the neutronics input files,
   //! while 'heat' sets temperature based on a thermal-fluids input (or restart) file.

src/coupled_driver.cpp

@@ -25,6 +25,9 @@ CoupledDriver::CoupledDriver(MPI_Comm comm, pugi::xml_node node)
   if (node.child("alpha_T"))
     alpha_T_ = node.child("alpha_T").text().as_double();
 
+  if (node.child("alpha_d"))
+    alpha_d_ = node.child("alpha_d").text().as_double();
+
   if (node.child("temperature_ic")) {
     auto s = std::string{node.child_value("temperature_ic")};
 
@@ -54,6 +57,8 @@ CoupledDriver::CoupledDriver(MPI_Comm comm, pugi::xml_node node)
   Expects(max_picard_iter_ >= 0);
   Expects(epsilon_ > 0);
   Expects(alpha_ > 0);
+  Expects(alpha_T_ > 0);
+  Expects(alpha_d_ > 0);
 }
 
 void CoupledDriver::execute()
@@ -207,7 +212,7 @@ void CoupledDriver::update_density()
     auto d = heat.density();
 
     if (heat.has_coupling_data())
-      densities_ = d;
+      densities_ = alpha_d_ * d + (1.0 - alpha_d_) * densities_prev_;
   }
 
   // Set density in the neutronics solver