28 #ifndef EWOMS_FV_BASE_AD_LOCAL_LINEARIZER_HH
29 #define EWOMS_FV_BASE_AD_LOCAL_LINEARIZER_HH
33 #include <opm/material/densead/Math.hpp>
34 #include <opm/common/Valgrind.hpp>
35 #include <opm/common/Unused.hpp>
37 #include <dune/istl/bvector.hh>
38 #include <dune/istl/matrix.hh>
40 #include <dune/common/fvector.hh>
41 #include <dune/common/fmatrix.hh>
45 template<
class TypeTag>
48 namespace Properties {
70 SET_PROP(AutoDiffLocalLinearizer, Evaluation)
78 typedef Opm::DenseAd::Evaluation<Scalar, numEq> type;
91 template<class TypeTag>
95 typedef typename GET_PROP_TYPE(TypeTag, LocalLinearizer) Implementation;
96 typedef typename
GET_PROP_TYPE(TypeTag, LocalResidual) LocalResidual;
100 typedef typename
GET_PROP_TYPE(TypeTag, PrimaryVariables) PrimaryVariables;
101 typedef typename
GET_PROP_TYPE(TypeTag, ElementContext) ElementContext;
104 typedef typename GridView::template Codim<0>::Entity Element;
108 typedef Dune::FieldVector<Scalar, numEq> ScalarVectorBlock;
109 typedef Dune::FieldMatrix<Scalar, numEq, numEq> ScalarMatrixBlock;
111 typedef Dune::BlockVector<ScalarVectorBlock> ScalarLocalBlockVector;
112 typedef Dune::Matrix<ScalarMatrixBlock> ScalarLocalBlockMatrix;
115 FvBaseAdLocalLinearizer()
116 : internalElemContext_(0)
121 FvBaseAdLocalLinearizer(
const FvBaseAdLocalLinearizer&) =
delete;
123 ~FvBaseAdLocalLinearizer()
124 {
delete internalElemContext_; }
142 simulatorPtr_ = &simulator;
143 delete internalElemContext_;
144 internalElemContext_ =
new ElementContext(simulator);
160 linearize(*internalElemContext_, element);
177 void linearize(ElementContext& elemCtx,
const Element& elem)
179 elemCtx.updateStencil(elem);
180 elemCtx.updateAllIntensiveQuantities();
183 model_().updatePVWeights(elemCtx);
190 unsigned numPrimaryDof = elemCtx.numPrimaryDof(0);
191 for (
unsigned focusDofIdx = 0; focusDofIdx < numPrimaryDof; focusDofIdx++) {
192 elemCtx.setFocusDofIndex(focusDofIdx);
193 elemCtx.updateAllExtensiveQuantities();
196 localResidual_.eval(elemCtx);
201 updateLocalLinearization_(elemCtx, focusDofIdx);
209 {
return localResidual_; }
215 {
return localResidual_; }
225 const ScalarMatrixBlock&
jacobian(
unsigned domainScvIdx,
unsigned rangeScvIdx)
const
226 {
return jacobian_[domainScvIdx][rangeScvIdx]; }
233 const ScalarVectorBlock&
residual(
unsigned dofIdx)
const
234 {
return residual_[dofIdx]; }
237 Implementation& asImp_()
238 {
return *
static_cast<Implementation*
>(
this); }
239 const Implementation& asImp_()
const
240 {
return *
static_cast<const Implementation*
>(
this); }
242 const Simulator& simulator_()
const
243 {
return *simulatorPtr_; }
244 const Problem& problem_()
const
245 {
return simulatorPtr_->problem(); }
246 const Model& model_()
const
247 {
return simulatorPtr_->model(); }
255 size_t numDof = elemCtx.numDof(0);
256 size_t numPrimaryDof = elemCtx.numPrimaryDof(0);
258 residual_.resize(numDof);
259 jacobian_.setSize(numDof, numPrimaryDof);
265 void reset_(
const ElementContext& elemCtx OPM_UNUSED)
276 unsigned focusDofIdx)
278 const auto& resid = localResidual_.residual();
280 for (
unsigned eqIdx = 0; eqIdx < numEq; eqIdx++)
281 residual_[focusDofIdx][eqIdx] = resid[focusDofIdx][eqIdx].value();
283 size_t numDof = elemCtx.numDof(0);
284 for (
unsigned dofIdx = 0; dofIdx < numDof; dofIdx++) {
285 for (
unsigned eqIdx = 0; eqIdx < numEq; eqIdx++) {
286 for (
unsigned pvIdx = 0; pvIdx < numEq; pvIdx++) {
291 jacobian_[dofIdx][focusDofIdx][eqIdx][pvIdx] = resid[dofIdx][eqIdx].derivative(pvIdx);
292 Opm::Valgrind::CheckDefined(jacobian_[dofIdx][focusDofIdx][eqIdx][pvIdx]);
301 ElementContext *internalElemContext_;
303 LocalResidual localResidual_;
305 ScalarLocalBlockVector residual_;
306 ScalarLocalBlockMatrix jacobian_;
void linearize(const Element &element)
Compute an element's local Jacobian matrix and evaluate its residual.
Definition: fvbaseadlocallinearizer.hh:158
const LocalResidual & localResidual() const
Return reference to the local residual.
Definition: fvbaseadlocallinearizer.hh:214
void linearize(ElementContext &elemCtx, const Element &elem)
Compute an element's local Jacobian matrix and evaluate its residual.
Definition: fvbaseadlocallinearizer.hh:177
void reset_(const ElementContext &elemCtx OPM_UNUSED)
Reset the all relevant internal attributes to 0.
Definition: fvbaseadlocallinearizer.hh:265
#define GET_PROP_VALUE(TypeTag, PropTagName)
Access the value attribute of a property for a type tag.
Definition: propertysystem.hh:469
const ScalarVectorBlock & residual(unsigned dofIdx) const
Returns the local residual of a sub-control volume.
Definition: fvbaseadlocallinearizer.hh:233
#define NEW_TYPE_TAG(...)
Define a new type tag.
Definition: propertysystem.hh:169
#define GET_PROP_TYPE(TypeTag, PropTagName)
Access the type attribute of a property for a type tag.
Definition: propertysystem.hh:486
void resize_(const ElementContext &elemCtx)
Resize all internal attributes to the size of the element.
Definition: fvbaseadlocallinearizer.hh:253
Declare the properties used by the infrastructure code of the finite volume discretizations.
const ScalarMatrixBlock & jacobian(unsigned domainScvIdx, unsigned rangeScvIdx) const
Returns the local Jacobian matrix of the residual of a sub-control volume.
Definition: fvbaseadlocallinearizer.hh:225
#define SET_PROP(EffTypeTagName, PropTagName)
Set a property for a specific type tag.
Definition: propertysystem.hh:297
Calculates the local residual and its Jacobian for a single element of the grid.
Definition: fvbaseadlocallinearizer.hh:46
void updateLocalLinearization_(const ElementContext &elemCtx, unsigned focusDofIdx)
Updates the current local Jacobian matrix with the partial derivatives of all equations for the degre...
Definition: fvbaseadlocallinearizer.hh:275
static void registerParameters()
Register all run-time parameters for the local jacobian.
Definition: fvbaseadlocallinearizer.hh:129
void init(Simulator &simulator)
Initialize the local Jacobian object.
Definition: fvbaseadlocallinearizer.hh:140
Manages the initializing and running of time dependent problems.
Definition: simulator.hh:75
#define NEW_PROP_TAG(PTagName)
Define a property tag.
Definition: propertysystem.hh:247
#define SET_TYPE_PROP(EffTypeTagName, PropTagName,...)
Set a property which defines a type.
Definition: propertysystem.hh:377
LocalResidual & localResidual()
Return reference to the local residual.
Definition: fvbaseadlocallinearizer.hh:208