dof_application.h

Go to the documentation of this file.

//---------------------------------------------------------------------------
//    $Id: dof_application.h 1348 2013-08-16 00:45:07Z secanell $
//
//    Copyright (C) 2006, 2007, 2008, 2009 by Guido Kanschat
//
//    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.
//
//---------------------------------------------------------------------------

#ifndef __deal2__appframe__dof_application_h
#define __deal2__appframe__dof_application_h

#include <appframe/base.h>
#include <appframe/residual.h>
#include <appframe/mesh_worker.h>
#include <boost/shared_ptr.hpp>
#include <lac/block_indices.h>
#include <grid/grid_out.h>
#include <grid/grid_tools.h>
#include <lac/constraint_matrix.h>
#include <fe/fe_values.h>
#include <numerics/data_out.h>
#include <numerics/solution_transfer.h>

namespace dealii
{
    template <int, int> class Triangulation;
    template <int, int> class Mapping;
    template <int, int> class FiniteElement;
    template <int, int> class DoFHandler;
}


namespace AppFrame
{
    template<int dim> class DoFApplication;
    
    struct VectorSelector
    {
        std::string name;
        unsigned int block;
        unsigned int components;
        
        dealii::UpdateFlags flags;
    };
    
    template <int dim>
    class LocalIntegrator
    {
    public:
        typedef typename MeshWorker::InfoObjects::IntegrationInfo<dim, FEValuesBase<dim, dim> > CellInfo;
        typedef typename MeshWorker::InfoObjects::IntegrationInfo<dim, FEFaceValuesBase<dim, dim> > FaceInfo;
        
        LocalIntegrator();
        virtual ~LocalIntegrator();
        
        virtual void cell(FEVector& cell_vector, const CellInfo& cell) = 0;
        
        virtual void bdry(FEVector& face_vector, const FaceInfo& face) = 0;
        
        virtual void face(FEVector& face_vector1, FEVector& face_vector2,
                          const FaceInfo& face1, const FaceInfo& face2) = 0;
                          
        virtual void post_loop(FEVector& dst, const FEVectors& srcs);
                          
        Quadrature<dim> cell_quadrature;
                          
        Quadrature<dim-1> bdry_quadrature;
        
        Quadrature<dim-1> face_quadrature;
        
        void set_gauss_quadrature(unsigned int n_cell_points,
                                  unsigned int n_bdry_points,
                                  unsigned int n_face_points);
    };
    
    template <int dim>
    class LocalResidual : public MeshWorker::Assembler::ResidualLocalBlocksToGlobalBlocks<double>,
    public MeshWorker::WorkerObjects::IntegrationWorker<dim>
    {
    private:
        SmartPointer<DoFApplication<dim> > app;
    public:
        void initialize(DoFApplication<dim>* a,
                        const MeshWorker::BlockInfo& bi,
                        FEVectors& m)
        {
            app = a;
            MeshWorker::Assembler::ResidualLocalBlocksToGlobalBlocks<double>::initialize(bi, m);
        }
        
        void cell(const typename DoFApplication<dim>::CellInfo& cell);
        
        void bdry(const typename DoFApplication<dim>::FaceInfo& face);
        
        void face(const typename DoFApplication<dim>::FaceInfo& face1,
                  const typename DoFApplication<dim>::FaceInfo& face2);
    };
    
    template <int dim>
    class LocalEstimate : public MeshWorker::Assembler::CellsAndFaces<double>,
    public MeshWorker::WorkerObjects::IntegrationWorker<dim>
    {
    private:
        SmartPointer<DoFApplication<dim> > app;
    public:
        void initialize(DoFApplication<dim>* a,
                        const MeshWorker::BlockInfo&,
                        FEVectors& m)
        {
            app = a;
            MeshWorker::Assembler::CellsAndFaces<double>::initialize(m, true);
        }
        
        void cell(const typename DoFApplication<dim>::CellInfo& cell);
        
        void bdry(const typename DoFApplication<dim>::FaceInfo& face);
        
        void face(const typename DoFApplication<dim>::FaceInfo& face1,
                  const typename DoFApplication<dim>::FaceInfo& face2);
    };
    
//     /**
//      * Integrate the residual using
//      * DoFApplication::integrate_1form(). Optionally, this class can use
//      * its own LocalIntegrator or use the standard residual of the
//      * DoFApplication.
//      *
//      * @author Guido Kanschat, 2008
//      */
//     template <int dim>
//     class DoFResidual : public Residual
//     {
//     public:
//         /**
//          * Constructor initializing
//          * this object to compute the
//          * residual using
//          * <tt>application</tt> and the
//          * provided <tt>local_residual</tt>.
//          */
//         DoFResidual(const DoFApplication<dim>& application,
//                     const LocalIntegrator<dim>& local_residual);
//         
//         virtual double operator() (FEVector& dst, const FEVectors& src) const;
//         
//     private:
//         /**
//          * The DoFApplication class
//          * providing the function
//          * DoFApplication::integrate_1form().
//          */
//         SmartPointer<const DoFApplication<dim> > app;
//         
//         /**
//          * The object integrating the
//          * local residuals.
//          */
//         LocalIntegrator<dim> local_residual;
//     };
//     
  
