← Back to branch summary

~ubuntu-branches/ubuntu/maverick/dolfin/maverick

« back to all changes in this revision

Viewing changes to demo/function/nonmatching-interpolation/cpp/P1.h

  • Committer: Bazaar Package Importer
  • Author(s): Johannes Ring
  • Date: 2009-10-12 14:13:18 UTC
  • mfrom: (1.1.1 upstream)
  • Revision ID: james.westby@ubuntu.com-20091012141318-hkbxl0sq555vqv7d
Tags: 0.9.4-1
* New upstream release. This version cleans up the design of the
  function class by adding a new abstraction for user-defined
  functions called Expression. A number of minor bugfixes and
  improvements have also been made.
* debian/watch: Update for new flat directory structure.
* Update debian/copyright.
* debian/rules: Use explicit paths to PETSc 3.0.0 and SLEPc 3.0.0.

Show diffs side-by-side

added added

removed removed

Lines of Context:
demo/function/nonmatching-interpolation/cpp/P1.h
 
1
// This code conforms with the UFC specification version 1.0
 
2
// and was automatically generated by FFC version 0.7.0.
 
3
//
 
4
// Warning: This code was generated with the option '-l dolfin'
 
5
// and contains DOLFIN-specific wrappers that depend on DOLFIN.
 
6
 
 
7
#ifndef __P1_H
 
8
#define __P1_H
 
9
 
 
10
#include <cmath>
 
11
#include <stdexcept>
 
12
#include <fstream>
 
13
#include <ufc.h>
 
14
    
 
15
/// This class defines the interface for a finite element.
 
16
 
 
17
class p1_0_finite_element_0: public ufc::finite_element
 
