Flow123d  master-f44eb46
assembly_convection.hh
Go to the documentation of this file.
1 /*!
2  *
3  * Copyright (C) 2015 Technical University of Liberec. All rights reserved.
4  *
5  * This program is free software; you can redistribute it and/or modify it under
6  * the terms of the GNU General Public License version 3 as published by the
7  * Free Software Foundation. (http://www.gnu.org/licenses/gpl-3.0.en.html)
8  *
9  * This program is distributed in the hope that it will be useful, but WITHOUT
10  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
11  * FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.
12  *
13  *
14  * @file assembly_convection.hh
15  * @brief
16  */
17 
18 #ifndef ASSEMBLY_CONVECTION_HH_
19 #define ASSEMBLY_CONVECTION_HH_
20 
23 #include "transport/transport.h"
24 #include "fem/fe_p.hh"
25 #include "fem/fe_values.hh"
27 #include "coupling/balance.hh"
29 
30 
31 /**
32  * Auxiliary container class for Finite element and related objects of given dimension.
33  */
34 template <unsigned int dim>
36 {
37 public:
40 
41  static constexpr const char * name() { return "MassAssemblyConvection"; }
42 
43  /// Constructor.
44  MassAssemblyConvection(EqFields *eq_fields, EqData *eq_data)
45  : AssemblyBase<dim>(0), eq_fields_(eq_fields), eq_data_(eq_data) {
49  }
50 
51  /// Destructor.
53 
54  /// Initialize auxiliary vectors and other data members
55  void initialize(ElementCacheMap *element_cache_map) {
56  this->element_cache_map_ = element_cache_map;
57  }
58 
59 
60  /// Assemble integral over element
61  inline void cell_integral(DHCellAccessor cell, unsigned int element_patch_idx)
62  {
63  ASSERT_EQ(cell.dim(), dim).error("Dimension of element mismatch!");
64 
65  ElementAccessor<3> elm = cell.elm();
66  // we have currently zero order P_Disc FE
67  ASSERT(cell.get_loc_dof_indices().size() == 1);
68  IntIdx local_p0_dof = cell.get_loc_dof_indices()[0];
69 
70  auto p = *( this->bulk_points(element_patch_idx).begin() );
71  for (unsigned int sbi=0; sbi<eq_data_->n_substances(); ++sbi)
72  eq_data_->balance_->add_mass_values(eq_data_->subst_idx[sbi], cell, {local_p0_dof},
73  {eq_fields_->cross_section(p)*eq_fields_->water_content(p)*elm.measure()}, 0);
74 
75  VecSetValue(eq_data_->mass_diag, eq_data_->dh_->get_local_to_global_map()[local_p0_dof],
76  eq_fields_->cross_section(p)*eq_fields_->water_content(p), INSERT_VALUES);
77  }
78 
79  /// Implements @p AssemblyBase::begin.
80  void begin() override
81  {
82  VecZeroEntries(eq_data_->mass_diag);
83  eq_data_->balance_->start_mass_assembly(eq_data_->subst_idx);
84  }
85 
86  /// Implements @p AssemblyBase::end.
87  void end() override
88  {
89  eq_data_->balance_->finish_mass_assembly(eq_data_->subst_idx);
90 
91  VecAssemblyBegin(eq_data_->mass_diag);
92  VecAssemblyEnd(eq_data_->mass_diag);
93 
94  eq_data_->is_mass_diag_changed = true;
95  }
96 
97  private:
98  shared_ptr<FiniteElement<dim>> fe_; ///< Finite element for the solution of the advection-diffusion equation.
99 
100  /// Data objects shared with TransportDG
103 
104  /// Sub field set contains fields used in calculation.
106 
107  template < template<IntDim...> class DimAssembly>
108  friend class GenericAssembly;
109 
110 };
111 
112 
113 
114 /**
115  * Auxiliary container class for Finite element and related objects of given dimension.
116  */
117 template <unsigned int dim>
119 {
120 public:
123 
124  static constexpr const char * name() { return "InitCondAssemblyConvection"; }
125 
126  /// Constructor.
128  : AssemblyBase<dim>(0), eq_fields_(eq_fields), eq_data_(eq_data) {
131  }
132 
133  /// Destructor.
135 
136  /// Initialize auxiliary vectors and other data members
137  void initialize(ElementCacheMap *element_cache_map) {
138  this->element_cache_map_ = element_cache_map;
139  for (unsigned int sbi=0; sbi<eq_data_->n_substances(); sbi++) {
140  vecs_.push_back(eq_fields_->conc_mobile_fe[sbi]->vec());
141  }
142  }
143 
144 
145  /// Assemble integral over element
146  inline void cell_integral(DHCellAccessor cell, unsigned int element_patch_idx)
147  {
148  ASSERT_EQ(cell.dim(), dim).error("Dimension of element mismatch!");
149 
150  LongIdx index = cell.local_idx();
151  auto p = *( this->bulk_points(element_patch_idx).begin() );
152 
153  for (unsigned int sbi=0; sbi<eq_data_->n_substances(); sbi++) {
154  vecs_[sbi].set( index, eq_fields_->init_conc[sbi](p) );
155  }
156  }
157 
158 private:
159  shared_ptr<FiniteElement<dim>> fe_; ///< Finite element for the solution of the advection-diffusion equation.
160 
161  /// Data objects shared with TransportDG
164 
165  std::vector<VectorMPI> vecs_; ///< Set of data vectors of conc_mobile_fe objects.
166 
167  /// Sub field set contains fields used in calculation.
169 
170  template < template<IntDim...> class DimAssembly>
171  friend class GenericAssembly;
172 
173 };
174 
175 
176 /**
177  * Auxiliary container class for Finite element and related objects of given dimension.
178  *
179  * Assembles concentration sources for each substance and set boundary conditions.
180  * note: the source of concentration is multiplied by time interval (gives the mass, not the flow like before)
181  */
182 template <unsigned int dim>
184 {
185 public:
188 
189  static constexpr const char * name() { return "ConcSourcesBdrAssemblyConvection"; }
190 
191  /// Constructor.
193  : AssemblyBase<dim>(0), eq_fields_(eq_fields), eq_data_(eq_data) {
200  this->used_fields_ += eq_fields_->bc_conc;
201  }
202 
203  /// Destructor.
205 
206  /// Initialize auxiliary vectors and other data members
207  void initialize(ElementCacheMap *element_cache_map) {
208  this->element_cache_map_ = element_cache_map;
209 
210  fe_ = std::make_shared< FE_P_disc<dim> >(0);
213  }
214 
215 
216  /// Assemble integral over element
217  inline void cell_integral(DHCellAccessor cell, unsigned int element_patch_idx)
218  {
219  ASSERT_EQ(cell.dim(), dim).error("Dimension of element mismatch!");
220  if (!eq_data_->sources_changed_) return;
221 
222  // read for all substances
223  double max_cfl=0;
224  double source, diag;
225  ElementAccessor<3> elm = cell.elm();
226  // we have currently zero order P_Disc FE
227  ASSERT(cell.get_loc_dof_indices().size() == 1);
228  IntIdx local_p0_dof = cell.get_loc_dof_indices()[0];
229 
230  auto p = *( this->bulk_points(element_patch_idx).begin() );
231  for (unsigned int sbi = 0; sbi < eq_data_->n_substances(); sbi++)
232  {
233  source = eq_fields_->cross_section(p) * (eq_fields_->sources_density[sbi](p)
234  + eq_fields_->sources_sigma[sbi](p) * eq_fields_->sources_conc[sbi](p));
235  // addition to RHS
236  eq_data_->corr_vec[sbi].set(local_p0_dof, source);
237  // addition to diagonal of the transport matrix
238  diag = eq_fields_->sources_sigma[sbi](p) * eq_fields_->cross_section(p);
239  eq_data_->tm_diag[sbi].set(local_p0_dof, - diag);
240 
241  // compute maximal cfl condition over all substances
242  max_cfl = std::max(max_cfl, fabs(diag));
243 
244  eq_data_->balance_->add_source_values(sbi, elm.region().bulk_idx(), {local_p0_dof},
245  {- eq_fields_->sources_sigma[sbi](p) * elm.measure() * eq_fields_->cross_section(p)},
246  {source * elm.measure()});
247  }
248 
249  eq_data_->cfl_source_.set(local_p0_dof, max_cfl);
250  }
251 
252  /// Assembles the fluxes on the boundary.
253  inline void boundary_side_integral(DHCellSide cell_side)
254  {
255  ASSERT_EQ(cell_side.dim(), dim).error("Dimension of element mismatch!");
256  if (!cell_side.cell().is_own()) return;
257 
258  // we have currently zero order P_Disc FE
259  ASSERT(cell_side.cell().get_loc_dof_indices().size() == 1);
260  IntIdx local_p0_dof = cell_side.cell().get_loc_dof_indices()[0];
261  LongIdx glob_p0_dof = eq_data_->dh_->get_local_to_global_map()[local_p0_dof];
262 
263  fe_values_side_.reinit(cell_side.side());
264 
265  unsigned int sbi;
266  auto p_side = *( this->boundary_points(cell_side).begin() );
267  auto p_bdr = p_side.point_bdr(cell_side.cond().element_accessor() );
268  double flux = eq_fields_->side_flux(p_side, fe_values_side_);
269  if (flux < 0.0) {
270  double aij = -(flux / cell_side.element().measure() );
271 
272  for (sbi=0; sbi<eq_data_->n_substances(); sbi++)
273  {
274  double value = eq_fields_->bc_conc[sbi](p_bdr);
275 
276  VecSetValue(eq_data_->bcvcorr[sbi], glob_p0_dof, value * aij, ADD_VALUES);
277 
278  // CAUTION: It seems that PETSc possibly optimize allocated space during assembly.
279  // So we have to add also values that may be non-zero in future due to changing velocity field.
280  eq_data_->balance_->add_flux_values(eq_data_->subst_idx[sbi], cell_side,
281  {local_p0_dof}, {0.0}, flux*value);
282  }
283  } else {
284  for (sbi=0; sbi<eq_data_->n_substances(); sbi++)
285  VecSetValue(eq_data_->bcvcorr[sbi], glob_p0_dof, 0, ADD_VALUES);
286 
287  for (sbi=0; sbi<eq_data_->n_substances(); sbi++)
288  eq_data_->balance_->add_flux_values(eq_data_->subst_idx[sbi], cell_side,
289  {local_p0_dof}, {flux}, 0.0);
290  }
291  }
292 
293  /// Implements @p AssemblyBase::begin.
294  void begin() override
295  {
300 
301  //TODO: would it be possible to check the change in data for chosen substance? (may be in multifields?)
302 
303  if (eq_data_->sources_changed_) eq_data_->balance_->start_source_assembly(eq_data_->subst_idx);
304  // Assembly bcvcorr vector
305  for(unsigned int sbi=0; sbi < eq_data_->n_substances(); sbi++) VecZeroEntries(eq_data_->bcvcorr[sbi]);
306 
307  eq_data_->balance_->start_flux_assembly(eq_data_->subst_idx);
308  }
309 
310  /// Implements @p AssemblyBase::end.
311  void end() override
312  {
313  eq_data_->balance_->finish_flux_assembly(eq_data_->subst_idx);
314  if (eq_data_->sources_changed_) eq_data_->balance_->finish_source_assembly(eq_data_->subst_idx);
315 
316  for (unsigned int sbi=0; sbi<eq_data_->n_substances(); sbi++) VecAssemblyBegin(eq_data_->bcvcorr[sbi]);
317  for (unsigned int sbi=0; sbi<eq_data_->n_substances(); sbi++) VecAssemblyEnd(eq_data_->bcvcorr[sbi]);
318 
319  // we are calling set_boundary_conditions() after next_time() and
320  // we are using data from t() before, so we need to set corresponding bc time
322  }
323 
324  private:
325  shared_ptr<FiniteElement<dim>> fe_; ///< Finite element for the solution of the advection-diffusion equation.
326 
327  /// Data objects shared with ConvectionTransport
330 
331  /// Sub field set contains fields used in calculation.
333 
334  FEValues<3> fe_values_side_; ///< FEValues of object (of P disc finite element type)
335 
336  template < template<IntDim...> class DimAssembly>
337  friend class GenericAssembly;
338 
339 };
340 
341 
342 /**
343  * Auxiliary container class for Finite element and related objects of given dimension.
344  *
345  *
346  * Assembly convection term part of the matrix and boundary matrix for application of boundary conditions.
347  *
348  * Discretization of the convection term use explicit time scheme and finite volumes with full upwinding.
349  * We count on with exchange between dimensions and mixing on edges where more then two elements connect (can happen for 2D and 1D elements in
350  * 3D embedding space)
351  *
352  * In order to get multiplication matrix for explicit transport one have to scale the convection part by the acctual time step and
353  * add time term, i. e. unit matrix (see. transport_matrix_step_mpi)
354  *
355  * Updates CFL time step constrain.
356  */
357 template <unsigned int dim>
359 {
360 public:
363 
364  static constexpr const char * name() { return "MatrixMpiAssemblyConvection"; }
365 
366  /// Constructor.
368  : AssemblyBase<dim>(0), eq_fields_(eq_fields), eq_data_(eq_data) {
371  }
372 
373  /// Destructor.
375 
376  /// Initialize auxiliary vectors and other data members
377  void initialize(ElementCacheMap *element_cache_map) {
378  this->element_cache_map_ = element_cache_map;
379 
380  fe_ = std::make_shared< FE_P_disc<dim> >(0);
384  for (unsigned int sid=0; sid<eq_data_->max_edg_sides; sid++)
385  {
386  fe_values_vec_[sid].initialize(*this->quad_low_, *fe_, u);
387  }
393  dof_indices_i_.resize(1);
394  dof_indices_j_.resize(1);
395  }
396 
397  /// Assembles the fluxes between sides of elements of the same dimension.
399  ASSERT_EQ(edge_side_range.begin()->element().dim(), dim).error("Dimension of element mismatch!");
400 
401  unsigned int sid=0, s1, s2, i_col;
402  edg_flux = 0.0;
403  for( DHCellSide edge_side : edge_side_range )
404  {
405  fe_values_vec_[sid].reinit(edge_side.side());
406  edge_side.cell().get_dof_indices(dof_indices_i_);
407  side_dofs_[sid] = dof_indices_i_[0];
408  elm_meassures_[sid] = edge_side.element().measure();
409  auto p = *( this->edge_points(edge_side).begin() );
411  if (side_flux_[sid] > 0.0) {
412  eq_data_->cfl_flow_.add_global(side_dofs_[sid], -(side_flux_[sid] / edge_side.element().measure()) );
413  edg_flux += side_flux_[sid];
414  }
415  ++sid;
416  }
417 
418  unsigned int n_edge_sides = edge_side_range.begin()->n_edge_sides();
419  if (n_edge_sides<2) return; // following loop has no effect if edge contains only 1 side
420  for( s1=0; s1<n_edge_sides; s1++ )
421  {
422  for( s2=0, i_col=0; s2<n_edge_sides; s2++ )
423  {
424  if (s2==s1) continue;
425 
426  if (side_flux_[s2] > 0.0 && side_flux_[s1] < 0.0)
427  row_values_[i_col] = -(side_flux_[s1] * side_flux_[s2] / ( edg_flux * elm_meassures_[s1] ) );
428  else row_values_[i_col] = 0;
429  all_elem_dofs_[i_col++] = side_dofs_[s2];
430  }
431  MatSetValues(eq_data_->tm, 1, &(side_dofs_[s1]), n_edge_sides-1, &(all_elem_dofs_[0]), &(row_values_[0]), INSERT_VALUES);
432  }
433  }
434 
435 
436  /// Assembles the fluxes between elements of different dimensions.
437  inline void dimjoin_intergral(DHCellAccessor cell_lower_dim, DHCellSide neighb_side) {
438  if (dim == 1) return;
439  ASSERT_EQ(cell_lower_dim.dim(), dim-1).error("Dimension of element mismatch!");
440 
441  auto p_high = *( this->coupling_points(neighb_side).begin() );
442  fe_values_side_.reinit(neighb_side.side());
443 
444  cell_lower_dim.get_dof_indices(dof_indices_i_);
445  neighb_side.cell().get_dof_indices(dof_indices_j_);
447 
448  // volume source - out-flow from higher dimension
449  if (flux > 0.0) aij = flux / cell_lower_dim.elm().measure();
450  else aij=0;
451  MatSetValue(eq_data_->tm, dof_indices_i_[0], dof_indices_j_[0], aij, INSERT_VALUES);
452  // out flow from higher dim. already accounted
453 
454  // volume drain - in-flow to higher dimension
455  if (flux < 0.0) {
456  eq_data_->cfl_flow_.add_global( dof_indices_i_[0], (flux / cell_lower_dim.elm().measure()) ); // diagonal drain
457  aij = (-flux) / neighb_side.element().measure();
458  } else aij=0;
459  MatSetValue(eq_data_->tm, dof_indices_j_[0], dof_indices_i_[0], aij, INSERT_VALUES);
460  }
461 
462 
463  /// Implements @p AssemblyBase::begin.
464  void begin() override
465  {
466  MatZeroEntries(eq_data_->tm);
468  }
469 
470  /// Implements @p AssemblyBase::end.
471  void end() override
472  {
473  for ( DHCellAccessor dh_cell : eq_data_->dh_->own_range() ) {
474  dh_cell.get_dof_indices(dof_indices_i_);
475  MatSetValue(eq_data_->tm, dof_indices_i_[0], dof_indices_i_[0], eq_data_->cfl_flow_.get(dh_cell.local_idx()), INSERT_VALUES);
476 
477  eq_data_->cfl_flow_.set(dh_cell.local_idx(), fabs(eq_data_->cfl_flow_.get(dh_cell.local_idx())) );
478  }
479 
480  MatAssemblyBegin(eq_data_->tm, MAT_FINAL_ASSEMBLY);
481  MatAssemblyEnd(eq_data_->tm, MAT_FINAL_ASSEMBLY);
482  VecAssemblyBegin(eq_data_->cfl_flow_.petsc_vec());
483  VecAssemblyEnd(eq_data_->cfl_flow_.petsc_vec());
484 
487  }
488 
489 private:
490  shared_ptr<FiniteElement<dim>> fe_; ///< Finite element for the solution of the advection-diffusion equation.
491 
492  /// Data objects shared with ConvectionTransport
495 
496  /// Sub field set contains fields used in calculation.
498 
499  FEValues<3> fe_values_side_; ///< FEValues of object (of P disc finite element type)
500  vector<FEValues<3>> fe_values_vec_; ///< Vector of FEValues of object (of P disc finite element types)
501  vector<LongIdx> dof_indices_i_, dof_indices_j_; ///< Global DOF indices.
502 
508  double aij;
509  double edg_flux, flux;
510 
511  template < template<IntDim...> class DimAssembly>
512  friend class GenericAssembly;
513 
514 };
515 
516 
517 #endif /* ASSEMBLY_CONVECTION_HH_ */
#define ASSERT(expr)
Definition: asserts.hh:351
#define ASSERT_EQ(a, b)
Definition of comparative assert macro (EQual) only for debug mode.
Definition: asserts.hh:333
Range< BulkPoint > bulk_points(unsigned int element_patch_idx) const
Return BulkPoint range of appropriate dimension.
Range< EdgePoint > edge_points(const DHCellSide &cell_side) const
Return EdgePoint range of appropriate dimension.
Quadrature * quad_low_
Quadrature used in assembling methods (dim-1).
Range< CouplingPoint > coupling_points(const DHCellSide &cell_side) const
Return CouplingPoint range of appropriate dimension.
int active_integrals_
Holds mask of active integrals.
Range< BoundaryPoint > boundary_points(const DHCellSide &cell_side) const
Return BoundaryPoint range of appropriate dimension.
ElementCacheMap * element_cache_map_
ElementCacheMap shared with GenericAssembly object.
ElementAccessor< 3 > element_accessor()
FEValues< 3 > fe_values_side_
FEValues of object (of P disc finite element type)
ConcSourcesBdrAssemblyConvection(EqFields *eq_fields, EqData *eq_data)
Constructor.
EqFields * eq_fields_
Data objects shared with ConvectionTransport.
void end() override
Implements AssemblyBase::end.
FieldSet used_fields_
Sub field set contains fields used in calculation.
void initialize(ElementCacheMap *element_cache_map)
Initialize auxiliary vectors and other data members.
ConvectionTransport::EqData EqData
void boundary_side_integral(DHCellSide cell_side)
Assembles the fluxes on the boundary.
shared_ptr< FiniteElement< dim > > fe_
Finite element for the solution of the advection-diffusion equation.
static constexpr const char * name()
void begin() override
Implements AssemblyBase::begin.
ConvectionTransport::EqFields EqFields
void cell_integral(DHCellAccessor cell, unsigned int element_patch_idx)
Assemble integral over element.
bool sources_changed_
Flag indicates that sources part of equation was changed during last time.
Definition: transport.h:161
bool is_convection_matrix_scaled
Flag indicates the state of object.
Definition: transport.h:177
Vec * bcvcorr
Boundary condition correction vector.
Definition: transport.h:168
VectorMPI cfl_source_
Parallel vector for source term contribution to CFL condition.
Definition: transport.h:164
std::shared_ptr< Balance > balance_
object for calculation and writing the mass balance to file, shared with EquationBase.
Definition: transport.h:147
VectorMPI cfl_flow_
Parallel vector for flow contribution to CFL condition.
Definition: transport.h:165
unsigned int max_edg_sides
Maximal number of edge sides (evaluate from dim 1,2,3)
Definition: transport.h:180
double transport_bc_time
Time of the last update of the boundary condition terms.
Definition: transport.h:169
unsigned int n_substances()
Returns number of transported substances.
Definition: transport.h:126
vector< VectorMPI > tm_diag
additions to PETSC transport matrix on the diagonal - from sources (for each substance)
Definition: transport.h:167
Mat tm
PETSc transport matrix.
Definition: transport.h:166
vector< VectorMPI > corr_vec
Definition: transport.h:163
vector< unsigned int > subst_idx
List of indices used to call balance methods for a set of quantities.
Definition: transport.h:156
std::shared_ptr< DOFHandlerMultiDim > dh_
Definition: transport.h:144
MultiField< 3, FieldValue< 3 >::Scalar > init_conc
Initial concentrations.
Definition: transport.h:106
FieldFEScalarVec conc_mobile_fe
Underlaying FieldFE for each substance of conc_mobile.
Definition: transport.h:112
double side_flux(SidePoint &side_p, FEValues< 3 > &fe_side_values)
Calculate flux on given side point.
Definition: transport.h:95
BCMultiField< 3, FieldValue< 3 >::Scalar > bc_conc
Definition: transport.h:103
Cell accessor allow iterate over DOF handler cells.
const ElementAccessor< 3 > elm() const
Return ElementAccessor to element of loc_ele_idx_.
bool is_own() const
Return true if accessor represents own element (false for ghost element)
LocDofVec get_loc_dof_indices() const
Returns the local indices of dofs associated to the cell on the local process.
unsigned int get_dof_indices(std::vector< LongIdx > &indices) const
Fill vector of the global indices of dofs associated to the cell.
unsigned int dim() const
Return dimension of element appropriate to cell.
unsigned int local_idx() const
Return local index to element (index of DOF handler).
Side accessor allows to iterate over sides of DOF handler cell.
Side side() const
Return Side of given cell and side_idx.
Boundary cond() const
const DHCellAccessor & cell() const
Return DHCellAccessor appropriate to the side.
unsigned int dim() const
Return dimension of element appropriate to the side.
ElementAccessor< 3 > element() const
double measure() const
Computes the measure of the element.
Region region() const
Definition: accessors.hh:198
Directing class of FieldValueCache.
void initialize(Quadrature &_quadrature, FiniteElement< DIM > &_fe, UpdateFlags _flags)
Initialize structures and calculates cell-independent data.
Definition: fe_values.cc:137
void reinit(const ElementAccessor< spacedim > &cell)
Update cell-dependent data (gradients, Jacobians etc.)
Definition: fe_values.cc:537
bool changed() const
Container for various descendants of FieldCommonBase.
Definition: field_set.hh:159
Generic class of assemblation.
void cell_integral(DHCellAccessor cell, unsigned int element_patch_idx)
Assemble integral over element.
InitCondAssemblyConvection(EqFields *eq_fields, EqData *eq_data)
Constructor.
ConvectionTransport::EqFields EqFields
shared_ptr< FiniteElement< dim > > fe_
Finite element for the solution of the advection-diffusion equation.
ConvectionTransport::EqData EqData
std::vector< VectorMPI > vecs_
Set of data vectors of conc_mobile_fe objects.
FieldSet used_fields_
Sub field set contains fields used in calculation.
EqFields * eq_fields_
Data objects shared with TransportDG.
void initialize(ElementCacheMap *element_cache_map)
Initialize auxiliary vectors and other data members.
static constexpr const char * name()
EqFields * eq_fields_
Data objects shared with TransportDG.
void cell_integral(DHCellAccessor cell, unsigned int element_patch_idx)
Assemble integral over element.
void end() override
Implements AssemblyBase::end.
void begin() override
Implements AssemblyBase::begin.
~MassAssemblyConvection()
Destructor.
FieldSet used_fields_
Sub field set contains fields used in calculation.
MassAssemblyConvection(EqFields *eq_fields, EqData *eq_data)
Constructor.
static constexpr const char * name()
ConvectionTransport::EqFields EqFields
void initialize(ElementCacheMap *element_cache_map)
Initialize auxiliary vectors and other data members.
shared_ptr< FiniteElement< dim > > fe_
Finite element for the solution of the advection-diffusion equation.
ConvectionTransport::EqData EqData
ConvectionTransport::EqFields EqFields
static constexpr const char * name()
void initialize(ElementCacheMap *element_cache_map)
Initialize auxiliary vectors and other data members.
std::vector< LongIdx > side_dofs_
void edge_integral(RangeConvert< DHEdgeSide, DHCellSide > edge_side_range)
Assembles the fluxes between sides of elements of the same dimension.
std::vector< LongIdx > all_elem_dofs_
MatrixMpiAssemblyConvection(EqFields *eq_fields, EqData *eq_data)
Constructor.
std::vector< double > elm_meassures_
ConvectionTransport::EqData EqData
std::vector< double > row_values_
void end() override
Implements AssemblyBase::end.
std::vector< double > side_flux_
EqFields * eq_fields_
Data objects shared with ConvectionTransport.
shared_ptr< FiniteElement< dim > > fe_
Finite element for the solution of the advection-diffusion equation.
void begin() override
Implements AssemblyBase::begin.
void dimjoin_intergral(DHCellAccessor cell_lower_dim, DHCellSide neighb_side)
Assembles the fluxes between elements of different dimensions.
vector< FEValues< 3 > > fe_values_vec_
Vector of FEValues of object (of P disc finite element types)
FEValues< 3 > fe_values_side_
FEValues of object (of P disc finite element type)
vector< LongIdx > dof_indices_j_
Global DOF indices.
FieldSet used_fields_
Sub field set contains fields used in calculation.
IterConvert< ObjectIn, ObjectOut > begin()
Iterator to begin item of range.
unsigned int bulk_idx() const
Returns index of the region in the bulk set.
Definition: region.hh:90
double t() const
double last_t() const
Field< 3, FieldValue< 3 >::Scalar > water_content
Water content - result of unsaturated water flow model or porosity.
MultiField< 3, FieldValue< 3 >::Scalar > sources_sigma
Concentration sources - Robin type, in_flux = sources_sigma * (sources_conc - mobile_conc)
MultiField< 3, FieldValue< 3 >::Scalar > sources_conc
MultiField< 3, FieldValue< 3 >::Scalar > sources_density
Concentration sources - density of substance source, only positive part is used.
Field< 3, FieldValue< 3 >::VectorFixed > flow_flux
Flow flux, can be result of water flow model.
Field< 3, FieldValue< 3 >::Scalar > cross_section
Pointer to DarcyFlow field cross_section.
void add_global(unsigned int pos, double val)
Add value to item on given global position.
Definition: vector_mpi.hh:132
double get(unsigned int pos) const
Return value on given position.
Definition: vector_mpi.hh:105
void set(unsigned int pos, double val)
Set value on given position.
Definition: vector_mpi.hh:111
void zero_entries()
Definition: vector_mpi.hh:82
Vec & petsc_vec()
Getter for PETSC vector of output data (e.g. can be used by scatters).
Definition: vector_mpi.hh:79
Definitions of basic Lagrangean finite elements with polynomial shape functions.
Class FEValues calculates finite element data on the actual cells such as shape function values,...
@ coupling
@ boundary
@ bulk
@ edge
int IntIdx
Definition: index_types.hh:25
int LongIdx
Define type that represents indices of large arrays (elements, nodes, dofs etc.)
Definition: index_types.hh:24
static constexpr bool value
Definition: json.hpp:87
unsigned int IntDim
A dimension index type.
Definition: mixed.hh:19
Definitions of particular quadrature rules on simplices.
UpdateFlags
Enum type UpdateFlags indicates which quantities are to be recomputed on each finite element cell.
Definition: update_flags.hh:68
@ update_values
Shape function values.
Definition: update_flags.hh:87
@ update_normal_vectors
Normal vectors.
@ update_side_JxW_values
Transformed quadrature weight for cell sides.
@ update_quadrature_points
Transformed quadrature points.