1 #ifndef DUNE_SPGRID_CAPABILITIES_HH
2 #define DUNE_SPGRID_CAPABILITIES_HH
3 
4 #if HAVE_MPI
5 #include <mpi.h>
6 #endif
7 
8 #include <dune/geometry/type.hh>
9 
10 #include <dune/grid/common/capabilities.hh>
11 
12 #include <dune/grid/extensions/superentityiterator.hh>
13 
14 #include <dune/grid/spgrid/declaration.hh>
15 
16 /** \file
17  *  \author Martin Nolte
18  *  \brief  capabilities for \ref Dune::SPGrid "SPGrid"
19  */
20 
21 namespace Dune
22 {
23 
24   // Capabilities
25   // ------------
26 
27   /** \brief namespace containing all capabilities */
28   namespace Capabilities
29   {
30 
31     /** \brief Do elements of a grid always have the same geometry type?
32      *
33      *  \tparam  Grid  grid for which the information is desired
34      */
35     template< class ct, int dim, template< int > class Ref, class Comm >
36     struct hasSingleGeometryType< SPGrid< ct, dim, Ref, Comm > >
37     {
38       /** \brief all elements in \ref Dune::SPGrid "SPGrid" have the same
39        *         geometry type */
40       static const bool v = true;
41       /** \brief \ref Dune::SPGrid "SPGrid" has only cube elements */
42       static const unsigned int topologyId = Impl::CubeTopology< dim >::type::id;
43     };
44 
45     /** \brief Is the grid Cartesian?
46      *
47      *  Cartesian grids satisfy the following properties:
48      *  - all geometries are affine
49      *  - The unit outer normal can be computed by the following code:
50      *  \code
51      *  FieldVector< ctype, dim > n( 0 );
52      *  n[ face / 2 ] = ctype( 2*(face % 2) - 1 );
53      *  \endcode
54      *  .
55      *
56      *  \tparam  Grid  grid for which the information is desired
57      */
58     template< class ct, int dim, template< int > class Ref, class Comm >
59     struct isCartesian< SPGrid< ct, dim, Ref, Comm > >
60     {
61       /** \brief \ref Dune::SPGrid "SPGrid" is a Cartesian grid */
62       static const bool v = true;
63     };
64 
65     /** \brief Does a grid implement entities of a codimension?
66      *
67      *  \tparam  Grid   grid for which the information is desired
68      *  \tparam  codim  codimension in question
69      */
70     template< class ct, int dim, template< int > class Ref, class Comm, int codim >
71     struct hasEntity< SPGrid< ct, dim, Ref, Comm >, codim >
72     {
73       /** \brief \ref Dune::SPGrid "SPGrid" implements entities for all
74        *         codimensions */
75       static const bool v = ((codim >= 0) && (codim <= dim));
76     };
77 
78 #if HAVE_MPI
79     /** \brief Can a parallel grid communicate on a given codimension?
80      *
81      *  \tparam  Grid   grid for which the information is desired
82      *  \tparam  codim  codimension in question
83      *
84      *  \note In order to communicate on a given codimension, the grid has to
85      *        implement entities for that codimension.
86      */
87     template< class ct, int dim, template< int > class Ref, int codim >
88     struct canCommunicate< SPGrid< ct, dim, Ref, MPI_Comm >, codim >
89     {
90       /** \brief \ref Dune::SPGrid "SPGrid" with MPI_Comm can communicate on
91        *         all codimensions */
92       static const bool v = ((codim >= 0) && (codim <= dim));
93     };
94 #endif // #if HAVE_MPI
95 
96     /** \brief Are all levels of a grid always conform?
97      *
98      *  \tparam  Grid  grid for which the information is desired
99      */
100     template< class ct, int dim, template< int > class Ref, class Comm >
101     struct isLevelwiseConforming< SPGrid< ct, dim, Ref, Comm > >
102     {
103       /** \brief All levels of a \ref Dune::SPGrid "SPGrid" are always conform */
104       static const bool v = true;
105     };
106 
107     /** \brief Is the leaf level of a grid always conform?
108      *
109      *  \tparam  Grid  grid for which the information is desired
110      */
111     template< class ct, int dim, template< int > class Ref, class Comm >
112     struct isLeafwiseConforming< SPGrid< ct, dim, Ref, Comm > >
113     {
114       /** \brief The leaf level of a \ref Dune::SPGrid "SPGrid" are always conform */
115       static const bool v = true;
116     };
117 
118     /** \brief Does a grid provide backup and restore facilities?
119      *
120      *  \tparam  Grid  grid for which the information is desired
121      */
122     template< class ct, int dim, template< int > class Ref, class Comm >
123     struct hasBackupRestoreFacilities< SPGrid< ct, dim, Ref, Comm > >
124     {
125       /** \brief \ref Dune::SPGrid "SPGrid" provides backup and restore facilities */
126       static const bool v = true;
127     };
128 
129     /** \brief Is a grid implementation thread safe?
130      *
131      *  \tparam  Grid  grid for which the information is desired
132      */
133     template< class ct, int dim, template< int > class Ref, class Comm >
134     struct threadSafe< SPGrid< ct, dim, Ref, Comm > >
135     {
136       /** \brief \ref Dune::SPGrid "SPGrid" is not thread safe */
137       static const bool v = false;
138     };
139 
140     /** \brief Is a grid implementation thread safe while not being modified?
141      *
142      *  \tparam  Grid  grid for which the information is desired
143      */
144     template< class ct, int dim, template< int > class Ref, class Comm >
145     struct viewThreadSafe< SPGrid< ct, dim, Ref, Comm > >
146     {
147       /** \brief \ref Dune::SPGrid "SPGrid" is not thread safe */
148       static const bool v = false;
149     };
150 
151 
152 
153     // non-standard capabilities (see dune-fem)
154     // ----------------------------------------
155 
156     template< class Grid >
157     struct hasHierarchicIndexSet;
158 
159     template< class ct, int dim, template< int > class Ref, class Comm >
160     struct hasHierarchicIndexSet< SPGrid< ct, dim, Ref, Comm > >
161     {
162       static const bool v = true;
163     };
164 
165 
166     template< class Grid >
167     struct supportsCallbackAdaptation;
168 
169     /** Does a grid implementation support callback adaptation?
170      *
171      *  \tparam  Grid  grid for which the information is desired
172      *
173      *  \note This is not a standard dune-grid capability.
174      */
175     template< class ct, int dim, template< int > class Ref, class Comm >
176     struct supportsCallbackAdaptation< SPGrid< ct, dim, Ref, Comm > >
177     {
178       /** \brief \ref Dune::SPGrid "SPGrid" supports callback adaptation */
179       static const bool v = true;
180     };
181 
182   } // namespace Capabilities
183 
184 
185 
186   // Extensions
187   // ----------
188 
189   namespace Extensions
190   {
191 
192     /** \brief Does a grid support superentity iterators of a codimension?
193      *
194      *  \tparam  Grid   grid for which the information is desired
195      *  \tparam  codim  codimension in question
196      */
197     template< class ct, int dim, template< int > class Ref, class Comm, int codim >
198     struct SuperEntityIterator< SPGrid< ct, dim, Ref, Comm >, codim >
199     {
200       /** \brief \ref Dune::SPGrid "SPGrid" supports superentity iterators for all
201        *         codimensions */
202       static const bool v = ((codim >= 0) && (codim <= dim));
203     };
204 
205   } // namespace Extensions
206 
207 } // namespace Dune
208 
209 #endif // #ifndef DUNE_SPGRID_CAPABILITIES_HH
210