I'm trying to implement a time-dependent solver that assembles the system
using the MeshWorker::loop tool. However, I cannot figure out how get the
values of the solution at the previous timestep at each quadrature point.
Currently, I store the old solution in a class called "AdvectionProblem"
such that---I've left out some code in the interest of brevity.
template<unsigned int dim>
class AdvectionProblem {
public:
...
void AssembleSystem()
private:
...
Vector<double> old_solution;
DoFHandler<dim> dofHandler;
};
Since we need the old_solution values to assemble the system, I create a
CellIntegrator class that MeshWorker::loop can use to integrate over each
cell (I also create similar objects for faces/boundaries):
template<unsigned int dim>
class IntegrateCell {
public:
IntegrateCell(dealii::Vector<double> const& old_solution);
virtual ~IntegrateCell() = default;
void operator()(dealii::MeshWorker::DoFInfo<dim>& dinfo,
dealii::MeshWorker::IntegrationInfo<dim>& info) const;
private:
/// Store a constant reference to the old solution
const dealii::Vector<double>& old_solution;
};
template<unsigned int dim>
IntegrateCell<dim>::IntegrateCell(Vector<double> const& old_solution) :
old_solution(old_solution) {}
template<unsigned int dim>
void IntegrateCell<dim>::operator()(MeshWorker::DoFInfo<dim>& dinfo,
MeshWorker::IntegrationInfo<dim>& info) const {
// get some useful objects
const FEValuesBase<dim>& feValues = info.fe_values();
FullMatrix<double>& localMatrix = dinfo.matrix(0).matrix;
Vector<double>& localVector = dinfo.vector(0).block(0);
const std::vector<double>& JxW = feValues.get_JxW_values();
// get the values of the old solution
std::vector<double> old_solution_values(feValues.n_quadrature_points);
feValues.get_function_values(old_solution, old_solution_values);
*// these are non-zero, IngtegrateCell is correctly getting the
old_solution std::cout << old_solution << std::endl;*
// loop over the quadrature points and compute the advection vector at
the current point
for( unsigned int point=0; point<feValues.n_quadrature_points; ++point ) {
*// these are always all zero*
* std::cout << "value: " << old_solution_values[point] << std::endl;*
}
}
I then try to assemble the system using the MeshWorker::loop tool:
template<unsigned int dim>
void AdvectionProblem<dim>::AssembleSystem() {
// an object that knows about data structures and local integration
MeshWorker::IntegrationInfoBox<dim> infoBox;
// initialize the quadrature formula
const unsigned int nGauss = dofHandler.get_fe().degree + 1;
infoBox.initialize_gauss_quadrature(nGauss, nGauss, nGauss);
// these are the values we need to integrate our system
infoBox.initialize_update_flags();
const UpdateFlags updateFlags = update_quadrature_points | update_values
| update_gradients | update_JxW_values;
infoBox.add_update_flags(updateFlags, true, true, true, true);
// initialize the finite element values
infoBox.initialize(fe, mapping);
// create an object that forwards local data to the assembler
MeshWorker::DoFInfo<dim> dofInfo(dofHandler);
// create teh assembler---tell it where to put local data
MeshWorker::Assembler::SystemSimple<SparseMatrix<double>, Vector<double>
> assembler;
assembler.initialize(systemMatrix, rhs);
IntegrateCell<dim> cellIntegration(old_solution);
*MeshWorker::loop<dim, dim, MeshWorker::DoFInfo<dim>,
MeshWorker::IntegrationInfoBox<dim> >(dofHandler.begin_active(),
dofHandler.end(), dofInfo, infoBox, cellIntegration, nullptr, nullptr,
assembler);*
}
*For some reason the feValues.get_function_values(old_solution,
old_solution_values); call in the CellIntegrator function always sets the
old_solution_values to zero. Does anyone see what I'm doing wrong?*
Thanks!
--
The deal.II project is located at http://www.dealii.org/
For mailing list/forum options, see
https://groups.google.com/d/forum/dealii?hl=en
---
You received this message because you are subscribed to the Google Groups
"deal.II User Group" group.
To unsubscribe from this group and stop receiving emails from it, send an email
to dealii+unsubscr...@googlegroups.com.
To view this discussion on the web visit
https://groups.google.com/d/msgid/dealii/aa9ce784-3b72-4398-aa8d-0f33e6235ae6%40googlegroups.com.
