8 #include "flow/mortar_assembly.hh"
20 quadrature_(*(eq_data->mesh))
25 for(
uint row_dim=0; row_dim<4; row_dim ++)
26 for(
uint col_dim=0; col_dim<4; col_dim ++) {
28 row_avg.fill(1.0 / (row_dim+1));
30 col_avg.fill(1.0 / (col_dim+1));
48 i_data.
dofs.resize(nsides);
55 for(
unsigned int i_side=0; i_side < nsides; i_side++ ) {
113 if (ele.
dim() > 2)
return;
115 if (isec_list.size() == 0)
return;
120 double m_sigma =
eq_fields_->sigma.value( ele_centre, ele);
121 double m_conductivity =
eq_fields_->conductivity.value( ele_centre, ele);
122 double m_crossection =
eq_fields_->cross_section.value( ele_centre, ele );
124 double master_sigma = 2*m_sigma*m_conductivity *
141 double cs_sqr_avg = 0.0;
142 double isec_sum = 0.0;
143 unsigned int slave_dim = 0;
145 for(; i < isec_list.size(); ++i) {
149 slave_dim = ele_slave.
dim();
150 if (slave_dim == ele.
dim())
break;
151 if (! non_zero)
continue;
156 cs_sqr_avg += cs*cs*isec_measure;
157 isec_sum += isec_measure;
161 if ( ! (slave_dim == 2 && ele.
dim() ==2 ) ) {
162 if( fabs(isec_sum - ele.
measure()) > 1E-5) {
164 for(
auto & isec : isec_list) {
180 master_sigma = master_sigma * (cs_sqr_avg / isec_sum)
189 if (i < isec_list.size()) {
192 for(; i < isec_list.size(); ++i) {
196 isec_sum += isec_measure;
202 master_sigma = 100 * m_conductivity/ m_crossection / isec_sum;
216 double s = -scale * row_ele.delta * col_ele.delta;
222 for(
uint i=0; i< row_ele.n_dirichlet; i++)
250 for(
uint icol=0; icol < n_cols; icol++ ) pressure[icol] =
eq_data_->full_solution.get(col_ele.
dofs[icol]);
253 for(
uint irow=0; irow < n_rows; irow++ )
eq_data_->full_solution.add( row_ele.
vel_dofs[irow], add_velocity[irow] );
261 for(
unsigned int i_side=0; i_side < ele->
n_sides(); i_side++ ) {
273 dofs[shift + i_side] = -
dofs[shift + i_side];
ElementAccessor< 3 > element_accessor()
Cell accessor allow iterate over DOF handler cells.
const ElementAccessor< 3 > elm() const
Return ElementAccessor to element of loc_ele_idx_.
unsigned int n_dofs() const
Return number of dofs on given cell.
LocDofVec get_loc_dof_indices() const
Returns the local indices of dofs associated to the cell on the local process.
double measure() const
Computes the measure of the element.
SideIter side(const unsigned int loc_index)
unsigned int idx() const
We need this method after replacing Region by RegionIdx, and movinf RegionDB instance into particular...
arma::vec::fixed< spacedim > centre() const
Computes the barycenter.
unsigned int n_sides() const
bool reinit(const IntersectionLocalBase *isec)
void set_solution_col(uint loc_col, double solution)
void eliminate_solution()
const arma::mat & get_matrix()
void set_rhs(arma::vec rhs)
const arma::vec & get_rhs()
void set_matrix(arma::mat matrix)
void reset()
Resets the matrix, RHS, dofs to zero and clears solution settings.
void set_solution_row(uint loc_row, double solution, double diag=0.0)
std::vector< std::vector< ILpair > > element_intersections_
AssemblyFieldsPtr eq_fields_
MixedMeshIntersections & mixed_mesh_
void fix_velocity_local(const IsecData &row_ele, const IsecData &col_ele)
void add_to_linsys(double scale)
vector< IsecData > isec_data_list
std::vector< arma::vec > col_average_
void pressure_diff(const DHCellAccessor &dh_cell, double delta)
P0_CouplingAssembler(AssemblyFieldsPtr eq_fields, AssemblyDataPtr eq_data)
IntersectionQuadratureP0 quadrature_
arma::mat & tensor_average(unsigned int row_dim, unsigned int col_dim)
void assembly(const DHCellAccessor &dh_cell)
void assembly(const DHCellAccessor &dh_cell)
vector< double > dirichlet
void add_sides(const DHCellAccessor &dh_cell, unsigned int shift, vector< int > &dofs, vector< double > &dirichlet)
bool is_boundary() const
Returns true for side on the boundary.
Wrappers for linear systems based on MPIAIJ and MATIS format.
#define WarningOut()
Macro defining 'warning' record of log.
std::shared_ptr< DarcyMH::EqFields > AssemblyFieldsPtr
std::shared_ptr< DarcyMH::EqData > AssemblyDataPtr
ArmaMat< double, N, M > mat
std::string format(CStringRef format_str, ArgList args)