<commit_before>/* $Id$ */
/* Author: <NAME>, Texas A&M University, 2008 */

/*    $Id$       */
/*    Version: $Name$                                          */
/*                                                                */
/*    Copyright (C) 2008 by the deal.II authors */
/*                                                                */
/*    This file is subject to QPL and may not be  distributed     */
/*    without copyright and license information. Please refer     */
/*    to the file deal.II/doc/license.html for the  text  and     */
/*    further information on this license.                        */


				 // @sect3{Include files}
                        
				 // As usual, we start by including 
				 // some well-known files:
#include <base/quadrature_lib.h>
#include <base/logstream.h>
#include <base/function.h>
#include <base/utilities.h>

#include <lac/block_vector.h>
#include <lac/full_matrix.h>
#include <lac/block_sparse_matrix.h>
#include <lac/solver_cg.h>
#include <lac/precondition.h>

#include <grid/tria.h>
#include <grid/grid_generator.h>
#include <grid/tria_accessor.h>
#include <grid/tria_iterator.h>
#include <grid/tria_boundary_lib.h>
#include <grid/grid_tools.h>
#include <grid/grid_refinement.h>

#include <dofs/dof_handler.h>
#include <dofs/dof_renumbering.h>
#include <dofs/dof_accessor.h>
#include <dofs/dof_tools.h>
#include <dofs/dof_constraints.h>

#include <fe/fe_q.h>
#include <fe/fe_system.h>
#include <fe/fe_values.h>
#include <fe/mapping_q1.h>

#include <numerics/vectors.h>
#include <numerics/matrices.h>
#include <numerics/data_out.h>
#include <numerics/error_estimator.h>

				 // Then we need to include the header file
				 // for the sparse direct solver UMFPACK:
#include <lac/sparse_direct.h>

				 // This includes the libary for the 
				 // incomplete LU factorization that will 
				 // be used as a preconditioner in 3D:
#include <lac/sparse_ilu.h>

				 // This is C++:
#include <fstream>
#include <sstream>

				 // As in all programs, the namespace dealii
				 // in included:
using namespace dealii;

				 // @sect3{Defining the inner preconditioner type}
                        
				 // As explained in the introduction, we are
				 // going to use different preconditioners for
				 // two and three space dimensions,
				 // respectively. We differentiate between
				 // them by the use of the spatial dimension
				 // as a template parameter. See step-4 for
				 // details on templates.  We are not going to
				 // create any preconditioner object here, all
				 // we do is to create class that holds a
				 // local typedef determining the
				 // preconditioner class so we can write our
				 // program in a dimension-independent way.
template <int dim>
struct InnerPreconditioner;

				 // In 2D, we are going to use a sparse direct
				 // solver as preconditioner:
template <>
struct InnerPreconditioner<2> 
{
    typedef SparseDirectUMFPACK type;
};

				 // And the ILU preconditioning in 3D, called
				 // by SparseILU:
template <>
struct InnerPreconditioner<3> 
{
    typedef SparseILU<double> type;
};


				 // @sect3{The <code>StokesProblem</code> class template}
                    
				 // This is an adaptation of step-20, so the
				 // main class and the data types are the same
				 // as used there. In this example we also use
				 // adaptive grid refinement, which is handled
				 // in complete analogy to step-6:
template <int dim>
class StokesProblem 
{
  public:
    StokesProblem (const unsigned int degree);
    void run ();
    
  private:
    void setup_dofs ();
    void assemble_system ();
    void solve ();
    void output_results (const unsigned int refinement_cycle) const;
    void refine_mesh ();
    
    const unsigned int   degree;
    
    Triangulation<dim>   triangulation;
    FESystem<dim>        fe;
    DoFHandler<dim>      dof_handler;

    ConstraintMatrix     hanging_node_constraints;
    
    BlockSparsityPattern      sparsity_pattern;
    BlockSparseMatrix<double> system_matrix;

    BlockVector<double> solution;
    BlockVector<double> system_rhs;

