36#ifndef OPENRS_BOUNDARYCONDITIONS_HEADER
37#define OPENRS_BOUNDARYCONDITIONS_HEADER
43#include <opm/common/ErrorMacros.hpp>
44#include <dune/common/fvector.hh>
58 enum BCType { Dirichlet, Neumann, Periodic };
64 template<
typename charT,
class traits>
65 void write(std::basic_ostream<charT,traits>& os)
const
67 os <<
"Type: " << type_ <<
" Value: " << value_;
73 : type_(Neumann), value_(0.0)
80 : type_(type), value_(value)
87 return type_ == Dirichlet;
93 return type_ == Neumann;
99 return type_ == Periodic;
108 template<
typename charT,
class traits,
typename T>
109 std::basic_ostream<charT,traits>&
110 operator<<(std::basic_ostream<charT,traits>& os,
126 :
BCBase<double>(Neumann, 0.0)
133 :
BCBase<double>(type, value)
148 OPM_ERROR_IF(!
isDirichlet(),
"Pressure boundary conditions are only valid for Dirichlet boundaries");
157 OPM_ERROR_IF(!
isNeumann(),
"Outflux boundary conditions are only valid for Neumann boundaries");
166 OPM_ERROR_IF(!
isPeriodic(),
"Pressure difference boundary conditions are only valid for periodic boundaries");
180 :
BCBase<double>(Dirichlet, 1.0)
187 :
BCBase<double>(type, value)
213 template <
int numComponents>
217 typedef Dune::FieldVector<double, numComponents> CompVec;
221 :
Base(
Base::Dirichlet, CompVec(0.0))
226 explicit SurfvolBC(Dune::FieldVector<double, numComponents> value)
252 : periodic_partner_bid_(num_different_boundary_ids, 0),
253 canonical_bid_(num_different_boundary_ids, 0)
257 void resize(
int new_size)
259 periodic_partner_bid_.resize(new_size, 0);
260 canonical_bid_.resize(new_size, 0);
265 return periodic_partner_bid_.empty() && canonical_bid_.empty();
270 periodic_partner_bid_.clear();
271 canonical_bid_.clear();
276 return periodic_partner_bid_.size();
279 void setPeriodicPartners(
int boundary_id_1,
int boundary_id_2)
281 assert(boundary_id_1 >= 0 && boundary_id_1 <
int(periodic_partner_bid_.size()));
282 assert(boundary_id_2 >= 0 && boundary_id_2 <
int(periodic_partner_bid_.size()));
283 assert(periodic_partner_bid_[boundary_id_1] == 0);
284 assert(periodic_partner_bid_[boundary_id_2] == 0);
285 periodic_partner_bid_[boundary_id_1] = boundary_id_2;
286 periodic_partner_bid_[boundary_id_2] = boundary_id_1;
289 int getPeriodicPartner(
int boundary_id)
const
291 assert(boundary_id >= 0 && boundary_id <
int(periodic_partner_bid_.size()));
292 return periodic_partner_bid_[boundary_id];
295 void setCanonicalBoundaryId(
int boundary_id,
int canonical_bid)
297 assert(boundary_id >= 0 && boundary_id <
int(canonical_bid_.size()));
298 canonical_bid_[boundary_id] = canonical_bid;
301 int getCanonicalBoundaryId(
int boundary_id)
const
303 assert(boundary_id >= 0 && boundary_id <
int(canonical_bid_.size()));
304 return canonical_bid_[boundary_id];
307 template<
typename charT,
class traits>
308 void write(std::basic_ostream<charT,traits>& os)
const
310 for (
int i = 0; i < int(periodic_partner_bid_.size()); ++i) {
311 os <<
"Partner of bid " << i <<
" is " << periodic_partner_bid_[i]
312 <<
" (canonical bid is " << canonical_bid_[i] <<
'\n';
317 std::vector<int> periodic_partner_bid_;
318 std::vector<int> canonical_bid_;
321 template<
typename charT,
class traits>
322 std::basic_ostream<charT,traits>&
323 operator<<(std::basic_ostream<charT,traits>& os,
331 template <
typename T>
341 template <
bool FC = false,
bool SC = false,
bool ZC = false,
int numComponents = 3>
343 private std::conditional<FC, std::vector<FlowBC>, DummyVec<FlowBC> >::type,
344 private std::conditional<SC, std::vector<SatBC>, DummyVec<SatBC> >::type,
345 private std::conditional<ZC, std::vector<SurfvolBC<numComponents> >, DummyVec<SurfvolBC<numComponents> > >::type
348 typedef typename std::conditional<FC, std::vector<FlowBC>,
DummyVec<FlowBC> >::type FlowConds;
349 typedef typename std::conditional<SC, std::vector<SatBC>,
DummyVec<SatBC> >::type SatConds;
351 const static bool HasFlowConds = FC;
352 const static bool HasSatConds = SC;
353 const static bool HasSurfvolConds = SC;
361 FlowConds(num_different_boundary_ids),
362 SatConds(num_different_boundary_ids),
363 SurfvolConds(num_different_boundary_ids)
367 void resize(
int new_size)
369 PeriodicConditionHandler::resize(new_size);
370 FlowConds::resize(new_size);
371 SatConds::resize(new_size);
372 SurfvolConds::resize(new_size);
377 return PeriodicConditionHandler::empty();
382 PeriodicConditionHandler::clear();
385 SurfvolConds::clear();
390 return PeriodicConditionHandler::size();
395 return FlowConds::operator[](i);
398 const FlowBC& flowCond(
int i)
const
400 return FlowConds::operator[](i);
403 template <
class BoundaryFace>
404 const FlowBC& flowCond(
const BoundaryFace& bf)
const
406 assert(bf.boundary());
407 return FlowConds::operator[](bf.boundaryId());
410 SatBC& satCond(
int i)
412 return SatConds::operator[](i);
415 const SatBC& satCond(
int i)
const
417 return SatConds::operator[](i);
420 template <
class BoundaryFace>
421 const SatBC& satCond(
const BoundaryFace& bf)
const
423 assert(bf.boundary());
424 return SatConds::operator[](bf.boundaryId());
429 return SurfvolConds::operator[](i);
434 return SurfvolConds::operator[](i);
437 template <
class BoundaryFace>
440 assert(bf.boundary());
441 return SurfvolConds::operator[](bf.boundaryId());
444 template<
typename charT,
class traits>
445 void write(std::basic_ostream<charT,traits>& os)
const
447 PeriodicConditionHandler::write(os);
451 template<
typename charT,
class traits,
bool F,
bool S>
452 std::basic_ostream<charT,traits>&
453 operator<<(std::basic_ostream<charT,traits>& os,
A class for building boundary conditions in a uniform way.
Definition BoundaryConditions.hpp:52
BCType
Enum for the allowed boundary condition types.
Definition BoundaryConditions.hpp:58
bool isPeriodic() const
Type query.
Definition BoundaryConditions.hpp:97
BCBase()
Default constructor, that makes a Neumann condition with value 0.0.
Definition BoundaryConditions.hpp:72
BCBase(BCType type, T value)
Constructor taking a type and value.
Definition BoundaryConditions.hpp:79
bool isDirichlet() const
Type query.
Definition BoundaryConditions.hpp:85
void write(std::basic_ostream< charT, traits > &os) const
Write type and value to an output stream.
Definition BoundaryConditions.hpp:65
bool isNeumann() const
Type query.
Definition BoundaryConditions.hpp:91
Definition BoundaryConditions.hpp:346
Definition BoundaryConditions.hpp:333
A class for representing a flow boundary condition.
Definition BoundaryConditions.hpp:122
double pressure() const
Query a Dirichlet condition.
Definition BoundaryConditions.hpp:145
double outflux() const
Query a Neumann condition.
Definition BoundaryConditions.hpp:154
FlowBC(BCType type, double value)
Constructor taking a type and value.
Definition BoundaryConditions.hpp:132
double pressureDifference() const
Query a Periodic condition.
Definition BoundaryConditions.hpp:163
FlowBC()
Default constructor, that makes a noflow condition (Neumann, value 0.0).
Definition BoundaryConditions.hpp:125
Definition BoundaryConditions.hpp:245
A class for representing a saturation boundary condition.
Definition BoundaryConditions.hpp:176
double saturation() const
Query a Dirichlet condition.
Definition BoundaryConditions.hpp:197
double saturationDifference() const
Query a Periodic condition.
Definition BoundaryConditions.hpp:204
SatBC()
Default constructor, that makes a Dirichlet condition with value 1.0.
Definition BoundaryConditions.hpp:179
SatBC(BCType type, double value)
Constructor taking a type and value.
Definition BoundaryConditions.hpp:186
A class for representing a surface volume boundary condition.
Definition BoundaryConditions.hpp:215
CompVec surfvol() const
Query a Dirichlet condition.
Definition BoundaryConditions.hpp:236
SurfvolBC()
Default constructor, that makes a Dirichlet condition with all zero component values.
Definition BoundaryConditions.hpp:220
bool isDirichlet() const
Forwarding the relevant type query.
Definition BoundaryConditions.hpp:85
SurfvolBC(Dune::FieldVector< double, numComponents > value)
Constructor taking a value.
Definition BoundaryConditions.hpp:226
Inverting small matrices.
Definition ImplicitAssembly.hpp:43
std::basic_ostream< charT, traits > & operator<<(std::basic_ostream< charT, traits > &os, const BCBase< T > &bc)
Stream insertion for BCBase.
Definition BoundaryConditions.hpp:110