/* $Id$ */
-/* Author: Guido Kanschat, University of Heidelberg, 2003 */
+/* Author: Guido Kanschat, University of Heidelberg, 2003 */
+/* Baerbel Janssen, University of Heidelberg, 2010 */
+/* Wolfgang Bangerth, Texas A&M University, 2010 */
/* $Id$ */
/* */
-/* Copyright (C) 2003, 2004, 2005, 2006, 2007, 2008, 2009 by the deal.II authors */
+/* Copyright (C) 2003, 2004, 2006, 2007, 2008, 2009, 2010 by the deal.II authors */
/* */
/* This file is subject to QPL and may not be distributed */
/* without copyright and license information. Please refer */
#include <base/quadrature_lib.h>
#include <base/function.h>
#include <base/logstream.h>
+#include <base/utilities.h>
+
+#include <lac/constraint_matrix.h>
#include <lac/vector.h>
#include <lac/full_matrix.h>
#include <lac/sparse_matrix.h>
#include <lac/solver_cg.h>
#include <lac/precondition.h>
+
#include <grid/tria.h>
#include <grid/tria_accessor.h>
#include <grid/tria_iterator.h>
-#include <grid/tria_boundary_lib.h>
#include <grid/grid_generator.h>
-#include <dofs/dof_handler.h>
+
#include <dofs/dof_accessor.h>
#include <dofs/dof_tools.h>
+
#include <fe/fe_q.h>
#include <fe/fe_values.h>
+
#include <numerics/vectors.h>
-#include <numerics/matrices.h>
#include <numerics/data_out.h>
- // These are the new include files
- // required for multi-level methods.
- // First, the file defining the
- // multigrid method itself.
+
+//These are the same include files
+//as in step-16 necessary for the
+//multi-level methods
#include <multigrid/multigrid.h>
- // The DoFHandler is replaced by an
- // MGDoFHandler which is defined
- // here.
#include <multigrid/mg_dof_handler.h>
#include <multigrid/mg_dof_accessor.h>
-
- // Then, we need some pre-made
- // transfer routines between grids.
#include <multigrid/mg_transfer.h>
- // And a file in which equivalents to the
- // DoFTools class are declared:
#include <multigrid/mg_tools.h>
-
#include <multigrid/mg_coarse.h>
#include <multigrid/mg_smoother.h>
#include <multigrid/mg_matrix.h>
// previous programs:
using namespace dealii;
- // This class is based on the same
- // class in step-5. Remark that we
- // replaced the DoFHandler by
- // MGDoFHandler. since this inherits
- // from DoFHandler, the new object
- // incorporates the old functionality
- // plus the new functions for degrees
- // of freedom on different
- // levels. Furthermore, we added
- // MGLevelObjects for sparsity
- // patterns and matrices.
+
+//This class is basically the same
+//class as in step-16. The only
+//difference is that here we solve Laplace's
+//problem on an adaptively refined grid.
template <int dim>
class LaplaceProblem
{
public:
- LaplaceProblem ();
+ LaplaceProblem (const unsigned int deg);
void run ();
-
+
private:
void setup_system ();
void assemble_system ();
- // We add this function for
- // assembling the multilevel
- // matrices.
void assemble_multigrid ();
void solve ();
+ void refine_local ();
+ void test_get_coarse_cell ();
void output_results (const unsigned int cycle) const;
Triangulation<dim> triangulation;
SparsityPattern sparsity_pattern;
SparseMatrix<double> system_matrix;
- // Here are the new objects for
- // handling level matrices: sparsity
- // patterns and matrices. We use number
- // type float to save memory. It's only
- // a preconditioner!
- MGLevelObject<SparsityPattern> mg_sparsity;
- MGLevelObject<SparseMatrix<float> > mg_matrices;
-
+ //This object holds the information f
+ //or the hanging nodes.
+ ConstraintMatrix constraints;
+
+ MGLevelObject<SparsityPattern> mg_sparsity;
+ MGLevelObject<SparseMatrix<double> > mg_matrices;
+
+ /* The matrices at the interface
+ * between two refinement levels,
+ * coupling coarse to fine.*/
+ MGLevelObject<SparseMatrix<double> > mg_interface_matrices_up;
+
Vector<double> solution;
Vector<double> system_rhs;
-};
+ const unsigned int degree;
+};
- // This function is as before.
template <int dim>
-LaplaceProblem<dim>::LaplaceProblem () :
- fe (1),
- mg_dof_handler (triangulation)
+LaplaceProblem<dim>::LaplaceProblem (const unsigned int deg) :
+ triangulation (Triangulation<dim>::limit_level_difference_at_vertices),
+ fe (deg),
+ mg_dof_handler (triangulation),
+ degree(deg)
{}
-
- // This is the function of step-5
- // augmented by the setup of the
- // multi-level matrices in the end.
template <int dim>
void LaplaceProblem<dim>::setup_system ()
{
// multilevel structure
deallog << "Number of degrees of freedom: "
<< mg_dof_handler.n_dofs();
+
for (unsigned int l=0;l<triangulation.n_levels();++l)
deallog << " " << 'L' << l << ": "
<< mg_dof_handler.n_dofs(l);
deallog << std::endl;
-
+
sparsity_pattern.reinit (mg_dof_handler.n_dofs(),
mg_dof_handler.n_dofs(),
mg_dof_handler.max_couplings_between_dofs());
- DoFTools::make_sparsity_pattern (static_cast<const DoFHandler<dim>&>(mg_dof_handler),
- sparsity_pattern);
- sparsity_pattern.compress();
-
- system_matrix.reinit (sparsity_pattern);
+ DoFTools::make_sparsity_pattern (
+ static_cast<const DoFHandler<dim>&>(mg_dof_handler),
+ sparsity_pattern);
solution.reinit (mg_dof_handler.n_dofs());
system_rhs.reinit (mg_dof_handler.n_dofs());
+ constraints.clear ();
+ DoFTools::make_hanging_node_constraints (mg_dof_handler, constraints);
+ VectorTools::interpolate_boundary_values (mg_dof_handler,
+ 0,
+ ZeroFunction<dim>(),
+ constraints);
+ constraints.close ();
+ constraints.condense (sparsity_pattern);
+ sparsity_pattern.compress();
+ system_matrix.reinit (sparsity_pattern);
+
// The multi-level objects are
// resized to hold matrices for
// every level. The coarse level is
// SparsityPattern before it can be
// destroyed.
const unsigned int nlevels = triangulation.n_levels();
+
+ mg_interface_matrices_up.resize(0, nlevels-1);
+ mg_interface_matrices_up.clear ();
mg_matrices.resize(0, nlevels-1);
+ mg_matrices.clear ();
mg_sparsity.resize(0, nlevels-1);
// Now, we have to build a matrix
mg_dof_handler.n_dofs(level),
mg_dof_handler.max_couplings_between_dofs());
MGTools::make_sparsity_pattern (mg_dof_handler, mg_sparsity[level], level);
+ CompressedSparsityPattern ci_sparsity;
+ if(level>0)
+ {
+ ci_sparsity.reinit(mg_dof_handler.n_dofs(level),
+ mg_dof_handler.n_dofs(level));
+ MGTools::make_sparsity_pattern(mg_dof_handler, ci_sparsity, level);
+ }
+ }
+
+//And the same for the mg matrices
+//for the interface. Note that there
+//is no such interface on the coarsest level
+ for(unsigned int level=0; level<nlevels; ++level)
+ {
mg_sparsity[level].compress();
mg_matrices[level].reinit(mg_sparsity[level]);
+ mg_interface_matrices_up[level].reinit(mg_sparsity[level]);
}
}
template <int dim>
void LaplaceProblem<dim>::assemble_system ()
{
- QGauss<dim> quadrature_formula(2);
+ QGauss<dim> quadrature_formula(1+degree);
FEValues<dim> fe_values (fe, quadrature_formula,
update_values | update_gradients |
- update_quadrature_points | update_JxW_values);
+ update_quadrature_points | update_JxW_values);
const unsigned int dofs_per_cell = fe.dofs_per_cell;
const unsigned int n_q_points = quadrature_formula.size();
std::vector<unsigned int> local_dof_indices (dofs_per_cell);
- typename DoFHandler<dim>::active_cell_iterator cell = mg_dof_handler.begin_active(),
- endc = mg_dof_handler.end();
+ typename DoFHandler<dim>::active_cell_iterator
+ cell = mg_dof_handler.begin_active(),
+ endc = mg_dof_handler.end();
for (; cell!=endc; ++cell)
{
cell_matrix = 0;
cell_rhs = 0;
- // As before, we want the FEValues
- // object to compute the quantities
- // which we told him to compute in
- // the constructor using the update
- // flags. Then, we loop over all
- // quadrature points and the local
- // matrix rows and columns for
- // computing the element
- // contribution. This is the same as
- // in step-4. For the right hand
- // side, we use a constant value of
- // 1.
fe_values.reinit (cell);
-
for (unsigned int q_point=0; q_point<n_q_points; ++q_point)
for (unsigned int i=0; i<dofs_per_cell; ++i)
{
cell_rhs(i) += (fe_values.shape_value(i,q_point)
* 1.0 * fe_values.JxW(q_point));
- };
+ }
cell->get_dof_indices (local_dof_indices);
- for (unsigned int i=0; i<dofs_per_cell; ++i)
- {
- for (unsigned int j=0; j<dofs_per_cell; ++j)
- system_matrix.add (local_dof_indices[i],
- local_dof_indices[j],
- cell_matrix(i,j));
-
- system_rhs(local_dof_indices[i]) += cell_rhs(i);
- };
- };
-
- // The Dirichlet boundary
- // conditions on the finest level
- // are handled as usual.
- std::map<unsigned int,double> boundary_values;
-
- VectorTools::interpolate_boundary_values (mg_dof_handler,
- 0,
- ZeroFunction<dim>(),
- boundary_values);
-
- MatrixTools::apply_boundary_values (boundary_values,
- system_matrix,
- solution,
- system_rhs);
+ constraints.distribute_local_to_global (cell_matrix, cell_rhs,
+ local_dof_indices,
+ system_matrix, system_rhs);
+ }
}
// the same as above. Only the loop
// goes over all existing cells now
// and the results must be entered
- // into the correct matrix.
+ // into the correct matrix. Here comes
+ // the difference to global refinement
+ // into play. We have to fill the interface
+ // matrices correctly.
// Since we only do multi-level
// preconditioning, no right-hand
template <int dim>
void LaplaceProblem<dim>::assemble_multigrid ()
{
- QGauss<dim> quadrature_formula(2);
+ QGauss<dim> quadrature_formula(1+degree);
FEValues<dim> fe_values (fe, quadrature_formula,
update_values | update_gradients |
- update_quadrature_points | update_JxW_values);
+ update_quadrature_points | update_JxW_values);
- const unsigned int dofs_per_cell = fe.dofs_per_cell;
- const unsigned int n_q_points = quadrature_formula.size();
+ const unsigned int dofs_per_cell = fe.dofs_per_cell;
+ const unsigned int n_q_points = quadrature_formula.size();
FullMatrix<double> cell_matrix (dofs_per_cell, dofs_per_cell);
std::vector<unsigned int> local_dof_indices (dofs_per_cell);
+ std::vector<std::vector<bool> > interface_dofs;
+ std::vector<std::vector<bool> > boundary_interface_dofs;
+ for (unsigned int level = 0; level<triangulation.n_levels(); ++level)
+ {
+ std::vector<bool> tmp (mg_dof_handler.n_dofs(level));
+ interface_dofs.push_back (tmp);
+ boundary_interface_dofs.push_back (tmp);
+ }
+ MGTools::extract_inner_interface_dofs (mg_dof_handler,
+ interface_dofs,
+ boundary_interface_dofs);
+
+ typename FunctionMap<dim>::type dirichlet_boundary;
+ ZeroFunction<dim> homogeneous_dirichlet_bc (1);
+ dirichlet_boundary[0] = &homogeneous_dirichlet_bc;
+
+ std::vector<std::set<unsigned int> > boundary_indices(triangulation.n_levels());
+ MGTools::make_boundary_list (mg_dof_handler, dirichlet_boundary,
+ boundary_indices);
+
+ std::vector<ConstraintMatrix> boundary_constraints (triangulation.n_levels());
+ std::vector<ConstraintMatrix> boundary_interface_constraints (triangulation.n_levels());
+ for (unsigned int level=0; level<triangulation.n_levels(); ++level)
+ {
+ boundary_constraints[level].add_lines (interface_dofs[level]);
+ boundary_constraints[level].add_lines (boundary_indices[level]);
+ boundary_constraints[level].close ();
+
+ boundary_interface_constraints[level]
+ .add_lines (boundary_interface_dofs[level]);
+ boundary_interface_constraints[level].close ();
+ }
+
typename MGDoFHandler<dim>::cell_iterator cell = mg_dof_handler.begin(),
endc = mg_dof_handler.end();
+
for (; cell!=endc; ++cell)
{
- // Remember the level of the
- // current cell.
- const unsigned int level = cell->level();
cell_matrix = 0;
// Compute the values specified
* fe_values.JxW(q_point);
}
-
// Oops! This is a tiny
// difference easily
// forgotten. The indices we
// 'mg' entered into the
// function call.
cell->get_mg_dof_indices (local_dof_indices);
+
+ const unsigned int level = cell->level();
+ boundary_constraints[level]
+ .distribute_local_to_global (cell_matrix,
+ local_dof_indices,
+ mg_matrices[level]);
+
for (unsigned int i=0; i<dofs_per_cell; ++i)
- {
- for (unsigned int j=0; j<dofs_per_cell; ++j)
- // And now add everything
- // to the matrix on the
- // right level.
- mg_matrices[level].add (local_dof_indices[i],
- local_dof_indices[j],
- cell_matrix(i,j));
- }
+ for (unsigned int j=0; j<dofs_per_cell; ++j)
+ if( !(interface_dofs[level][local_dof_indices[i]]==true &&
+ interface_dofs[level][local_dof_indices[j]]==false))
+ cell_matrix(i,j) = 0;
+
+ boundary_interface_constraints[level]
+ .distribute_local_to_global (cell_matrix,
+ local_dof_indices,
+ mg_interface_matrices_up[level]);
}
-
- // Now we have to eliminate the
- // boundary nodes on all
- // levels. This is done exactly in
- // the same fashion as on the
- // finest level. Therefore, we need
- // a function map first. On the
- // other hand, since we use
- // multigrid only as a
- // preconditioner, we always use
- // homogeneous boundary conditions
- // and the ZeroFunction will be
- // sufficient in all cases; using a
- // different function here would
- // not hurt on the other hand,
- // since the values are ignored.
-
- typename FunctionMap<dim>::type dirichlet_boundary;
- ZeroFunction<dim> homogeneous_dirichlet_bc (1);
- dirichlet_boundary[0] = &homogeneous_dirichlet_bc;
-
- // Next we generate the set of
- // dof indices to be eliminated on
- // each level.
- std::vector<std::set<unsigned int> > boundary_indices(triangulation.n_levels());
- MGTools::make_boundary_list (mg_dof_handler,
- dirichlet_boundary,
- boundary_indices);
-
- // And finally we eliminate these
- // degrees of freedom from every
- // matrix in the multilevel hierarchy.
- const unsigned int nlevels = triangulation.n_levels();
- for (unsigned int level=0;level<nlevels;++level)
- {
- MGTools::apply_boundary_values(boundary_indices[level],
- mg_matrices[level],
- true);
- }
}
-
-
- // The solution process again looks
- // mostly like in the previous
- // examples. However, we will now use
- // a preconditioned conjugate
- // gradient algorithm. It is not very
- // difficult to make this change:
template <int dim>
void LaplaceProblem<dim>::solve ()
{
// the transfer of functions
// between different grid
// levels.
- MGTransferPrebuilt<Vector<double> > mg_transfer;
+ MGTransferPrebuilt<Vector<double> > mg_transfer(constraints);
mg_transfer.build_matrices(mg_dof_handler);
// Next, we need a coarse grid
// direct solver, even if its
// implementation is not very
// clever.
- FullMatrix<float> coarse_matrix;
+ FullMatrix<double> coarse_matrix;
coarse_matrix.copy_from (mg_matrices[0]);
- MGCoarseGridHouseholder<float, Vector<double> > mg_coarse;
+ MGCoarseGridHouseholder<double, Vector<double> > mg_coarse;
mg_coarse.initialize(coarse_matrix);
// The final ingredient for the
// here. Names are getting quite
// long here, so we help with
// typedefs.
- typedef PreconditionSOR<SparseMatrix<float> > RELAXATION;
- MGSmootherRelaxation<SparseMatrix<float>, RELAXATION, Vector<double> >
+ typedef PreconditionSOR<SparseMatrix<double> > RELAXATION;
+// typedef PreconditionJacobi<SparseMatrix<double> > RELAXATION;
+// typedef SparseILU<double> RELAXATION;
+ MGSmootherRelaxation<SparseMatrix<double>, RELAXATION, Vector<double> >
mg_smoother(vector_memory);
// Initialize the smoother with our
- // level matrices and the (required)
+ // level matrices and the required,
// additional data for the
- // relaxation method with default
+ // relaxaton method with default
// values.
- RELAXATION::AdditionalData smoother_data;
+ RELAXATION::AdditionalData smoother_data;//(0, 9,false);
mg_smoother.initialize(mg_matrices, smoother_data);
-
+
// Do two smoothing steps per level
- mg_smoother.set_steps(2);
+ mg_smoother.set_steps(1);
// Since the SOR method is not
// symmetric, but we use conjugate
// gradient iteration below, here
// symmetric operator even for
// nonsymmetric smoothers.
mg_smoother.set_symmetric(true);
+ mg_smoother.set_variable(true);
+
// We must wrap our matrices in an
// object having the required
// multiplication functions.
- MGMatrix<SparseMatrix<float>, Vector<double> >
+ MGMatrix<SparseMatrix<double>, Vector<double> >
mg_matrix(&mg_matrices);
+ //do the same for the interface matrices
+ MGMatrix<SparseMatrix<double>, Vector<double> >
+ mg_interface_up(&mg_interface_matrices_up);
+ MGMatrix<SparseMatrix<double>, Vector<double> >
+ mg_interface_down(&mg_interface_matrices_up);
// Now, we are ready to set up the
// V-cycle operator and the
// multilevel preconditioner.
mg_transfer,
mg_smoother,
mg_smoother);
+ mg.set_edge_matrices(mg_interface_down, mg_interface_up);
+
PreconditionMG<dim, Vector<double>,
MGTransferPrebuilt<Vector<double> > >
- preconditioner(mg_dof_handler, mg, mg_transfer);
-
- // Finally, create the solver
- // object and solve the system
- SolverControl solver_control (1000, 1e-12);
- SolverCG<> cg (solver_control);
+ preconditioner(mg_dof_handler, mg, mg_transfer);
-
- cg.solve (system_matrix, solution, system_rhs,
- preconditioner);
+ // Finally, create the solver
+ // object and solve the system
+ReductionControl solver_control (100, 1.e-20, 1.e-10, true, true);
+SolverCG<> cg (solver_control);
+
+solution = 0;
+
+cg.solve (system_matrix, solution, system_rhs,
+ preconditioner);
+constraints.distribute (solution);
+
+std::cout << " " << solver_control.last_step()
+<< " CG iterations needed to obtain convergence."
+<< std::endl;
}
+template <int dim>
+void LaplaceProblem<dim>::refine_local ()
+{
+ bool cell_refined = false;
+ for (typename Triangulation<dim>::active_cell_iterator
+ cell = triangulation.begin_active();
+ cell != triangulation.end(); ++cell)
+ {
+ if(false)
+ {
+ for (unsigned int vertex=0;
+ vertex < GeometryInfo<dim>::vertices_per_cell;
+ ++vertex)
+ {
+ if(cell->vertex(vertex)[dim-1]==0 && cell->vertex(vertex)[0]==0 && cell->vertex(vertex)[dim-2]==0)
+ {
+ cell->set_refine_flag ();
+ cell_refined = true;
+ break;
+ }
+ }
+ }
+ else if(true) //Kreis
+ {
+ for (unsigned int vertex=0;
+ vertex < GeometryInfo<dim>::vertices_per_cell;
+ ++vertex)
+ {
+ const Point<dim> p = cell->vertex(vertex);
+ const Point<dim> origin = (dim == 2 ?
+ Point<dim>(0,0) :
+ Point<dim>(0,0,0));
+ const double dist = p.distance(origin);
+ if(dist<0.25/M_PI)
+ {
+ cell->set_refine_flag ();
+ cell_refined = true;
+ break;
+ }
+ }
+ }
+ else if(false) //linke Diagonale
+ {
+ for (unsigned int vertex=0;
+ vertex < GeometryInfo<dim>::vertices_per_cell;
+ ++vertex)
+ {
+ const Point<dim> p = cell->vertex(vertex);
+ if(p[0]==p[1])
+ {
+ cell->set_refine_flag ();
+ cell_refined = true;
+ break;
+ }
+ }
+ }
+ else if(false) //inneres Quadrat
+ {
+ for (unsigned int vertex=0;
+ vertex < GeometryInfo<dim>::vertices_per_cell;
+ ++vertex)
+ {
+ const Point<dim> p = cell->vertex(vertex);
+ const double dist = std::max(std::fabs(p[0]),std::fabs(p[1]));
+ if(dist<0.5)
+ {
+ cell->set_refine_flag ();
+ cell_refined = true;
+ break;
+ }
+ }
+ }
+ else if(false) //Raute
+ {
+ for (unsigned int vertex=0;
+ vertex < GeometryInfo<dim>::vertices_per_cell;
+ ++vertex)
+ {
+ const Point<dim> p = cell->vertex(vertex);
+ const double dist = std::fabs(p[0])+std::fabs(p[1]);
+ if(dist<0.5)
+ {
+ cell->set_refine_flag ();
+ cell_refined = true;
+ break;
+ }
+ }
+ }
+ else //erster Quadrant
+ {
+ const Point<dim> p = cell->center();
+ bool positive = p(0) > 0;
+ if (dim>1 && p(1) <= 0)
+ positive = false;
+ if (dim>2 && p(2) <= 0)
+ positive = false;
+ if (positive)
+ {
+ cell->set_refine_flag();
+ cell_refined = true;
+ }
+ }
+ }
+ //Wenn nichts verfeinert wurde bisher, global verfeinern!
+ if(!cell_refined)
+ for (typename Triangulation<dim>::active_cell_iterator
+ cell = triangulation.begin_active();
+ cell != triangulation.end(); ++cell)
+ cell->set_refine_flag();
+
+ triangulation.execute_coarsening_and_refinement ();
+}
- // Here is the data output, which is
- // a simplified version of step-5. We
- // do a standard gnuplot output for
- // each grid produced in the
- // refinement process.
template <int dim>
void LaplaceProblem<dim>::output_results (const unsigned int cycle) const
{
- // Construct and initialize a DataOut object
DataOut<dim> data_out;
data_out.attach_dof_handler (mg_dof_handler);
data_out.add_data_vector (solution, "solution");
data_out.build_patches ();
- // The following block generates
- // the file name and opens the
- // file:
std::ostringstream filename;
filename << "solution-"
<< cycle
- << ".gnuplot";
+ << ".vtk";
std::ofstream output (filename.str().c_str());
- data_out.write_gnuplot (output);
+ data_out.write_vtk (output);
}
-
-
template <int dim>
void LaplaceProblem<dim>::run ()
{
- for (unsigned int cycle=0; cycle<6; ++cycle)
+ for (unsigned int cycle=0; cycle<9; ++cycle)
{
- deallog << "Cycle " << cycle << std::endl;
-
if (cycle == 0)
{
- // Generate a simple hyperball grid.
- GridGenerator::hyper_ball(triangulation);
- static const HyperBallBoundary<dim> boundary;
- triangulation.set_boundary (0, boundary);
+ GridGenerator::hyper_cube(triangulation, -1, 1);
+ triangulation.refine_global (1);
}
- triangulation.refine_global (1);
+ refine_local ();
+
+ std::cout << "Cycle " << cycle
+ << " with " << triangulation.n_active_cells()
+ << " cells."
+ << std::endl;
+
setup_system ();
assemble_system ();
assemble_multigrid ();
};
}
-
- // The main function looks mostly
- // like the one in the previous
- // example, so we won't comment on it
- // further.
+
int main ()
{
- LaplaceProblem<2> laplace_problem_2d;
- laplace_problem_2d.run ();
-
+ try
+ {
+ deallog.depth_console (0);
+
+ LaplaceProblem<2> laplace_problem(1);
+ laplace_problem.run ();
+ }
+ catch (std::exception &exc)
+ {
+ std::cerr << std::endl << std::endl
+ << "----------------------------------------------------"
+ << std::endl;
+ std::cerr << "Exception on processing: " << std::endl
+ << exc.what() << std::endl
+ << "Aborting!" << std::endl
+ << "----------------------------------------------------"
+ << std::endl;
+
+ return 1;
+ }
+ catch (...)
+ {
+ std::cerr << std::endl << std::endl
+ << "----------------------------------------------------"
+ << std::endl;
+ std::cerr << "Unknown exception!" << std::endl
+ << "Aborting!" << std::endl
+ << "----------------------------------------------------"
+ << std::endl;
+ return 1;
+ }
+
return 0;
}