www.mooseframework.org
FluxLimitedTVDAdvection.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://www.mooseframework.org
3 //*
4 //* All rights reserved, see COPYRIGHT for full restrictions
5 //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 //*
7 //* Licensed under LGPL 2.1, please see LICENSE for details
8 //* https://www.gnu.org/licenses/lgpl-2.1.html
9 
11 #include "Assembly.h"
12 
14 
15 template <>
16 InputParameters
18 {
19  InputParameters params = validParams<Kernel>();
20  params.addClassDescription("Conservative form of $\\nabla \\cdot \\vec{v} u$ (advection), using "
21  "the Flux Limited TVD scheme invented by Kuzmin and Turek");
22  params.addRequiredParam<UserObjectName>("advective_flux_calculator",
23  "AdvectiveFluxCalculator UserObject");
24  return params;
25 }
26 
27 FluxLimitedTVDAdvection::FluxLimitedTVDAdvection(const InputParameters & parameters)
28  : Kernel(parameters),
29  _fluo(getUserObject<AdvectiveFluxCalculatorBase>("advective_flux_calculator"))
30 {
31 }
32 
33 Real
35 {
36  mooseError("FluxLimitedTVDAdvection::computeQpResidual() called\n");
37  return 0.0;
38 }
39 
40 void
42 {
43  prepareVectorTag(_assembly, _var.number());
44  precalculateResidual();
45 
46  // get the residual contributions from _fluo
47  for (unsigned i = 0; i < _current_elem->n_nodes(); ++i)
48  {
49  const dof_id_type node_id_i = _current_elem->node_id(i);
50  _local_re(i) = _fluo.getFluxOut(node_id_i) / _fluo.getValence(node_id_i);
51  }
52 
53  accumulateTaggedLocalResidual();
54 
55  if (_has_save_in)
56  {
57  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
58  for (unsigned int i = 0; i < _save_in.size(); i++)
59  _save_in[i]->sys().solution().add_vector(_local_re, _save_in[i]->dofIndices());
60  }
61 }
62 
63 void
65 {
66  prepareMatrixTag(_assembly, _var.number(), _var.number());
67  precalculateJacobian();
68 
69  // Run through the nodes of this element using "i", getting the Jacobian contributions
70  // d(residual_i)/du(node_j) for all nodes j that can have a nonzero Jacobian contribution. Some
71  // of these node_j will live in this element, but some will live in other elements connected with
72  // node "i", and some will live in the next layer of nodes (eg, in 1D residual_3 could have
73  // contributions from node1, node2, node3, node4 and node5).
74  for (unsigned i = 0; i < _current_elem->n_nodes(); ++i)
75  {
76  // global id of node "i"
77  const dof_id_type node_id_i = _current_elem->node_id(i);
78  // dof number of _var on node "i"
79  std::vector<dof_id_type> idof_indices(
80  1, _current_elem->node_ref(i).dof_number(_sys.number(), _var.number(), 0));
81  // number of times node "i" is encountered in a sweep over elements
82  const unsigned valence = _fluo.getValence(node_id_i);
83 
84  // retrieve the derivative information from _fluo
85  const std::map<dof_id_type, Real> derivs = _fluo.getdFluxOutdu(node_id_i);
86 
87  // now build up the dof numbers of all the "j" nodes and the derivative matrix
88  // d(residual_i)/d(u_j)
89  std::vector<dof_id_type> jdof_indices(derivs.size());
90  DenseMatrix<Number> deriv_matrix(1, derivs.size());
91  unsigned j = 0;
92  for (const auto & node_j_deriv : derivs)
93  {
94  // global id of j:
95  const dof_id_type node_id_j = node_j_deriv.first;
96  // dof at node j:
97  jdof_indices[j] =
98  _mesh.getMesh().node_ref(node_id_j).dof_number(_sys.number(), _var.number(), 0);
99  // derivative must be divided by valence, otherwise the loop over elements will multiple-count
100  deriv_matrix(0, j) = node_j_deriv.second / valence;
101  j++;
102  }
103  // Add the result to the system's Jacobian matrix
104  _assembly.cacheJacobianBlock(deriv_matrix, idof_indices, jdof_indices, _var.scalingFactor());
105  }
106 }
unsigned getValence(dof_id_type node_i) const
Returns the valence of the global node i Valence is the number of times the node is encountered in a ...
Real getFluxOut(dof_id_type node_i) const
Returns the flux out of lobal node id.
FluxLimitedTVDAdvection(const InputParameters &parameters)
virtual Real computeQpResidual() override
const AdvectiveFluxCalculatorBase & _fluo
The user object that computes Kuzmin and Turek&#39;s K_ij, R+ and R-, etc quantities. ...
virtual void computeResidual() override
Advection of the variable with velocity set in the AdvectiveFluxCalculator.
Base class to compute Advective fluxes.
const std::map< dof_id_type, Real > & getdFluxOutdu(dof_id_type node_i) const
Returns r where r[j] = d(flux out of global node i)/du(global node j) used in Jacobian computations...
registerMooseObject("PorousFlowApp", FluxLimitedTVDAdvection)
virtual void computeJacobian() override
InputParameters validParams< FluxLimitedTVDAdvection >()