24#ifndef DUMUX_FVTRANSPORT2P2C_MULTIPHYSICS_HH
25#define DUMUX_FVTRANSPORT2P2C_MULTIPHYSICS_HH
54template<
class TypeTag>
69 dim = GridView::dimension
73 wPhaseIdx = Indices::wPhaseIdx, nPhaseIdx = Indices::nPhaseIdx,
74 wCompIdx = Indices::wPhaseIdx, nCompIdx = Indices::nPhaseIdx
77 using PhaseVector = Dune::FieldVector<Scalar, 2>;
87 virtual void update(
const Scalar t, Scalar& dt, TransportSolutionType& updateVec,
bool impet =
false);
118template<
class TypeTag>
124 unsigned int size = problem().gridView().size(0);
125 if (this->enableLocalTimeStepping())
127 if (this->timeStepData_.size() != size)
128 this->timeStepData_.resize(size);
131 this->impet_ = impet;
132 this->averagedFaces_ = 0.;
135 updateVec.resize(getPropValue<TypeTag, Properties::NumComponents>());
136 updateVec[wCompIdx].resize(problem().gridView().size(0));
137 updateVec[nCompIdx].resize(problem().gridView().size(0));
138 updateVec[wCompIdx] = 0;
139 updateVec[nCompIdx] = 0;
142 int restrictingCell = -1;
144 PhaseVector entries(0.), timestepFlux(0.);
146 for (
const auto& element : elements(problem().gridView()))
149 int globalIdxI = problem().variables().index(element);
150 CellData& cellDataI = problem().variables().cellData(globalIdxI);
152 if (impet || cellDataI.subdomain() == 2)
155 double sumfactorin = 0;
156 double sumfactorout = 0;
159 for (
const auto& intersection : intersections(problem().gridView(), element))
161 int indexInInside = intersection.indexInInside();
164 if (intersection.neighbor())
165 this->getFlux(entries, timestepFlux, intersection, cellDataI);
168 if (intersection.boundary())
169 this->getFluxOnBoundary(entries, timestepFlux, intersection, cellDataI);
171 if (this->enableLocalTimeStepping())
173 LocalTimesteppingData& localData = this->timeStepData_[globalIdxI];
175 if (localData.faceTargetDt[indexInInside] < this->accumulatedDt_ + this->dtThreshold_)
177 localData.faceFluxes[indexInInside] = entries;
183 updateVec[wCompIdx][globalIdxI] += entries[wCompIdx];
184 updateVec[nCompIdx][globalIdxI] += entries[nCompIdx];
188 sumfactorin += timestepFlux[0];
189 sumfactorout += timestepFlux[1];
193 if (this->enableLocalTimeStepping())
195 LocalTimesteppingData& localData = this->timeStepData_[globalIdxI];
196 for (
int i=0; i < 2*dim; i++)
198 updateVec[wCompIdx][globalIdxI] += localData.faceFluxes[i][wCompIdx];
199 updateVec[nCompIdx][globalIdxI] += localData.faceFluxes[i][nCompIdx];
204 PrimaryVariables q(NAN);
205 problem().source(q, element);
206 updateVec[wCompIdx][globalIdxI] += q[Indices::contiWEqIdx];
207 updateVec[nCompIdx][globalIdxI] += q[Indices::contiNEqIdx];
211 sumfactorin = max(sumfactorin,sumfactorout)
212 / problem().spatialParams().porosity(element);
215 if (this->enableLocalTimeStepping())
217 this->timeStepData_[globalIdxI].dt = 1./sumfactorin;
218 if ( 1./sumfactorin < dt)
221 restrictingCell= globalIdxI;
226 if ( 1./sumfactorin < dt)
229 restrictingCell= globalIdxI;
238 using ElementMapper =
typename SolutionTypes::ElementMapper;
240 for (
int i = 0; i < updateVec.size(); i++)
242 DataHandle dataHandle(problem().variables().elementMapper(), updateVec[i]);
243 problem().gridView().template communicate<DataHandle>(dataHandle,
244 Dune::InteriorBorder_All_Interface,
245 Dune::ForwardCommunication);
247 dt = problem().gridView().comm().min(dt);
252 Dune::dinfo <<
"Timestep restricted by CellIdx " << restrictingCell <<
253 " leads to dt = "<<dt * getParam<Scalar>(
"Impet.CFLFactor")<< std::endl;
254 if(this->averagedFaces_ != 0)
255 Dune::dinfo <<
" Averageing done for " << this->averagedFaces_ <<
" faces. "<< std::endl;
Contains a class to exchange entries of a vector.
Finite volume discretization of the component transport equation.
typename Properties::Detail::GetPropImpl< TypeTag, Property >::type GetProp
get the type of a property (equivalent to old macro GET_PROP(...))
Definition: propertysystem.hh:140
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:149
A data handle class to exchange entries of a vector.
Definition: vectorcommdatahandle.hh:78
Compositional transport step in a Finite Volume discretization.
Definition: fvtransport.hh:62
Problem & problem_
Definition: fvtransport.hh:302
Definition: fvtransport.hh:114
Compositional transport step in a finite volume discretization.
Definition: fvtransportmultiphysics.hh:56
virtual ~FVTransport2P2CMultiPhysics()
Definition: fvtransportmultiphysics.hh:96
FVTransport2P2CMultiPhysics(Problem &problem)
Constructs a FVTransport2P2CMultiPhysics object.
Definition: fvtransportmultiphysics.hh:93
virtual void update(const Scalar t, Scalar &dt, TransportSolutionType &updateVec, bool impet=false)
Calculate the update vector and determine timestep size.
Definition: fvtransportmultiphysics.hh:119