18
{
 
19
public:
 
20
 
 
21
  /// Constructor
 
22
  p1_0_finite_element_0() : ufc::finite_element()
 
23
  {
 
24
    // Do nothing
 
25
  }
 
26
 
 
27
  /// Destructor
 
28
  virtual ~p1_0_finite_element_0()
 
29
  {
 
30
    // Do nothing
 
31
  }
 
32
 
 
33
  /// Return a string identifying the finite element
 
34
  virtual const char* signature() const
 
35
  {
 
36
    return "FiniteElement('Lagrange', 'triangle', 1)";
 
37
  }
 
38
 
 
39
  /// Return the cell shape
 
40
  virtual ufc::shape cell_shape() const
 
41
  {
 
42
    return ufc::triangle;
 
43
  }
 
44
 
 
45
  /// Return the dimension of the finite element function space
 
46
  virtual unsigned int space_dimension() const
 
47
  {
 
48
    return 3;
 
49
  }
 
50
 
 
51
  /// Return the rank of the value space
 
52
  virtual unsigned int value_rank() const
 
53
  {
 
54
    return 0;
 
55
  }
 
56
 
 
57
  /// Return the dimension of the value space for axis i
 
58
  virtual unsigned int value_dimension(unsigned int i) const
 
59
  {
 
60
    return 1;
 
61
  }
 
62
 
 
63
  /// Evaluate basis function i at given point in cell
 
64
  virtual void evaluate_basis(unsigned int i,
 
65
                              double* values,
 
66
                              const double* coordinates,
 
67
                              const ufc::cell& c) const
 
68
  {
 
69
    // Extract vertex coordinates
 
70
    const double * const * element_coordinates = c.coordinates;
 
71
    
 
72
    // Compute Jacobian of affine map from reference cell
 
73
    const double J_00 = element_coordinates[1][0] - element_coordinates[0][0];
 
74
    const double J_01 = element_coordinates[2][0] - element_coordinates[0][0];
 
75
    const double J_10 = element_coordinates[1][1] - element_coordinates[0][1];
 
76
    const double J_11 = element_coordinates[2][1] - element_coordinates[0][1];
 
77
    
 
78
    // Compute determinant of Jacobian
 
79
    const double detJ = J_00*J_11 - J_01*J_10;
 
80
    
 
81
    // Compute inverse of Jacobian
 
82
    
 
83
    // Get coordinates and map to the reference (UFC) element
 
84
    double x = (element_coordinates[0][1]*element_coordinates[2][0] -\
 
85
                element_coordinates[0][0]*element_coordinates[2][1] +\
 
86
                J_11*coordinates[0] - J_01*coordinates[1]) / detJ;
 
87
    double y = (element_coordinates[1][1]*element_coordinates[0][0] -\
 
88
                element_coordinates[1][0]*element_coordinates[0][1] -\
 
89
                J_10*coordinates[0] + J_00*coordinates[1]) / detJ;
 
90
    
 
91
    // Map coordinates to the reference square
 
92
    if (std::abs(y - 1.0) < 1e-14)
 
93
      x = -1.0;
 
94
    else
 
95
      x = 2.0 *x/(1.0 - y) - 1.0;
 
96
    y = 2.0*y - 1.0;
 
97
    
 
98
    // Reset values
 
99
    *values = 0;
 
100
    
 
101
    // Map degree of freedom to element degree of freedom
 
102
    const unsigned int dof = i;
 
103
    
 
104
    // Generate scalings
 
105
    const double scalings_y_0 = 1;
 
106
    const double scalings_y_1 = scalings_y_0*(0.5 - 0.5*y);
 
107
    
 
108
    // Compute psitilde_a
 
109
    const double psitilde_a_0 = 1;
 
110
    const double psitilde_a_1 = x;
 
111
    
 
112
    // Compute psitilde_bs
 
113
    const double psitilde_bs_0_0 = 1;
 
114
    const double psitilde_bs_0_1 = 1.5*y + 0.5;
 
115
    const double psitilde_bs_1_0 = 1;
 
116
    
 
117
    // Compute basisvalues
 
118
    const double basisvalue0 = 0.707106781186548*psitilde_a_0*scalings_y_0*psitilde_bs_0_0;
 
119
    const double basisvalue1 = 1.73205080756888*psitilde_a_1*scalings_y_1*psitilde_bs_1_0;
 
120
    const double basisvalue2 = psitilde_a_0*scalings_y_0*psitilde_bs_0_1;
 
121
    
 
122
    // Table(s) of coefficients
 
123
    static const double coefficients0[3][3] = \
 
124
    {{0.471404520791032, -0.288675134594813, -0.166666666666667},
 
125
    {0.471404520791032, 0.288675134594813, -0.166666666666667},
 
126
    {0.471404520791032, 0, 0.333333333333333}};
 
127
    
 
128
    // Extract relevant coefficients
 
129
    const double coeff0_0 = coefficients0[dof][0];
 
130
    const double coeff0_1 = coefficients0[dof][1];
 
131
    const double coeff0_2 = coefficients0[dof][2];
 
132
    
 
133
    // Compute value(s)
 
134
    *values = coeff0_0*basisvalue0 + coeff0_1*basisvalue1 + coeff0_2*basisvalue2;
 
135
  }
 
136
 
 
137
  /// Evaluate all basis functions at given point in cell
 
138
  virtual void evaluate_basis_all(double* values,
 
139
                                  const double* coordinates,
 
140
                                  const ufc::cell& c) const
 
141
  {
 
142
    throw std::runtime_error("The vectorised version of evaluate_basis() is not yet implemented.");
 
143
  }
 
144
 
 
145
  /// Evaluate order n derivatives of basis function i at given point in cell
 
146
  virtual void evaluate_basis_derivatives(unsigned int i,
 
147
                                          unsigned int n,
 
148
                                          double* values,
 
149
                                          const double* coordinates,
 
150
                                          const ufc::cell& c) const
 
151
  {
 
152
    // Extract vertex coordinates
 
153
    const double * const * element_coordinates = c.coordinates;
 
154
    
 
155
    // Compute Jacobian of affine map from reference cell
 
156
    const double J_00 = element_coordinates[1][0] - element_coordinates[0][0];
 
157
    const double J_01 = element_coordinates[2][0] - element_coordinates[0][0];
 
158
    const double J_10 = element_coordinates[1][1] - element_coordinates[0][1];
 
159
    const double J_11 = element_coordinates[2][1] - element_coordinates[0][1];
 
160
    
 
161
    // Compute determinant of Jacobian
 
162
    const double detJ = J_00*J_11 - J_01*J_10;
 
163
    
 
164
    // Compute inverse of Jacobian
 
165
    
 
166
    // Get coordinates and map to the reference (UFC) element
 
167
    double x = (element_coordinates[0][1]*element_coordinates[2][0] -\
 
168
                element_coordinates[0][0]*element_coordinates[2][1] +\
 
169
                J_11*coordinates[0] - J_01*coordinates[1]) / detJ;
 
170
    double y = (element_coordinates[1][1]*element_coordinates[0][0] -\
 
171
                element_coordinates[1][0]*element_coordinates[0][1] -\
 
172
                J_10*coordinates[0] + J_00*coordinates[1]) / detJ;
 
173
    
 
174
    // Map coordinates to the reference square
 
175
    if (std::abs(y - 1.0) < 1e-14)
 
176
      x = -1.0;
 
177
    else
 
178
      x = 2.0 *x/(1.0 - y) - 1.0;
 
179
    y = 2.0*y - 1.0;
 
180
    
 
181
    // Compute number of derivatives
 
182
    unsigned int num_derivatives = 1;
 
183
    
 
184
    for (unsigned int j = 0; j < n; j++)
 
185
      num_derivatives *= 2;
 
186
    
 
187
    
 
188
    // Declare pointer to two dimensional array that holds combinations of derivatives and initialise
 
189
    unsigned int **combinations = new unsigned int *[num_derivatives];
 
190
    
 
191
    for (unsigned int j = 0; j < num_derivatives; j++)
 
192
    {
 
193
      combinations[j] = new unsigned int [n];
 
194
      for (unsigned int k = 0; k < n; k++)
 
195
        combinations[j][k] = 0;
 
196
    }
 
197
    
 
198
    // Generate combinations of derivatives
 
199
    for (unsigned int row = 1; row < num_derivatives; row++)
 
200
    {
 
201
      for (unsigned int num = 0; num < row; num++)
 
202
      {
 
203
        for (unsigned int col = n-1; col+1 > 0; col--)
 
204
        {
 
205
          if (combinations[row][col] + 1 > 1)
 
206
            combinations[row][col] = 0;
 
207
          else
 
208
          {
 
209
            combinations[row][col] += 1;
 
210
            break;
 
211
          }
 
212
        }
 
213
      }
 
214
    }
 
215
    
 
216
    // Compute inverse of Jacobian
 
217
    const double Jinv[2][2] =  {{J_11 / detJ, -J_01 / detJ}, {-J_10 / detJ, J_00 / detJ}};
 
218
    
 
219
    // Declare transformation matrix
 
220
    // Declare pointer to two dimensional array and initialise
 
221
    double **transform = new double *[num_derivatives];
 
222
    
 
223
    for (unsigned int j = 0; j < num_derivatives; j++)
 
224
    {
 
225
      transform[j] = new double [num_derivatives];
 
226
      for (unsigned int k = 0; k < num_derivatives; k++)
 
227
        transform[j][k] = 1;
 
228
    }
 
229
    
 
230
    // Construct transformation matrix
 
231
    for (unsigned int row = 0; row < num_derivatives; row++)
 
232
    {
 
233
      for (unsigned int col = 0; col < num_derivatives; col++)
 
234
      {
 
235
        for (unsigned int k = 0; k < n; k++)
 
236
          transform[row][col] *= Jinv[combinations[col][k]][combinations[row][k]];
 
237
      }
 
238
    }
 
239
    
 
240
    // Reset values
 
241
    for (unsigned int j = 0; j < 1*num_derivatives; j++)
 
242
      values[j] = 0;
 
243
    
 
244
    // Map degree of freedom to element degree of freedom
 
245
    const unsigned int dof = i;
 
246
    
 
247
    // Generate scalings
 
248
    const double scalings_y_0 = 1;
 
249
    const double scalings_y_1 = scalings_y_0*(0.5 - 0.5*y);
 
250
    
 
251
    // Compute psitilde_a
 
252
    const double psitilde_a_0 = 1;
 
253
    const double psitilde_a_1 = x;
 
254
    
 
255
    // Compute psitilde_bs
 
256
    const double psitilde_bs_0_0 = 1;
 
257
    const double psitilde_bs_0_1 = 1.5*y + 0.5;
 
258
    const double psitilde_bs_1_0 = 1;
 
259
    
 
260
    // Compute basisvalues
 
261
    const double basisvalue0 = 0.707106781186548*psitilde_a_0*scalings_y_0*psitilde_bs_0_0;
 
262
    const double basisvalue1 = 1.73205080756888*psitilde_a_1*scalings_y_1*psitilde_bs_1_0;
 
263
    const double basisvalue2 = psitilde_a_0*scalings_y_0*psitilde_bs_0_1;
 
264
    
 
265
    // Table(s) of coefficients
 
266
    static const double coefficients0[3][3] = \
 
267
    {{0.471404520791032, -0.288675134594813, -0.166666666666667},
 
268
    {0.471404520791032, 0.288675134594813, -0.166666666666667},
 
269
    {0.471404520791032, 0, 0.333333333333333}};
 
270
    
 
271
    // Interesting (new) part
 
272
    // Tables of derivatives of the polynomial base (transpose)
 
273
    static const double dmats0[3][3] = \
 
274
    {{0, 0, 0},
 
275
    {4.89897948556636, 0, 0},
 
276
    {0, 0, 0}};
 
277
    
 
278
    static const double dmats1[3][3] = \
 
279
    {{0, 0, 0},
 
280
    {2.44948974278318, 0, 0},
 
281
    {4.24264068711928, 0, 0}};
 
282
    
 
283
    // Compute reference derivatives
 
284
    // Declare pointer to array of derivatives on FIAT element
 
285
    double *derivatives = new double [num_derivatives];
 
286
    
 
287
    // Declare coefficients
 
288
    double coeff0_0 = 0;
 
289
    double coeff0_1 = 0;
 
290
    double coeff0_2 = 0;
 
291
    
 
292
    // Declare new coefficients
 
293
    double new_coeff0_0 = 0;
 
294
    double new_coeff0_1 = 0;
 
295
    double new_coeff0_2 = 0;
 
296
    
 
297
    // Loop possible derivatives
 
298
    for (unsigned int deriv_num = 0; deriv_num < num_derivatives; deriv_num++)
 
299
    {
 
300
      // Get values from coefficients array
 
301
      new_coeff0_0 = coefficients0[dof][0];
 
302
      new_coeff0_1 = coefficients0[dof][1];
 
303
      new_coeff0_2 = coefficients0[dof][2];
 
304
    
 
305
      // Loop derivative order
 
306
      for (unsigned int j = 0; j < n; j++)
 
307
      {
 
308
        // Update old coefficients
 
309
        coeff0_0 = new_coeff0_0;
 
310
        coeff0_1 = new_coeff0_1;
 
311
        coeff0_2 = new_coeff0_2;
 
312
    
 
313
        if(combinations[deriv_num][j] == 0)
 
314
        {
 
315
          new_coeff0_0 = coeff0_0*dmats0[0][0] + coeff0_1*dmats0[1][0] + coeff0_2*dmats0[2][0];
 
316
          new_coeff0_1 = coeff0_0*dmats0[0][1] + coeff0_1*dmats0[1][1] + coeff0_2*dmats0[2][1];
 
317
          new_coeff0_2 = coeff0_0*dmats0[0][2] + coeff0_1*dmats0[1][2] + coeff0_2*dmats0[2][2];
 
318
        }
 
319
        if(combinations[deriv_num][j] == 1)
 
320
        {
 
321
          new_coeff0_0 = coeff0_0*dmats1[0][0] + coeff0_1*dmats1[1][0] + coeff0_2*dmats1[2][0];
 
322
          new_coeff0_1 = coeff0_0*dmats1[0][1] + coeff0_1*dmats1[1][1] + coeff0_2*dmats1[2][1];
 
323
          new_coeff0_2 = coeff0_0*dmats1[0][2] + coeff0_1*dmats1[1][2] + coeff0_2*dmats1[2][2];
 
324
        }
 
325
    
 
326
      }
 
327
      // Compute derivatives on reference element as dot product of coefficients and basisvalues
 
328
      derivatives[deriv_num] = new_coeff0_0*basisvalue0 + new_coeff0_1*basisvalue1 + new_coeff0_2*basisvalue2;
 
329
    }
 
330
    
 
331
    // Transform derivatives back to physical element
 
332
    for (unsigned int row = 0; row < num_derivatives; row++)
 
333
    {
 
334
      for (unsigned int col = 0; col < num_derivatives; col++)
 
335
      {
 
336
        values[row] += transform[row][col]*derivatives[col];
 
337
      }
 
338
    }
 
339
    // Delete pointer to array of derivatives on FIAT element
 
340
    delete [] derivatives;
 
341
    
 
342
    // Delete pointer to array of combinations of derivatives and transform
 
343
    for (unsigned int row = 0; row < num_derivatives; row++)
 
344
    {
 
345
      delete [] combinations[row];
 
346
      delete [] transform[row];
 
347
    }
 
348
    
 
349
    delete [] combinations;
 
350
    delete [] transform;
 
351
  }
 
352
 
 
353
  /// Evaluate order n derivatives of all basis functions at given point in cell
 
354
  virtual void evaluate_basis_derivatives_all(unsigned int n,
 
355
                                              double* values,
 
356
                                              const double* coordinates,
 
357
                                              const ufc::cell& c) const
 
358
  {
 
359
    throw std::runtime_error("The vectorised version of evaluate_basis_derivatives() is not yet implemented.");
 
360
  }
 
361
 
 
362
  /// Evaluate linear functional for dof i on the function f
 
363
  virtual double evaluate_dof(unsigned int i,
 
364
                              const ufc::function& f,
 
365
                              const ufc::cell& c) const
 
366
  {
 
367
    // The reference points, direction and weights:
 
368
    static const double X[3][1][2] = {{{0, 0}}, {{1, 0}}, {{0, 1}}};
 
369
    static const double W[3][1] = {{1}, {1}, {1}};
 
370
    static const double D[3][1][1] = {{{1}}, {{1}}, {{1}}};
 
371
    
 
372
    const double * const * x = c.coordinates;
 
373
    double result = 0.0;
 
374
    // Iterate over the points:
 
375
    // Evaluate basis functions for affine mapping
 
376
    const double w0 = 1.0 - X[i][0][0] - X[i][0][1];
 
377
    const double w1 = X[i][0][0];
 
378
    const double w2 = X[i][0][1];
 
379
    
 
380
    // Compute affine mapping y = F(X)
 
381
    double y[2];
 
382
    y[0] = w0*x[0][0] + w1*x[1][0] + w2*x[2][0];
 
383
    y[1] = w0*x[0][1] + w1*x[1][1] + w2*x[2][1];
 
384
    
 
385
    // Evaluate function at physical points
 
386
    double values[1];
 
387
    f.evaluate(values, y, c);
 
388
    
 
389
    // Map function values using appropriate mapping
 
390
    // Affine map: Do nothing
 
391
    
 
392
    // Note that we do not map the weights (yet).
 
393
    
 
394
    // Take directional components
 
395
    for(int k = 0; k < 1; k++)
 
396
      result += values[k]*D[i][0][k];
 
397
    // Multiply by weights
 
398
    result *= W[i][0];
 
399
    
 
400
    return result;
 
401
  }
 
402
 
 
403
  /// Evaluate linear functionals for all dofs on the function f
 
404
  virtual void evaluate_dofs(double* values,
 
405
                             const ufc::function& f,
 
406
                             const ufc::cell& c) const
 
407
  {
 
408
    throw std::runtime_error("Not implemented (introduced in UFC v1.1).");
 
409
  }
 
410
 
 
411
  /// Interpolate vertex values from dof values
 
412
  virtual void interpolate_vertex_values(double* vertex_values,
 
413
                                         const double* dof_values,
 
414
                                         const ufc::cell& c) const
 
415
  {
 
416
    // Evaluate at vertices and use affine mapping
 
417
    vertex_values[0] = dof_values[0];
 
418
    vertex_values[1] = dof_values[1];
 
419
    vertex_values[2] = dof_values[2];
 
420
  }
 
421
 
 
422
  /// Return the number of sub elements (for a mixed element)
 
423
  virtual unsigned int num_sub_elements() const
 
424
  {
 
425
    return 1;
 
426
  }
 
427
 
 
428
  /// Create a new finite element for sub element i (for a mixed element)
 
429
  virtual ufc::finite_element* create_sub_element(unsigned int i) const
 
430
  {
 
431
    return new p1_0_finite_element_0();
 
432
  }
 
433
 
 
434
};
 
