← Back to branch summary

~njansson/dolfin/hpc

« back to all changes in this revision

Viewing changes to src/test/tmp/L2Norm.h

  • Committer: Anders Logg
  • Date: 2007-01-10 09:04:44 UTC
  • mfrom: (1689.1.221 trunk)
  • Revision ID: logg@simula.no-20070110090444-ecyux3n1qnei4i8h
Remove old head

Show diffs side-by-side

added added

removed removed

Lines of Context:
src/test/tmp/L2Norm.h
 
1
// Automatically generated by FFC, the FEniCS Form Compiler, version 0.3.5.
 
2
// For further information, go to http://www/fenics.org/ffc/.
 
3
// Licensed under the GNU GPL Version 2.
 
4
 
 
5
#ifndef __L2NORM_H
 
6
#define __L2NORM_H
 
7
 
 
8
#include <dolfin/Mesh.h>
 
9
#include <dolfin/Cell.h>
 
10
#include <dolfin/Point.h>
 
11
#include <dolfin/AffineMap.h>
 
12
#include <dolfin/FiniteElement.h>
 
13
#include <dolfin/FiniteElementSpec.h>
 
14
#include <dolfin/BilinearForm.h>
 
15
#include <dolfin/LinearForm.h>
 
16
#include <dolfin/Functional.h>
 
17
#include <dolfin/FEM.h>
 
