From d9480c73bd92a53f9ba28f5f62604678eb3c2f38 Mon Sep 17 00:00:00 2001 From: frohne Date: Wed, 19 Oct 2011 16:25:35 +0000 Subject: [PATCH] Using Trilinos vectors and sparse-matrices git-svn-id: https://svn.dealii.org/trunk@24644 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/examples/step-41/step-41.cc | 54 +++++++++++++++++++---------- 1 file changed, 36 insertions(+), 18 deletions(-) diff --git a/deal.II/examples/step-41/step-41.cc b/deal.II/examples/step-41/step-41.cc index 6ae88cbfd3..885346b30e 100644 --- a/deal.II/examples/step-41/step-41.cc +++ b/deal.II/examples/step-41/step-41.cc @@ -38,6 +38,11 @@ #include #include +#include +#include +#include +#include + #include #include #include @@ -89,7 +94,7 @@ class Step4 void assemble_system (); void projection_active_set (); void solve (); - void output_results (Vector vector_to_plot, const std::string& title) const; + void output_results (TrilinosWrappers::Vector vector_to_plot, const std::string& title) const; Triangulation triangulation; FE_Q fe; @@ -98,14 +103,14 @@ class Step4 ConstraintMatrix constraints; SparsityPattern sparsity_pattern; - SparseMatrix system_matrix; - SparseMatrix system_matrix_complete; + TrilinosWrappers::SparseMatrix system_matrix; + TrilinosWrappers::SparseMatrix system_matrix_complete; - Vector solution; - Vector system_rhs; - Vector system_rhs_complete; - Vector resid_vector; - Vector active_set; + TrilinosWrappers::Vector solution; + TrilinosWrappers::Vector system_rhs; + TrilinosWrappers::Vector system_rhs_complete; + TrilinosWrappers::Vector resid_vector; + TrilinosWrappers::Vector active_set; std::map boundary_values; }; @@ -233,7 +238,8 @@ template double RightHandSide::value (const Point &p, const unsigned int /*component*/) const { - double return_value = 0; +// double return_value = -2.0*p.square () - 2.0; + double return_value = -10; // for (unsigned int i=0; i double Obstacle::value (const Point &p, const unsigned int /*component*/) const { - return 2.0*p.square() - 0.5; +// return 2.0*p.square() - 0.5; + + double return_value = 0; + if (p (0) < -0.5) + return_value = -0.2; + else if (p (0) >= -0.5 && p (0) < 0.0) + return_value = -0.4; + else if (p (0) >= 0.0 && p (0) < 0.5) + return_value = -0.6; + else + return_value = -0.8; + + return return_value; } @@ -348,7 +366,7 @@ template void Step4::make_grid () { GridGenerator::hyper_cube (triangulation, -1, 1); - triangulation.refine_global (6); + triangulation.refine_global (7); std::cout << " Number of active cells: " << triangulation.n_active_cells() @@ -470,7 +488,7 @@ void Step4::assemble_system () const unsigned int n_q_points = quadrature_formula.size(); FullMatrix cell_matrix (dofs_per_cell, dofs_per_cell); - Vector cell_rhs (dofs_per_cell); + TrilinosWrappers::Vector cell_rhs (dofs_per_cell); std::vector local_dof_indices (dofs_per_cell); @@ -663,7 +681,7 @@ void Step4::projection_active_set () Point point (cell->vertex (v)[0], cell->vertex (v)[1]); double obstacle_value = obstacle.value (point); - if (solution (index_x) >= obstacle_value && resid_vector (index_x) <= 0) + if (solution (index_x) <= obstacle_value && resid_vector (index_x) >= 0) { constraints.add_line (index_x); constraints.set_inhomogeneity (index_x, obstacle_value); @@ -692,9 +710,9 @@ template void Step4::solve () { ReductionControl reduction_control (100, 1e-12, 1e-2); - SolverCG<> solver (reduction_control); - SolverBicgstab<> solver_bicgstab (reduction_control); - PreconditionSSOR > precondition; + SolverCG solver (reduction_control); +// SolverBicgstab<> solver_bicgstab (reduction_control); + TrilinosWrappers::PreconditionSSOR precondition; precondition.initialize (system_matrix, 1.2); solver.solve (system_matrix, solution, system_rhs, precondition); @@ -736,7 +754,7 @@ void Step4::solve () // than 2 or 3, but we neglect this here for // the sake of brevity). template -void Step4::output_results (Vector vector_to_plot, const std::string& title) const +void Step4::output_results (TrilinosWrappers::Vector vector_to_plot, const std::string& title) const { DataOut data_out; @@ -822,7 +840,7 @@ void Step4::run () system_matrix_complete.vmult_add (resid_vector, solution); for (unsigned int k = 0; k 0) resid_vector (k) = 0; std::ostringstream filename_residuum; -- 2.39.5