435
 
 
436
/// This class defines the interface for a local-to-global mapping of
 
437
/// degrees of freedom (dofs).
 
438
 
 
439
class p1_0_dof_map_0: public ufc::dof_map
 
440
{
 
441
private:
 
442
 
 
443
  unsigned int __global_dimension;
 
444
 
 
445
public:
 
446
 
 
447
  /// Constructor
 
448
  p1_0_dof_map_0() : ufc::dof_map()
 
449
  {
 
450
    __global_dimension = 0;
 
451
  }
 
452
 
 
453
  /// Destructor
 
454
  virtual ~p1_0_dof_map_0()
 
455
  {
 
456
    // Do nothing
 
457
  }
 
458
 
 
459
  /// Return a string identifying the dof map
 
460
  virtual const char* signature() const
 
461
  {
 
462
    return "FFC dof map for FiniteElement('Lagrange', 'triangle', 1)";
 
463
  }
 
464
 
 
465
  /// Return true iff mesh entities of topological dimension d are needed
 
466
  virtual bool needs_mesh_entities(unsigned int d) const
 
467
  {
 
468
    switch ( d )
 
469
    {
 
470
    case 0:
 
471
      return true;
 
472
      break;
 
473
    case 1:
 
474
      return false;
 
475
      break;
 
476
    case 2:
 
477
      return false;
 
478
      break;
 
479
    }
 
480
    return false;
 
481
  }
 
482
 
 
483
  /// Initialize dof map for mesh (return true iff init_cell() is needed)
 
484
  virtual bool init_mesh(const ufc::mesh& m)
 
485
  {
 
486
    __global_dimension = m.num_entities[0];
 
487
    return false;
 
488
  }
 
489
 
 
490
  /// Initialize dof map for given cell
 
491
  virtual void init_cell(const ufc::mesh& m,
 
492
                         const ufc::cell& c)
 
493
  {
 
494
    // Do nothing
 
495
  }
 
496
 
 
497
  /// Finish initialization of dof map for cells
 
498
  virtual void init_cell_finalize()
 
499
  {
 
500
    // Do nothing
 
501
  }
 
502
 
 
503
  /// Return the dimension of the global finite element function space
 
504
  virtual unsigned int global_dimension() const
 
505
  {
 
506
    return __global_dimension;
 
507
  }
 
508
 
 
509
  /// Return the dimension of the local finite element function space for a cell
 
510
  virtual unsigned int local_dimension(const ufc::cell& c) const
 
511
  {
 
512
    return 3;
 
513
  }
 
514
 
 
515
  /// Return the maximum dimension of the local finite element function space
 
516
  virtual unsigned int max_local_dimension() const
 
517
  {
 
518
    return 3;
 
519
  }
 
520
 
 
521
  // Return the geometric dimension of the coordinates this dof map provides
 
522
  virtual unsigned int geometric_dimension() const
 
523
  {
 
524
    return 2;
 
525
  }
 
526
 
 
527
  /// Return the number of dofs on each cell facet
 
528
  virtual unsigned int num_facet_dofs() const
 
529
  {
 
530
    return 2;
 
531
  }
 
532
 
 
533
  /// Return the number of dofs associated with each cell entity of dimension d
 
534
  virtual unsigned int num_entity_dofs(unsigned int d) const
 
535
  {
 
536
    throw std::runtime_error("Not implemented (introduced in UFC v1.1).");
 
537
  }
 
538
 
 
539
  /// Tabulate the local-to-global mapping of dofs on a cell
 
540
  virtual void tabulate_dofs(unsigned int* dofs,
 
541
                             const ufc::mesh& m,
 
542
                             const ufc::cell& c) const
 
543
  {
 
544
    dofs[0] = c.entity_indices[0][0];
 
545
    dofs[1] = c.entity_indices[0][1];
 
546
    dofs[2] = c.entity_indices[0][2];
 
547
  }
 
548
 
 
549
  /// Tabulate the local-to-local mapping from facet dofs to cell dofs
 
550
  virtual void tabulate_facet_dofs(unsigned int* dofs,
 
551
                                   unsigned int facet) const
 
552
  {
 
553
    switch ( facet )
 
554
    {
 
555
    case 0:
 
556
      dofs[0] = 1;
 
557
      dofs[1] = 2;
 
558
      break;
 
559
    case 1:
 
560
      dofs[0] = 0;
 
561
      dofs[1] = 2;
 
562
      break;
 
563
    case 2:
 
564
      dofs[0] = 0;
 
565
      dofs[1] = 1;
 
566
      break;
 
567
    }
 
568
  }
 
569
 
 
570
  /// Tabulate the local-to-local mapping of dofs on entity (d, i)
 
571
  virtual void tabulate_entity_dofs(unsigned int* dofs,
 
572
                                    unsigned int d, unsigned int i) const
 
573
  {
 
574
    throw std::runtime_error("Not implemented (introduced in UFC v1.1).");
 
575
  }
 
576
 
 
577
  /// Tabulate the coordinates of all dofs on a cell
 
578
  virtual void tabulate_coordinates(double** coordinates,
 
579
                                    const ufc::cell& c) const
 
580
  {
 
581
    const double * const * x = c.coordinates;
 
582
    coordinates[0][0] = x[0][0];
 
583
    coordinates[0][1] = x[0][1];
 
584
    coordinates[1][0] = x[1][0];
 
585
    coordinates[1][1] = x[1][1];
 
586
    coordinates[2][0] = x[2][0];
 
587
    coordinates[2][1] = x[2][1];
 
588
  }
 
589
 
 
590
  /// Return the number of sub dof maps (for a mixed element)
 
591
  virtual unsigned int num_sub_dof_maps() const
 
592
  {
 
593
    return 1;
 
594
  }
 
595
 
 
596
  /// Create a new dof_map for sub dof map i (for a mixed element)
 
597
  virtual ufc::dof_map* create_sub_dof_map(unsigned int i) const
 
598
  {
 
599
    return new p1_0_dof_map_0();
 
600
  }
 
601
 
 
602
};
 
