This file is indexed.

/usr/include/deal.II/fe/fe_nothing.h is in libdeal.ii-dev 6.3.1-1.1.

This file is owned by root:root, with mode 0o644.

The actual contents of the file can be viewed below.

  1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
//---------------------------------------------------------------------------
//    $Id: fe_nothing.h 19771 2009-10-08 18:41:41Z bangerth $
//    Version: $Name$
//
//    Copyright (C) 2009 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.
//
//---------------------------------------------------------------------------
#ifndef __deal2__fe_nothing_h
#define __deal2__fe_nothing_h

#include <base/config.h>
#include <fe/fe.h>

DEAL_II_NAMESPACE_OPEN


/*!@addtogroup fe */
/*@{*/

/**
 * Definition of a finite element with zero degrees of freedom.  This
 * class is useful (in the context of an hp method) to represent empty
 * cells in the triangulation on which no degrees of freedom should
 * be allocated.  Thus a triangulation may be divided into two regions:
 * an active region where normal elements are used, and an inactive
 * region where FE_Nothing elements are used.  The hp::DoFHandler will
 * therefore assign no degrees of freedom to the FE_Nothing cells, and
 * this subregion is therefore implicitly deleted from the computation.
 *
 * Note that some care must be taken that the resulting mesh topology
 * continues to make sense when FE_Nothing elements are introduced.
 * This is particularly true when dealing with hanging node constraints,
 * because the library makes some basic assumptions about the nature
 * of those constraints.  The following geometries are acceptable:
 *
 * +---------+----+----+
 * |         | 0  |    |
 * |    1    +----+----+
 * |         | 0  |    |
 * +---------+----+----+
 *
 * +---------+----+----+
 * |         | 1  |    |
 * |    0    +----+----+
 * |         | 1  |    |
 * +---------+----+----+
 *
 * Here, 0 denotes an FE_Nothing cell, and 1 denotes some other
 * element type.  The library has no difficulty computing the necessary
 * hanging node constraints in these cases (i.e. no constraint).  
 * However, the following geometry is NOT acceptable (at least 
 * in the current implementation):
 *
 * +---------+----+----+
 * |         | 0  |    |
 * |    1    +----+----+
 * |         | 1  |    |
 * +---------+----+----+
 * 
 * The distinction lies in the mixed nature of the child faces,
 * a case we have not implemented as of yet.
 *
 * @author Joshua White
 */
template <int dim>
class FE_Nothing : public FiniteElement<dim>
{
  public:

                                    /**
                                      * Constructor. Argument denotes the
                                      * number of components to give this
                                      * finite element (default = 1).  
                                      */
    FE_Nothing (unsigned int n_components = 1);
    
                                     /**
                                      * A sort of virtual copy
                                      * constructor. Some places in
                                      * the library, for example the
                                      * constructors of FESystem as
                                      * well as the hp::FECollection
                                      * class, need to make copied of
                                      * finite elements without
                                      * knowing their exact type. They
                                      * do so through this function.
                                      */
    virtual
    FiniteElement<dim> *
    clone() const;

                                     /**
                                      * Return a string that uniquely
                                      * identifies a finite
                                      * element. In this case it is
                                      * <code>FE_Nothing@<dim@></code>.
                                      */
    virtual
    std::string
    get_name() const;

                                     /**
                                      * Number of base elements in a
                                      * mixed discretization. In this case
                                      * we only have one.
                                      */
    virtual
    unsigned int
    n_base_elements () const;
    
                                     /**
                                      * Access to base element
                                      * objects. Since the element is
                                      * scalar, then
                                      * <code>base_element(0)</code> is
                                      * @p this.
                                      */
    virtual
    const FiniteElement<dim> &
    base_element (const unsigned int index) const;

                                     /**
                                      * This index denotes how often
                                      * the base element @p index is
                                      * used in a composed element. If
                                      * the element is scalar, then
                                      * the result is always equal to
                                      * one. See the documentation for
                                      * the n_base_elements()
                                      * function for more details.
                                      */
    virtual
    unsigned int
    element_multiplicity (const unsigned int index) const;

                                     /**
                                      * Determine the values a finite
                                      * element should compute on
                                      * initialization of data for
                                      * FEValues.
                                      *
                                      * Given a set of flags
                                      * indicating what quantities are
                                      * requested from a FEValues
                                      * object, update_once() and
                                      * update_each() compute which
                                      * values must really be
                                      * computed. Then, the
                                      * <tt>fill_*_values</tt> functions
                                      * are called with the result of
                                      * these.
                                      *
                                      * In this case, since the element 
                                      * has zero degrees of freedom and
                                      * no information can be computed on
                                      * it, this function simply returns
                                      * the default (empty) set of update
                                      * flags.
                                      */

    virtual
    UpdateFlags
    update_once (const UpdateFlags flags) const;

                                     /**
                                      * Complementary function for
                                      * update_once().
                                      *
                                      * While update_once() returns
                                      * the values to be computed on
                                      * the unit cell for yielding the
                                      * required data, this function
                                      * determines the values that
                                      * must be recomputed on each
                                      * cell.
                                      *
                                      * Refer to update_once() for
                                      * more details.
                                      */
    virtual
    UpdateFlags
    update_each (const UpdateFlags flags) const;

                                     /**
                                      * Return the value of the
                                      * @p ith shape function at the
                                      * point @p p. @p p is a point
                                      * on the reference element. Because the
                                      * current element has no degrees of freedom,
                                      * this function should obviously not be
                                      * called in practice.  All this function
                                      * really does, therefore, is trigger an
                                      * exception.
                                      */
    virtual
    double
    shape_value (const unsigned int i, const Point<dim> &p) const;

