Fixing RK-2 (refs #1929)

The scheme was using F^{n} in stage 2 instead of F^{1/2} + adding a new test - re-golding 2D rk2 test (b/c of the bug fix) r22416
idaholab · Feb 14, 2014 · d1f3bab · d1f3bab
1 parent 84b727a
commit d1f3bab
Show file tree

Hide file tree

Showing 7 changed files with 120 additions and 19 deletions.
diff --git a/framework/include/timeintegrators/RungeKutta2.h b/framework/include/timeintegrators/RungeKutta2.h
@@ -39,9 +39,6 @@ class RungeKutta2 : public TimeIntegrator
   virtual void postStep(NumericVector<Number> & residual);
 
 protected:
-  NumericVector<Number> & _sln_half;
-  NumericVector<Number> & _Re_half;
-
   unsigned int _stage;
 };
 

diff --git a/framework/src/timeintegrators/RungeKutta2.C b/framework/src/timeintegrators/RungeKutta2.C
@@ -26,8 +26,6 @@ InputParameters validParams<RungeKutta2>()
 
 RungeKutta2::RungeKutta2(const std::string & name, InputParameters parameters) :
     TimeIntegrator(name, parameters),
-    _sln_half(_nl.addVector("sln_half", true, GHOSTED)),
-    _Re_half(_nl.addVector("Re_half", true, GHOSTED)),
     _stage(0)
 {
   _fe_problem.setConstJacobian(true);
@@ -46,11 +44,21 @@ void
 RungeKutta2::computeTimeDerivatives()
 {
   _u_dot  = *_nl.currentSolution();
-  _u_dot -= _nl.solutionOld();
-  _u_dot *= 1. / _dt;
-  _u_dot.close();
+  if (_stage == 1)
+  {
+    _u_dot -= _nl.solutionOld();
+    _u_dot *= 2. / _dt;
+
+    _du_dot_du = 2. / _dt;
+  }
+  else
+  {
+    _u_dot -= _nl.solutionOlder();
+    _u_dot *= 1. / _dt;
 
-  _du_dot_du = 1. / _dt;
+    _du_dot_du = 1. / _dt;
+  }
+  _u_dot.close();
   _du_dot_du.close();
 }
 
@@ -62,29 +70,23 @@ RungeKutta2::solve()
   Real time_half = (time + time_old) / 2.;
 
   _stage = 1;
-  // make sure that time derivative contribution is zero in the first pre-solve step
-  _u_dot.zero();
-  _u_dot.close();
-  _du_dot_du.zero();
-  _du_dot_du.close();
-
   _fe_problem.time() = time_half;
-  _nl.sys().solve();
+  _fe_problem.solve();
+
+  _fe_problem.copyOldSolutions();
 
   // ---------------------------------
   std::cout << " 2. stage" << std::endl;
   _stage = 2;
   _fe_problem.time() = time;
   _fe_problem.timeOld() = time_half;
-  _nl.sys().solve();
+  _fe_problem.solve();
 }
 
 void
 RungeKutta2::postStep(NumericVector<Number> & residual)
 {
   residual += _Re_non_time;
-  if (_stage == 1)
-    residual *= 0.5;
   residual += _Re_time;
   residual.close();
 }
diff --git a/test/tests/executioners/time_integration/gold/rk2-1d-linear_out.e b/test/tests/executioners/time_integration/gold/rk2-1d-linear_out.e
diff --git a/test/tests/executioners/time_integration/gold/rk2-2d-quadratic_out.e b/test/tests/executioners/time_integration/gold/rk2-2d-quadratic_out.e
diff --git a/test/tests/executioners/time_integration/rk2-1d-linear.i b/test/tests/executioners/time_integration/rk2-1d-linear.i
@@ -0,0 +1,95 @@
+[Mesh]
+  type = GeneratedMesh
+  dim = 1
+  xmin = -1
+  xmax = 1
+  nx = 20
+  elem_type = EDGE2
+[]
+
+[Functions]
+  [./ic]
+    type = ParsedFunction
+    value = 0
+  [../]
+
+  [./forcing_fn]
+    type = ParsedFunction
+    value = x
+  [../]
+
+  [./exact_fn]
+    type = ParsedFunction
+    value = t*x
+  [../]
+[]
+
+[Variables]
+  [./u]
+    order = FIRST
+    family = LAGRANGE
+
+    [./InitialCondition]
+      type = FunctionIC
+      function = ic
+    [../]
+  [../]
+[]
+
+[Kernels]
+  [./ie]
+    type = TimeDerivative
+    variable = u
+    implicit = true
+  [../]
+
+  [./diff]
+    type = Diffusion
+    variable = u
+    implicit = false
+  [../]
+
+  [./ffn]
+    type = UserForcingFunction
+    variable = u
+    function = forcing_fn
+    implicit = false
+  [../]
+[]
+
+[BCs]
+  [./all]
+    type = FunctionDirichletBC
+    variable = u
+    boundary = '0 1'
+    function = exact_fn
+  [../]
+[]
+
+[Postprocessors]
+  [./l2_err]
+    type = ElementL2Error
+    variable = u
+    function = exact_fn
+  [../]
+  
+  [./res]
+    type = Residual
+  [../]
+[]
+
+[Executioner]
+  type = Transient
+  scheme = 'rk-2'
+
+  start_time = 0.0
+  num_steps = 10
+  dt = 0.001
+  l_tol = 1e-15
+[]
+
+[Output]
+  output_initial = true
+  interval = 1
+  exodus = true
+[]
diff --git a/test/tests/executioners/time_integration/rk2-2d-quadratic.i b/test/tests/executioners/time_integration/rk2-2d-quadratic.i
@@ -86,6 +86,7 @@
   start_time = 0.0
   num_steps = 10
   dt = 0.0001
+  l_tol = 1e-8
 []
 
 [Output]

diff --git a/test/tests/executioners/time_integration/tests b/test/tests/executioners/time_integration/tests
@@ -106,6 +106,12 @@
   [../]
 
   # RK-2
+  [./rk2-1d-linear]
+    type = 'Exodiff'
+    input = 'rk2-1d-linear.i'
+    exodiff = 'rk2-1d-linear_out.e'
+  [../]
+
   [./rk2-2d-quadratic]
     type = 'Exodiff'
     input = 'rk2-2d-quadratic.i'