				     // This one is new: We shall use a
				     // so-called shared pointer structure to
				     // access the preconditioner. Shared
				     // pointers are essentially just a
				     // convenient form of pointers. Several
				     // shared pointers can point to the same
				     // object (just like regular pointers),
				     // but when the last shared pointer
				     // object to point to a preconditioner
				     // object is deleted (for example if a
				     // shared pointer object goes out of
				     // scope, if the class of which it is a
				     // member is destroyed, or if the pointer
				     // is assigned a different preconditioner
				     // object) then the preconditioner object
				     // pointed to is also destroyed. This
				     // ensures that we don't have to manually
				     // track in how many places a
				     // preconditioner object is still
				     // referenced, it can never create a
				     // memory leak, and can never produce a
				     // dangling pointer to an already
				     // destroyed object:
    boost::shared_ptr<typename InnerPreconditioner<dim>::type> A_preconditioner;
};

				 // @sect3{%Boundary values and right hand side}

				 // As in step-20 and most other
				 // example programs, the next task is
				 // to define the data for the PDE:
				 // For the Stokes problem, we are
				 // going to use natural boundary
				 // values on parts of the boundary
				 // (i.e. homogenous Neumann-type) for
				 // which we won't have to do anything
				 // special (the homogeneity implies
				 // that the corresponding terms in
				 // the weak form are simply zero),
				 // and boundary conditions on the
				 // velocity (Dirichlet-type) on the
				 // rest of the boundary, as described
				 // in the introduction.
				 //
				 // In order to enforce the Dirichlet
				 // boundary values on the velocity,
				 // we will use the
				 // VectorTools::interpolate_boundary_values
				 // function as usual which requires
				 // us to write a function object with
				 // as many components as the finite
				 // element has. In other words, we
				 // have to define the function on the
				 // $(u,p)$-space, but we are going to
				 // filter out the pressure component
				 // when interpolating the boundary
				 // values.
                    <commit_msg>Removed % char at one point.<commit_after>/* $Id$ */
/* Author: <NAME>, Texas A&M University, 2008 */

/*    $Id$       */
/*    Version: $Name$                                          */
/*                                                                */
/*    Copyright (C) 2008 by the deal.II authors */
/*                                                                */
/*    This file is subject to QPL and may not be  distributed     */
/*    without copyright and license information. Please refer     */
/*    to the file deal.II/doc/license.html for the  text  and     */
/*    further information on this license.                        */


				 // @sect3{Include files}
                        
				 // As usual, we start by including 
				 // some well-known files:
#include <base/quadrature_lib.h>
#include <base/logstream.h>
#include <base/function.h>
#include <base/utilities.h>

#include <lac/block_vector.h>
#include <lac/full_matrix.h>
#include <lac/block_sparse_matrix.h>
#include <lac/solver_cg.h>
#include <lac/precondition.h>

#include <grid/tria.h>
#include <grid/grid_generator.h>
#include <grid/tria_accessor.h>
#include <grid/tria_iterator.h>
#include <grid/tria_boundary_lib.h>
#include <grid/grid_tools.h>
#include <grid/grid_refinement.h>

#include <dofs/dof_handler.h>
#include <dofs/dof_renumbering.h>
#include <dofs/dof_accessor.h>
#include <dofs/dof_tools.h>
#include <dofs/dof_constraints.h>

#include <fe/fe_q.h>
#include <fe/fe_system.h>
#include <fe/fe_values.h>
#include <fe/mapping_q1.h>

#include <numerics/vectors.h>
#include <numerics/matrices.h>
#include <numerics/data_out.h>
#include <numerics/error_estimator.h>

				 // Then we need to include the header file
				 // for the sparse direct solver UMFPACK:
#include <lac/sparse_direct.h>

				 // This includes the libary for the 
				 // incomplete LU factorization that will 
				 // be used as a preconditioner in 3D:
#include <lac/sparse_ilu.h>

				 // This is C++:
#include <fstream>
#include <sstream>

				 // As in all programs, the namespace dealii
				 // in included:
