double estimated_error;
};
+ public:
+
struct ErrorOnCell {
double part[8];
ErrorOnCell ();
};
-struct CellwiseError
+ struct CellwiseError
{
CellwiseError (const unsigned int n_errors);
std::vector<ErrorOnCell> errors;
ErrorOnCell* next_free_slot;
};
+ protected:
+
Vector<float> estimated_error_per_cell;
FullMatrix<double> embedding_matrix;
FullMatrix<double> interpolation_matrix;
const unsigned int) const {
const double pi = 3.1415926536;
if (p(0)==0)
- return sin(pi*get_time()/0.4)*sin(pi*get_time()/0.4);
+ return sin(pi*this->get_time()/0.4)*sin(pi*this->get_time()/0.4);
else
return 0;
};
const unsigned int) const {
const double pi = 3.1415926536;
if (p(0)==0)
- return 2*pi/0.4*sin(pi*get_time()/0.4)*cos(pi*get_time()/0.4);
+ return 2*pi/0.4*sin(pi*this->get_time()/0.4)*cos(pi*this->get_time()/0.4);
else
return 0;
};
virtual double value (const Point<dim> &p,
const unsigned int) const {
const double pi = 3.1415926536;
- if ((get_time()<0.2) && (p(0)==0))
- return sin(pi*get_time()/0.2)*sin(pi*get_time()/0.2);
+ if ((this->get_time()<0.2) && (p(0)==0))
+ return sin(pi*this->get_time()/0.2)*sin(pi*this->get_time()/0.2);
else
return 0;
};
virtual double value (const Point<dim> &p,
const unsigned int) const {
const double pi = 3.1415926536;
- if ((get_time()<0.2) && (p(0)==0))
- return 2*pi/0.2*sin(pi*get_time()/0.2)*cos(pi*get_time()/0.2);
+ if ((this->get_time()<0.2) && (p(0)==0))
+ return 2*pi/0.2*sin(pi*this->get_time()/0.2)*cos(pi*this->get_time()/0.2);
else
return 0;
};
const unsigned int) const {
const double pi = 3.1415926536;
if ((0.4 <= p(1)) && (p(1) <= 0.6) && (p(0) <= 0.5))
- return (p(1)-0.4)*(0.6-p(1)) * sin(pi*get_time()/0.2);
+ return (p(1)-0.4)*(0.6-p(1)) * sin(pi*this->get_time()/0.2);
else
return 0;
};
const unsigned int) const {
const double pi = 3.1415926536;
if ((0.4 <= p(1)) && (p(1) <= 0.6) && (p(0) <= 0.5))
- return pi/0.2*(p(1)-0.4)*(0.6-p(1)) * cos(pi*get_time()/0.2);
+ return pi/0.2*(p(1)-0.4)*(0.6-p(1)) * cos(pi*this->get_time()/0.2);
else
return 0;
};
const double period = 60;
- if ((get_time()>=period) || (r>=a))
+ if ((this->get_time()>=period) || (r>=a))
return 0;
- const double s = cos(r/a*pi/2)*sin(pi*get_time()/period);
+ const double s = cos(r/a*pi/2)*sin(pi*this->get_time()/period);
return s*s;
};
};
const double a = 5000000;
const double period = 60;
- if ((get_time()>=period) || (r>=a))
+ if ((this->get_time()>=period) || (r>=a))
return 0;
else
return (2*pi/period*cos(r/a*pi/2)*cos(r/a*pi/2)*
- sin(pi*get_time()/period)*cos(pi*get_time()/period));
+ sin(pi*this->get_time()/period)*cos(pi*this->get_time()/period));
};
};
template <int dim>
TimeStepBase_Wave<dim>::TimeStepBase_Wave ():
TimeStepBase_Tria<dim> (),
- parameters (*static_cast<WaveParameters<dim>*>(0))
+ parameters (parameters)
{}
case TimeStepBase::dual_problem:
{
Assert (((this->next_action == TimeStepBase::primal_problem) &&
- (static_cast<const TimeStep_Wave<dim>*>(&get_timestep_primal())
+ (static_cast<const TimeStep_Wave<dim>*>(&this->get_timestep_primal())
== this))
||
((this->next_action == TimeStepBase::dual_problem) &&
- (static_cast<const TimeStep_Wave<dim>*>(&get_timestep_dual())
+ (static_cast<const TimeStep_Wave<dim>*>(&this->get_timestep_dual())
== this)),
ExcInternalError());
case TimeStepBase::postprocess:
{
this->sweep_info->get_timers().postprocessing.start();
- std::ifstream tmp_in(tmp_filename_base(branch_signature()).c_str());
+ std::ifstream tmp_in(this->tmp_filename_base(branch_signature()).c_str());
u.block_read (tmp_in);
v.block_read (tmp_in);
tmp_in.close ();
Assert (u.size() != 0, ExcInternalError());
Assert (v.size() != 0, ExcInternalError());
- std::ofstream tmp_out(tmp_filename_base(branch_signature()).c_str());
+ std::ofstream tmp_out(this->tmp_filename_base(branch_signature()).c_str());
u.block_write (tmp_out);
v.block_write (tmp_out);
tmp_out.close ();
template <int dim>
void TimeStep_Wave<dim>::end_sweep ()
{
- std::string tmp_filename = tmp_filename_base(branch_signature());
+ std::string tmp_filename = this->tmp_filename_base(branch_signature());
remove (tmp_filename.c_str());
}
system_matrix.copy_from (this->mass_matrix);
this->constraints.condense (static_cast<SparseMatrix<double>&>(system_matrix));
const unsigned int
- solver_steps1 = solve (system_matrix, this->u, tmp_u_bar),
- solver_steps2 = solve (system_matrix, this->v, tmp_v_bar);
+ solver_steps1 = this->solve (system_matrix, this->u, tmp_u_bar),
+ solver_steps2 = this->solve (system_matrix, this->v, tmp_v_bar);
this->statistic_data = typename TimeStep_Wave<dim>::StatisticData (this->tria->n_active_cells(),
this->dof_handler->n_dofs(),
solver_steps1, solver_steps2,
- compute_energy ());
+ this->compute_energy ());
}
else
this->statistic_data = typename TimeStep_Wave<dim>::StatisticData (this->tria->n_active_cells(),
this->sweep_info->get_data().dual_dofs += this->dof_handler->n_dofs() * 2;
- const double time_step = get_forward_timestep ();
+ const double time_step = this->get_forward_timestep ();
Vector<double> right_hand_side1 (this->dof_handler->n_dofs());
Vector<double> right_hand_side2 (this->dof_handler->n_dofs());
old_u.reinit (this->dof_handler->n_dofs());
old_v.reinit (this->dof_handler->n_dofs());
- transfer_old_solutions (old_u, old_v);
+ this->transfer_old_solutions (old_u, old_v);
};
assemble_vectors (right_hand_side1, right_hand_side2);
right_hand_side1);
};
- const unsigned int solver_steps1 = solve (system_matrix, this->v, right_hand_side1);
+ const unsigned int solver_steps1 = this->solve (system_matrix, this->v, right_hand_side1);
system_matrix.copy_from (this->mass_matrix);
this->constraints.condense (static_cast<SparseMatrix<double>&>(system_matrix));
this->u -= old_u;
};
- const unsigned int solver_steps2 = solve (system_matrix, this->u, right_hand_side2);
+ const unsigned int solver_steps2 = this->solve (system_matrix, this->u, right_hand_side2);
this->statistic_data = typename TimeStep_Wave<dim>::StatisticData (this->tria->n_active_cells(),
this->dof_handler->n_dofs(),
solver_steps1,
solver_steps2,
- compute_energy ());
+ this->compute_energy ());
deallog << "." << std::endl;
}
{
Assert (this->system_sparsity.empty(), ExcInternalError());
- create_matrices ();
+ this->create_matrices ();
};
this->sweep_info->get_timers().dual_problem.stop();
}
this->parameters.dual_functional->reset (*this);
this->parameters.dual_functional->compute_functionals (dual1, dual2);
- const double timestep = get_forward_timestep();
+ const double timestep = this->get_forward_timestep();
right_hand_side1.add (timestep, dual2);
right_hand_side1.add (this->parameters.theta * timestep * timestep, dual1);
= static_cast<const TimeStepBase_Wave<dim>*>(this->next_timestep)->get_timestep_dual();
const unsigned int dofs_per_cell = this->fe.dofs_per_cell;
- const double time_step = get_forward_timestep();
+ const double time_step = this->get_forward_timestep();
if (!old_cell->has_children() && !new_cell->has_children())
{
= static_cast<const TimeStepBase_Wave<dim>*>(this->next_timestep)->get_timestep_dual();
const unsigned int dofs_per_cell = this->fe.dofs_per_cell;
- const double time_step = get_forward_timestep();
+ const double time_step = this->get_forward_timestep();
FullMatrix<double> cell_matrix (dofs_per_cell, dofs_per_cell);
unsigned int level_difference = 1;
const unsigned int dofs_per_cell = this->fe.dofs_per_cell;
- const double time_step = get_forward_timestep();
+ const double time_step = this->get_forward_timestep();
FullMatrix<double> cell_matrix(dofs_per_cell, dofs_per_cell);
Vector<double> local_old_dof_values_u (dofs_per_cell);
Assert (estimated_error_per_cell.size()!=0,
ExcInternalError());
- std::ofstream tmp_out(tmp_filename_base(branch_signature()).c_str());
+ std::ofstream tmp_out(this->tmp_filename_base(branch_signature()).c_str());
estimated_error_per_cell.block_write (tmp_out);
tmp_out.close ();
void
TimeStep_ErrorEstimation<dim>::get_error_indicators (Vector<float> &indicators) const
{
- std::ifstream in (tmp_filename_base(branch_signature()).c_str());
+ std::ifstream in (this->tmp_filename_base(branch_signature()).c_str());
indicators.block_read (in);
}
neumann_boundary[1] = &homogeneous_neumann_bc;
const TimeStep_Wave<dim> &target = (which_variables==0 ?
- static_cast<const TimeStep_Wave<dim>&>(get_timestep_primal()) :
- static_cast<const TimeStep_Wave<dim>&>(get_timestep_dual ()));
+ static_cast<const TimeStep_Wave<dim>&>(this->get_timestep_primal()) :
+ static_cast<const TimeStep_Wave<dim>&>(this->get_timestep_dual ()));
KellyErrorEstimator<dim>::estimate (*target.dof_handler,
target.quadrature_face,
void TimeStep_ErrorEstimation<dim>::estimate_error_dual () {
CellwiseError cellwise_error (this->tria->n_active_cells());
- const TimeStep_Primal<dim> &primal_problem = get_timestep_primal(),
+ const TimeStep_Primal<dim> &primal_problem = this->get_timestep_primal(),
&primal_problem_old = static_cast<const TimeStepBase_Wave<dim>*>
(this->previous_timestep)->get_timestep_primal();
- const TimeStep_Dual<dim> &dual_problem = get_timestep_dual(),
+ const TimeStep_Dual<dim> &dual_problem = this->get_timestep_dual(),
&dual_problem_old = static_cast<const TimeStepBase_Wave<dim>*>
(this->previous_timestep)->get_timestep_dual();
};
-const TimeStep_Primal<dim> &primal_problem = get_timestep_primal(),
+const TimeStep_Primal<dim> &primal_problem = this->get_timestep_primal(),
&primal_problem_old = static_cast<const TimeStepBase_Wave<dim>*>
(this->previous_timestep)->get_timestep_primal();
- const TimeStep_Dual<dim> &dual_problem = get_timestep_dual(),
+ const TimeStep_Dual<dim> &dual_problem = this->get_timestep_dual(),
&dual_problem_old = static_cast<const TimeStepBase_Wave<dim>*>
(this->previous_timestep)->get_timestep_dual();
- const FiniteElement<dim> &primal_fe = get_timestep_primal().fe,
- &dual_fe = get_timestep_dual().fe;
+ const FiniteElement<dim> &primal_fe = this->get_timestep_primal().fe,
+ &dual_fe = this->get_timestep_dual().fe;
const unsigned int dofs_per_cell_primal = primal_fe.dofs_per_cell,
dofs_per_cell_dual = dual_fe.dofs_per_cell;
const Vector<double> &local_v_bar_old,
CellwiseError &cellwise_error,
FEValues<dim> &fe_values) const {
- const TimeStep_Primal<dim> &primal_problem = get_timestep_primal();
- const TimeStep_Dual<dim> &dual_problem = get_timestep_dual();
+ const TimeStep_Primal<dim> &primal_problem = this->get_timestep_primal();
+ const TimeStep_Dual<dim> &dual_problem = this->get_timestep_dual();
- const FiniteElement<dim> &dual_fe = get_timestep_dual().fe;
+ const FiniteElement<dim> &dual_fe = this->get_timestep_dual().fe;
const unsigned int dofs_per_cell_dual = dual_fe.dofs_per_cell;
ErrorOnCell error_on_cell;
- const unsigned int dofs_per_cell = get_timestep_dual().fe.dofs_per_cell;
+ const unsigned int dofs_per_cell = this->get_timestep_dual().fe.dofs_per_cell;
Vector<double> tmp1(dofs_per_cell);
Vector<double> tmp2(dofs_per_cell);
tmp1 = local_v;
tmp1 += local_v_old;
- error_on_cell.part[2] = -(get_backward_timestep() / 4 *
+ error_on_cell.part[2] = -(this->get_backward_timestep() / 4 *
mass_matrix.matrix_scalar_product (tmp1, tmp2));
tmp2 = local_u_bar;
tmp2 -= local_u_bar_old;
- error_on_cell.part[3] = -(get_backward_timestep() / 12 *
+ error_on_cell.part[3] = -(this->get_backward_timestep() / 12 *
mass_matrix.matrix_scalar_product (tmp1, tmp2));
tmp1 = local_u;
tmp1 += local_u_old;
- error_on_cell.part[4] = (get_backward_timestep() / 4 *
+ error_on_cell.part[4] = (this->get_backward_timestep() / 4 *
laplace_matrix.matrix_scalar_product (tmp1, tmp2));
tmp2 = local_v_bar;
tmp2 -= local_v_bar_old;
- error_on_cell.part[5] = (get_backward_timestep() / 12 *
+ error_on_cell.part[5] = (this->get_backward_timestep() / 12 *
laplace_matrix.matrix_scalar_product (tmp1, tmp2));
template <int dim>
void TimeStep_ErrorEstimation<dim>::make_interpolation_matrices () {
- const FiniteElement<dim> &primal_fe = get_timestep_primal().fe,
- &dual_fe = get_timestep_dual().fe;
+ const FiniteElement<dim> &primal_fe = this->get_timestep_primal().fe,
+ &dual_fe = this->get_timestep_dual().fe;
embedding_matrix.reinit (dual_fe.dofs_per_cell,
primal_fe.dofs_per_cell);
case TimeStepBase_Tria<dim>::grid_refinement:
if (sleep_level == 1)
- save_refine_flags ();
+ this->save_refine_flags ();
break;
default:
if ((this->sweep_no < this->parameters.number_of_sweeps-1) ||
(this->parameters.refinement_strategy == WaveParameters<dim>::dual_estimator))
- estimate_error ();
+ this->estimate_error ();
this->sweep_info->get_timers().postprocessing.start();
for (typename std::list<EvaluationBase<dim>*>::const_iterator i = this->parameters.eval_list.begin();
i != this->parameters.eval_list.end(); ++i)
{
- (*i)->reset_timelevel (get_timestep_primal());
+ (*i)->reset_timelevel (this->get_timestep_primal());
statistic_data.evaluation_results.push_back ((*i)->evaluate());
};
typename DataOut<dim>::OutputFormat output_format
= DataOut<dim>::parse_output_format (this->parameters.output_format);
- out.attach_dof_handler (*get_timestep_primal().dof_handler);
- out.add_data_vector (get_timestep_primal().u, "u");
- out.add_data_vector (get_timestep_primal().v, "v");
+ out.attach_dof_handler (*this->get_timestep_primal().dof_handler);
+ out.add_data_vector (this->get_timestep_primal().u, "u");
+ out.add_data_vector (this->get_timestep_primal().v, "v");
Vector<double> u_bar, v_bar;
&&
(this->sweep_no >= this->parameters.initial_energy_estimator_sweeps))
{
- u_bar.reinit (get_timestep_primal().u.size());
- v_bar.reinit (get_timestep_primal().u.size());
+ u_bar.reinit (this->get_timestep_primal().u.size());
+ v_bar.reinit (this->get_timestep_primal().u.size());
if (this->parameters.primal_fe == this->parameters.dual_fe)
{
- u_bar = get_timestep_dual().u;
- v_bar = get_timestep_dual().v;
+ u_bar = this->get_timestep_dual().u;
+ v_bar = this->get_timestep_dual().v;
}
else
interpolate_dual_solution (u_bar, v_bar);
}
else
{
- estimated_error.reinit (get_timestep_primal().dof_handler->n_dofs());
- DoFTools::distribute_cell_to_dof_vector (*get_timestep_primal().dof_handler,
+ estimated_error.reinit (this->get_timestep_primal().dof_handler->n_dofs());
+ DoFTools::distribute_cell_to_dof_vector (*this->get_timestep_primal().dof_handler,
this->estimated_error_per_cell,
estimated_error);
};
this->sweep_data->data_out_stack->new_parameter_value (this->time,
(this->timestep_no == 0 ?
0 :
- get_backward_timestep() *
+ this->get_backward_timestep() *
this->parameters.write_stacked_interval));
- this->sweep_data->data_out_stack->attach_dof_handler (*get_timestep_primal().dof_handler);
- this->sweep_data->data_out_stack->add_data_vector (get_timestep_primal().u, "u");
- this->sweep_data->data_out_stack->add_data_vector (get_timestep_primal().v, "v");
+ this->sweep_data->data_out_stack->attach_dof_handler (*this->get_timestep_primal().dof_handler);
+ this->sweep_data->data_out_stack->add_data_vector (this->get_timestep_primal().u, "u");
+ this->sweep_data->data_out_stack->add_data_vector (this->get_timestep_primal().v, "v");
if ((this->parameters.refinement_strategy == WaveParameters<dim>::dual_estimator)
&&
{
if (this->parameters.primal_fe == this->parameters.dual_fe)
{
- this->sweep_data->data_out_stack->add_data_vector (get_timestep_dual().u, "dual_u");
- this->sweep_data->data_out_stack->add_data_vector (get_timestep_dual().v, "dual_v");
+ this->sweep_data->data_out_stack->add_data_vector (this->get_timestep_dual().u, "dual_u");
+ this->sweep_data->data_out_stack->add_data_vector (this->get_timestep_dual().v, "dual_v");
}
else
{
- Vector<double> u_bar(get_timestep_primal().dof_handler->n_dofs());
- Vector<double> v_bar(get_timestep_primal().dof_handler->n_dofs());
+ Vector<double> u_bar(this->get_timestep_primal().dof_handler->n_dofs());
+ Vector<double> v_bar(this->get_timestep_primal().dof_handler->n_dofs());
interpolate_dual_solution (u_bar, v_bar);
template <int dim>
void TimeStep_Postprocess<dim>::end_sweep ()
{
- std::string tmp_filename = tmp_filename_base(branch_signature());
+ std::string tmp_filename = this->tmp_filename_base(branch_signature());
remove (tmp_filename.c_str());
}
template <int dim>
void TimeStep_Postprocess<dim>::interpolate_dual_solution (Vector<double> &interpolated_u_bar,
Vector<double> &interpolated_v_bar) const {
- const unsigned int n_primal_dofs = get_timestep_primal().dof_handler->n_dofs();
+ const unsigned int n_primal_dofs = this->get_timestep_primal().dof_handler->n_dofs();
interpolated_u_bar.reinit (n_primal_dofs);
interpolated_v_bar.reinit (n_primal_dofs);
- const TimeStep_Wave<dim> &target = get_timestep_dual ();
+ const TimeStep_Wave<dim> &target = this->get_timestep_dual ();
typename DoFHandler<dim>::active_cell_iterator primal_cell, dual_cell, endc;
- primal_cell = get_timestep_primal().dof_handler->begin_active();
- endc = get_timestep_primal().dof_handler->end();
+ primal_cell = this->get_timestep_primal().dof_handler->begin_active();
+ endc = this->get_timestep_primal().dof_handler->end();
dual_cell = target.dof_handler->begin_active();
for (; primal_cell != endc; ++primal_cell, ++dual_cell)
this->sweep_info->get_data().primal_dofs += this->dof_handler->n_dofs() * 2;
- const double time_step = get_backward_timestep ();
+ const double time_step = this->get_backward_timestep ();
Vector<double> right_hand_side1 (this->dof_handler->n_dofs());
Vector<double> right_hand_side2 (this->dof_handler->n_dofs());
old_u.reinit (this->dof_handler->n_dofs());
old_v.reinit (this->dof_handler->n_dofs());
- transfer_old_solutions (old_u, old_v);
+ this->transfer_old_solutions (old_u, old_v);
};
-assemble_vectors (right_hand_side1, right_hand_side2);
+ assemble_vectors (right_hand_side1, right_hand_side2);
UserMatrix system_matrix (this->system_sparsity, this->parameters.preconditioning);
system_matrix.copy_from (this->mass_matrix);
right_hand_side1);
};
- const unsigned int solver_steps1 = solve (system_matrix, this->u, right_hand_side1);
+ const unsigned int solver_steps1 = this->solve (system_matrix, this->u, right_hand_side1);
system_matrix.copy_from (this->mass_matrix);
this->constraints.condense (static_cast<SparseMatrix<double>&>(system_matrix));
this->v -= old_v;
};
- const unsigned int solver_steps2 = solve (system_matrix, this->v, right_hand_side2);
+ const unsigned int solver_steps2 = this->solve (system_matrix, this->v, right_hand_side2);
this->statistic_data = typename TimeStep_Wave<dim>::StatisticData (this->tria->n_active_cells(),
this->dof_handler->n_dofs(),
solver_steps1,
solver_steps2,
- compute_energy ());
+ this->compute_energy ());
deallog << "." << std::endl;
}
{
Assert (this->system_sparsity.empty(), ExcInternalError());
- create_matrices ();
+ this->create_matrices ();
};
this->sweep_info->get_timers().primal_problem.stop();
}
= static_cast<const TimeStepBase_Wave<dim>*>(this->previous_timestep)->get_timestep_primal();
const unsigned int dofs_per_cell = this->fe.dofs_per_cell;
- const double time_step = get_backward_timestep();
+ const double time_step = this->get_backward_timestep();
if (!old_cell->has_children() && !new_cell->has_children())
{
= static_cast<const TimeStepBase_Wave<dim>*>(this->previous_timestep)->get_timestep_primal();
const unsigned int dofs_per_cell = this->fe.dofs_per_cell;
- const double time_step = get_backward_timestep();
+ const double time_step = this->get_backward_timestep();
FullMatrix<double> cell_matrix (dofs_per_cell, dofs_per_cell);
unsigned int level_difference = 1;
const unsigned int dofs_per_cell = this->fe.dofs_per_cell;
- const double time_step = get_backward_timestep();
+ const double time_step = this->get_backward_timestep();
FullMatrix<double> cell_matrix(dofs_per_cell, dofs_per_cell);
Vector<double> local_old_dof_values_u (dofs_per_cell);