Make potential_reduction_length optional

Author: davidrudlstorfer

src/beaminteraction/src/potential/4C_beaminteraction_potential_input.cpp

@@ -73,11 +73,10 @@ void BeamPotential::set_valid_parameters(std::map<std::string, Core::IO::InputSp
                                              "slave be assigned to beam elements?",
                                  .default_value = MasterSlaveChoice::smaller_eleGID_is_slave}));
 
-  // TODO add None as default value
-  beampotential.specs.emplace_back(parameter<double>("POTENTIAL_REDUCTION_LENGTH",
-      {.description = "Within this length of the master beam end point the potential is smoothly "
-                      "reduced to one half to account for infinitely long master beam surrogates.",
-          .default_value = -1.0}));
+  beampotential.specs.emplace_back(parameter<std::optional<double>>("POTENTIAL_REDUCTION_LENGTH",
+      {.description =
+              "Within this length of the master beam end point the potential is smoothly "
+              "reduced to one half to account for infinitely long master beam surrogates."}));
 
   beampotential.move_into_collection(list);
 

src/beaminteraction/src/potential/4C_beaminteraction_potential_pair_beam_to_beam.cpp

@@ -1068,7 +1068,7 @@ void BeamInteraction::BeamToBeamPotentialPair<numnodes, numnodalvalues, T>::
   }
 
   if (params()->strategy() == BeamPotential::Strategy::single_length_specific_small_separations &&
-      params()->potential_reduction_length() != -1.0)
+      params()->potential_reduction_length().has_value())
   {
     FOUR_C_THROW(
         "The potential reduction strategy is currently not implemented for the beam interaction "
@@ -1084,20 +1084,20 @@ void BeamInteraction::BeamToBeamPotentialPair<numnodes, numnodalvalues, T>::
   const std::optional<double> cutoff_radius = params()->cutoff_radius();
 
   // get potential reduction length
-  const double potential_reduction_length = params()->potential_reduction_length();
+  const std::optional<double> potential_reduction_length = params()->potential_reduction_length();
 
   // get length from current master beam element to beam edge
   double length_prior_left = 0.0;
   double length_prior_right = 0.0;
 
-  if (potential_reduction_length > 0.0)
+  if (potential_reduction_length.has_value())
   {
     std::tie(length_prior_left, length_prior_right) =
         params()->ele_gid_prior_length_map_.at(element2()->id());
 
     if ((length_prior_left >= 0.0) && (length_prior_right >= 0.0) &&
-        ((ele2length_ - 2 * potential_reduction_length + length_prior_left + length_prior_right) <
-            0.0))
+        ((ele2length_ - 2 * potential_reduction_length.value() + length_prior_left +
+             length_prior_right) < 0.0))
     {
       FOUR_C_THROW(
           "ERROR: Master beam endpoint reduction strategy would interfere on current master beam "
@@ -3029,8 +3029,8 @@ bool BeamInteraction::BeamToBeamPotentialPair<numnodes, numnodalvalues,
     Core::LinAlg::Matrix<1, 3, T> const& xi_master_partial_r_slave,
     Core::LinAlg::Matrix<1, 3, T> const& xi_master_partial_r_master,
     Core::LinAlg::Matrix<1, 3, T> const& xi_master_partial_r_xi_master,
-    double potential_reduction_length, double length_prior_right, double length_prior_left,
-    T& interaction_potential_GP,
+    std::optional<double> potential_reduction_length, double length_prior_right,
+    double length_prior_left, T& interaction_potential_GP,
     Core::LinAlg::Matrix<1, numnodes * numnodalvalues, double> const& N_i_slave,
     Core::LinAlg::Matrix<1, numnodes * numnodalvalues, double> const& N_i_xi_slave,
     Core::LinAlg::Matrix<1, numnodes * numnodalvalues, T> const& N_i_master,
@@ -3227,35 +3227,37 @@ bool BeamInteraction::BeamToBeamPotentialPair<numnodes, numnodalvalues,
   }
 
   // determine potential reduction factor for master beam endpoint strategy
-  if (potential_reduction_length > 0.0)
+  if (potential_reduction_length.has_value())
   {
     T left_length_to_edge = length_prior_left + ele2length_ * 0.5 * (1 + xi_master);
     T right_length_to_edge = length_prior_right + ele2length_ * 0.5 * (1 - xi_master);
     T length_to_edge = -1.0;
     bool right_node = false;
 
     // determine potential reduction factor (master beam can only consist of one element!)
-    if (left_length_to_edge < potential_reduction_length)
+    if (left_length_to_edge < potential_reduction_length.value())
     {
       potential_reduction_factor_GP =
-          0.5 - 0.5 * std::cos(M_PI * left_length_to_edge / potential_reduction_length);
+          0.5 - 0.5 * std::cos(M_PI * left_length_to_edge / potential_reduction_length.value());
       length_to_edge = left_length_to_edge;
     }
-    else if (right_length_to_edge < potential_reduction_length)
+    else if (right_length_to_edge < potential_reduction_length.value())
     {
       potential_reduction_factor_GP =
-          0.5 - 0.5 * std::cos(M_PI * right_length_to_edge / potential_reduction_length);
+          0.5 - 0.5 * std::cos(M_PI * right_length_to_edge / potential_reduction_length.value());
       length_to_edge = right_length_to_edge;
       right_node = true;
     }
 
-    if ((length_to_edge < potential_reduction_length) && (length_to_edge != -1.0))
+    if ((length_to_edge < potential_reduction_length.value()) && (length_to_edge != -1.0))
     {
-      pot_red_fac_deriv_l_edge = 0.5 * M_PI / potential_reduction_length *
-                                 std::sin(M_PI * length_to_edge / potential_reduction_length);
-      pot_red_fac_2ndderiv_l_edge = 0.5 * M_PI * M_PI /
-                                    (potential_reduction_length * potential_reduction_length) *
-                                    std::cos(M_PI * length_to_edge / potential_reduction_length);
+      pot_red_fac_deriv_l_edge =
+          0.5 * M_PI / potential_reduction_length.value() *
+          std::sin(M_PI * length_to_edge / potential_reduction_length.value());
+      pot_red_fac_2ndderiv_l_edge =
+          0.5 * M_PI * M_PI /
+          (potential_reduction_length.value() * potential_reduction_length.value()) *
+          std::cos(M_PI * length_to_edge / potential_reduction_length.value());
     }
 
     l_edge_deriv_xi_master = 0.5 * ele2length_;

src/beaminteraction/src/potential/4C_beaminteraction_potential_pair_beam_to_beam.hpp

@@ -284,8 +284,8 @@ namespace BeamInteraction
         Core::LinAlg::Matrix<1, 3, T> const& xi_master_partial_r_slave,
         Core::LinAlg::Matrix<1, 3, T> const& xi_master_partial_r_master,
         Core::LinAlg::Matrix<1, 3, T> const& xi_master_partial_r_xi_master,
-        double potential_reduction_length, double length_prior_right, double length_prior_left,
-        T& interaction_potential_GP,
+        std::optional<double> potential_reduction_length, double length_prior_right,
+        double length_prior_left, T& interaction_potential_GP,
         Core::LinAlg::Matrix<1, numnodes * numnodalvalues, double> const& N_i_slave,
         Core::LinAlg::Matrix<1, numnodes * numnodalvalues, double> const& N_i_xi_slave,
         Core::LinAlg::Matrix<1, numnodes * numnodalvalues, T> const& N_i_master,

src/beaminteraction/src/potential/4C_beaminteraction_potential_params.cpp

@@ -36,7 +36,7 @@ BeamInteraction::BeamPotentialParams::BeamPotentialParams()
       choice_master_slave_(BeamPotential::MasterSlaveChoice::vague),
       visualization_output_(false),
       params_runtime_visualization_output_btb_potential_(nullptr),
-      potential_reduction_length_(0.0)
+      potential_reduction_length_(std::nullopt)
 {
   // empty constructor
 }
@@ -121,7 +121,7 @@ void BeamInteraction::BeamPotentialParams::init(const double restart_time)
   /****************************************************************************/
   cutoff_radius_ = beam_potential_params_list.get<std::optional<double>>("CUTOFF_RADIUS");
 
-  if (cutoff_radius_.has_value() && cutoff_radius_.value() <= 0.0)
+  if (cutoff_radius_.has_value() and cutoff_radius_.value() <= 0.0)
     FOUR_C_THROW("Invalid cutoff radius! Must be positive value or null to deactivate.");
 
   /****************************************************************************/
@@ -189,10 +189,11 @@ void BeamInteraction::BeamPotentialParams::init(const double restart_time)
 
   /****************************************************************************/
   potential_reduction_length_ =
-      beam_potential_params_list.get<double>("POTENTIAL_REDUCTION_LENGTH");
+      beam_potential_params_list.get<std::optional<double>>("POTENTIAL_REDUCTION_LENGTH");
 
-  if (potential_reduction_length_ != -1.0 and potential_reduction_length_ <= 0.0)
-    FOUR_C_THROW("Invalid potential reduction length! Must be positive value or -1 to deactivate.");
+  if (potential_reduction_length_.has_value() and potential_reduction_length_.value() <= 0.0)
+    FOUR_C_THROW(
+        "Invalid potential reduction length! Must be positive value or none to deactivate.");
 
   isinit_ = true;
 }

src/beaminteraction/src/potential/4C_beaminteraction_potential_params.hpp

@@ -131,7 +131,7 @@ namespace BeamInteraction
       return params_runtime_visualization_output_btb_potential_;
     }
 
-    inline double potential_reduction_length() const
+    inline std::optional<double> potential_reduction_length() const
     {
       throw_error_if_not_init_and_setup();
       return potential_reduction_length_;
@@ -191,7 +191,7 @@ namespace BeamInteraction
     //! within this length starting from the master beam end point the potential is smoothly
     //! reduced to zero to account for infinitely long master beam surrogates and enable an
     //! axial pull-off force.
-    double potential_reduction_length_;
+    std::optional<double> potential_reduction_length_;
   };
 
 }  // namespace BeamInteraction

src/beaminteraction/src/potential/4C_beaminteraction_potential_submodel_evaluator.cpp

@@ -143,7 +143,7 @@ void BeamInteraction::SUBMODELEVALUATOR::BeamPotential::post_setup()
 {
   check_init_setup();
 
-  if (beam_potential_params().potential_reduction_length() != -1.0)
+  if (beam_potential_params().potential_reduction_length().has_value())
     setup_potential_reduction_strategy();
 
   nearby_elements_map_.clear();
@@ -599,7 +599,7 @@ void BeamInteraction::SUBMODELEVALUATOR::BeamPotential::post_read_restart()
 {
   check_init_setup();
 
-  if (beam_potential_params().potential_reduction_length() != -1.0)
+  if (beam_potential_params().potential_reduction_length().has_value())
     setup_potential_reduction_strategy();
 
   nearby_elements_map_.clear();