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/material/common/Valgrind.hpp>
36#include <dune/istl/bvector.hh>
37#include <dune/istl/matrix.hh>
39#include <dune/common/fvector.hh>
40#include <dune/common/fmatrix.hh>
44template<
class TypeTag>
45class FvBaseAdLocalLinearizer;
48namespace Opm::Properties {
57template<
class TypeTag>
62template<
class TypeTag>
71 using type = DenseAd::Evaluation<Scalar, numEq>;
86template<
class TypeTag>
99 using Element =
typename GridView::template Codim<0>::Entity;
103 using ScalarVectorBlock = Dune::FieldVector<Scalar, numEq>;
107 using ScalarLocalBlockVector = Dune::BlockVector<ScalarVectorBlock>;
108 using ScalarLocalBlockMatrix = Dune::Matrix<ScalarMatrixBlock>;
112 : internalElemContext_(0)
120 {
delete internalElemContext_; }
136 void init(Simulator& simulator)
138 simulatorPtr_ = &simulator;
139 delete internalElemContext_;
140 internalElemContext_ =
new ElementContext(simulator);
156 linearize(*internalElemContext_, element);
176 elemCtx.updateAllIntensiveQuantities();
179 model_().updatePVWeights(
elemCtx);
185 unsigned numPrimaryDof =
elemCtx.numPrimaryDof(0);
188 elemCtx.updateAllExtensiveQuantities();
204 {
return localResidual_; }
210 {
return localResidual_; }
229 {
return residual_[
dofIdx]; }
232 Implementation& asImp_()
233 {
return *
static_cast<Implementation*
>(
this); }
234 const Implementation& asImp_()
const
235 {
return *
static_cast<const Implementation*
>(
this); }
237 const Simulator& simulator_()
const
238 {
return *simulatorPtr_; }
239 const Problem& problem_()
const
240 {
return simulatorPtr_->problem(); }
241 const Model& model_()
const
242 {
return simulatorPtr_->model(); }
250 size_t numDof =
elemCtx.numDof(0);
251 size_t numPrimaryDof =
elemCtx.numPrimaryDof(0);
253 residual_.resize(numDof);
254 if (jacobian_.N() != numDof || jacobian_.M() != numPrimaryDof)
255 jacobian_.setSize(numDof, numPrimaryDof);
274 const auto&
resid = localResidual_.residual();
279 size_t numDof =
elemCtx.numDof(0);
297 ElementContext *internalElemContext_;
299 LocalResidual localResidual_;
301 ScalarLocalBlockVector residual_;
302 ScalarLocalBlockMatrix jacobian_;
Calculates the local residual and its Jacobian for a single element of the grid.
Definition fvbaseadlocallinearizer.hh:88
const ScalarVectorBlock & residual(unsigned dofIdx) const
Returns the local residual of a sub-control volume.
Definition fvbaseadlocallinearizer.hh:228
void linearize(const Element &element)
Compute an element's local Jacobian matrix and evaluate its residual.
Definition fvbaseadlocallinearizer.hh:154
LocalResidual & localResidual()
Return reference to the local residual.
Definition fvbaseadlocallinearizer.hh:203
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:271
void reset_(const ElementContext &)
Reset the all relevant internal attributes to 0.
Definition fvbaseadlocallinearizer.hh:261
void resize_(const ElementContext &elemCtx)
Resize all internal attributes to the size of the element.
Definition fvbaseadlocallinearizer.hh:248
const LocalResidual & localResidual() const
Return reference to the local residual.
Definition fvbaseadlocallinearizer.hh:209
void init(Simulator &simulator)
Initialize the local Jacobian object.
Definition fvbaseadlocallinearizer.hh:136
static void registerParameters()
Register all run-time parameters for the local jacobian.
Definition fvbaseadlocallinearizer.hh:125
void linearize(ElementContext &elemCtx, const Element &elem)
Compute an element's local Jacobian matrix and evaluate its residual.
Definition fvbaseadlocallinearizer.hh:173
const ScalarMatrixBlock & jacobian(unsigned domainScvIdx, unsigned rangeScvIdx) const
Returns the local Jacobian matrix of the residual of a sub-control volume.
Definition fvbaseadlocallinearizer.hh:220
Manages the initializing and running of time dependent problems.
Definition simulator.hh:92
Declare the properties used by the infrastructure code of the finite volume discretizations.
This file contains a set of helper functions used by VFPProd / VFPInj.
Definition blackoilboundaryratevector.hh:37
constexpr auto getPropValue()
get the value data member of a property
Definition propertysystem.hh:242
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type::type GetPropType
get the type alias defined in the property (equivalent to old macro GET_PROP_TYPE(....
Definition propertysystem.hh:235
Representation of a function evaluation and all necessary derivatives with regard to the intensive qu...
Definition fvbaseproperties.hh:66
The type of the local linearizer.
Definition fvbaseproperties.hh:97
Definition fvbaseadlocallinearizer.hh:53