    template <int dim>
    class DoFApplication
    :
    public ApplicationBase
    {
    public:
        typedef typename MeshWorker::InfoObjects::IntegrationInfo<dim, FEValuesBase<dim> > CellInfo;
        typedef typename MeshWorker::InfoObjects::IntegrationInfo<dim, FEFaceValuesBase<dim> > FaceInfo;
        DoFApplication(boost::shared_ptr<ApplicationData> data = boost::shared_ptr<ApplicationData>(),
                       bool multigrid = false);
        
        DoFApplication(bool multigrid);
        
        DoFApplication(DoFApplication<dim>&,
                       bool triangulation_only,
                       bool multigrid = false);
        
        ~DoFApplication();
        
        virtual void declare_parameters(ParameterHandler& param);
        
        void _initialize (ParameterHandler& param);
        
        void add_vector_for_transfer(FEVector*);
        void delete_vector_for_transfer();
        virtual void remesh_dofs();
        
        virtual void initialize (ParameterHandler& param);
        virtual void init_vector (FEVector& dst) const;
        
        virtual void remesh();
        virtual double evaluate(const FEVectors&);
        virtual double estimate(const FEVectors& src);
        
        virtual double residual(FEVector& dst, const FEVectors& src, bool apply_boundaries = true);
        
        virtual void residual_constraints(FEVector& dst) const;
        

        
        virtual void grid_out(const std::string& basename);
        
        virtual void data_out(const std::string& basename,
                              const FEVectors&   src);
        
        virtual void data_out(const std::string&              basename,
                              const FEVectors&                src,
                              const std::vector<std::string>& solution_names);
        
        virtual void data_out(const std::string&              basename,
                              
                              const FEVector&                 solution,
                              const std::vector<std::string>& solution_names,
                              
                              const FEVector&                 postprocessing       = FEVector(),
                              const std::vector<std::string>& postprocessing_names = std::vector<std::string>());
        
        std::vector< DataComponentInterpretation::DataComponentInterpretation > solution_interpretations;
        
        std::vector< DataComponentInterpretation::DataComponentInterpretation > postprocessing_interpretations;
        
        bool output_materials_and_levels;
        
        Vector<double> output_materials;
        
        Vector<double> output_levels;
        
               
        unsigned int memory_consumption() const;
        
        virtual double cell_estimate(const CellInfo&);
        
        virtual double bdry_estimate(const FaceInfo&);
        
        virtual double face_estimate(const FaceInfo&,
                                     const FaceInfo&);
        
        virtual void cell_residual(FEVector& cell_vector,
                                   const CellInfo& cell);
        
        virtual void bdry_residual(FEVector& face_vector,
                                   const FaceInfo& face);
        
        virtual void face_residual(FEVector& face_vector1,
                                   FEVector& face_vector2,
                                   const FaceInfo& face1,
                                   const FaceInfo& face2);
        template<class BOX, class WORKER>
        void integrate_1form(BOX& box, WORKER& local_forms, FEVector& dst, const FEVectors& src);
        
        void integrate_functional(LocalEstimate<dim>& local_forms,
                                  FEVectors& dst,
                                  const FEVectors& src);
        
        void store_triangulation(Triangulation<dim>& new_tr)
        {
            new_tr.copy_triangulation(*this->tr);
        }
        
        void transfer_solution_to_coarse_mesh(Triangulation<dim>& tr_coarse,
                                              AppFrame::FEVector& coarse_solution,
                                              AppFrame::FEVector& refined_solution);
        
        void read_init_solution(AppFrame::FEVector& dst, bool &good_solution) const;
        
        unsigned int verbosity;
        
    protected:
        void constrain_boundary(FEVector& v, bool homogeneous) const;
        
        std::map<unsigned int, double> boundary_constraints;
        
        ConstraintMatrix hanging_node_constraints; 
        
        boost::shared_ptr<Mapping<dim> > mapping;
        
        unsigned int mapping_degree;    
        
        boost::shared_ptr<Triangulation<dim> > tr;
        
        boost::shared_ptr<FiniteElement<dim> > element;
        
        boost::shared_ptr<DoFHandler<dim> > dof;
        
        boost::shared_ptr<MGDoFHandler<dim> > mg_dof;
        
        MeshWorker::BlockInfo block_info;
        
        //      Vector<double> cell_errors;
        Vector<float> cell_errors;
        Vector<float> face_errors;
        
        GridOut g_out;
        
        DataOut<dim, DoFHandler<dim> > d_out;
        
        std::vector<DataComponentInterpretation::DataComponentInterpretation>
        data_interpretation;
        
        std::string refinement;
        
        unsigned int initial_refinement;
        
        double refinement_threshold;
        
        double coarsening_threshold;
        
        bool sort_cuthill;
        
        Point<dim> sort_direction;
        
        std::vector<FEVector*> transfer_vectors;
        Quadrature<dim> quadrature_residual_cell;
        
        Quadrature<dim-1> quadrature_residual_bdry;
        
        Quadrature<dim-1> quadrature_residual_face;
        
        bool boundary_fluxes;
        
        bool interior_fluxes;
        
    public:
        
        template <class INFO, typename TYPE>
        void fill_local_data(const INFO& info,
                             const std::vector<VectorSelector>& data_vector,
                             std::vector<std::vector<std::vector<TYPE> > >& data) const;
                             
    private:
        virtual void sort_dofs();
        
        virtual void distribute_face_to_cell_errors();
        
        virtual double global_from_local_errors() const;
        friend class LocalResidual<dim>;
        friend class LocalEstimate<dim>;
    };   
}

#endif