603
 
 
604
/// This class defines the interface for the assembly of the global
 
605
/// tensor corresponding to a form with r + n arguments, that is, a
 
606
/// mapping
 
607
///
 
608
///     a : V1 x V2 x ... Vr x W1 x W2 x ... x Wn -> R
 
609
///
 
610
/// with arguments v1, v2, ..., vr, w1, w2, ..., wn. The rank r
 
611
/// global tensor A is defined by
 
612
///
 
613
///     A = a(V1, V2, ..., Vr, w1, w2, ..., wn),
 
614
///
 
615
/// where each argument Vj represents the application to the
 
616
/// sequence of basis functions of Vj and w1, w2, ..., wn are given
 
617
/// fixed functions (coefficients).
 
618
 
 
619
class p1_form_0: public ufc::form
 
620
{
 
621
public:
 
622
 
 
623
  /// Constructor
 
624
  p1_form_0() : ufc::form()
 
625
  {
 
626
    // Do nothing
 
627
  }
 
628
 
 
629
  /// Destructor
 
630
  virtual ~p1_form_0()
 
631
  {
 
632
    // Do nothing
 
633
  }
 
634
 
 
635
  /// Return a string identifying the form
 
636
  virtual const char* signature() const
 
637
  {
 
638
    return "None";
 
639
  }
 
640
 
 
641
  /// Return the rank of the global tensor (r)
 
642
  virtual unsigned int rank() const
 
643
  {
 
644
    return -1;
 
645
  }
 
646
 
 
647
  /// Return the number of coefficients (n)
 
648
  virtual unsigned int num_coefficients() const
 
649
  {
 
650
    return 0;
 
651
  }
 
652
 
 
653
  /// Return the number of cell integrals
 
654
  virtual unsigned int num_cell_integrals() const
 
655
  {
 
656
    return 0;
 
657
  }
 
658
 
 
659
  /// Return the number of exterior facet integrals
 
660
  virtual unsigned int num_exterior_facet_integrals() const
 
661
  {
 
662
    return 0;
 
663
  }
 
664
 
 
665
  /// Return the number of interior facet integrals
 
666
  virtual unsigned int num_interior_facet_integrals() const
 
667
  {
 
668
    return 0;
 
669
  }
 
670
 
 
671
  /// Create a new finite element for argument function i
 
672
  virtual ufc::finite_element* create_finite_element(unsigned int i) const
 
673
  {
 
674
    return 0;
 
675
  }
 
676
 
 
677
  /// Create a new dof map for argument function i
 
678
  virtual ufc::dof_map* create_dof_map(unsigned int i) const
 
679
  {
 
680
    return 0;
 
681
  }
 
682
 
 
683
  /// Create a new cell integral on sub domain i
 
684
  virtual ufc::cell_integral* create_cell_integral(unsigned int i) const
 
685
  {
 
686
    return 0;
 
687
  }
 
688
 
 
689
  /// Create a new exterior facet integral on sub domain i
 
690
  virtual ufc::exterior_facet_integral* create_exterior_facet_integral(unsigned int i) const
 
691
  {
 
692
    return 0;
 
693
  }
 
694
 
 
695
  /// Create a new interior facet integral on sub domain i
 
696
  virtual ufc::interior_facet_integral* create_interior_facet_integral(unsigned int i) const
 
697
  {
 
698
    return 0;
 
699
  }
 
700
 
 
701
};
 