using namespace dealii;

				 // @sect3{Defining the inner preconditioner type}
                        
				 // As explained in the introduction, we are
				 // going to use different preconditioners for
				 // two and three space dimensions,
				 // respectively. We differentiate between
				 // them by the use of the spatial dimension
				 // as a template parameter. See step-4 for
				 // details on templates.  We are not going to
				 // create any preconditioner object here, all
				 // we do is to create class that holds a
				 // local typedef determining the
				 // preconditioner class so we can write our
				 // program in a dimension-independent way.
template <int dim>
struct InnerPreconditioner;

				 // In 2D, we are going to use a sparse direct
				 // solver as preconditioner:
template <>
struct InnerPreconditioner<2> 
{
    typedef SparseDirectUMFPACK type;
};

				 // And the ILU preconditioning in 3D, called
				 // by SparseILU:
template <>
struct InnerPreconditioner<3> 
{
    typedef SparseILU<double> type;
};


				 // @sect3{The <code>StokesProblem</code> class template}
                    
				 // This is an adaptation of step-20, so the
				 // main class and the data types are the same
				 // as used there. In this example we also use
				 // adaptive grid refinement, which is handled
				 // in complete analogy to step-6:
template <int dim>
class StokesProblem 
{
  public:
    StokesProblem (const unsigned int degree);
    void run ();
    
  private:
    void setup_dofs ();
    void assemble_system ();
    void solve ();
    void output_results (const unsigned int refinement_cycle) const;
    void refine_mesh ();
    
    const unsigned int   degree;
    
    Triangulation<dim>   triangulation;
    FESystem<dim>        fe;
    DoFHandler<dim>      dof_handler;

    ConstraintMatrix     hanging_node_constraints;
    
    BlockSparsityPattern      sparsity_pattern;
    BlockSparseMatrix<double> system_matrix;

    BlockVector<double> solution;
    BlockVector<double> system_rhs;

				     // This one is new: We shall use a
				     // so-called shared pointer structure to
				     // access the preconditioner. Shared
				     // pointers are essentially just a
				     // convenient form of pointers. Several
				     // shared pointers can point to the same
				     // object (just like regular pointers),
				     // but when the last shared pointer
				     // object to point to a preconditioner
				     // object is deleted (for example if a
				     // shared pointer object goes out of
				     // scope, if the class of which it is a
				     // member is destroyed, or if the pointer
				     // is assigned a different preconditioner
				     // object) then the preconditioner object
				     // pointed to is also destroyed. This
				     // ensures that we don't have to manually
				     // track in how many places a
				     // preconditioner object is still
				     // referenced, it can never create a
				     // memory leak, and can never produce a
				     // dangling pointer to an already
				     // destroyed object:
    boost::shared_ptr<typename InnerPreconditioner<dim>::type> A_preconditioner;
};

				 // @sect3{Boundary values and right hand side}

				 // As in step-20 and most other
				 // example programs, the next task is
				 // to define the data for the PDE:
				 // For the Stokes problem, we are
				 // going to use natural boundary
				 // values on parts of the boundary
				 // (i.e. homogenous Neumann-type) for
				 // which we won't have to do anything
				 // special (the homogeneity implies
				 // that the corresponding terms in
				 // the weak form are simply zero),
				 // and boundary conditions on the
				 // velocity (Dirichlet-type) on the
				 // rest of the boundary, as described
				 // in the introduction.
				 //
				 // In order to enforce the Dirichlet
				 // boundary values on the velocity,
				 // we will use the
				 // VectorTools::interpolate_boundary_values
				 // function as usual which requires
				 // us to write a function object with
				 // as many components as the finite
				 // element has. In other words, we
				 // have to define the function on the
				 // $(u,p)$-space, but we are going to
				 // filter out the pressure component
				 // when interpolating the boundary
				 // values.
                    
{'fraction_non_alphanumeric': 0.0896343352959405, 'fraction_numerical': 0.0037186241090796405, 'mean_word_length': 2.094941261088468, 'pattern_counts': {'":': 0, '<': 105, '<?xml version=': 0, '>': 105, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 10, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '6941554', 'n_tokens_mistral': 4274, 'n_tokens_neox': 3400, 'n_words': 1340}