www.mooseframework.org
FDKernel.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 "FDKernel.h"
16 #include "Assembly.h"
17 #include "MooseVariable.h"
18 #include "Problem.h"
19 #include "SubProblem.h"
20 #include "SystemBase.h"
21 
22 #include "libmesh/threads.h"
23 #include "libmesh/quadrature.h"
24 
25 template <>
28 {
30  return params;
31 }
32 
33 FDKernel::FDKernel(const InputParameters & parameters) : Kernel(parameters)
34 {
35  _scale = 1.490116119384766e-08; // HACK: sqrt of the machine epsilon for double precision
36 #ifdef LIBMESH_HAVE_PETSC
37  _scale = PETSC_SQRT_MACHINE_EPSILON;
38 #endif
39 }
40 
41 DenseVector<Number>
42 FDKernel::perturbedResidual(unsigned int varnum,
43  unsigned int perturbationj,
44  Real perturbation_scale,
45  Real & perturbation)
46 {
47  DenseVector<Number> re;
48  re.resize(_var.dofIndices().size());
49  re.zero();
50 
51  MooseVariable & var = _sys.getVariable(_tid, varnum);
52  var.computePerturbedElemValues(perturbationj, perturbation_scale, perturbation);
54  for (_i = 0; _i < _test.size(); _i++)
55  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
56  re(_i) += _JxW[_qp] * _coord[_qp] * computeQpResidual();
58  return re;
59 }
60 
61 void
63 {
65 }
66 
67 void
68 FDKernel::computeOffDiagJacobian(unsigned int jvar_index)
69 {
70 
71  DenseMatrix<Number> & ke = _assembly.jacobianBlock(_var.number(), jvar_index);
72  DenseMatrix<Number> local_ke;
73  local_ke.resize(ke.m(), ke.n());
74  local_ke.zero();
75 
76  // FIXME: pull out the already computed element residual instead of recomputing it
77  Real h;
78  DenseVector<Number> re = perturbedResidual(_var.number(), 0, 0.0, h);
79  for (_j = 0; _j < _phi.size(); _j++)
80  {
81  DenseVector<Number> p_re = perturbedResidual(jvar_index, _j, _scale, h);
82  for (_i = 0; _i < _test.size(); _i++)
83  local_ke(_i, _j) = (p_re(_i) - re(_i)) / h;
84  }
85  ke += local_ke;
86 
87  if (jvar_index == _var.number())
88  {
89  _local_ke = local_ke;
91  {
92  unsigned int rows = ke.m();
93  DenseVector<Number> diag(rows);
94  for (unsigned int i = 0; i < rows; i++)
95  diag(i) = _local_ke(i, i);
96 
97  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
98  for (const auto & var : _diag_save_in)
99  var->sys().solution().add_vector(diag, var->dofIndices());
100  }
101  }
102 }
103 
104 void
106 {
107  // FIXME: implement me.
108  /*
109  DenseMatrix<Number> & ke = _assembly.jacobianBlock(_var.number(), jvar);
110  MooseVariableScalar & jv = _sys.getScalarVariable(_tid, jvar);
111 
112  for (_i = 0; _i < _test.size(); _i++)
113  for (_j = 0; _j < jv.order(); _j++)
114  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
115  ke(_i, _j) += _JxW[_qp] * _coord[_qp] * computeQpOffDiagJacobian(jvar);
116  */
117 }
void restoreUnperturbedElemValues()
Restore the values the variable had before a call to computePerturbedElemValues().
virtual Real computeQpResidual()=0
Compute this Kernel&#39;s contribution to the residual at the current quadrature point.
std::vector< MooseVariable * > _diag_save_in
Definition: KernelBase.h:176
virtual DenseVector< Number > perturbedResidual(unsigned int ivar, unsigned int i, Real perturbation_scale, Real &perturbation)
Computes the residual when the current state of j-th variable at element node i is perturbed by pertu...
Definition: FDKernel.C:42
Class for stuff related to variables.
Definition: MooseVariable.h:43
const MooseArray< Real > & _JxW
The current quadrature point weight value.
Definition: KernelBase.h:140
DenseMatrix< Number > _local_ke
Holds residual entries as they are accumulated by this Kernel.
Definition: KernelBase.h:167
const VariableTestValue & _test
the current test function
Definition: KernelBase.h:152
SystemBase & _sys
Reference to the EquationSystem object.
Definition: KernelBase.h:111
THREAD_ID _tid
The thread ID for this kernel.
Definition: KernelBase.h:114
QBase *& _qrule
active quadrature rule
Definition: KernelBase.h:137
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
Real _scale
Definition: FDKernel.h:52
void computePerturbedElemValues(unsigned i, Real scale, Real &h)
Compute values at interior quadrature points when this variable&#39;s elem dof i is perturbed by h in the...
bool _has_diag_save_in
The aux variables to save the diagonal Jacobian contributions to.
Definition: KernelBase.h:175
Assembly & _assembly
Reference to this Kernel&#39;s assembly object.
Definition: KernelBase.h:117
virtual MooseVariable & getVariable(THREAD_ID tid, const std::string &var_name)
Gets a reference to a variable of with specified name.
Definition: SystemBase.C:103
unsigned int _i
current index for the test function
Definition: KernelBase.h:146
InputParameters validParams< FDKernel >()
Definition: FDKernel.C:27
virtual void computeOffDiagJacobianScalar(unsigned int jvar)
Computes jacobian block with respect to a scalar variable.
Definition: FDKernel.C:105
const MooseArray< Real > & _coord
The scaling factor to convert from cartesian to another coordinate system (e.g rz, spherical, etc.)
Definition: KernelBase.h:143
const VariablePhiValue & _phi
the current shape functions
Definition: KernelBase.h:158
virtual void computeOffDiagJacobian(unsigned int jvar)
Computes d-residual / d-jvar...
Definition: FDKernel.C:68
unsigned int _j
current index for the shape function
Definition: KernelBase.h:149
InputParameters validParams< Kernel >()
Definition: Kernel.C:30
unsigned int number() const
Get variable number coming from libMesh.
std::vector< dof_id_type > & dofIndices()
unsigned int size() const
The number of elements that can currently be stored in the array.
Definition: MooseArray.h:250
Definition: Kernel.h:25
virtual void computeJacobian()
Compute this Kernel&#39;s contribution to the diagonal Jacobian entries.
Definition: FDKernel.C:62
DenseMatrix< Number > & jacobianBlock(unsigned int ivar, unsigned int jvar)
Definition: Assembly.C:887
virtual void precalculateResidual()
Following methods are used for Kernels that need to perform a per-element calculation.
Definition: Kernel.C:138
MooseVariable & _var
Reference to this Kernel&#39;s MooseVariable object.
Definition: KernelBase.h:120
FDKernel(const InputParameters &parameters)
Definition: FDKernel.C:33
unsigned int _qp
The current quadrature point index.
Definition: KernelBase.h:131