plot_properties.hh

#ifndef PLOT_STRESSES_HH
#define PLOT_STRESSES_HH
 
#include "composite_properties_lop.hh"
 
namespace Dune{
    namespace Composites{
 
        template<class MODEL, typename V, typename GV1, typename GFS, typename MBE>
            void plotProperties(MODEL& model, const GV1& gv1, const GFS& gfs,  MBE& mbe){
                using Dune::PDELab::Backend::native;
 
                typedef Dune::PDELab::istl::VectorBackend<Dune::PDELab::istl::Blocking::none,1> Scalar_VectorBackend;
                typedef double RF;
                std::string vtk_output = config.get<std::string>("ModelName","TestModel");
                Dune::Timer timer;
 
                using GV = Dune::PDELab::AllEntitySet<GV1>; 
                auto gv = GV(gv1);
 
                // Initially we consider only a DG stress plot
 
                //Function spaces for Stresses
                typedef typename GV::Grid::ctype e_ctype;
                typedef Dune::PDELab::QkDGLocalFiniteElementMap<e_ctype,double,0,3> FEM_stress;
                FEM_stress fem_stress;
 
                typedef Dune::PDELab::GridFunctionSpace<GV, FEM_stress, Dune::PDELab::NoConstraints, Scalar_VectorBackend> SCALAR_P1GFS;
                SCALAR_P1GFS materialType(gv,fem_stress); materialType.name("materialType");
                SCALAR_P1GFS orientation(gv,fem_stress); orientation.name("orientation");
 
                typedef Dune::PDELab::istl::VectorBackend<Dune::PDELab::istl::Blocking::fixed,2> VectorBackend;
 
                typedef Dune::PDELab::CompositeGridFunctionSpace <VectorBackend,Dune::PDELab::EntityBlockedOrderingTag, 
                        SCALAR_P1GFS, SCALAR_P1GFS> P1GFS;
 
                P1GFS p1gfs(materialType,orientation);
 
                typedef Dune::PDELab::EmptyTransformation NoTrafo;
 
                Dune::PDELab::composite_properties<GV,MODEL> lopProperties(gv,model);
 
                typedef Dune::PDELab::GridOperator<P1GFS,P1GFS,Dune::PDELab::composite_properties<GV,MODEL>,MBE,RF,RF,RF,NoTrafo,NoTrafo> GO;
                GO go(p1gfs,p1gfs,lopProperties,mbe);
 
                typedef typename GO::Traits::Range V1;
                V1 properties(p1gfs,0.0), dummy(p1gfs,0.0); 
                go.residual(dummy,properties);
 
                Dune::SubsamplingVTKWriter<GV1> vtkwriter(gv1,0);
                Dune::PDELab::addSolutionToVTKWriter(vtkwriter,p1gfs,properties);
                vtkwriter.write(model.vtk_properties_output,Dune::VTK::appendedraw);
            }
 
    }
}
 
#endif