                                     /**
                                      * Fill the fields of
                                      * FEValues. This function
                                      * performs all the operations
                                      * needed to compute the data of an
                                      * FEValues object.
                                      *
                                      * In the current case, this function
                                      * returns no meaningful information, 
                                      * since the element has no degrees of
                                      * freedom.
                                      */
    virtual
    void
    fill_fe_values (const Mapping<dim> & mapping,
		    const typename Triangulation<dim>::cell_iterator & cell,
		    const Quadrature<dim> & quadrature,
		    typename Mapping<dim>::InternalDataBase & mapping_data,
		    typename Mapping<dim>::InternalDataBase & fedata,
		    FEValuesData<dim,dim> & data,
		    CellSimilarity::Similarity & cell_similarity) const;

                                     /**
                                      * Fill the fields of
                                      * FEFaceValues. This function
                                      * performs all the operations
                                      * needed to compute the data of an
                                      * FEFaceValues object.
                                      *
                                      * In the current case, this function
                                      * returns no meaningful information, 
                                      * since the element has no degrees of
                                      * freedom.
                                      */
    virtual
    void
    fill_fe_face_values (const Mapping<dim> & mapping,
			 const typename Triangulation<dim> :: cell_iterator & cell,
			 const unsigned int face,
			 const Quadrature<dim-1> & quadrature,
			 typename Mapping<dim> :: InternalDataBase & mapping_data,
			 typename Mapping<dim> :: InternalDataBase & fedata,
			 FEValuesData<dim,dim> & data) const;

                                     /**
                                      * Fill the fields of
                                      * FESubFaceValues. This function
                                      * performs all the operations
                                      * needed to compute the data of an
                                      * FESubFaceValues object.
                                      *
                                      * In the current case, this function
                                      * returns no meaningful information, 
                                      * since the element has no degrees of
                                      * freedom.
                                      */
    virtual
    void
    fill_fe_subface_values (const Mapping<dim> & mapping,
			    const typename Triangulation<dim>::cell_iterator & cell,
			    const unsigned int face,
			    const unsigned int subface,
			    const Quadrature<dim-1> & quadrature,
			    typename Mapping<dim>::InternalDataBase & mapping_data,
			    typename Mapping<dim>::InternalDataBase & fedata,
			    FEValuesData<dim,dim> & data) const;

                                     /**
                                      * Prepare internal data
                                      * structures and fill in values
                                      * independent of the
                                      * cell. Returns a pointer to an
                                      * object of which the caller of
                                      * this function then has to
                                      * assume ownership (which
                                      * includes destruction when it
                                      * is no more needed).
                                      *
                                      * In the current case, this function 
                                      * just returns a default pointer, since
                                      * no meaningful data exists for this 
                                      * element.
                                      */
    virtual
    typename Mapping<dim>::InternalDataBase *
    get_data (const UpdateFlags     update_flags,
	      const Mapping<dim>    & mapping,
	      const Quadrature<dim> & quadrature) const;
        
                                     /**
                                      * Return whether this element dominates
                                      * the one given as argument when they
                                      * meet at a common face,
                                      * whether it is the other way around,
                                      * whether neither dominates, or if
                                      * either could dominate.
                                      *
                                      * For a definition of domination, see
                                      * FiniteElementBase::Domination and in
                                      * particular the @ref hp_paper "hp paper".
                                      *
                                      * In the current case, this element
                                      * is always assumed to dominate, unless 
                                      * it is also of type FE_Nothing().  In
                                      * that situation, either element can
                                      * dominate.
                                      */
    virtual
    FiniteElementDomination::Domination
    compare_for_face_domination (const FiniteElement<dim> & fe_other) const;
    
    
    
    virtual
    std::vector<std::pair<unsigned int, unsigned int> >
    hp_vertex_dof_identities (const FiniteElement<dim> &fe_other) const;
    
    virtual
    std::vector<std::pair<unsigned int, unsigned int> >
    hp_line_dof_identities (const FiniteElement<dim> &fe_other) const;
    
    virtual
    std::vector<std::pair<unsigned int, unsigned int> >
    hp_quad_dof_identities (const FiniteElement<dim> &fe_other) const;
    
    virtual
    bool
    hp_constraints_are_implemented () const;

                                      /**
                                      * Return the matrix
                                      * interpolating from a face of
                                      * of one element to the face of
                                      * the neighboring element. 
                                      * The size of the matrix is
                                      * then <tt>source.#dofs_per_face</tt> times
                                      * <tt>this->#dofs_per_face</tt>.
                                      *
                                      * Since the current finite element has no
                                      * degrees of freedom, the interpolation
                                      * matrix is necessarily empty.
                                      */

    virtual
    void
    get_face_interpolation_matrix (const FiniteElement<dim> &source_fe,
                                   FullMatrix<double>       &interpolation_matrix) const;
                                   
                                   
                                     /**
                                      * Return the matrix
                                      * interpolating from a face of
                                      * of one element to the subface of
                                      * the neighboring element. 
                                      * The size of the matrix is
                                      * then <tt>source.#dofs_per_face</tt> times
                                      * <tt>this->#dofs_per_face</tt>.
                                      *
                                      * Since the current finite element has no
                                      * degrees of freedom, the interpolation
                                      * matrix is necessarily empty.
                                      */
                                      
    virtual
    void
    get_subface_interpolation_matrix (const FiniteElement<dim> & source_fe,
                                      const unsigned int index,
                                      FullMatrix<double>  &interpolation_matrix) const;

    
};


/*@}*/

DEAL_II_NAMESPACE_CLOSE

#endif