702
 
 
703
// DOLFIN wrappers
 
704
 
 
705
// Standard library includes
 
706
#include <string>
 
707
 
 
708
// DOLFIN includes
 
709
#include <dolfin/common/NoDeleter.h>
 
710
#include <dolfin/fem/FiniteElement.h>
 
711
#include <dolfin/fem/DofMap.h>
 
712
#include <dolfin/fem/Form.h>
 
713
#include <dolfin/function/FunctionSpace.h>
 
714
#include <dolfin/function/GenericFunction.h>
 
715
#include <dolfin/function/CoefficientAssigner.h>
 
716
 
 
717
namespace P1
 
718
{
 
719
 
 
720
class Form_0_FunctionSpace_0: public dolfin::FunctionSpace
 
721
{
 
722
public:
 
723
 
 
724
  Form_0_FunctionSpace_0(const dolfin::Mesh& mesh):
 
725
      dolfin::FunctionSpace(dolfin::reference_to_no_delete_pointer(mesh),
 
726
                            boost::shared_ptr<const dolfin::FiniteElement>(new dolfin::FiniteElement(boost::shared_ptr<ufc::finite_element>(new p1_0_finite_element_0()))),
 
727
                            boost::shared_ptr<const dolfin::DofMap>(new dolfin::DofMap(boost::shared_ptr<ufc::dof_map>(new p1_0_dof_map_0()), dolfin::reference_to_no_delete_pointer(mesh))))
 
728
  {
 
729
    // Do nothing
 
730
  }
 
731
 
 
732
  Form_0_FunctionSpace_0(dolfin::Mesh& mesh):
 
733
    dolfin::FunctionSpace(dolfin::reference_to_no_delete_pointer(mesh),
 
734
                          boost::shared_ptr<const dolfin::FiniteElement>(new dolfin::FiniteElement(boost::shared_ptr<ufc::finite_element>(new p1_0_finite_element_0()))),
 
735
                          boost::shared_ptr<const dolfin::DofMap>(new dolfin::DofMap(boost::shared_ptr<ufc::dof_map>(new p1_0_dof_map_0()), dolfin::reference_to_no_delete_pointer(mesh))))
 
736
  {
 
737
    // Do nothing
 
738
  }
 
739
 
 
740
  Form_0_FunctionSpace_0(boost::shared_ptr<dolfin::Mesh> mesh):
 
741
      dolfin::FunctionSpace(mesh,
 
742
                            boost::shared_ptr<const dolfin::FiniteElement>(new dolfin::FiniteElement(boost::shared_ptr<ufc::finite_element>(new p1_0_finite_element_0()))),
 
743
                            boost::shared_ptr<const dolfin::DofMap>(new dolfin::DofMap(boost::shared_ptr<ufc::dof_map>(new p1_0_dof_map_0()), mesh)))
 
744
  {
 
745
      // Do nothing
 
746
  }
 
747
 
 
748
  Form_0_FunctionSpace_0(boost::shared_ptr<const dolfin::Mesh> mesh):
 
749
      dolfin::FunctionSpace(mesh,
 
750
                            boost::shared_ptr<const dolfin::FiniteElement>(new dolfin::FiniteElement(boost::shared_ptr<ufc::finite_element>(new p1_0_finite_element_0()))),
 
751
                            boost::shared_ptr<const dolfin::DofMap>(new dolfin::DofMap(boost::shared_ptr<ufc::dof_map>(new p1_0_dof_map_0()), mesh)))
 
752
  {
 
753
      // Do nothing
 
754
  }
 
755
 
 
756
 
 
757
  ~Form_0_FunctionSpace_0()
 
758
  {
 
759
  }
 
760
 
 
761
};
 
762
 
 
763
class Form_0: public dolfin::Form
 
764
{
 
765
public:
 
766
 
 
767
  // Constructor
 
768
  Form_0(const dolfin::FunctionSpace& V0):
 
769
    dolfin::Form(1, 0)
 
770
  {
 
771
    _function_spaces[0] = reference_to_no_delete_pointer(V0);
 
772
 
 
773
    _ufc_form = boost::shared_ptr<const ufc::form>(new p1_form_0());
 
774
  }
 
775
 
 
776
  // Constructor
 
777
  Form_0(boost::shared_ptr<const dolfin::FunctionSpace> V0):
 
778
    dolfin::Form(1, 0)
 
779
  {
 
780
    _function_spaces[0] = V0;
 
781
 
 
782
    _ufc_form = boost::shared_ptr<const ufc::form>(new p1_form_0());
 
783
  }
 
784
 
 
785
  // Destructor
 
786
  ~Form_0()
 
787
  {}
 
788
 
 
789
  /// Return the number of the coefficient with this name
 
790
  virtual dolfin::uint coefficient_number(const std::string& name) const
 
791
  {
 
792
 
 
793
    dolfin::error("No coefficients.");
 
794
    return 0;
 
795
  }
 
796
 
 
797
  /// Return the name of the coefficient with this number
 
798
  virtual std::string coefficient_name(dolfin::uint i) const
 
799
  {
 
800
 
 
801
    dolfin::error("No coefficients.");
 
802
    return "unnamed";
 
803
  }
 
804
 
 
805
  // Typedefs
 
806
  typedef Form_0_FunctionSpace_0 TestSpace;
 
807
 
 
808
  // Coefficients
 
809
};
 
810
 
 
811
// Class typedefs
 
812
typedef Form_0 LinearForm;
 
813
typedef Form_0::TestSpace FunctionSpace;
 
814
 
 
815
} // namespace P1
 
816
 
 
817
#endif