Making shear and P wave speed user defined quantities for seismic for…

…ce and non reflecting BC idaholab#183

include/bcs/NonReflectingBC.h

@@ -74,11 +74,11 @@ class NonReflectingBC : public IntegratedBC
   /// Density of the soil
   const MaterialProperty<Real> & _density;
 
-  /// Shear wave speed of the soil
-  const MaterialProperty<Real> & _shear_wave_speed;
-
   /// P wave speed of the soil
-  const MaterialProperty<Real> & _P_wave_speed;
+  const Real & _p_wave_speed;
+
+  /// Shear wave speed of the soil
+  const Real & _shear_wave_speed;
 };
 
 #endif // NONREFLECTINGBC_H

include/bcs/SeismicForce.h

@@ -49,11 +49,11 @@ class SeismicForce : public IntegratedBC
   /// density of the soil
   const MaterialProperty<Real> & _density;
 
-  /// shear wave speed of the soil
-  const MaterialProperty<Real> & _shear_wave_speed;
-
   /// P wave speed of the soil
-  const MaterialProperty<Real> & _P_wave_speed;
+  const Real & _p_wave_speed;
+
+  /// shear wave speed of the soil
+  const Real & _shear_wave_speed;
 };
 
 #endif // SEISMICFORCE_H

src/actions/SeismicForceAction.C

@@ -48,6 +48,10 @@ validParams<SeismicForceAction>()
                                              "input_component.");
   params.addParam<Real>(
       "alpha", 0.0, "The alpha parameter required for HHT time integration scheme.");
+  params.addRequiredRangeCheckedParam<Real>(
+      "p_wave_speed", "p_wave_speed>0.0", "P-wave speed of the material.");
+  params.addRequiredRangeCheckedParam<Real>(
+      "shear_wave_speed", "shear_wave_speed>0.0", "shear wave speed of the material.");
   return params;
 }
 
@@ -97,6 +101,8 @@ SeismicForceAction::act()
   InputParameters params = _factory.getValidParams("SeismicForce");
   params.set<std::vector<BoundaryName>>("boundary") = _boundary;
   params.set<Real>("alpha") = getParam<Real>("alpha");
+  params.set<Real>("p_wave_speed") = getParam<Real>("p_wave_speed");
+  params.set<Real>("shear_wave_speed") = getParam<Real>("shear_wave_speed");
 
   // Loop through each component and create SeismicForce object
   std::string prefix = "mastodon_seismic_force";

src/bcs/NonReflectingBC.C

@@ -52,6 +52,10 @@ NonReflectingBC::commonParameters()
   params.addRequiredParam<Real>("beta", "The beta parameter for Newmark time integration.");
   params.addRequiredParam<Real>("gamma", "The gamma parameter for Newmark time integration.");
   params.addParam<Real>("alpha", 0.0, "The alpha parameter for HHT time integration.");
+  params.addRequiredRangeCheckedParam<Real>(
+      "p_wave_speed", "p_wave_speed>0.0", "P-wave speed of the material.");
+  params.addRequiredRangeCheckedParam<Real>(
+      "shear_wave_speed", "shear_wave_speed>0.0", "shear wave speed of the material.");
   return params;
 }
 
@@ -68,8 +72,8 @@ NonReflectingBC::NonReflectingBC(const InputParameters & parameters)
     _gamma(getParam<Real>("gamma")),
     _alpha(getParam<Real>("alpha")),
     _density(getMaterialProperty<Real>("density")),
-    _shear_wave_speed(getMaterialProperty<Real>("shear_wave_speed")),
-    _P_wave_speed(getMaterialProperty<Real>("P_wave_speed"))
+    _p_wave_speed(getParam<Real>("p_wave_speed")),
+    _shear_wave_speed(getParam<Real>("shear_wave_speed"))
 {
 
   // Error checking on variable vectors
@@ -122,17 +126,17 @@ NonReflectingBC::computeQpResidual()
   // residual is test[i][_qp] *( density* V_p * normal component of velocity +
   // density * V_s* tangential component of velocity)
   return _test[_i][_qp] * _density[_qp] *
-         (_P_wave_speed[_qp] * normal_vel * _normals[_qp](_component) +
-          _shear_wave_speed[_qp] * (vel[_component] - normal_vel * _normals[_qp](_component)));
+         (_p_wave_speed * normal_vel * _normals[_qp](_component) +
+          _shear_wave_speed * (vel[_component] - normal_vel * _normals[_qp](_component)));
 }
 
 Real
 NonReflectingBC::computeQpJacobian()
 {
   return _test[_i][_qp] * _density[_qp] *
-         ((_P_wave_speed[_qp] - _shear_wave_speed[_qp]) * _normals[_qp](_component) *
+         ((_p_wave_speed - _shear_wave_speed) * _normals[_qp](_component) *
               _normals[_qp](_component) +
-          _shear_wave_speed[_qp]) *
+          _shear_wave_speed) *
          (1. + _alpha) * _gamma / _beta / _dt * _phi[_j][_qp];
 }
 