18
 
 
19
namespace dolfin { namespace L2Norm {
 
20
 
 
21
/// This class contains the form to be evaluated, including
 
22
/// contributions from the interior and boundary of the domain.
 
23
 
 
24
class LinearForm : public dolfin::LinearForm
 
25
{
 
26
public:
 
27
 
 
28
  class TestElement;
 
29
 
 
30
  class FunctionElement_0;
 
31
 
 
32
  class FunctionElement_1;
 
33
 
 
34
  LinearForm(Function& w0, Function& w1);
 
35
  
 
36
 
 
37
  bool interior_contribution() const;
 
38
 
 
39
  void eval(real block[], const AffineMap& map, real det) const;
 
40
 
 
41
  bool boundary_contribution() const;
 
42
 
 
43
  void eval(real block[], const AffineMap& map, real det, unsigned int facet) const;
 
44
 
 
45
  bool interior_boundary_contribution() const;
 
46
 
 
47
  void eval(real block[], const AffineMap& map0, const AffineMap& map1, real det, unsigned int facet0, unsigned int facet1, unsigned int alignment) const;
 
48
 
 
49
};
 
50
 
 
51
class LinearForm::TestElement : public dolfin::FiniteElement
 
52
{
 
53
public:
 
54
 
 
55
  TestElement() : dolfin::FiniteElement(), tensordims(0), subelements(0)
 
56
  {
 
57
    // Element is scalar, don't need to initialize tensordims
 
58
 
 
59
    // Element is simple, don't need to initialize subelements
 
60
  }
 
61
 
 
62
  ~TestElement()
 
63
  {
 
64
    if ( tensordims ) delete [] tensordims;
 
65
    if ( subelements )
 
66
    {
 
67
      for (unsigned int i = 0; i < elementdim(); i++)
 
68
        delete subelements[i];
 
69
      delete [] subelements;
 
70
    }
 
71
  }
 
72
 
 
73
  inline unsigned int spacedim() const
 
74
  {
 
75
    return 1;
 
76
  }
 
77
 
 
78
  inline unsigned int shapedim() const
 
79
  {
 
80
    return 2;
 
81
  }
 
82
 
 
83
  inline unsigned int tensordim(unsigned int i) const
 
84
  {
 
85
    dolfin_error("Element is scalar.");
 
86
    return 0;
 
87
  }
 
88
 
 
89
  inline unsigned int elementdim() const
 
90
  {
 
91
    return 1;
 
92
  }
 
93
 
 
94
  inline unsigned int rank() const
 
95
  {
 
96
    return 0;
 
97
  }
 
98
 
 
99
  void nodemap(int nodes[], const Cell& cell, const Mesh& mesh) const
 
100
  {
 
101
    nodes[0] = cell.index();
 
102
  }
 
103
 
 
104
  void pointmap(Point points[], unsigned int components[], const AffineMap& map) const
 
105
  {
 
106
    points[0] = map(3.333333333333334e-01, 3.333333333333334e-01);
 
107
    components[0] = 0;
 
108
  }
 
109
 
 
110
  void vertexeval(uint vertex_nodes[], unsigned int vertex, const Mesh& mesh) const
 
111
  {
 
112
    // FIXME: Temporary fix for Lagrange elements
 
113
    vertex_nodes[0] = vertex;
 
114
  }
 
115
 
 
116
  const FiniteElement& operator[] (unsigned int i) const
 
117
  {
 
118
    return *this;
 
119
  }
 
120
 
 
121
  FiniteElement& operator[] (unsigned int i)
 
122
  {
 
123
    return *this;
 
124
  }
 
125
 
 
126
  FiniteElementSpec spec() const
 
127
  {
 
128
    FiniteElementSpec s("Discontinuous Lagrange", "triangle", 0);
 
129
    return s;
 
130
  }
 
131
  
 
132
private:
 
133
 
 
134
  unsigned int* tensordims;
 
135
  FiniteElement** subelements;
 
136
 
 
137
};
 
138
 
 
139
class LinearForm::FunctionElement_0 : public dolfin::FiniteElement
 
140
{
 
141
public:
 
142
 
 
143
  FunctionElement_0() : dolfin::FiniteElement(), tensordims(0), subelements(0)
 
144
  {
 
145
    // Element is scalar, don't need to initialize tensordims
 
146
 
 
147
    // Element is simple, don't need to initialize subelements
 
148
  }
 
149
 
 
150
  ~FunctionElement_0()
 
151
  {
 
152
    if ( tensordims ) delete [] tensordims;
 
153
    if ( subelements )
 
154
    {
 
155
      for (unsigned int i = 0; i < elementdim(); i++)
 
156
        delete subelements[i];
 
157
      delete [] subelements;
 
158
    }
 
159
  }
 
160
 
 
161
  inline unsigned int spacedim() const
 
162
  {
 
163
    return 3;
 
164
  }
 
165
 
 
166
  inline unsigned int shapedim() const
 
167
  {
 
168
    return 2;
 
169
  }
 
170
 
 
171
  inline unsigned int tensordim(unsigned int i) const
 
172
  {
 
173
    dolfin_error("Element is scalar.");
 
174
    return 0;
 
175
  }
 
176
 
 
177
  inline unsigned int elementdim() const
 
178
  {
 
179
    return 1;
 
180
  }
 
181
 
 
182
  inline unsigned int rank() const
 
183
  {
 
184
    return 0;
 
185
  }
 
186
 
 
187
  void nodemap(int nodes[], const Cell& cell, const Mesh& mesh) const
 
188
  {
 
189
    nodes[0] = cell.entities(0)[0];
 
190
    nodes[1] = cell.entities(0)[1];
 
191
    nodes[2] = cell.entities(0)[2];
 
192
  }
 
193
 
 
194
  void pointmap(Point points[], unsigned int components[], const AffineMap& map) const
 
195
  {
 
196
    points[0] = map(0.000000000000000e+00, 0.000000000000000e+00);
 
197
    points[1] = map(1.000000000000000e+00, 0.000000000000000e+00);
 
198
    points[2] = map(0.000000000000000e+00, 1.000000000000000e+00);
 
199
    components[0] = 0;
 
200
    components[1] = 0;
 
201
    components[2] = 0;
 
202
  }
 
203
 
 
204
  void vertexeval(uint vertex_nodes[], unsigned int vertex, const Mesh& mesh) const
 
205
  {
 
206
    // FIXME: Temporary fix for Lagrange elements
 
207
    vertex_nodes[0] = vertex;
 
208
  }
 
209
 
 
210
  const FiniteElement& operator[] (unsigned int i) const
 
211
  {
 
212
    return *this;
 
213
  }
 
214
 
 
215
  FiniteElement& operator[] (unsigned int i)
 
216
  {
 
217
    return *this;
 
218
  }
 
219
 
 
220
  FiniteElementSpec spec() const
 
221
  {
 
222
    FiniteElementSpec s("Lagrange", "triangle", 1);
 
223
    return s;
 
224
  }
 
225
  
 
226
private:
 
227
 
 
228
  unsigned int* tensordims;
 
229
  FiniteElement** subelements;
 
230
 
 
231
};
 
232
 
 
233
class LinearForm::FunctionElement_1 : public dolfin::FiniteElement
 
234
{
 
235
public:
 
236
 
 
237
  FunctionElement_1() : dolfin::FiniteElement(), tensordims(0), subelements(0)
 
238
  {
 
239
    // Element is scalar, don't need to initialize tensordims
 
240
 
 
241
    // Element is simple, don't need to initialize subelements
 
242
  }
 
243
 
 
244
  ~FunctionElement_1()
 
245
  {
 
246
    if ( tensordims ) delete [] tensordims;
 
247
    if ( subelements )
 
248
    {
 
249
      for (unsigned int i = 0; i < elementdim(); i++)
 
250
        delete subelements[i];
 
251
      delete [] subelements;
 
252
    }
 
253
  }
 
254
 
 
255
  inline unsigned int spacedim() const
 
256
  {
 
257
    return 3;
 
258
  }
 
259
 
 
260
  inline unsigned int shapedim() const
 
261
  {
 
262
    return 2;
 
263
  }
 
264
 
 
265
  inline unsigned int tensordim(unsigned int i) const
 
266
  {
 
267
    dolfin_error("Element is scalar.");
 
268
    return 0;
 
269
  }
 
270
 
 
271
  inline unsigned int elementdim() const
 
272
  {
 
273
    return 1;
 
274
  }
 
275
 
 
276
  inline unsigned int rank() const
 
277
  {
 
278
    return 0;
 
279
  }
 
280
 
 
281
  void nodemap(int nodes[], const Cell& cell, const Mesh& mesh) const
 
282
  {
 
283
    nodes[0] = cell.entities(0)[0];
 
284
    nodes[1] = cell.entities(0)[1];
 
285
    nodes[2] = cell.entities(0)[2];
 
286
  }
 
287
 
 
288
  void pointmap(Point points[], unsigned int components[], const AffineMap& map) const
 
289
  {
 
290
    points[0] = map(0.000000000000000e+00, 0.000000000000000e+00);
 
291
    points[1] = map(1.000000000000000e+00, 0.000000000000000e+00);
 
292
    points[2] = map(0.000000000000000e+00, 1.000000000000000e+00);
 
293
    components[0] = 0;
 
294
    components[1] = 0;
 
295
    components[2] = 0;
 
296
  }
 
297
 
 
298
  void vertexeval(uint vertex_nodes[], unsigned int vertex, const Mesh& mesh) const
 
299
  {
 
300
    // FIXME: Temporary fix for Lagrange elements
 
301
    vertex_nodes[0] = vertex;
 
302
  }
 
303
 
 
304
  const FiniteElement& operator[] (unsigned int i) const
 
305
  {
 
306
    return *this;
 
307
  }
 
308
 
 
309
  FiniteElement& operator[] (unsigned int i)
 
310
  {
 
311
    return *this;
 
312
  }
 
313
 
 
314
  FiniteElementSpec spec() const
 
315
  {
 
316
    FiniteElementSpec s("Lagrange", "triangle", 1);
 
317
    return s;
 
318
  }
 
319
  
 
320
private:
 
321
 
 
322
  unsigned int* tensordims;
 
323
  FiniteElement** subelements;
 
324
 
 
325
};
 
326
 
 
327
LinearForm::LinearForm(Function& w0, Function& w1) : dolfin::LinearForm(2)
 
328
{
 
329
  // Create finite element for test space
 
330
  _test = new TestElement();
 
331
 
 
332
  // Add functions
 
333
  initFunction(0, w0, new FunctionElement_0());
 
334
  initFunction(1, w1, new FunctionElement_1());
 
335
}
 
336
 
 
337
// Contribution from the interior
 
338
bool LinearForm::interior_contribution() const { return true; }
 
339
 
 
340
void LinearForm::eval(real block[], const AffineMap& map, real det) const
 
341
{
 
342
  // Compute coefficients
 
343
  const real c0_0 = c[0][0];
 
344
  const real c0_1 = c[0][1];
 
345
  const real c0_2 = c[0][2];
 
346
  const real c1_0 = c[1][0];
 
347
  const real c1_1 = c[1][1];
 
348
  const real c1_2 = c[1][2];
 
349
 
 
350
  // Compute geometry tensors
 
351
  const real G0_0_0 = det*c0_0*c0_0 + det*c1_0*c1_0;
 
352
  const real G0_0_1 = det*c0_0*c0_1 + det*c1_0*c1_1;
 
353
  const real G0_0_2 = det*c0_0*c0_2 + det*c1_0*c1_2;
 
354
  const real G0_1_0 = det*c0_1*c0_0 + det*c1_1*c1_0;
 
355
  const real G0_1_1 = det*c0_1*c0_1 + det*c1_1*c1_1;
 
356
  const real G0_1_2 = det*c0_1*c0_2 + det*c1_1*c1_2;
 
357
  const real G0_2_0 = det*c0_2*c0_0 + det*c1_2*c1_0;
 
358
  const real G0_2_1 = det*c0_2*c0_1 + det*c1_2*c1_1;
 
359
  const real G0_2_2 = det*c0_2*c0_2 + det*c1_2*c1_2;
 
360
  const real G1_0_0 = det*c0_0*c1_0 + det*c1_0*c0_0;
 
361
  const real G1_0_1 = det*c0_0*c1_1 + det*c1_0*c0_1;
 
362
  const real G1_0_2 = det*c0_0*c1_2 + det*c1_0*c0_2;
 
363
  const real G1_1_0 = det*c0_1*c1_0 + det*c1_1*c0_0;
 
364
  const real G1_1_1 = det*c0_1*c1_1 + det*c1_1*c0_1;
 
365
  const real G1_1_2 = det*c0_1*c1_2 + det*c1_1*c0_2;
 
366
  const real G1_2_0 = det*c0_2*c1_0 + det*c1_2*c0_0;
 
367
  const real G1_2_1 = det*c0_2*c1_1 + det*c1_2*c0_1;
 
368
  const real G1_2_2 = det*c0_2*c1_2 + det*c1_2*c0_2;
 
369
 
 
370
  // Compute element tensor
 
371
  block[0] = 8.333333333333318e-02*G0_0_0 + 4.166666666666659e-02*G0_0_1 + 4.166666666666658e-02*G0_0_2 + 4.166666666666659e-02*G0_1_0 + 8.333333333333318e-02*G0_1_1 + 4.166666666666659e-02*G0_1_2 + 4.166666666666658e-02*G0_2_0 + 4.166666666666659e-02*G0_2_1 + 8.333333333333316e-02*G0_2_2 - 8.333333333333318e-02*G1_0_0 - 4.166666666666659e-02*G1_0_1 - 4.166666666666658e-02*G1_0_2 - 4.166666666666659e-02*G1_1_0 - 8.333333333333318e-02*G1_1_1 - 4.166666666666659e-02*G1_1_2 - 4.166666666666658e-02*G1_2_0 - 4.166666666666659e-02*G1_2_1 - 8.333333333333316e-02*G1_2_2;
 
372
}
 
373
 
 
374
// No contribution from the boundary
 
375
bool LinearForm::boundary_contribution() const { return false; }
 
376
 
 
377
void LinearForm::eval(real block[], const AffineMap& map, real det, unsigned int facet) const {}
 
378
 
 
379
// No contribution from interior boundaries
 
380
bool LinearForm::interior_boundary_contribution() const { return false; }
 
381
 
 
382
void LinearForm::eval(real block[], const AffineMap& map0, const AffineMap& map1, real det, unsigned int facet0, unsigned int facet1, unsigned int alignment) const {}
 
383
 
 
384
} }
 
385
 
 
386
#endif