18 #ifndef ASSEMBLY_FLOW_OUTPUT_HH_
19 #define ASSEMBLY_FLOW_OUTPUT_HH_
53 template <
unsigned int dim>
60 static constexpr
const char *
name() {
return "L2DifferenceAssembly"; }
81 fe_p0_ = std::make_shared< FE_P_disc<dim> >(0);
92 ASSERT_EQ(cell.
dim(), dim).error(
"Dimension of element mismatch!");
99 for (
unsigned int li = 0; li < ele->
n_sides(); li++) {
106 double velocity_diff=0, divergence_diff=0, pressure_diff=0, diff;
110 double mean_x_squared=0;
111 for(
unsigned int i_node=0; i_node < ele->
n_nodes(); i_node++ )
112 for(
unsigned int j_node=0; j_node < ele->
n_nodes(); j_node++ )
114 mean_x_squared += (i_node == j_node ? 2.0 : 1.0) / ( 6 * dim )
115 * arma::dot( *ele.
node(i_node), *ele.
node(j_node));
118 unsigned int i_point=0, oposite_node;
119 for (
auto p : this->
bulk_points(element_patch_idx) )
131 for(
unsigned int i_shape=0; i_shape < ele->
n_sides(); i_shape++) {
138 + arma::dot( ele.
centre(), *ele.
node(oposite_node) )
148 for(
unsigned int i_shape=0; i_shape < ele->
n_sides(); i_shape++) {
159 for(
unsigned int i_shape=0; i_shape < ele->
n_sides(); i_shape++) diff +=
fluxes_[ i_shape ];
171 velocity_data.set( idx, sqrt(velocity_diff) );
178 pressure_data.set( idx, sqrt(pressure_diff) );
182 div_data.set( idx, sqrt(divergence_diff) );
192 for(
unsigned int j=0; j<3; j++){
205 for(
unsigned int j=0; j<3; j++){
212 os <<
"l2 norm output\n\n"
249 template <
template<
IntDim...>
class DimAssembly>
258 template <
unsigned int dim>
265 static constexpr
const char *
name() {
return "OutputInternalFlowAssembly"; }
286 ASSERT_EQ(cell.
dim(), dim).error(
"Dimension of element mismatch!");
291 std::stringstream ss;
296 auto p = *( this->
bulk_points(element_patch_idx).begin() );
298 for (
unsigned int i = 0; i < 3; i++)
307 for (
unsigned int i = 0; i < ele->
n_sides(); i++) {
310 uint old_opp_node = new_to_old_node[new_opp_node];
311 uint old_iside = ele->
n_sides() - old_opp_node - 1;
312 old_to_new_side[old_iside] = i;
317 for (
unsigned int i = 0; i < ele->
n_sides(); i++) {
318 uint new_lid = ele->
n_sides() + 1 + old_to_new_side[i];
322 for (
unsigned int i = 0; i < ele->
n_sides(); i++) {
323 uint new_iside = old_to_new_side[i];
328 string line = ss.str();
335 DebugOut() <<
"Internal flow data output\n";
345 eq_data_->
raw_output_file <<
"// fields:\n//ele_id ele_presure flux_in_barycenter[3] n_sides side_pressures[n] side_fluxes[n]\n";
350 for (
unsigned int i_elem=0; i_elem<
eq_data_->
flow_data_->dh_->n_own_cells(); ++i_elem) {
368 template <
template<
IntDim...>
class DimAssembly>
#define ASSERT_EQ(a, b)
Definition of comparative assert macro (EQual) only for debug mode.
Range< BulkPoint > bulk_points(unsigned int element_patch_idx) const
Return BulkPoint range of appropriate dimension.
Quadrature * quad_
Quadrature used in assembling methods.
int active_integrals_
Holds mask of active integrals.
ElementCacheMap * element_cache_map_
ElementCacheMap shared with GenericAssembly object.
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.
DHCellAccessor cell_with_other_dh(const DOFHandlerMultiDim *dh) const
Create new accessor with same local idx and given DOF handler. Actual and given DOF handler must be c...
unsigned int dim() const
Return dimension of element appropriate to cell.
unsigned int elm_idx() const
Return serial idx to element of loc_ele_idx_.
Output specific field stuff.
std::shared_ptr< FieldFE< 3, FieldValue< 3 >::Scalar > > div_diff_ptr
std::shared_ptr< SubDOFHandlerMultiDim > dh_
std::shared_ptr< FieldFE< 3, FieldValue< 3 >::Scalar > > vel_diff_ptr
std::vector< int > velocity_mask
std::shared_ptr< FieldFE< 3, FieldValue< 3 >::Scalar > > pressure_diff_ptr
std::shared_ptr< DarcyLMH::EqData > flow_data_
std::shared_ptr< DarcyLMH::EqData > flow_data_
ofstream raw_output_file
Raw data output file.
TimeGovernor * time_
Time is shared with flow equation.
std::vector< std::string > raw_output_strings_
Output lines of cells.
Field< 3, FieldValue< 3 >::Scalar > ref_divergence
Field< 3, FieldValue< 3 >::Scalar > cross_section
Field< 3, FieldValue< 3 >::Scalar > conductivity
Field< 3, FieldValue< 3 >::VectorFixed > field_ele_velocity
Field< 3, FieldValue< 3 >::VectorFixed > ref_velocity
Precompute l2 difference outputs.
Field< 3, FieldValue< 3 >::Scalar > ref_pressure
auto & orig_nodes_order() const
unsigned int input_id() const
Return the element ID in the input mesh. Should be only used for special output.
NodeAccessor< 3 > node(unsigned int ni) const
double measure() const
Computes the measure of the element.
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.
Directing class of FieldValueCache.
unsigned int n_sides() const
unsigned int n_nodes() const
double JxW(const unsigned int point_no)
Return the product of Jacobian determinant and the quadrature weight at given quadrature point.
const FEValuesViews::Vector< spacedim > & vector_view(unsigned int i) const
Accessor to vector values of multicomponent FE.
void initialize(Quadrature &_quadrature, FiniteElement< DIM > &_fe, UpdateFlags _flags)
Initialize structures and calculates cell-independent data.
void reinit(const ElementAccessor< spacedim > &cell)
Update cell-dependent data (gradients, Jacobians etc.)
arma::vec::fixed< spacedim > point(const unsigned int point_no)
Return coordinates of the quadrature point in the actual cell system.
Raviart-Thomas element of order 0.
Container for various descendants of FieldCommonBase.
Dedicated class for storing path to input and output files.
Generic class of assemblation.
void end() override
Implements AssemblyBase::end.
L2DifferenceAssembly(EqFields *eq_fields, EqData *eq_data)
Constructor.
arma::vec3 flux_in_q_point_
Precomputed flux in quadrature point.
FieldSet used_fields_
Sub field set contains fields used in calculation.
DarcyLMH::EqFields EqFields
std::vector< double > fluxes_
Precomputed fluxes on element sides.
DarcyFlowMHOutput::DiffEqData EqData
void initialize(ElementCacheMap *element_cache_map)
Initialize auxiliary vectors and other data members.
void begin() override
Implements AssemblyBase::begin.
double cross_
Field results.
virtual ~L2DifferenceAssembly()
Destructor.
static constexpr const char * name()
void cell_integral(DHCellAccessor cell, unsigned int element_patch_idx)
Assemble integral over element.
FE_RT0< dim > fe_rt_
FEValues for velocity.
EqFields * eq_fields_
Data objects shared with Flow equation.
std::shared_ptr< FE_P_disc< dim > > fe_p0_
arma::vec3 q_point_
Local coords of quadrature point.
arma::vec3 ref_flux_
Field result.
FieldSet used_fields_
Sub field set contains fields used in calculation.
void cell_integral(DHCellAccessor cell, unsigned int element_patch_idx)
Assemble integral over element.
DarcyFlowMHOutput::RawOutputEqData EqData
void initialize(ElementCacheMap *element_cache_map)
Initialize auxiliary vectors and other data members.
arma::vec3 flux_in_center_
static constexpr const char * name()
void begin() override
Implements AssemblyBase::begin.
OutputInternalFlowAssembly(EqFields *eq_fields, EqData *eq_data)
Constructor.
virtual ~OutputInternalFlowAssembly()
Destructor.
DarcyLMH::EqFields EqFields
void end() override
Implements AssemblyBase::end.
EqFields * eq_fields_
Data objects shared with Flow equation.
static unsigned int oposite_node(unsigned int sid)
Lumped mixed-hybrid model of linear Darcy flow, possibly unsteady.
Output class for darcy_flow_mh model.
Definitions of basic Lagrangean finite elements with polynomial shape functions.
Definitions of Raviart-Thomas finite elements.
Class FESystem for compound finite elements.
Class FEValues calculates finite element data on the actual cells such as shape function values,...
arma::Col< IntIdx > LocDofVec
#define DebugOut()
Macro defining 'debug' record of log.
unsigned int IntDim
A dimension index type.
std::string format(CStringRef format_str, ArgList args)
Definitions of particular quadrature rules on simplices.
UpdateFlags
Enum type UpdateFlags indicates which quantities are to be recomputed on each finite element cell.
@ update_values
Shape function values.
@ update_JxW_values
Transformed quadrature weights.
@ update_quadrature_points
Transformed quadrature points.