@@ -150,7 +154,7 @@ NonReflectingBC::computeQpOffDiagJacobian(unsigned int jvar)
     }
 
   if (active)
-    return _test[_i][_qp] * _density[_qp] * (_P_wave_speed[_qp] - _shear_wave_speed[_qp]) *
+    return _test[_i][_qp] * _density[_qp] * (_p_wave_speed - _shear_wave_speed) *
            _normals[_qp](_component) * _normals[_qp](coupled_component) * (1. + _alpha) * _gamma /
            _beta / _dt * _phi[_j][_qp];
 

src/bcs/SeismicForce.C

@@ -33,6 +33,10 @@ validParams<SeismicForce>()
   params.addParam<FunctionName>("velocity_function",
                                 "The function that describes the input ground velocity.");
   params.addCoupledVar("velocity", "The variable that describes the input velocity.");
+  params.addRequiredRangeCheckedParam<Real>(
+      "p_wave_speed", "p_wave_speed>0.0", "P-wave speed of the material.");
+  params.addRequiredRangeCheckedParam<Real>(
+      "shear_wave_speed", "shear_wave_speed>0.0", "shear wave speed of the material.");
   params.addParam<Real>(
       "alpha", 0.0, "The alpha parameter required for HHT time integration scheme.");
   params.set<bool>("use_displaced_mesh") = true;
@@ -46,8 +50,8 @@ SeismicForce::SeismicForce(const InputParameters & parameters)
     _factor(getParam<Real>("factor")),
     _alpha(getParam<Real>("alpha")),
     _density(getMaterialProperty<Real>("density")),
-    _shear_wave_speed(getMaterialProperty<Real>("shear_wave_speed")),
-    _P_wave_speed(getMaterialProperty<Real>("P_wave_speed"))
+    _p_wave_speed(getParam<Real>("p_wave_speed")),
+    _shear_wave_speed(getParam<Real>("shear_wave_speed"))
 {
   if (_component >= _mesh.dimension())
     mooseError("Invalid value for 'component' (",
@@ -98,6 +102,6 @@ SeismicForce::computeQpResidual()
     tangential_vel[i] = vel[i] - normal_vel * _normals[_qp](i);
 
   return _factor * _test[_i][_qp] * _density[_qp] *
-         (_P_wave_speed[_qp] * normal_vel * _normals[_qp](_component) +
-          _shear_wave_speed[_qp] * tangential_vel[_component]);
+         (_p_wave_speed * normal_vel * _normals[_qp](_component) +
+          _shear_wave_speed * tangential_vel[_component]);
 }

test/tests/bcs/nonreflecting_bc/non_reflecting_bc_error.i

@@ -23,6 +23,8 @@
    boundary = 'left'
    beta = 1
    gamma = 1
+   shear_wave_speed = 1.0
+   p_wave_speed = 3.5
  [../]
 []
 

test/tests/bcs/nonreflecting_bc/non_reflecting_bc_test.i

@@ -161,6 +161,8 @@
       beta = 0.25
       gamma = 0.5
       boundary = 'back'
+      shear_wave_speed = 10.0
+      p_wave_speed = 16.329931618554521
     [../]
   [../]
   [./Periodic]

test/tests/bcs/seismic_force/seismic_force_action_error.i

@@ -25,6 +25,8 @@
       displacements = 'disp_x disp_y'
       input_components = '0 1'
       boundary = 'front'
+      shear_wave_speed = 1.0
+      p_wave_speed = 3.5
     [../]
   [../]
 []

test/tests/bcs/seismic_force/seismic_force_error.i

@@ -24,6 +24,8 @@
     type = SeismicForce
     variable = 'disp_x'
     boundary = 'left'
+    shear_wave_speed = 1.0
+    p_wave_speed = 3.5
   [../]
 []
 

test/tests/bcs/seismic_force/seismic_force_test_x_y.i

@@ -151,6 +151,8 @@
       input_components = '0 1'
       boundary = 'front'
       velocity_functions = 'x_vel y_vel' # input velocity functions
+      shear_wave_speed = 1.0
+      p_wave_speed = 3.5 
     [../]
   [../]
   [./back_x]

test/tests/bcs/seismic_force/seismic_force_test_z.i

@@ -196,6 +196,8 @@
        input_components = '2'
        boundary = 'front'
        velocity_functions = 'z_vel'
+       shear_wave_speed = 3718.204815052883833
+       p_wave_speed = 7381.527628151251754 
    [../]
   [../]
   [./NonReflectingBC]
@@ -206,6 +208,8 @@
         beta = 0.25
         gamma = 0.5
         boundary = 'back'
+        shear_wave_speed = 3718.204815052883833
+        p_wave_speed = 7381.527628151251754 
     [../]
   [../]
   [./Periodic]