1 #include <config.h>
2
3 #include <sstream>
4 #include <string>
5
6 #include <dune/grid/io/file/dgfparser/dgfparser.hh>
7
8 #include <dune/fem/function/adaptivefunction.hh>
9 #include <dune/fem/function/common/gridfunctionadapter.hh>
10 #include <dune/fem/gridpart/leafgridpart.hh>
11 #include <dune/fem/io/file/vtkio.hh>
12 #include <dune/fem/misc/mpimanager.hh>
13 #include <dune/fem/space/lagrange.hh>
14 #include <dune/fem/space/discontinuousgalerkin.hh>
15 #include <dune/fem/space/combinedspace.hh>
16 #include <dune/fem/space/common/interpolate.hh>
17 #include <dune/fem/test/exactsolution.hh>
18
19
20 // dgfUnitCube
21 // -----------
22
dgfUnitCube(int dimWorld,int cells)23 inline static std::string dgfUnitCube ( int dimWorld, int cells )
24 {
25 std::string dgf = "DGF\nINTERVAL\n";
26 for( int i = 0; i < dimWorld; ++i )
27 dgf += " 0";
28 dgf += "\n";
29 for( int i = 0; i < dimWorld; ++i )
30 dgf += " 1";
31 dgf += "\n";
32 for( int i = 0; i < dimWorld; ++i )
33 dgf += (" " + std::to_string( cells ));
34 dgf += "\n#\n";
35 return dgf;
36 }
37
38
39
40 // FunctionSpace
41 // -------------
42
43 template< class GridPart, int dimRange >
44 using FunctionSpace = Dune::Fem::FunctionSpace< typename GridPart::ctype, double, GridPart::dimensionworld, dimRange >;
45
46
47
48 // main
49 // ----
50
main(int argc,char ** argv)51 int main ( int argc, char **argv )
52 {
53 Dune::Fem::MPIManager::initialize( argc, argv );
54
55 // construct unit cube
56 typedef typename Dune::GridSelector::GridType GridType;
57 std::istringstream dgf( dgfUnitCube( GridType::dimensionworld, 4 ) );
58 Dune::GridPtr< GridType > grid( dgf );
59
60 // create leaf grid part
61 typedef Dune::Fem::LeafGridPart< GridType > GridPartType;
62 GridPartType gridPart( *grid );
63
64 using FunctionSpaceType = FunctionSpace< GridPartType, 1 > ;
65
66 // construct FE-DG space
67 typedef Dune::Fem::LagrangeDiscreteFunctionSpace< FunctionSpaceType, GridPartType, 1 > ContinuousSpaceType;
68 typedef Dune::Fem::DiscontinuousGalerkinSpace< FunctionSpaceType, GridPartType, 0 > DGSpaceType;
69 typedef Dune::Fem::EnrichedDiscreteFunctionSpace< DGSpaceType, ContinuousSpaceType > DiscreteFunctionSpaceType;
70 DiscreteFunctionSpaceType dfSpace( gridPart );
71
72 // interpolate a function
73 Dune::Fem::ExactSolution< FunctionSpaceType > uExact;
74 const auto uGridExact = gridFunctionAdapter( "exact solution", uExact, gridPart, 3 );
75 Dune::Fem::AdaptiveDiscreteFunction< DiscreteFunctionSpaceType > u( "solution", dfSpace );
76
77 std::cout << "Grid elements " << grid->size( 0 ) << std::endl;
78 std::cout << "Grid vertices " << grid->size( GridPartType::dimension ) << std::endl;
79 std::cout << "DOFs " << dfSpace.size() << std::endl;
80 interpolate( uGridExact, u );
81
82 Dune::Fem::GlobalRefine::apply( *grid, 2 );
83
84 // output analytical function and interolation to vtk file
85 Dune::Fem::VTKIO< GridPartType > vtkIO( gridPart, Dune::VTK::nonconforming );
86 vtkIO.addVertexData( uGridExact );
87 vtkIO.addVertexData( u );
88 vtkIO.write( "test-tuplespace" );
89
90 return 0;
91 }
92