#include <numerics/matrices.h>
#include <numerics/data_out.h>
#include <numerics/derivative_approximation.h>
+#include <numerics/solution_transfer.h>
#include <fstream>
#include <sstream>
double get_maximal_velocity () const;
void solve ();
void output_results () const;
+ void refine_mesh ();
const unsigned int degree;
if (setup_matrices == true)
{
+ system_matrix.clear ();
+
sparsity_pattern.reinit (3,3);
sparsity_pattern.block(0,0).reinit (n_u, n_u, n_couplings);
sparsity_pattern.block(1,0).reinit (n_p, n_u, n_couplings);
<< " CG iterations for temperature."
<< std::endl;
}
-
-
- old_solution = solution;
}
+template <int dim>
+void
+BoussinesqFlowProblem<dim>::refine_mesh ()
+{
+ return;
+
+ Vector<float> estimated_error_per_cell (triangulation.n_active_cells());
+
+//TODO do this better
+ DerivativeApproximation::approximate_gradient (dof_handler,
+ old_solution,
+ estimated_error_per_cell,
+ dim+1);
+
+ GridRefinement::refine_and_coarsen_fixed_fraction (triangulation,
+ estimated_error_per_cell,
+ 0.3, 0.03,
+ triangulation.n_active_cells());
+
+ SolutionTransfer<dim, double> soltrans(dof_handler);
+
+ triangulation.prepare_coarsening_and_refinement();
+
+ Vector<double> x_old_solution (dof_handler.n_dofs());
+ x_old_solution = old_solution;
+
+ soltrans.prepare_for_coarsening_and_refinement(x_old_solution);
+
+ triangulation.execute_coarsening_and_refinement ();
+ setup_dofs (true);
+
+ Vector<double> tmp (dof_handler.n_dofs());
+ soltrans.interpolate(x_old_solution, tmp);
+
+ rebuild_preconditioner = true;
+
+ old_solution = tmp;
+}
+
+
+
template <int dim>
double
BoussinesqFlowProblem<dim>::get_maximal_velocity () const
}
- for (unsigned int pre_refinement=0; pre_refinement<5-dim; ++pre_refinement)
+ for (unsigned int pre_refinement=0; pre_refinement<3-dim; ++pre_refinement)
{
setup_dofs(false);
Vector<float> estimated_error_per_cell (triangulation.n_active_cells());
+//TODO do this better
DerivativeApproximation::approximate_gradient (dof_handler,
old_solution,
estimated_error_per_cell,
time += time_step;
++timestep_number;
+
+ old_solution = solution;
std::cout << std::endl;
+
+ if (timestep_number % 1 == 0)
+ refine_mesh ();
}
while (time <= 500);
}