www.mooseframework.org
AdamsPredictor.C
Go to the documentation of this file.
1 /****************************************************************/
2 /* DO NOT MODIFY THIS HEADER */
3 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
4 /* */
5 /* (c) 2010 Battelle Energy Alliance, LLC */
6 /* ALL RIGHTS RESERVED */
7 /* */
8 /* Prepared by Battelle Energy Alliance, LLC */
9 /* Under Contract No. DE-AC07-05ID14517 */
10 /* With the U. S. Department of Energy */
11 /* */
12 /* See COPYRIGHT for full restrictions */
13 /****************************************************************/
14 
15 #include "AdamsPredictor.h"
16 #include "NonlinearSystem.h"
17 
18 #include "libmesh/numeric_vector.h"
19 
20 template <>
23 {
25  params.addParam<int>("order", 2, "The maximum reachable order of the Adams-Bashforth Predictor");
26  return params;
27 }
28 
30  : Predictor(parameters),
31  _order(getParam<int>("order")),
32  _current_old_solution(_nl.addVector("AB2_current_old_solution", true, GHOSTED)),
33  _older_solution(_nl.addVector("AB2_older_solution", true, GHOSTED)),
34  _oldest_solution(_nl.addVector("AB2_rejected_solution", true, GHOSTED)),
35  _tmp_previous_solution(_nl.addVector("tmp_previous_solution", true, GHOSTED)),
36  _tmp_residual_old(_nl.addVector("tmp_residual_old", true, GHOSTED)),
37  _tmp_third_vector(_nl.addVector("tmp_third_vector", true, GHOSTED)),
38  _t_step_old(declareRestartableData<int>("t_step_old", 0)),
39  _dt_older(declareRestartableData<Real>("dt_older", 0)),
40  _dtstorage(declareRestartableData<Real>("dtstorage", 0))
41 {
42 }
43 
44 void
46 {
47  // if the time step number hasn't changed then do nothing
48  if (_t_step == _t_step_old)
49  return;
50 
51  // Otherwise move back the previous old solution and copy the current old solution,
52  // This will probably not work with DT2, but I don't need to get it to work with dt2.
54 
56  // Set older solution to hold the previous old solution
58  // Set current old solution to hold what it says.
60  // Same thing for dt
63 }
64 
65 bool
67 {
69  return false;
70 
71  // AB2 can only be applied if there are enough old solutions
72  // AB1 could potentially be used for the time step prior?
73  // It would be possible to do VSVO Adams, Kevin has the info
74  // Doing so requires a time stack of some sort....
75  if (_dt == 0 || _dt_old == 0 || _dt_older == 0 || _t_step < 2)
76  return false;
77  else
78  return true;
79 }
80 
81 void
82 AdamsPredictor::apply(NumericVector<Number> & sln)
83 {
84  // localize current solution to working vec
85  sln.localize(_solution_predictor);
86  // NumericVector<Number> & vector1 = _tmp_previous_solution;
87  NumericVector<Number> & vector2 = _tmp_residual_old;
88  NumericVector<Number> & vector3 = _tmp_third_vector;
89 
90  Real commonpart = _dt / _dt_old;
91  Real firstpart = (1 + .5 * commonpart);
92  Real secondpart = .5 * _dt / _dt_older;
93 
94  _older_solution.localize(vector2);
95  _oldest_solution.localize(vector3);
96 
97  _solution_predictor *= 1 + commonpart * firstpart;
98  vector2 *= -1. * commonpart * (firstpart + secondpart);
99  vector3 *= commonpart * secondpart;
100 
101  _solution_predictor += vector2;
102  _solution_predictor += vector3;
103 
104  _solution_predictor.localize(sln);
105 }
AdamsPredictor(const InputParameters &parameters)
NumericVector< Number > & _tmp_third_vector
virtual void apply(NumericVector< Number > &sln) override
InputParameters validParams< AdamsPredictor >()
NumericVector< Number > & _current_old_solution
Base class for predictors.
Definition: Predictor.h:39
virtual NumericVector< Number > & solutionOld()=0
NonlinearSystemBase & _nl
Definition: Predictor.h:54
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
Real & _dt_old
Definition: Predictor.h:58
NumericVector< Number > & _solution_predictor
Definition: Predictor.h:62
InputParameters validParams< Predictor >()
Definition: Predictor.C:24
Real & _dt
Definition: Predictor.h:57
virtual void timestepSetup() override
NumericVector< Number > & _oldest_solution
void addParam(const std::string &name, const S &value, const std::string &doc_string)
These methods add an option parameter and a documentation string to the InputParameters object...
NumericVector< Number > & _tmp_residual_old
int & _t_step
Definition: Predictor.h:56
NumericVector< Number > & _older_solution
virtual bool shouldApply()
Definition: Predictor.C:69
virtual bool shouldApply() override