From: bangerth Date: Tue, 13 May 2008 02:52:08 +0000 (+0000) Subject: Simplify a bit. X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=aa1e48a93959b872e732b5c02288d5472e4b9dc0;p=dealii-svn.git Simplify a bit. git-svn-id: https://svn.dealii.org/trunk@16087 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/examples/step-33/step-33.cc b/deal.II/examples/step-33/step-33.cc index e02d932901..4878af89b8 100644 --- a/deal.II/examples/step-33/step-33.cc +++ b/deal.II/examples/step-33/step-33.cc @@ -306,7 +306,6 @@ class ConsLaw void load_parameters(const char *); private: - void build_fe(); void setup_system (); void initialize_system (); void assemble_system (double &res_norm); @@ -323,7 +322,7 @@ class ConsLaw const MappingQ1 mapping; - FESystem *fe_ptr; + FESystem fe; DoFHandler dof_handler; @@ -348,10 +347,8 @@ class ConsLaw public: - void assemble_cell_term(const FEValues& fe_v, - std::vector &dofs, - unsigned int cell_no - ); + void assemble_cell_term (const FEValues &fe_v, + const std::vector &dofs); void assemble_face_term( int face_no, @@ -490,11 +487,8 @@ void ConsLaw::initialize() { // to the right hand side, and adding in the Jacobian // contributions. template -void ConsLaw::assemble_cell_term( - const FEValues &fe_v, - std::vector &dofs, - unsigned int /*cell_no*/ -) +void ConsLaw::assemble_cell_term (const FEValues &fe_v, + const std::vector &dofs) { unsigned int dofs_per_cell = fe_v.dofs_per_cell; unsigned int n_q_points = fe_v.n_quadrature_points; @@ -650,7 +644,7 @@ void ConsLaw::assemble_cell_term( Matrix->SumIntoGlobalValues(dofs[i], dofs_per_cell, values, - reinterpret_cast(&dofs[0])); + reinterpret_cast(const_cast(&dofs[0]))); // Add minus the residual to the right hand side. right_hand_side(dofs[i]) -= F_i.val(); @@ -881,7 +875,6 @@ void ConsLaw::assemble_face_term( template void ConsLaw::assemble_system (double &res_norm) { - FESystem &fe = *fe_ptr; const unsigned int dofs_per_cell = dof_handler.get_fe().dofs_per_cell; // We track the dofs on this cell and (if necessary) @@ -976,8 +969,7 @@ void ConsLaw::assemble_system (double &res_norm) cell_diameter = cell->diameter(); assemble_cell_term(fe_v, - dofs, - cell_no); + dofs); // We use the DG style loop through faces // to determine if we need to apply a @@ -1120,7 +1112,7 @@ template ConsLaw::ConsLaw () : mapping (), - fe_ptr(NULL), + fe (FE_Q(1), EulerEquations::n_components), dof_handler (triangulation), quadrature (2), face_quadrature (2), @@ -1134,19 +1126,12 @@ ConsLaw::ConsLaw () theta(0.5) {} - // At one time this example could work for both DG and - // continuous finite elements. The choice was made here. -template -void ConsLaw::build_fe() { - fe_ptr = new FESystem(FE_Q(1), EulerEquations::n_components); -} // Bye bye Conservation law. template ConsLaw::~ConsLaw () { dof_handler.clear (); - delete fe_ptr; } // @sect3{Initialize System} @@ -1162,7 +1147,7 @@ void ConsLaw::initialize_system () // First we need to distribute the // DoFs. dof_handler.clear(); - dof_handler.distribute_dofs (*fe_ptr); + dof_handler.distribute_dofs (fe); // Size all of the fields. solution.reinit (dof_handler.n_dofs()); @@ -1192,7 +1177,7 @@ void ConsLaw::setup_system () sparsity_pattern.reinit (dof_handler.n_dofs(), dof_handler.n_dofs(), (GeometryInfo::faces_per_cell - *GeometryInfo::subfaces_per_face+1)*fe_ptr->dofs_per_cell); + *GeometryInfo::subfaces_per_face+1)*fe.dofs_per_cell); // Since the continuous sparsity pattern is // a subset of the DG one, and since we need @@ -1369,7 +1354,7 @@ void ConsLaw::postprocess() { QGauss quadrature_formula(4); - const std::vector > &us = fe_ptr->base_element(0).get_unit_support_points(); + const std::vector > &us = fe.base_element(0).get_unit_support_points(); Quadrature unit_support(us); @@ -1378,10 +1363,10 @@ void ConsLaw::postprocess() { int n_uq_points = unit_support.n_quadrature_points; FEValues fe_v ( - mapping, *fe_ptr, quadrature_formula, update_flags); + mapping, fe, quadrature_formula, update_flags); FEValues fe_v_unit ( - mapping, *fe_ptr, unit_support, update_flags1); + mapping, fe, unit_support, update_flags1); std::vector > U(n_uq_points, Vector(EulerEquations::n_components)); @@ -1451,7 +1436,7 @@ void ConsLaw::estimate() { FEValues fe_v ( - mapping, *fe_ptr, quadrature_formula, update_flags); + mapping, fe, quadrature_formula, update_flags); std::vector > U(n_q_points, Vector(EulerEquations::n_components)); @@ -1522,7 +1507,7 @@ void ConsLaw::refine_grid () triangulation.execute_coarsening_and_refinement (); dof_handler.clear(); - dof_handler.distribute_dofs (*fe_ptr); + dof_handler.distribute_dofs (fe); { Vector new_solution(1); @@ -1926,8 +1911,6 @@ void ConsLaw::run () grid_in.read_ucd(input_file); input_file.close(); - build_fe(); - unsigned int nstep = 0; // Initialize fields and matrices.