DoFHandler<dim> dof_handler;
ConstraintMatrix constraints;
- MatrixFree<dim,double> system_mf_storage;
+ std_cxx11::shared_ptr<MatrixFree<dim,double> > system_mf_storage;
typedef LaplaceOperator<dim,degree_finite_element,double> SystemMatrixType;
SystemMatrixType system_matrix;
#endif
fe (degree_finite_element),
dof_handler (triangulation),
+ system_mf_storage(new MatrixFree<dim,double>()),
pcout (std::cout, Utilities::MPI::this_mpi_process(MPI_COMM_WORLD) == 0),
// The LaplaceProblem class holds an additional output stream that
// collects detailed timings about the setup phase. This stream, called
MatrixFree<dim,double>::AdditionalData::none;
additional_data.mapping_update_flags = (update_gradients | update_JxW_values |
update_quadrature_points);
- system_mf_storage.reinit (dof_handler, constraints, QGauss<1>(fe.degree+1),
- additional_data);
+ system_mf_storage->reinit (dof_handler, constraints, QGauss<1>(fe.degree+1),
+ additional_data);
}
system_matrix.initialize (system_mf_storage);
mg_mf_storage[level].reinit(dof_handler, level_constraints,
QGauss<1>(fe.degree+1), additional_data);
- mg_matrices[level].initialize(mg_mf_storage[level], mg_constrained_dofs, level);
+ mg_matrices[level].initialize(
+ std_cxx11::make_shared<MatrixFree<dim,float>>(mg_mf_storage[level]),
+ mg_constrained_dofs, level);
mg_matrices[level].evaluate_coefficient(Coefficient<dim>());
}
setup_time += time.wall_time();
Timer time;
system_rhs = 0;
- FEEvaluation<dim,degree_finite_element> phi(system_mf_storage);
- for (unsigned int cell=0; cell<system_mf_storage.n_macro_cells(); ++cell)
+ FEEvaluation<dim,degree_finite_element> phi(*system_mf_storage);
+ for (unsigned int cell=0; cell<system_mf_storage->n_macro_cells(); ++cell)
{
phi.reinit(cell);
for (unsigned int q=0; q<phi.n_q_points; ++q)
/**
* Initialize operator on fine scale.
*/
- void initialize (const MatrixFree<dim,Number> &data);
+ void initialize (std_cxx11::shared_ptr<const MatrixFree<dim,Number> > data);
/**
* Initialize operator on a level @p level.
*/
- void initialize (const MatrixFree<dim,Number> &data,
+ void initialize (std_cxx11::shared_ptr<const MatrixFree<dim,Number> > data,
const MGConstrainedDoFs &mg_constrained_dofs,
const unsigned int level);
/**
* MatrixFree object to be used with this operator.
*/
- SmartPointer<const MatrixFree<dim,Number>, Base<dim,Number> > data;
+ std_cxx11::shared_ptr<const MatrixFree<dim,Number> > data;
/**
* A shared pointer to a diagonal matrix that stores the inverse of
template <int dim, typename Number>
void
Base<dim,Number>::
- initialize (const MatrixFree<dim,Number> &data_)
+ initialize (std_cxx11::shared_ptr<const MatrixFree<dim,Number> > data_)
{
- data = SmartPointer<const MatrixFree<dim,Number>, Base<dim,Number> >(&data_,typeid(*this).name());
+ data = data_;
edge_constrained_indices.clear();
have_interface_matrices = false;
}
template <int dim, typename Number>
void
Base<dim,Number>::
- initialize (const MatrixFree<dim,Number> &data_,
+ initialize (std_cxx11::shared_ptr<const MatrixFree<dim,Number> > data_,
const MGConstrainedDoFs &mg_constrained_dofs,
const unsigned int level)
{
AssertThrow (level != numbers::invalid_unsigned_int,
ExcMessage("level is not set"));
- if (data_.n_macro_cells() > 0)
+ if (data_->n_macro_cells() > 0)
{
AssertDimension(static_cast<int>(level),
- data_.get_cell_iterator(0,0)->level());
+ data_->get_cell_iterator(0,0)->level());
}
- data = SmartPointer<const MatrixFree<dim,Number>, Base<dim,Number> >(&data_,typeid(*this).name());
+ data = data_;
// setup edge_constrained indices
std::vector<types::global_dof_index> interface_indices;
if (locally_owned.is_element(interface_indices[i]))
edge_constrained_indices.push_back(locally_owned.index_within_set(interface_indices[i]));
have_interface_matrices = Utilities::MPI::max((unsigned int)edge_constrained_indices.size(),
- data_.get_vector_partitioner()->get_mpi_communicator()) > 0;
+ data->get_vector_partitioner()->get_mpi_communicator()) > 0;
}
* which stores quadrature point data.
*/
template <int dim, typename VectorType>
- void project (const MatrixFree<dim,typename VectorType::value_type> &data,
+ void project (std_cxx11::shared_ptr<const MatrixFree<dim,typename VectorType::value_type> > data,
const ConstraintMatrix &constraints,
const unsigned int n_q_points_1d,
const std_cxx11::function< VectorizedArray<typename VectorType::value_type> (const unsigned int, const unsigned int)> func,
* Same as above but for <code>n_q_points_1d = matrix_free.get_dof_handler().get_fe().degree+1</code>.
*/
template <int dim, typename VectorType>
- void project (const MatrixFree<dim,typename VectorType::value_type> &data,
+ void project (std_cxx11::shared_ptr<const MatrixFree<dim,typename VectorType::value_type> > data,
const ConstraintMatrix &constraints,
const std_cxx11::function< VectorizedArray<typename VectorType::value_type> (const unsigned int, const unsigned int)> func,
VectorType &vec_result);
additional_data.tasks_parallel_scheme =
MatrixFree<dim,number>::AdditionalData::partition_color;
additional_data.mapping_update_flags = (update_values | update_JxW_values);
- MatrixFree<dim, number> matrix_free;
- matrix_free.reinit (mapping, dof, constraints,
- QGauss<1>(fe_degree+2), additional_data);
+ std_cxx11::shared_ptr<MatrixFree<dim, number> > matrix_free(
+ new MatrixFree<dim, number> ());
+ matrix_free->reinit (mapping, dof, constraints,
+ QGauss<1>(fe_degree+2), additional_data);
typedef MatrixFreeOperators::MassOperator<dim, fe_degree, fe_degree+2, components, number> MatrixType;
MatrixType mass_matrix;
mass_matrix.initialize(matrix_free);
typedef LinearAlgebra::distributed::Vector<number> LocalVectorType;
LocalVectorType vec, rhs, inhomogeneities;
- matrix_free.initialize_dof_vector(vec);
- matrix_free.initialize_dof_vector(rhs);
- matrix_free.initialize_dof_vector(inhomogeneities);
+ matrix_free->initialize_dof_vector(vec);
+ matrix_free->initialize_dof_vector(rhs);
+ matrix_free->initialize_dof_vector(inhomogeneities);
constraints.distribute(inhomogeneities);
inhomogeneities*=-1.;
additional_data.tasks_parallel_scheme =
MatrixFree<dim,Number>::AdditionalData::partition_color;
additional_data.mapping_update_flags = (update_values | update_JxW_values);
- MatrixFree<dim, Number> matrix_free;
- matrix_free.reinit (mapping, dof, constraints,
- QGauss<1>(fe_degree+2), additional_data);
+ std_cxx11::shared_ptr<MatrixFree<dim, Number> > matrix_free(
+ new MatrixFree<dim, Number>());
+ matrix_free->reinit (mapping, dof, constraints,
+ QGauss<1>(fe_degree+2), additional_data);
typedef MatrixFreeOperators::MassOperator<dim, fe_degree, fe_degree+2, 1, Number> MatrixType;
MatrixType mass_matrix;
mass_matrix.initialize(matrix_free);
typedef LinearAlgebra::distributed::Vector<Number> LocalVectorType;
LocalVectorType vec, rhs, inhomogeneities;
- matrix_free.initialize_dof_vector(vec);
- matrix_free.initialize_dof_vector(rhs);
- matrix_free.initialize_dof_vector(inhomogeneities);
+ matrix_free->initialize_dof_vector(vec);
+ matrix_free->initialize_dof_vector(rhs);
+ matrix_free->initialize_dof_vector(inhomogeneities);
constraints.distribute(inhomogeneities);
inhomogeneities*=-1.;
template <int dim, typename VectorType, int spacedim, int fe_degree, int n_q_points_1d>
- void project_parallel (const MatrixFree<dim,typename VectorType::value_type> &matrix_free,
+ void project_parallel (std_cxx11::shared_ptr<const MatrixFree<dim,typename VectorType::value_type> > matrix_free,
const ConstraintMatrix &constraints,
const std_cxx11::function< VectorizedArray<typename VectorType::value_type> (const unsigned int, const unsigned int)> func,
VectorType &vec_result)
{
- const DoFHandler<dim,spacedim> &dof = matrix_free.get_dof_handler();
+ const DoFHandler<dim,spacedim> &dof = matrix_free->get_dof_handler();
typedef typename VectorType::value_type Number;
Assert (dof.get_fe().n_components() == 1,
typedef LinearAlgebra::distributed::Vector<Number> LocalVectorType;
LocalVectorType vec, rhs, inhomogeneities;
- matrix_free.initialize_dof_vector(vec);
- matrix_free.initialize_dof_vector(rhs);
- matrix_free.initialize_dof_vector(inhomogeneities);
+ matrix_free->initialize_dof_vector(vec);
+ matrix_free->initialize_dof_vector(rhs);
+ matrix_free->initialize_dof_vector(inhomogeneities);
constraints.distribute(inhomogeneities);
inhomogeneities*=-1.;
// assemble right hand side:
{
- FEEvaluation<dim,fe_degree,n_q_points_1d,1,Number> fe_eval(matrix_free);
- const unsigned int n_cells = matrix_free.n_macro_cells();
+ FEEvaluation<dim,fe_degree,n_q_points_1d,1,Number> fe_eval(*matrix_free);
+ const unsigned int n_cells = matrix_free->n_macro_cells();
const unsigned int n_q_points = fe_eval.n_q_points;
for (unsigned int cell=0; cell<n_cells; ++cell)
template <int dim, typename VectorType>
- void project (const MatrixFree<dim,typename VectorType::value_type> &matrix_free,
+ void project (std_cxx11::shared_ptr<const MatrixFree<dim,typename VectorType::value_type> > matrix_free,
const ConstraintMatrix &constraints,
const unsigned int n_q_points_1d,
const std_cxx11::function< VectorizedArray<typename VectorType::value_type> (const unsigned int, const unsigned int)> func,
VectorType &vec_result)
{
(void) n_q_points_1d;
- const unsigned int fe_degree = matrix_free.get_dof_handler().get_fe().degree;
+ const unsigned int fe_degree = matrix_free->get_dof_handler().get_fe().degree;
Assert (fe_degree+1 == n_q_points_1d,
ExcNotImplemented());
template <int dim, typename VectorType>
- void project (const MatrixFree<dim,typename VectorType::value_type> &matrix_free,
+ void project (std_cxx11::shared_ptr<const MatrixFree<dim,typename VectorType::value_type> > matrix_free,
const ConstraintMatrix &constraints,
const std_cxx11::function< VectorizedArray<typename VectorType::value_type> (const unsigned int, const unsigned int)> func,
VectorType &vec_result)
{
project (matrix_free,
constraints,
- matrix_free.get_dof_handler().get_fe().degree+1,
+ matrix_free->get_dof_handler().get_fe().degree+1,
func,
vec_result);
}
template
void project<deal_II_dimension, VEC>(
- const MatrixFree<deal_II_dimension,VEC::value_type> &matrix_free,
+ std_cxx11::shared_ptr<const MatrixFree<deal_II_dimension,VEC::value_type> > matrix_free,
const ConstraintMatrix &constraints,
const std_cxx11::function< VectorizedArray<VEC::value_type> (const unsigned int, const unsigned int)>,
VEC &);
template
void project<deal_II_dimension, VEC>(
- const MatrixFree<deal_II_dimension,VEC::value_type> &matrix_free,
+ std_cxx11::shared_ptr<const MatrixFree<deal_II_dimension,VEC::value_type> > matrix_free,
const ConstraintMatrix &constraints,
const unsigned int,
const std_cxx11::function< VectorizedArray<VEC::value_type> (const unsigned int, const unsigned int)>,
//std::cout << "Number of degrees of freedom: " << dof.n_dofs() << std::endl;
//std::cout << "Number of constraints: " << constraints.n_constraints() << std::endl;
- MatrixFree<dim,number> mf_data;
+ std_cxx11::shared_ptr<MatrixFree<dim,number> > mf_data(new MatrixFree<dim,number> ());
{
const QGauss<1> quad (fe_degree+1);
typename MatrixFree<dim,number>::AdditionalData data;
data.tasks_parallel_scheme =
MatrixFree<dim,number>::AdditionalData::none;
data.tasks_block_size = 7;
- mf_data.reinit (dof, constraints, quad, data);
+ mf_data->reinit (dof, constraints, quad, data);
}
MatrixFreeOperators::LaplaceOperator<dim,fe_degree,fe_degree+1, 1, number> mf;
mf.initialize(mf_data);
mf.compute_diagonal();
LinearAlgebra::distributed::Vector<number> in, out, ref;
- mf_data.initialize_dof_vector (in);
+ mf_data->initialize_dof_vector (in);
out.reinit (in);
ref.reinit (in);
//std::cout << "Number of degrees of freedom: " << dof.n_dofs() << std::endl;
//std::cout << "Number of constraints: " << constraints.n_constraints() << std::endl;
- MatrixFree<dim,number> mf_data;
+ std_cxx11::shared_ptr<MatrixFree<dim,number> > mf_data(new MatrixFree<dim,number> ());
{
const QGauss<1> quad (fe_degree+1);
typename MatrixFree<dim,number>::AdditionalData data;
MatrixFree<dim,number>::AdditionalData::none;
data.tasks_block_size = 7;
data.mapping_update_flags = update_quadrature_points | update_gradients | update_JxW_values;
- mf_data.reinit (dof, constraints, quad, data);
+ mf_data->reinit (dof, constraints, quad, data);
}
std_cxx11::shared_ptr<Table<2, VectorizedArray<number> > > coefficient;
coefficient = std_cxx11::make_shared<Table<2, VectorizedArray<number> > >();
{
- FEEvaluation<dim,fe_degree,fe_degree+1,1,number> fe_eval(mf_data);
+ FEEvaluation<dim,fe_degree,fe_degree+1,1,number> fe_eval(*mf_data);
- const unsigned int n_cells = mf_data.n_macro_cells();
+ const unsigned int n_cells = mf_data->n_macro_cells();
const unsigned int n_q_points = fe_eval.n_q_points;
coefficient->reinit(n_cells, n_q_points);
mf.set_coefficient(coefficient);
mf.compute_diagonal();
LinearAlgebra::distributed::Vector<number> in, out, ref;
- mf_data.initialize_dof_vector (in);
+ mf_data->initialize_dof_vector (in);
out.reinit (in);
ref.reinit (in);
//std::cout << "Number of degrees of freedom: " << dof.n_dofs() << std::endl;
//std::cout << "Number of constraints: " << constraints.n_constraints() << std::endl;
- MatrixFree<dim,number> mf_data;
+ std_cxx11::shared_ptr<MatrixFree<dim,number> > mf_data(new MatrixFree<dim,number> ());
{
const QGauss<1> quad (fe_degree+2);
typename MatrixFree<dim,number>::AdditionalData data;
data.tasks_parallel_scheme =
MatrixFree<dim,number>::AdditionalData::none;
data.tasks_block_size = 7;
- mf_data.reinit (dof, constraints, quad, data);
+ mf_data->reinit (dof, constraints, quad, data);
}
MatrixFreeOperators::MassOperator<dim,fe_degree, fe_degree+2, 1, number> mf;
mf.initialize(mf_data);
mf.compute_diagonal();
LinearAlgebra::distributed::Vector<number> in, out, ref;
- mf_data.initialize_dof_vector (in);
+ mf_data->initialize_dof_vector (in);
out.reinit (in);
ref.reinit (in);
//std::cout << "Number of degrees of freedom: " << dof.n_dofs() << std::endl;
//std::cout << "Number of constraints: " << constraints.n_constraints() << std::endl;
- MatrixFree<dim,number> mf_data;
+ std_cxx11::shared_ptr<MatrixFree<dim,number> > mf_data(new MatrixFree<dim,number> ());
{
const QGauss<1> quad (fe_degree+2);
typename MatrixFree<dim,number>::AdditionalData data;
data.tasks_parallel_scheme =
MatrixFree<dim,number>::AdditionalData::none;
data.tasks_block_size = 7;
- mf_data.reinit (dof, constraints, quad, data);
+ mf_data->reinit (dof, constraints, quad, data);
}
MatrixFreeOperators::MassOperator<dim,fe_degree, fe_degree+2, 1, number> mf;
= mf.get_matrix_diagonal_inverse()->get_vector();
LinearAlgebra::distributed::Vector<number> in, out, ref;
- mf_data.initialize_dof_vector (in);
+ mf_data->initialize_dof_vector (in);
out.reinit (in);
ref.reinit (in);
//std::cout << "Number of degrees of freedom: " << dof.n_dofs() << std::endl;
//std::cout << "Number of constraints: " << constraints.n_constraints() << std::endl;
- MatrixFree<dim,number> mf_data;
+ std_cxx11::shared_ptr<MatrixFree<dim,number> > mf_data(new MatrixFree<dim,number> ());
{
const QGauss<1> quad (fe_degree+2);
typename MatrixFree<dim,number>::AdditionalData data;
data.tasks_parallel_scheme =
MatrixFree<dim,number>::AdditionalData::none;
data.tasks_block_size = 7;
- mf_data.reinit (dof, constraints, quad, data);
+ mf_data->reinit (dof, constraints, quad, data);
}
MatrixFreeOperators::MassOperator<dim,fe_degree, fe_degree+2, 1, number> mf;
= mf.get_matrix_diagonal_inverse()->get_vector();
LinearAlgebra::distributed::Vector<number> in, out, ref;
- mf_data.initialize_dof_vector (in);
+ mf_data->initialize_dof_vector (in);
out.reinit (in);
ref.reinit (in);
MappingQ<dim> mapping(fe_degree+1);
LaplaceOperator<dim,fe_degree,n_q_points_1d,1,number> fine_matrix;
- MatrixFree<dim,number> fine_level_data;
+ std_cxx11::shared_ptr<MatrixFree<dim,number> > fine_level_data(new MatrixFree<dim,number> ());
typename MatrixFree<dim,number>::AdditionalData fine_level_additional_data;
fine_level_additional_data.tasks_parallel_scheme = MatrixFree<dim,number>::AdditionalData::none;
- fine_level_data.reinit (mapping, dof, constraints, QGauss<1>(n_q_points_1d),
- fine_level_additional_data);
+ fine_level_data->reinit (mapping, dof, constraints, QGauss<1>(n_q_points_1d),
+ fine_level_additional_data);
fine_matrix.initialize(fine_level_data);
fine_matrix.compute_diagonal();
mg_level_data[level].reinit (mapping, dof, level_constraints, QGauss<1>(n_q_points_1d),
mg_additional_data);
- mg_matrices[level].initialize(mg_level_data[level],
+ mg_matrices[level].initialize(std_cxx11::make_shared<MatrixFree<dim,number> >(
+ mg_level_data[level]),
mg_constrained_dofs,
level);
mg_matrices[level].compute_diagonal();
MappingQ<dim> mapping(fe_degree+1);
LaplaceOperator<dim,fe_degree,n_q_points_1d,1,number> fine_matrix;
- MatrixFree<dim,number> fine_level_data;
+ std_cxx11::shared_ptr<MatrixFree<dim,number> > fine_level_data(new MatrixFree<dim,number> ());
typename MatrixFree<dim,number>::AdditionalData fine_level_additional_data;
fine_level_additional_data.tasks_parallel_scheme = MatrixFree<dim,number>::AdditionalData::none;
fine_level_additional_data.tasks_block_size = 3;
- fine_level_data.reinit (mapping, dof, constraints, QGauss<1>(n_q_points_1d),
- fine_level_additional_data);
+ fine_level_data->reinit (mapping, dof, constraints, QGauss<1>(n_q_points_1d),
+ fine_level_additional_data);
fine_matrix.initialize(fine_level_data);
fine_matrix.compute_diagonal();
mg_level_data[level].reinit (mapping, dof, level_constraints, QGauss<1>(n_q_points_1d),
mg_additional_data);
- mg_matrices[level].initialize(mg_level_data[level],
+ mg_matrices[level].initialize(std_cxx11::make_shared<MatrixFree<dim,number> >(
+ mg_level_data[level]),
mg_constrained_dofs,
level);
mg_matrices[level].compute_diagonal();
MappingQ<dim> mapping(fe_degree+1);
LaplaceOperator<dim,fe_degree,n_q_points_1d,1,number> fine_matrix;
- MatrixFree<dim,number> fine_level_data;
+ std_cxx11::shared_ptr<MatrixFree<dim,number> > fine_level_data(new MatrixFree<dim,number> ());
typename MatrixFree<dim,number>::AdditionalData fine_level_additional_data;
fine_level_additional_data.tasks_parallel_scheme = MatrixFree<dim,number>::AdditionalData::none;
fine_level_additional_data.tasks_block_size = 3;
- fine_level_data.reinit (mapping, dof, constraints, QGauss<1>(n_q_points_1d),
- fine_level_additional_data);
+ fine_level_data->reinit (mapping, dof, constraints, QGauss<1>(n_q_points_1d),
+ fine_level_additional_data);
fine_matrix.initialize(fine_level_data);
fine_matrix.compute_diagonal();
mg_level_data[level].reinit (mapping, dof, level_constraints, QGauss<1>(n_q_points_1d),
mg_additional_data);
- mg_matrices[level].initialize(mg_level_data[level],
+ mg_matrices[level].initialize(std_cxx11::make_shared<MatrixFree<dim,number> > (
+ mg_level_data[level]),
mg_constrained_dofs,
level);
mg_matrices[level].compute_diagonal();