From ddec5491cc56b09fcab51b5a775ec32af5cde5be Mon Sep 17 00:00:00 2001 From: bangerth Date: Sun, 21 Sep 2008 21:39:19 +0000 Subject: [PATCH] Remove the pressure boundary values -- they are only ever zero, and don't provide much insight anyway. git-svn-id: https://svn.dealii.org/trunk@16886 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/examples/step-31/step-31.cc | 76 +++-------------------------- deal.II/examples/step-32/step-32.cc | 47 ------------------ 2 files changed, 8 insertions(+), 115 deletions(-) diff --git a/deal.II/examples/step-31/step-31.cc b/deal.II/examples/step-31/step-31.cc index f0ff402482..02c1f70567 100644 --- a/deal.II/examples/step-31/step-31.cc +++ b/deal.II/examples/step-31/step-31.cc @@ -90,13 +90,14 @@ using namespace dealii; // step-22. Regarding the details, though, // there are some differences. - // The first - // thing is that we don't set any boundary - // conditions on the velocity, as is - // explained in the introduction. So - // what is left are two conditions for - // pressure p and temperature - // T. + // The first thing is that we don't set any + // non-homogenous boundary conditions on the + // velocity, since as is explained in the + // introduction we will use no-flux + // conditions + // $\mathbf{n}\cdot\mathbf{u}=0$. So what is + // left are two conditions for pressure + // p and temperature T. // Secondly, we set an initial // condition for all problem variables, @@ -115,29 +116,6 @@ namespace EquationData const double Rayleigh_number = 10; - template - class PressureBoundaryValues : public Function - { - public: - PressureBoundaryValues () : Function(1) {} - - virtual double value (const Point &p, - const unsigned int component = 0) const; - }; - - - template - double - PressureBoundaryValues::value (const Point &/*p*/, - const unsigned int /*component*/) const - { - return 0; - } - - - - - // @sect4{Initial values} template class TemperatureInitialValues : public Function @@ -1278,7 +1256,6 @@ void BoussinesqFlowProblem::assemble_stokes_system () // individual blocks (velocity, // pressure, temperature) from // the total FE system. - const EquationData::PressureBoundaryValues pressure_boundary_values; std::vector boundary_values (n_face_q_points); std::vector old_temperature_values(n_q_points); @@ -1377,43 +1354,6 @@ void BoussinesqFlowProblem::assemble_stokes_system () stokes_fe_values.JxW(q); } - - // Next follows the assembly - // of the face terms, result - // from Neumann boundary - // conditions. Since these - // terms only enter the right - // hand side vector and not - // the matrix, there is no - // substantial benefit from - // extracting the data - // before using it, so - // we remain in the lines - // of step-20 at this point. - for (unsigned int face_no=0; - face_no::faces_per_cell; - ++face_no) - if (cell->at_boundary(face_no)) - { - stokes_fe_face_values.reinit (cell, face_no); - - pressure_boundary_values - .value_list (stokes_fe_face_values.get_quadrature_points(), - boundary_values); - - for (unsigned int q=0; q - phi_i_u = stokes_fe_face_values[velocities].value (i, q); - - local_rhs(i) += -(phi_i_u * - stokes_fe_face_values.normal_vector(q) * - boundary_values[q] * - stokes_fe_face_values.JxW(q)); - } - } - // The last step in the loop // over all cells is to // enter the local contributions diff --git a/deal.II/examples/step-32/step-32.cc b/deal.II/examples/step-32/step-32.cc index 1fd4eeada7..83cdeb17dd 100644 --- a/deal.II/examples/step-32/step-32.cc +++ b/deal.II/examples/step-32/step-32.cc @@ -85,29 +85,6 @@ namespace EquationData const double Rayleigh_number = 10; - template - class PressureBoundaryValues : public Function - { - public: - PressureBoundaryValues () : Function(1) {} - - virtual double value (const Point &p, - const unsigned int component = 0) const; - }; - - - template - double - PressureBoundaryValues::value (const Point &/*p*/, - const unsigned int /*component*/) const - { - return 0; - } - - - - - // @sect4{Initial values} template class TemperatureInitialValues : public Function @@ -856,7 +833,6 @@ void BoussinesqFlowProblem::assemble_stokes_system () std::vector local_dof_indices (dofs_per_cell); - const EquationData::PressureBoundaryValues pressure_boundary_values; std::vector boundary_values (n_face_q_points); std::vector old_temperature_values(n_q_points); @@ -919,29 +895,6 @@ void BoussinesqFlowProblem::assemble_stokes_system () gravity * phi_u[i] * old_temperature)* stokes_fe_values.JxW(q); } - for (unsigned int face_no=0; - face_no::faces_per_cell; - ++face_no) - if (cell->at_boundary(face_no)) - { - stokes_fe_face_values.reinit (cell, face_no); - - pressure_boundary_values - .value_list (stokes_fe_face_values.get_quadrature_points(), - boundary_values); - - for (unsigned int q=0; q - phi_i_u = stokes_fe_face_values[velocities].value (i, q); - - local_rhs(i) += -(phi_i_u * - stokes_fe_face_values.normal_vector(q) * - boundary_values[q] * - stokes_fe_face_values.JxW(q)); - } - } cell->get_dof_indices (local_dof_indices); -- 2.39.5