~ubuntu-branches/ubuntu/raring/ffc/raring
» Revision
8
: test/regression/references/NeumannProblem.h
« back to all changes in this revision
Viewing changes to test/regression/references/NeumannProblem.h
- browse files at revision 8
- compare with another revision
- download diff
- download tarball
- view history from revision 8
- Committer: Bazaar Package Importer
- Author(s): Johannes Ring
- Date: 2010-07-01 19:54:32 UTC
- mfrom: (1.1.5 upstream)
- Revision ID: james.westby@ubuntu.com-20100701195432-xz3gw5nprdj79jcb
Tags: 0.9.3-1
* New upstream release.
* debian/control:
- Minor fix in Vcs fields.
- Bump Standards-Version to 3.9.0 (no changes needed).
- Update version for python-ufc, python-fiat, and python-ufl in
Depends field.
* Switch to dpkg-source 3.0 (quilt) format.
* Update debian/copyright and debian/copyright_hints.
* debian/control:
- Minor fix in Vcs fields.
- Bump Standards-Version to 3.9.0 (no changes needed).
- Update version for python-ufc, python-fiat, and python-ufl in
Depends field.
* Switch to dpkg-source 3.0 (quilt) format.
* Update debian/copyright and debian/copyright_hints.
- files added:
- bench/MassH1_2D_1.ufl
- bench/results
- debian/source
- files removed:
- bench/NavierStokes.form
- files modified:
- ChangeLog
added
removed
test/regression/references/NeumannProblem.h
1
// This code conforms with the UFC specification version 1.4
2
// and was automatically generated by FFC version 0.9.2.
1
// This code conforms with the UFC specification version 1.4.1
2
// and was automatically generated by FFC version 0.9.3.
3
3
//
4
4
// This code was generated with the following parameters:
5
5
//
105
105
106
106
// Reset values.
107
107
*values = 0.00000000;
108
109
// Map degree of freedom to element degree of freedom
110
const unsigned int dof = i;
111
112
// Array of basisvalues.
113
double basisvalues[3] = {0.00000000, 0.00000000, 0.00000000};
114
115
// Declare helper variables.
116
unsigned int rr = 0;
117
unsigned int ss = 0;
118
double tmp0 = (1.00000000 + Y + 2.00000000*X)/2.00000000;
119
120
// Compute basisvalues.
121
basisvalues[0] = 1.00000000;
122
basisvalues[1] = tmp0;
123
for (unsigned int r = 0; r < 1; r++)
124
{
125
rr = (r + 1)*(r + 1 + 1)/2 + 1;
126
ss = r*(r + 1)/2;
127
basisvalues[rr] = basisvalues[ss]*(0.50000000 + r + Y*(1.50000000 + r));
128
}// end loop over 'r'
129
for (unsigned int r = 0; r < 2; r++)
130
{
131
for (unsigned int s = 0; s < 2 - r; s++)
132
{
133
rr = (r + s)*(r + s + 1)/2 + s;
134
basisvalues[rr] *= std::sqrt((0.50000000 + r)*(1.00000000 + r + s));
135
}// end loop over 's'
136
}// end loop over 'r'
137
138
// Table(s) of coefficients.
139
static const double coefficients0[3][3] = \
140
{{0.47140452, -0.28867513, -0.16666667},
141
{0.47140452, 0.28867513, -0.16666667},
142
{0.47140452, 0.00000000, 0.33333333}};
143
144
// Compute value(s).
145
for (unsigned int r = 0; r < 3; r++)
146
{
147
*values += coefficients0[dof][r]*basisvalues[r];
148
}// end loop over 'r'
108
switch (i)
109
{
110
case 0:
111
{
112
113
// Array of basisvalues.
114
double basisvalues[3] = {0.00000000, 0.00000000, 0.00000000};
115
116
// Declare helper variables.
117
unsigned int rr = 0;
118
unsigned int ss = 0;
119
double tmp0 = (1.00000000 + Y + 2.00000000*X)/2.00000000;
120
121
// Compute basisvalues.
122
basisvalues[0] = 1.00000000;
123
basisvalues[1] = tmp0;
124
for (unsigned int r = 0; r < 1; r++)
125
{
126
rr = (r + 1)*(r + 1 + 1)/2 + 1;
127
ss = r*(r + 1)/2;
128
basisvalues[rr] = basisvalues[ss]*(0.50000000 + r + Y*(1.50000000 + r));
129
}// end loop over 'r'
130
for (unsigned int r = 0; r < 2; r++)
131
{
132
for (unsigned int s = 0; s < 2 - r; s++)
133
{
134
rr = (r + s)*(r + s + 1)/2 + s;
135
basisvalues[rr] *= std::sqrt((0.50000000 + r)*(1.00000000 + r + s));
136
}// end loop over 's'
137
}// end loop over 'r'
138
139
// Table(s) of coefficients.
140
static const double coefficients0[3] = \
141
{0.47140452, -0.28867513, -0.16666667};
142
143
// Compute value(s).
144
for (unsigned int r = 0; r < 3; r++)
145
{
146
*values += coefficients0[r]*basisvalues[r];
147
}// end loop over 'r'
148
break;
149
}
150
case 1:
151
{
152
153
// Array of basisvalues.
154
double basisvalues[3] = {0.00000000, 0.00000000, 0.00000000};
155
156
// Declare helper variables.
157
unsigned int rr = 0;
158
unsigned int ss = 0;
159
double tmp0 = (1.00000000 + Y + 2.00000000*X)/2.00000000;
160
161
// Compute basisvalues.
162
basisvalues[0] = 1.00000000;
163
basisvalues[1] = tmp0;
164
for (unsigned int r = 0; r < 1; r++)
165
{
166
rr = (r + 1)*(r + 1 + 1)/2 + 1;
167
ss = r*(r + 1)/2;
168
basisvalues[rr] = basisvalues[ss]*(0.50000000 + r + Y*(1.50000000 + r));
169
}// end loop over 'r'
170
for (unsigned int r = 0; r < 2; r++)
171
{
172
for (unsigned int s = 0; s < 2 - r; s++)
173
{
174
rr = (r + s)*(r + s + 1)/2 + s;
175
basisvalues[rr] *= std::sqrt((0.50000000 + r)*(1.00000000 + r + s));
176
}// end loop over 's'
177
}// end loop over 'r'
178
179
// Table(s) of coefficients.
180
static const double coefficients0[3] = \
181
{0.47140452, 0.28867513, -0.16666667};
182
183
// Compute value(s).
184
for (unsigned int r = 0; r < 3; r++)
185
{
186
*values += coefficients0[r]*basisvalues[r];
187
}// end loop over 'r'
188
break;
189
}
190
case 2:
191
{
192
193
// Array of basisvalues.
194
double basisvalues[3] = {0.00000000, 0.00000000, 0.00000000};
195
196
// Declare helper variables.
197
unsigned int rr = 0;
198
unsigned int ss = 0;
199
double tmp0 = (1.00000000 + Y + 2.00000000*X)/2.00000000;
200
201
// Compute basisvalues.
202
basisvalues[0] = 1.00000000;
203
basisvalues[1] = tmp0;
204
for (unsigned int r = 0; r < 1; r++)
205
{
206
rr = (r + 1)*(r + 1 + 1)/2 + 1;
207
ss = r*(r + 1)/2;
208
basisvalues[rr] = basisvalues[ss]*(0.50000000 + r + Y*(1.50000000 + r));
209
}// end loop over 'r'
210
for (unsigned int r = 0; r < 2; r++)
211
{
212
for (unsigned int s = 0; s < 2 - r; s++)
213
{
214
rr = (r + s)*(r + s + 1)/2 + s;
215
basisvalues[rr] *= std::sqrt((0.50000000 + r)*(1.00000000 + r + s));
216
}// end loop over 's'
217
}// end loop over 'r'
218
219
// Table(s) of coefficients.
220
static const double coefficients0[3] = \
221
{0.47140452, 0.00000000, 0.33333333};
222
223
// Compute value(s).
224
for (unsigned int r = 0; r < 3; r++)
225
{
226
*values += coefficients0[r]*basisvalues[r];
227
}// end loop over 'r'
228
break;
229
}
230
}
231
149
232
}
150
233
151
234
/// Evaluate all basis functions at given point in cell
261
344
values[r] = 0.00000000;
262
345
}// end loop over 'r'
263
346
264
// Map degree of freedom to element degree of freedom
265
const unsigned int dof = i;
266
267
// Array of basisvalues.
268
double basisvalues[3] = {0.00000000, 0.00000000, 0.00000000};
269
270
// Declare helper variables.
271
unsigned int rr = 0;
272
unsigned int ss = 0;
273
double tmp0 = (1.00000000 + Y + 2.00000000*X)/2.00000000;
274
275
// Compute basisvalues.
276
basisvalues[0] = 1.00000000;
277
basisvalues[1] = tmp0;
278
for (unsigned int r = 0; r < 1; r++)
279
{
280
rr = (r + 1)*(r + 1 + 1)/2 + 1;
281
ss = r*(r + 1)/2;
282
basisvalues[rr] = basisvalues[ss]*(0.50000000 + r + Y*(1.50000000 + r));
283
}// end loop over 'r'
284
for (unsigned int r = 0; r < 2; r++)
285
{
286
for (unsigned int s = 0; s < 2 - r; s++)
287
{
288
rr = (r + s)*(r + s + 1)/2 + s;
289
basisvalues[rr] *= std::sqrt((0.50000000 + r)*(1.00000000 + r + s));
290
}// end loop over 's'
291
}// end loop over 'r'
292
293
// Table(s) of coefficients.
294
static const double coefficients0[3][3] = \
295
{{0.47140452, -0.28867513, -0.16666667},
296
{0.47140452, 0.28867513, -0.16666667},
297
{0.47140452, 0.00000000, 0.33333333}};
298
299
// Tables of derivatives of the polynomial base (transpose).
300
static const double dmats0[3][3] = \
301
{{0.00000000, 0.00000000, 0.00000000},
302
{4.89897949, 0.00000000, 0.00000000},
303
{0.00000000, 0.00000000, 0.00000000}};
304
305
static const double dmats1[3][3] = \
306
{{0.00000000, 0.00000000, 0.00000000},
307
{2.44948974, 0.00000000, 0.00000000},
308
{4.24264069, 0.00000000, 0.00000000}};
309
310
// Compute reference derivatives.
311
// Declare pointer to array of derivatives on FIAT element.
312
double *derivatives = new double[num_derivatives];
313
for (unsigned int r = 0; r < num_derivatives; r++)
314
{
315
derivatives[r] = 0.00000000;
316
}// end loop over 'r'
317
318
// Declare derivative matrix (of polynomial basis).
319
double dmats[3][3] = \
320
{{1.00000000, 0.00000000, 0.00000000},
321
{0.00000000, 1.00000000, 0.00000000},
322
{0.00000000, 0.00000000, 1.00000000}};
323
324
// Declare (auxiliary) derivative matrix (of polynomial basis).
325
double dmats_old[3][3] = \
326
{{1.00000000, 0.00000000, 0.00000000},
327
{0.00000000, 1.00000000, 0.00000000},
328
{0.00000000, 0.00000000, 1.00000000}};
329
330
// Loop possible derivatives.
331
for (unsigned int r = 0; r < num_derivatives; r++)
332
{
333
// Resetting dmats values to compute next derivative.
334
for (unsigned int t = 0; t < 3; t++)
335
{
336
for (unsigned int u = 0; u < 3; u++)
337
{
338
dmats[t][u] = 0.00000000;
339
if (t == u)
340
{
341
dmats[t][u] = 1.00000000;
342
}
343
344
}// end loop over 'u'
345
}// end loop over 't'
346
347
// Looping derivative order to generate dmats.
348
for (unsigned int s = 0; s < n; s++)
349
{
350
// Updating dmats_old with new values and resetting dmats.
351
for (unsigned int t = 0; t < 3; t++)
352
{
353
for (unsigned int u = 0; u < 3; u++)
354
{
355
dmats_old[t][u] = dmats[t][u];
356
dmats[t][u] = 0.00000000;
357
}// end loop over 'u'
358
}// end loop over 't'
359
360
// Update dmats using an inner product.
361
if (combinations[r][s] == 0)
362
{
363
for (unsigned int t = 0; t < 3; t++)
364
{
365
for (unsigned int u = 0; u < 3; u++)
366
{
367
for (unsigned int tu = 0; tu < 3; tu++)
368
{
369
dmats[t][u] += dmats0[t][tu]*dmats_old[tu][u];
370
}// end loop over 'tu'
371
}// end loop over 'u'
372
}// end loop over 't'
373
}
374
375
if (combinations[r][s] == 1)
376
{
377
for (unsigned int t = 0; t < 3; t++)
378
{
379
for (unsigned int u = 0; u < 3; u++)
380
{
381
for (unsigned int tu = 0; tu < 3; tu++)
382
{
383
dmats[t][u] += dmats1[t][tu]*dmats_old[tu][u];
384
}// end loop over 'tu'
385
}// end loop over 'u'
386
}// end loop over 't'
387
}
388
389
}// end loop over 's'
390
for (unsigned int s = 0; s < 3; s++)
391
{
392
for (unsigned int t = 0; t < 3; t++)
393
{
394
derivatives[r] += coefficients0[dof][s]*dmats[s][t]*basisvalues[t];
395
}// end loop over 't'
396
}// end loop over 's'
397
}// end loop over 'r'
398
399
// Transform derivatives back to physical element
400
for (unsigned int r = 0; r < num_derivatives; r++)
401
{
402
for (unsigned int s = 0; s < num_derivatives; s++)
403
{
404
values[r] += transform[r][s]*derivatives[s];
405
}// end loop over 's'
406
}// end loop over 'r'
407
408
// Delete pointer to array of derivatives on FIAT element
409
delete [] derivatives;
410
411
// Delete pointer to array of combinations of derivatives and transform
412
for (unsigned int r = 0; r < num_derivatives; r++)
413
{
414
delete [] combinations[r];
415
}// end loop over 'r'
416
delete [] combinations;
417
for (unsigned int r = 0; r < num_derivatives; r++)
418
{
419
delete [] transform[r];
420
}// end loop over 'r'
421
delete [] transform;
347
switch (i)
348
{
349
case 0:
350
{
351
352
// Array of basisvalues.
353
double basisvalues[3] = {0.00000000, 0.00000000, 0.00000000};
354
355
// Declare helper variables.
356
unsigned int rr = 0;
357
unsigned int ss = 0;
358
double tmp0 = (1.00000000 + Y + 2.00000000*X)/2.00000000;
359
360
// Compute basisvalues.
361
basisvalues[0] = 1.00000000;
362
basisvalues[1] = tmp0;
363
for (unsigned int r = 0; r < 1; r++)
364
{
365
rr = (r + 1)*(r + 1 + 1)/2 + 1;
366
ss = r*(r + 1)/2;
367
basisvalues[rr] = basisvalues[ss]*(0.50000000 + r + Y*(1.50000000 + r));
368
}// end loop over 'r'
369
for (unsigned int r = 0; r < 2; r++)
370
{
371
for (unsigned int s = 0; s < 2 - r; s++)
372
{
373
rr = (r + s)*(r + s + 1)/2 + s;
374
basisvalues[rr] *= std::sqrt((0.50000000 + r)*(1.00000000 + r + s));
375
}// end loop over 's'
376
}// end loop over 'r'
377
378
// Table(s) of coefficients.
379
static const double coefficients0[3] = \
380
{0.47140452, -0.28867513, -0.16666667};
381
382
// Tables of derivatives of the polynomial base (transpose).
383
static const double dmats0[3][3] = \
384
{{0.00000000, 0.00000000, 0.00000000},
385
{4.89897949, 0.00000000, 0.00000000},
386
{0.00000000, 0.00000000, 0.00000000}};
387
388
static const double dmats1[3][3] = \
389
{{0.00000000, 0.00000000, 0.00000000},
390
{2.44948974, 0.00000000, 0.00000000},
391
{4.24264069, 0.00000000, 0.00000000}};
392
393
// Compute reference derivatives.
394
// Declare pointer to array of derivatives on FIAT element.
395
double *derivatives = new double[num_derivatives];
396
for (unsigned int r = 0; r < num_derivatives; r++)
397
{
398
derivatives[r] = 0.00000000;
399
}// end loop over 'r'
400
401
// Declare derivative matrix (of polynomial basis).
402
double dmats[3][3] = \
403
{{1.00000000, 0.00000000, 0.00000000},
404
{0.00000000, 1.00000000, 0.00000000},
405
{0.00000000, 0.00000000, 1.00000000}};
406
407
// Declare (auxiliary) derivative matrix (of polynomial basis).
408
double dmats_old[3][3] = \
409
{{1.00000000, 0.00000000, 0.00000000},
410
{0.00000000, 1.00000000, 0.00000000},
411
{0.00000000, 0.00000000, 1.00000000}};
412
413
// Loop possible derivatives.
414
for (unsigned int r = 0; r < num_derivatives; r++)
415
{
416
// Resetting dmats values to compute next derivative.
417
for (unsigned int t = 0; t < 3; t++)
418
{
419
for (unsigned int u = 0; u < 3; u++)
420
{
421
dmats[t][u] = 0.00000000;
422
if (t == u)
423
{
424
dmats[t][u] = 1.00000000;
425
}
426
427
}// end loop over 'u'
428
}// end loop over 't'
429
430
// Looping derivative order to generate dmats.
431
for (unsigned int s = 0; s < n; s++)
432
{
433
// Updating dmats_old with new values and resetting dmats.
434
for (unsigned int t = 0; t < 3; t++)
435
{
436
for (unsigned int u = 0; u < 3; u++)
437
{
438
dmats_old[t][u] = dmats[t][u];
439
dmats[t][u] = 0.00000000;
440
}// end loop over 'u'
441
}// end loop over 't'
442
443
// Update dmats using an inner product.
444
if (combinations[r][s] == 0)
445
{
446
for (unsigned int t = 0; t < 3; t++)
447
{
448
for (unsigned int u = 0; u < 3; u++)
449
{
450
for (unsigned int tu = 0; tu < 3; tu++)
451
{
452
dmats[t][u] += dmats0[t][tu]*dmats_old[tu][u];
453
}// end loop over 'tu'
454
}// end loop over 'u'
455
}// end loop over 't'
456
}
457
458
if (combinations[r][s] == 1)
459
{
460
for (unsigned int t = 0; t < 3; t++)
461
{
462
for (unsigned int u = 0; u < 3; u++)
463
{
464
for (unsigned int tu = 0; tu < 3; tu++)
465
{
466
dmats[t][u] += dmats1[t][tu]*dmats_old[tu][u];
467
}// end loop over 'tu'
468
}// end loop over 'u'
469
}// end loop over 't'
470
}
471
472
}// end loop over 's'
473
for (unsigned int s = 0; s < 3; s++)
474
{
475
for (unsigned int t = 0; t < 3; t++)
476
{
477
derivatives[r] += coefficients0[s]*dmats[s][t]*basisvalues[t];
478
}// end loop over 't'
479
}// end loop over 's'
480
}// end loop over 'r'
481
482
// Transform derivatives back to physical element
483
for (unsigned int r = 0; r < num_derivatives; r++)
484
{
485
for (unsigned int s = 0; s < num_derivatives; s++)
486
{
487
values[r] += transform[r][s]*derivatives[s];
488
}// end loop over 's'
489
}// end loop over 'r'
490
491
// Delete pointer to array of derivatives on FIAT element
492
delete [] derivatives;
493
494
// Delete pointer to array of combinations of derivatives and transform
495
for (unsigned int r = 0; r < num_derivatives; r++)
496
{
497
delete [] combinations[r];
498
}// end loop over 'r'
499
delete [] combinations;
500
for (unsigned int r = 0; r < num_derivatives; r++)
501
{
502
delete [] transform[r];
503
}// end loop over 'r'
504
delete [] transform;
505
break;
506
}
507
case 1:
508
{
509
510
// Array of basisvalues.
511
double basisvalues[3] = {0.00000000, 0.00000000, 0.00000000};
512
513
// Declare helper variables.
514
unsigned int rr = 0;
515
unsigned int ss = 0;
516
double tmp0 = (1.00000000 + Y + 2.00000000*X)/2.00000000;
517
518
// Compute basisvalues.
519
basisvalues[0] = 1.00000000;
520
basisvalues[1] = tmp0;
521
for (unsigned int r = 0; r < 1; r++)
522
{
523
rr = (r + 1)*(r + 1 + 1)/2 + 1;
524
ss = r*(r + 1)/2;
525
basisvalues[rr] = basisvalues[ss]*(0.50000000 + r + Y*(1.50000000 + r));
526
}// end loop over 'r'
527
for (unsigned int r = 0; r < 2; r++)
528
{
529
for (unsigned int s = 0; s < 2 - r; s++)
530
{
531
rr = (r + s)*(r + s + 1)/2 + s;
532
basisvalues[rr] *= std::sqrt((0.50000000 + r)*(1.00000000 + r + s));
533
}// end loop over 's'
534
}// end loop over 'r'
535
536
// Table(s) of coefficients.
537
static const double coefficients0[3] = \
538
{0.47140452, 0.28867513, -0.16666667};
539
540
// Tables of derivatives of the polynomial base (transpose).
541
static const double dmats0[3][3] = \
542
{{0.00000000, 0.00000000, 0.00000000},
543
{4.89897949, 0.00000000, 0.00000000},
544
{0.00000000, 0.00000000, 0.00000000}};
545
546
static const double dmats1[3][3] = \
547
{{0.00000000, 0.00000000, 0.00000000},
548
{2.44948974, 0.00000000, 0.00000000},
549
{4.24264069, 0.00000000, 0.00000000}};
550
551
// Compute reference derivatives.
552
// Declare pointer to array of derivatives on FIAT element.
553
double *derivatives = new double[num_derivatives];
554
for (unsigned int r = 0; r < num_derivatives; r++)
555
{
556
derivatives[r] = 0.00000000;
557
}// end loop over 'r'
558
559
// Declare derivative matrix (of polynomial basis).
560
double dmats[3][3] = \
561
{{1.00000000, 0.00000000, 0.00000000},
562
{0.00000000, 1.00000000, 0.00000000},
563
{0.00000000, 0.00000000, 1.00000000}};
564
565
// Declare (auxiliary) derivative matrix (of polynomial basis).
566
double dmats_old[3][3] = \
567
{{1.00000000, 0.00000000, 0.00000000},
568
{0.00000000, 1.00000000, 0.00000000},
569
{0.00000000, 0.00000000, 1.00000000}};
570
571
// Loop possible derivatives.
572
for (unsigned int r = 0; r < num_derivatives; r++)
573
{
574
// Resetting dmats values to compute next derivative.
575
for (unsigned int t = 0; t < 3; t++)
576
{
577
for (unsigned int u = 0; u < 3; u++)
578
{
579
dmats[t][u] = 0.00000000;
580
if (t == u)
581
{
582
dmats[t][u] = 1.00000000;
583
}
584
585
}// end loop over 'u'
586
}// end loop over 't'
587
588
// Looping derivative order to generate dmats.
589
for (unsigned int s = 0; s < n; s++)
590
{
591
// Updating dmats_old with new values and resetting dmats.
592
for (unsigned int t = 0; t < 3; t++)
593
{
594
for (unsigned int u = 0; u < 3; u++)
595
{
596
dmats_old[t][u] = dmats[t][u];
597
dmats[t][u] = 0.00000000;
598
}// end loop over 'u'
599
}// end loop over 't'
600
601
// Update dmats using an inner product.
602
if (combinations[r][s] == 0)
603
{
604
for (unsigned int t = 0; t < 3; t++)
605
{
606
for (unsigned int u = 0; u < 3; u++)
607
{
608
for (unsigned int tu = 0; tu < 3; tu++)
609
{
610
dmats[t][u] += dmats0[t][tu]*dmats_old[tu][u];
611
}// end loop over 'tu'
612
}// end loop over 'u'
613
}// end loop over 't'
614
}
615
616
if (combinations[r][s] == 1)
617
{
618
for (unsigned int t = 0; t < 3; t++)
619
{
620
for (unsigned int u = 0; u < 3; u++)
621
{
622
for (unsigned int tu = 0; tu < 3; tu++)
623
{
624
dmats[t][u] += dmats1[t][tu]*dmats_old[tu][u];
625
}// end loop over 'tu'
626
}// end loop over 'u'
627
}// end loop over 't'
628
}
629
630
}// end loop over 's'
631
for (unsigned int s = 0; s < 3; s++)
632
{
633
for (unsigned int t = 0; t < 3; t++)
634
{
635
derivatives[r] += coefficients0[s]*dmats[s][t]*basisvalues[t];
636
}// end loop over 't'
637
}// end loop over 's'
638
}// end loop over 'r'
639
640
// Transform derivatives back to physical element
641
for (unsigned int r = 0; r < num_derivatives; r++)
642
{
643
for (unsigned int s = 0; s < num_derivatives; s++)
644
{
645
values[r] += transform[r][s]*derivatives[s];
646
}// end loop over 's'
647
}// end loop over 'r'
648
649
// Delete pointer to array of derivatives on FIAT element
650
delete [] derivatives;
651
652
// Delete pointer to array of combinations of derivatives and transform
653
for (unsigned int r = 0; r < num_derivatives; r++)
654
{
655
delete [] combinations[r];
656
}// end loop over 'r'
657
delete [] combinations;
658
for (unsigned int r = 0; r < num_derivatives; r++)
659
{
660
delete [] transform[r];
661
}// end loop over 'r'
662
delete [] transform;
663
break;
664
}
665
case 2:
666
{
667
668
// Array of basisvalues.
669
double basisvalues[3] = {0.00000000, 0.00000000, 0.00000000};
670
671
// Declare helper variables.
672
unsigned int rr = 0;
673
unsigned int ss = 0;
674
double tmp0 = (1.00000000 + Y + 2.00000000*X)/2.00000000;
675
676
// Compute basisvalues.
677
basisvalues[0] = 1.00000000;
678
basisvalues[1] = tmp0;
679
for (unsigned int r = 0; r < 1; r++)
680
{
681
rr = (r + 1)*(r + 1 + 1)/2 + 1;
682
ss = r*(r + 1)/2;
683
basisvalues[rr] = basisvalues[ss]*(0.50000000 + r + Y*(1.50000000 + r));
684
}// end loop over 'r'
685
for (unsigned int r = 0; r < 2; r++)
686
{
687
for (unsigned int s = 0; s < 2 - r; s++)
688
{
689
rr = (r + s)*(r + s + 1)/2 + s;
690
basisvalues[rr] *= std::sqrt((0.50000000 + r)*(1.00000000 + r + s));
691
}// end loop over 's'
692
}// end loop over 'r'
693
694
// Table(s) of coefficients.
695
static const double coefficients0[3] = \
696
{0.47140452, 0.00000000, 0.33333333};
697
698
// Tables of derivatives of the polynomial base (transpose).
699
static const double dmats0[3][3] = \
700
{{0.00000000, 0.00000000, 0.00000000},
701
{4.89897949, 0.00000000, 0.00000000},
702
{0.00000000, 0.00000000, 0.00000000}};
703
704
static const double dmats1[3][3] = \
705
{{0.00000000, 0.00000000, 0.00000000},
706
{2.44948974, 0.00000000, 0.00000000},
707
{4.24264069, 0.00000000, 0.00000000}};
708
709
// Compute reference derivatives.
710
// Declare pointer to array of derivatives on FIAT element.
711
double *derivatives = new double[num_derivatives];
712
for (unsigned int r = 0; r < num_derivatives; r++)
713
{
714
derivatives[r] = 0.00000000;
715
}// end loop over 'r'
716
717
// Declare derivative matrix (of polynomial basis).
718
double dmats[3][3] = \
719
{{1.00000000, 0.00000000, 0.00000000},
720
{0.00000000, 1.00000000, 0.00000000},
721
{0.00000000, 0.00000000, 1.00000000}};
722
723
// Declare (auxiliary) derivative matrix (of polynomial basis).
724
double dmats_old[3][3] = \
725
{{1.00000000, 0.00000000, 0.00000000},
726
{0.00000000, 1.00000000, 0.00000000},
727
{0.00000000, 0.00000000, 1.00000000}};
728
729
// Loop possible derivatives.
730
for (unsigned int r = 0; r < num_derivatives; r++)
731
{
732
// Resetting dmats values to compute next derivative.
733
for (unsigned int t = 0; t < 3; t++)
734
{
735
for (unsigned int u = 0; u < 3; u++)
736
{
737
dmats[t][u] = 0.00000000;
738
if (t == u)
739
{
740
dmats[t][u] = 1.00000000;
741
}
742
743
}// end loop over 'u'
744
}// end loop over 't'
745
746
// Looping derivative order to generate dmats.
747
for (unsigned int s = 0; s < n; s++)
748
{
749
// Updating dmats_old with new values and resetting dmats.
750
for (unsigned int t = 0; t < 3; t++)
751
{
752
for (unsigned int u = 0; u < 3; u++)
753
{
754
dmats_old[t][u] = dmats[t][u];
755
dmats[t][u] = 0.00000000;
756
}// end loop over 'u'
757
}// end loop over 't'
758
759
// Update dmats using an inner product.
760
if (combinations[r][s] == 0)
761
{
762
for (unsigned int t = 0; t < 3; t++)
763
{
764
for (unsigned int u = 0; u < 3; u++)
765
{
766
for (unsigned int tu = 0; tu < 3; tu++)
767
{
768
dmats[t][u] += dmats0[t][tu]*dmats_old[tu][u];
769
}// end loop over 'tu'
770
}// end loop over 'u'
771
}// end loop over 't'
772
}
773
774
if (combinations[r][s] == 1)
775
{
776
for (unsigned int t = 0; t < 3; t++)
777
{
778
for (unsigned int u = 0; u < 3; u++)
779
{
780
for (unsigned int tu = 0; tu < 3; tu++)
781
{
782
dmats[t][u] += dmats1[t][tu]*dmats_old[tu][u];
783
}// end loop over 'tu'
784
}// end loop over 'u'
785
}// end loop over 't'
786
}
787
788
}// end loop over 's'
789
for (unsigned int s = 0; s < 3; s++)
790
{
791
for (unsigned int t = 0; t < 3; t++)
792
{
793
derivatives[r] += coefficients0[s]*dmats[s][t]*basisvalues[t];
794
}// end loop over 't'
795
}// end loop over 's'
796
}// end loop over 'r'
797
798
// Transform derivatives back to physical element
799
for (unsigned int r = 0; r < num_derivatives; r++)
800
{
801
for (unsigned int s = 0; s < num_derivatives; s++)
802
{
803
values[r] += transform[r][s]*derivatives[s];
804
}// end loop over 's'
805
}// end loop over 'r'
806
807
// Delete pointer to array of derivatives on FIAT element
808
delete [] derivatives;
809
810
// Delete pointer to array of combinations of derivatives and transform
811
for (unsigned int r = 0; r < num_derivatives; r++)
812
{
813
delete [] combinations[r];
814
}// end loop over 'r'
815
delete [] combinations;
816
for (unsigned int r = 0; r < num_derivatives; r++)
817
{
818
delete [] transform[r];
819
}// end loop over 'r'
820
delete [] transform;
821
break;
822
}
823
}
824
422
825
}
423
826
424
827
/// Evaluate order n derivatives of all basis functions at given point in cell
635
1038
// Reset values.
636
1039
values[0] = 0.00000000;
637
1040
values[1] = 0.00000000;
638
if (0 <= i && i <= 2)
639
{
640
// Map degree of freedom to element degree of freedom
641
const unsigned int dof = i;
642
643
// Array of basisvalues.
644
double basisvalues[3] = {0.00000000, 0.00000000, 0.00000000};
645
646
// Declare helper variables.
647
unsigned int rr = 0;
648
unsigned int ss = 0;
649
double tmp0 = (1.00000000 + Y + 2.00000000*X)/2.00000000;
650
651
// Compute basisvalues.
652
basisvalues[0] = 1.00000000;
653
basisvalues[1] = tmp0;
654
for (unsigned int r = 0; r < 1; r++)
655
{
656
rr = (r + 1)*(r + 1 + 1)/2 + 1;
657
ss = r*(r + 1)/2;
658
basisvalues[rr] = basisvalues[ss]*(0.50000000 + r + Y*(1.50000000 + r));
659
}// end loop over 'r'
660
for (unsigned int r = 0; r < 2; r++)
661
{
662
for (unsigned int s = 0; s < 2 - r; s++)
663
{
664
rr = (r + s)*(r + s + 1)/2 + s;
665
basisvalues[rr] *= std::sqrt((0.50000000 + r)*(1.00000000 + r + s));
666
}// end loop over 's'
667
}// end loop over 'r'
668
669
// Table(s) of coefficients.
670
static const double coefficients0[3][3] = \
671
{{0.47140452, -0.28867513, -0.16666667},
672
{0.47140452, 0.28867513, -0.16666667},
673
{0.47140452, 0.00000000, 0.33333333}};
674
675
// Compute value(s).
676
for (unsigned int r = 0; r < 3; r++)
677
{
678
values[0] += coefficients0[dof][r]*basisvalues[r];
679
}// end loop over 'r'
680
}
681
682
if (3 <= i && i <= 5)
683
{
684
// Map degree of freedom to element degree of freedom
685
const unsigned int dof = i - 3;
686
687
// Array of basisvalues.
688
double basisvalues[3] = {0.00000000, 0.00000000, 0.00000000};
689
690
// Declare helper variables.
691
unsigned int rr = 0;
692
unsigned int ss = 0;
693
double tmp0 = (1.00000000 + Y + 2.00000000*X)/2.00000000;
694
695
// Compute basisvalues.
696
basisvalues[0] = 1.00000000;
697
basisvalues[1] = tmp0;
698
for (unsigned int r = 0; r < 1; r++)
699
{
700
rr = (r + 1)*(r + 1 + 1)/2 + 1;
701
ss = r*(r + 1)/2;
702
basisvalues[rr] = basisvalues[ss]*(0.50000000 + r + Y*(1.50000000 + r));
703
}// end loop over 'r'
704
for (unsigned int r = 0; r < 2; r++)
705
{
706
for (unsigned int s = 0; s < 2 - r; s++)
707
{
708
rr = (r + s)*(r + s + 1)/2 + s;
709
basisvalues[rr] *= std::sqrt((0.50000000 + r)*(1.00000000 + r + s));
710
}// end loop over 's'
711
}// end loop over 'r'
712
713
// Table(s) of coefficients.
714
static const double coefficients0[3][3] = \
715
{{0.47140452, -0.28867513, -0.16666667},
716
{0.47140452, 0.28867513, -0.16666667},
717
{0.47140452, 0.00000000, 0.33333333}};
718
719
// Compute value(s).
720
for (unsigned int r = 0; r < 3; r++)
721
{
722
values[1] += coefficients0[dof][r]*basisvalues[r];
723
}// end loop over 'r'
1041
switch (i)
1042
{
1043
case 0:
1044
{
1045
1046
// Array of basisvalues.
1047
double basisvalues[3] = {0.00000000, 0.00000000, 0.00000000};
1048
1049
// Declare helper variables.
1050
unsigned int rr = 0;
1051
unsigned int ss = 0;
1052
double tmp0 = (1.00000000 + Y + 2.00000000*X)/2.00000000;
1053
1054
// Compute basisvalues.
1055
basisvalues[0] = 1.00000000;
1056
basisvalues[1] = tmp0;
1057
for (unsigned int r = 0; r < 1; r++)
1058
{
1059
rr = (r + 1)*(r + 1 + 1)/2 + 1;
1060
ss = r*(r + 1)/2;
1061
basisvalues[rr] = basisvalues[ss]*(0.50000000 + r + Y*(1.50000000 + r));
1062
}// end loop over 'r'
1063
for (unsigned int r = 0; r < 2; r++)
1064
{
1065
for (unsigned int s = 0; s < 2 - r; s++)
1066
{
1067
rr = (r + s)*(r + s + 1)/2 + s;
1068
basisvalues[rr] *= std::sqrt((0.50000000 + r)*(1.00000000 + r + s));
1069
}// end loop over 's'
1070
}// end loop over 'r'
1071
1072
// Table(s) of coefficients.
1073
static const double coefficients0[3] = \
1074
{0.47140452, -0.28867513, -0.16666667};
1075
1076
// Compute value(s).
1077
for (unsigned int r = 0; r < 3; r++)
1078
{
1079
values[0] += coefficients0[r]*basisvalues[r];
1080
}// end loop over 'r'
1081
break;
1082
}
1083
case 1:
1084
{
1085
1086
// Array of basisvalues.
1087
double basisvalues[3] = {0.00000000, 0.00000000, 0.00000000};
1088
1089
// Declare helper variables.
1090
unsigned int rr = 0;
1091
unsigned int ss = 0;
1092
double tmp0 = (1.00000000 + Y + 2.00000000*X)/2.00000000;
1093
1094
// Compute basisvalues.
1095
basisvalues[0] = 1.00000000;
1096
basisvalues[1] = tmp0;
1097
for (unsigned int r = 0; r < 1; r++)
1098
{
1099
rr = (r + 1)*(r + 1 + 1)/2 + 1;
1100
ss = r*(r + 1)/2;
1101
basisvalues[rr] = basisvalues[ss]*(0.50000000 + r + Y*(1.50000000 + r));
1102
}// end loop over 'r'
1103
for (unsigned int r = 0; r < 2; r++)
1104
{
1105
for (unsigned int s = 0; s < 2 - r; s++)
1106
{
1107
rr = (r + s)*(r + s + 1)/2 + s;
1108
basisvalues[rr] *= std::sqrt((0.50000000 + r)*(1.00000000 + r + s));
1109
}// end loop over 's'
1110
}// end loop over 'r'
1111
1112
// Table(s) of coefficients.
1113
static const double coefficients0[3] = \
1114
{0.47140452, 0.28867513, -0.16666667};
1115
1116
// Compute value(s).
1117
for (unsigned int r = 0; r < 3; r++)
1118
{
1119
values[0] += coefficients0[r]*basisvalues[r];
1120
}// end loop over 'r'
1121
break;
1122
}
1123
case 2:
1124
{
1125
1126
// Array of basisvalues.
1127
double basisvalues[3] = {0.00000000, 0.00000000, 0.00000000};
1128
1129
// Declare helper variables.
1130
unsigned int rr = 0;
1131
unsigned int ss = 0;
1132
double tmp0 = (1.00000000 + Y + 2.00000000*X)/2.00000000;
1133
1134
// Compute basisvalues.
1135
basisvalues[0] = 1.00000000;
1136
basisvalues[1] = tmp0;
1137
for (unsigned int r = 0; r < 1; r++)
1138
{
1139
rr = (r + 1)*(r + 1 + 1)/2 + 1;
1140
ss = r*(r + 1)/2;
1141
basisvalues[rr] = basisvalues[ss]*(0.50000000 + r + Y*(1.50000000 + r));
1142
}// end loop over 'r'
1143
for (unsigned int r = 0; r < 2; r++)
1144
{
1145
for (unsigned int s = 0; s < 2 - r; s++)
1146
{
1147
rr = (r + s)*(r + s + 1)/2 + s;
1148
basisvalues[rr] *= std::sqrt((0.50000000 + r)*(1.00000000 + r + s));
1149
}// end loop over 's'
1150
}// end loop over 'r'
1151
1152
// Table(s) of coefficients.
1153
static const double coefficients0[3] = \
1154
{0.47140452, 0.00000000, 0.33333333};
1155
1156
// Compute value(s).
1157
for (unsigned int r = 0; r < 3; r++)
1158
{
1159
values[0] += coefficients0[r]*basisvalues[r];
1160
}// end loop over 'r'
1161
break;
1162
}
1163
case 3:
1164
{
1165
1166
// Array of basisvalues.
1167
double basisvalues[3] = {0.00000000, 0.00000000, 0.00000000};
1168
1169
// Declare helper variables.
1170
unsigned int rr = 0;
1171
unsigned int ss = 0;
1172
double tmp0 = (1.00000000 + Y + 2.00000000*X)/2.00000000;
1173
1174
// Compute basisvalues.
1175
basisvalues[0] = 1.00000000;
1176
basisvalues[1] = tmp0;
1177
for (unsigned int r = 0; r < 1; r++)
1178
{
1179
rr = (r + 1)*(r + 1 + 1)/2 + 1;
1180
ss = r*(r + 1)/2;
1181
basisvalues[rr] = basisvalues[ss]*(0.50000000 + r + Y*(1.50000000 + r));
1182
}// end loop over 'r'
1183
for (unsigned int r = 0; r < 2; r++)
1184
{
1185
for (unsigned int s = 0; s < 2 - r; s++)
1186
{
1187
rr = (r + s)*(r + s + 1)/2 + s;
1188
basisvalues[rr] *= std::sqrt((0.50000000 + r)*(1.00000000 + r + s));
1189
}// end loop over 's'
1190
}// end loop over 'r'
1191
1192
// Table(s) of coefficients.
1193
static const double coefficients0[3] = \
1194
{0.47140452, -0.28867513, -0.16666667};
1195
1196
// Compute value(s).
1197
for (unsigned int r = 0; r < 3; r++)
1198
{
1199
values[1] += coefficients0[r]*basisvalues[r];
1200
}// end loop over 'r'
1201
break;
1202
}
1203
case 4:
1204
{
1205
1206
// Array of basisvalues.
1207
double basisvalues[3] = {0.00000000, 0.00000000, 0.00000000};
1208
1209
// Declare helper variables.
1210
unsigned int rr = 0;
1211
unsigned int ss = 0;
1212
double tmp0 = (1.00000000 + Y + 2.00000000*X)/2.00000000;
1213
1214
// Compute basisvalues.
1215
basisvalues[0] = 1.00000000;
1216
basisvalues[1] = tmp0;
1217
for (unsigned int r = 0; r < 1; r++)
1218
{
1219
rr = (r + 1)*(r + 1 + 1)/2 + 1;
1220
ss = r*(r + 1)/2;
1221
basisvalues[rr] = basisvalues[ss]*(0.50000000 + r + Y*(1.50000000 + r));
1222
}// end loop over 'r'
1223
for (unsigned int r = 0; r < 2; r++)
1224
{
1225
for (unsigned int s = 0; s < 2 - r; s++)
1226
{
1227
rr = (r + s)*(r + s + 1)/2 + s;
1228
basisvalues[rr] *= std::sqrt((0.50000000 + r)*(1.00000000 + r + s));
1229
}// end loop over 's'
1230
}// end loop over 'r'
1231
1232
// Table(s) of coefficients.
1233
static const double coefficients0[3] = \
1234
{0.47140452, 0.28867513, -0.16666667};
1235
1236
// Compute value(s).
1237
for (unsigned int r = 0; r < 3; r++)
1238
{
1239
values[1] += coefficients0[r]*basisvalues[r];
1240
}// end loop over 'r'
1241
break;
1242
}
1243
case 5:
1244
{
1245
1246
// Array of basisvalues.
1247
double basisvalues[3] = {0.00000000, 0.00000000, 0.00000000};
1248
1249
// Declare helper variables.
1250
unsigned int rr = 0;
1251
unsigned int ss = 0;
1252
double tmp0 = (1.00000000 + Y + 2.00000000*X)/2.00000000;
1253
1254
// Compute basisvalues.
1255
basisvalues[0] = 1.00000000;
1256
basisvalues[1] = tmp0;
1257
for (unsigned int r = 0; r < 1; r++)
1258
{
1259
rr = (r + 1)*(r + 1 + 1)/2 + 1;
1260
ss = r*(r + 1)/2;
1261
basisvalues[rr] = basisvalues[ss]*(0.50000000 + r + Y*(1.50000000 + r));
1262
}// end loop over 'r'
1263
for (unsigned int r = 0; r < 2; r++)
1264
{
1265
for (unsigned int s = 0; s < 2 - r; s++)
1266
{
1267
rr = (r + s)*(r + s + 1)/2 + s;
1268
basisvalues[rr] *= std::sqrt((0.50000000 + r)*(1.00000000 + r + s));
1269
}// end loop over 's'
1270
}// end loop over 'r'
1271
1272
// Table(s) of coefficients.
1273
static const double coefficients0[3] = \
1274
{0.47140452, 0.00000000, 0.33333333};
1275
1276
// Compute value(s).
1277
for (unsigned int r = 0; r < 3; r++)
1278
{
1279
values[1] += coefficients0[r]*basisvalues[r];
1280
}// end loop over 'r'
1281
break;
1282
}
724
1283
}
725
1284
726
1285
}
841
1400
values[r] = 0.00000000;
842
1401
}// end loop over 'r'
843
1402
844
if (0 <= i && i <= 2)
845
{
846
// Map degree of freedom to element degree of freedom
847
const unsigned int dof = i;
848
849
// Array of basisvalues.
850
double basisvalues[3] = {0.00000000, 0.00000000, 0.00000000};
851
852
// Declare helper variables.
853
unsigned int rr = 0;
854
unsigned int ss = 0;
855
double tmp0 = (1.00000000 + Y + 2.00000000*X)/2.00000000;
856
857
// Compute basisvalues.
858
basisvalues[0] = 1.00000000;
859
basisvalues[1] = tmp0;
860
for (unsigned int r = 0; r < 1; r++)
861
{
862
rr = (r + 1)*(r + 1 + 1)/2 + 1;
863
ss = r*(r + 1)/2;
864
basisvalues[rr] = basisvalues[ss]*(0.50000000 + r + Y*(1.50000000 + r));
865
}// end loop over 'r'
866
for (unsigned int r = 0; r < 2; r++)
867
{
868
for (unsigned int s = 0; s < 2 - r; s++)
869
{
870
rr = (r + s)*(r + s + 1)/2 + s;
871
basisvalues[rr] *= std::sqrt((0.50000000 + r)*(1.00000000 + r + s));
872
}// end loop over 's'
873
}// end loop over 'r'
874
875
// Table(s) of coefficients.
876
static const double coefficients0[3][3] = \
877
{{0.47140452, -0.28867513, -0.16666667},
878
{0.47140452, 0.28867513, -0.16666667},
879
{0.47140452, 0.00000000, 0.33333333}};
880
881
// Tables of derivatives of the polynomial base (transpose).
882
static const double dmats0[3][3] = \
883
{{0.00000000, 0.00000000, 0.00000000},
884
{4.89897949, 0.00000000, 0.00000000},
885
{0.00000000, 0.00000000, 0.00000000}};
886
887
static const double dmats1[3][3] = \
888
{{0.00000000, 0.00000000, 0.00000000},
889
{2.44948974, 0.00000000, 0.00000000},
890
{4.24264069, 0.00000000, 0.00000000}};
891
892
// Compute reference derivatives.
893
// Declare pointer to array of derivatives on FIAT element.
894
double *derivatives = new double[num_derivatives];
895
for (unsigned int r = 0; r < num_derivatives; r++)
896
{
897
derivatives[r] = 0.00000000;
898
}// end loop over 'r'
899
900
// Declare derivative matrix (of polynomial basis).
901
double dmats[3][3] = \
902
{{1.00000000, 0.00000000, 0.00000000},
903
{0.00000000, 1.00000000, 0.00000000},
904
{0.00000000, 0.00000000, 1.00000000}};
905
906
// Declare (auxiliary) derivative matrix (of polynomial basis).
907
double dmats_old[3][3] = \
908
{{1.00000000, 0.00000000, 0.00000000},
909
{0.00000000, 1.00000000, 0.00000000},
910
{0.00000000, 0.00000000, 1.00000000}};
911
912
// Loop possible derivatives.
913
for (unsigned int r = 0; r < num_derivatives; r++)
914
{
915
// Resetting dmats values to compute next derivative.
916
for (unsigned int t = 0; t < 3; t++)
917
{
918
for (unsigned int u = 0; u < 3; u++)
919
{
920
dmats[t][u] = 0.00000000;
921
if (t == u)
922
{
923
dmats[t][u] = 1.00000000;
924
}
925
926
}// end loop over 'u'
927
}// end loop over 't'
928
929
// Looping derivative order to generate dmats.
930
for (unsigned int s = 0; s < n; s++)
931
{
932
// Updating dmats_old with new values and resetting dmats.
933
for (unsigned int t = 0; t < 3; t++)
934
{
935
for (unsigned int u = 0; u < 3; u++)
936
{
937
dmats_old[t][u] = dmats[t][u];
938
dmats[t][u] = 0.00000000;
939
}// end loop over 'u'
940
}// end loop over 't'
941
942
// Update dmats using an inner product.
943
if (combinations[r][s] == 0)
944
{
945
for (unsigned int t = 0; t < 3; t++)
946
{
947
for (unsigned int u = 0; u < 3; u++)
948
{
949
for (unsigned int tu = 0; tu < 3; tu++)
950
{
951
dmats[t][u] += dmats0[t][tu]*dmats_old[tu][u];
952
}// end loop over 'tu'
953
}// end loop over 'u'
954
}// end loop over 't'
955
}
956
957
if (combinations[r][s] == 1)
958
{
959
for (unsigned int t = 0; t < 3; t++)
960
{
961
for (unsigned int u = 0; u < 3; u++)
962
{
963
for (unsigned int tu = 0; tu < 3; tu++)
964
{
965
dmats[t][u] += dmats1[t][tu]*dmats_old[tu][u];
966
}// end loop over 'tu'
967
}// end loop over 'u'
968
}// end loop over 't'
969
}
970
971
}// end loop over 's'
972
for (unsigned int s = 0; s < 3; s++)
973
{
974
for (unsigned int t = 0; t < 3; t++)
975
{
976
derivatives[r] += coefficients0[dof][s]*dmats[s][t]*basisvalues[t];
977
}// end loop over 't'
978
}// end loop over 's'
979
}// end loop over 'r'
980
981
// Transform derivatives back to physical element
982
for (unsigned int r = 0; r < num_derivatives; r++)
983
{
984
for (unsigned int s = 0; s < num_derivatives; s++)
985
{
986
values[r] += transform[r][s]*derivatives[s];
987
}// end loop over 's'
988
}// end loop over 'r'
989
990
// Delete pointer to array of derivatives on FIAT element
991
delete [] derivatives;
992
993
// Delete pointer to array of combinations of derivatives and transform
994
for (unsigned int r = 0; r < num_derivatives; r++)
995
{
996
delete [] combinations[r];
997
}// end loop over 'r'
998
delete [] combinations;
999
for (unsigned int r = 0; r < num_derivatives; r++)
1000
{
1001
delete [] transform[r];
1002
}// end loop over 'r'
1003
delete [] transform;
1004
}
1005
1006
if (3 <= i && i <= 5)
1007
{
1008
// Map degree of freedom to element degree of freedom
1009
const unsigned int dof = i - 3;
1010
1011
// Array of basisvalues.
1012
double basisvalues[3] = {0.00000000, 0.00000000, 0.00000000};
1013
1014
// Declare helper variables.
1015
unsigned int rr = 0;
1016
unsigned int ss = 0;
1017
double tmp0 = (1.00000000 + Y + 2.00000000*X)/2.00000000;
1018
1019
// Compute basisvalues.
1020
basisvalues[0] = 1.00000000;
1021
basisvalues[1] = tmp0;
1022
for (unsigned int r = 0; r < 1; r++)
1023
{
1024
rr = (r + 1)*(r + 1 + 1)/2 + 1;
1025
ss = r*(r + 1)/2;
1026
basisvalues[rr] = basisvalues[ss]*(0.50000000 + r + Y*(1.50000000 + r));
1027
}// end loop over 'r'
1028
for (unsigned int r = 0; r < 2; r++)
1029
{
1030
for (unsigned int s = 0; s < 2 - r; s++)
1031
{
1032
rr = (r + s)*(r + s + 1)/2 + s;
1033
basisvalues[rr] *= std::sqrt((0.50000000 + r)*(1.00000000 + r + s));
1034
}// end loop over 's'
1035
}// end loop over 'r'
1036
1037
// Table(s) of coefficients.
1038
static const double coefficients0[3][3] = \
1039
{{0.47140452, -0.28867513, -0.16666667},
1040
{0.47140452, 0.28867513, -0.16666667},
1041
{0.47140452, 0.00000000, 0.33333333}};
1042
1043
// Tables of derivatives of the polynomial base (transpose).
1044
static const double dmats0[3][3] = \
1045
{{0.00000000, 0.00000000, 0.00000000},
1046
{4.89897949, 0.00000000, 0.00000000},
1047
{0.00000000, 0.00000000, 0.00000000}};
1048
1049
static const double dmats1[3][3] = \
1050
{{0.00000000, 0.00000000, 0.00000000},
1051
{2.44948974, 0.00000000, 0.00000000},
1052
{4.24264069, 0.00000000, 0.00000000}};
1053
1054
// Compute reference derivatives.
1055
// Declare pointer to array of derivatives on FIAT element.
1056
double *derivatives = new double[num_derivatives];
1057
for (unsigned int r = 0; r < num_derivatives; r++)
1058
{
1059
derivatives[r] = 0.00000000;
1060
}// end loop over 'r'
1061
1062
// Declare derivative matrix (of polynomial basis).
1063
double dmats[3][3] = \
1064
{{1.00000000, 0.00000000, 0.00000000},
1065
{0.00000000, 1.00000000, 0.00000000},
1066
{0.00000000, 0.00000000, 1.00000000}};
1067
1068
// Declare (auxiliary) derivative matrix (of polynomial basis).
1069
double dmats_old[3][3] = \
1070
{{1.00000000, 0.00000000, 0.00000000},
1071
{0.00000000, 1.00000000, 0.00000000},
1072
{0.00000000, 0.00000000, 1.00000000}};
1073
1074
// Loop possible derivatives.
1075
for (unsigned int r = 0; r < num_derivatives; r++)
1076
{
1077
// Resetting dmats values to compute next derivative.
1078
for (unsigned int t = 0; t < 3; t++)
1079
{
1080
for (unsigned int u = 0; u < 3; u++)
1081
{
1082
dmats[t][u] = 0.00000000;
1083
if (t == u)
1084
{
1085
dmats[t][u] = 1.00000000;
1086
}
1087
1088
}// end loop over 'u'
1089
}// end loop over 't'
1090
1091
// Looping derivative order to generate dmats.
1092
for (unsigned int s = 0; s < n; s++)
1093
{
1094
// Updating dmats_old with new values and resetting dmats.
1095
for (unsigned int t = 0; t < 3; t++)
1096
{
1097
for (unsigned int u = 0; u < 3; u++)
1098
{
1099
dmats_old[t][u] = dmats[t][u];
1100
dmats[t][u] = 0.00000000;
1101
}// end loop over 'u'
1102
}// end loop over 't'
1103
1104
// Update dmats using an inner product.
1105
if (combinations[r][s] == 0)
1106
{
1107
for (unsigned int t = 0; t < 3; t++)
1108
{
1109
for (unsigned int u = 0; u < 3; u++)
1110
{
1111
for (unsigned int tu = 0; tu < 3; tu++)
1112
{
1113
dmats[t][u] += dmats0[t][tu]*dmats_old[tu][u];
1114
}// end loop over 'tu'
1115
}// end loop over 'u'
1116
}// end loop over 't'
1117
}
1118
1119
if (combinations[r][s] == 1)
1120
{
1121
for (unsigned int t = 0; t < 3; t++)
1122
{
1123
for (unsigned int u = 0; u < 3; u++)
1124
{
1125
for (unsigned int tu = 0; tu < 3; tu++)
1126
{
1127
dmats[t][u] += dmats1[t][tu]*dmats_old[tu][u];
1128
}// end loop over 'tu'
1129
}// end loop over 'u'
1130
}// end loop over 't'
1131
}
1132
1133
}// end loop over 's'
1134
for (unsigned int s = 0; s < 3; s++)
1135
{
1136
for (unsigned int t = 0; t < 3; t++)
1137
{
1138
derivatives[r] += coefficients0[dof][s]*dmats[s][t]*basisvalues[t];
1139
}// end loop over 't'
1140
}// end loop over 's'
1141
}// end loop over 'r'
1142
1143
// Transform derivatives back to physical element
1144
for (unsigned int r = 0; r < num_derivatives; r++)
1145
{
1146
for (unsigned int s = 0; s < num_derivatives; s++)
1147
{
1148
values[num_derivatives + r] += transform[r][s]*derivatives[s];
1149
}// end loop over 's'
1150
}// end loop over 'r'
1151
1152
// Delete pointer to array of derivatives on FIAT element
1153
delete [] derivatives;
1154
1155
// Delete pointer to array of combinations of derivatives and transform
1156
for (unsigned int r = 0; r < num_derivatives; r++)
1157
{
1158
delete [] combinations[r];
1159
}// end loop over 'r'
1160
delete [] combinations;
1161
for (unsigned int r = 0; r < num_derivatives; r++)
1162
{
1163
delete [] transform[r];
1164
}// end loop over 'r'
1165
delete [] transform;
1403
switch (i)
1404
{
1405
case 0:
1406
{
1407
1408
// Array of basisvalues.
1409
double basisvalues[3] = {0.00000000, 0.00000000, 0.00000000};
1410
1411
// Declare helper variables.
1412
unsigned int rr = 0;
1413
unsigned int ss = 0;
1414
double tmp0 = (1.00000000 + Y + 2.00000000*X)/2.00000000;
1415
1416
// Compute basisvalues.
1417
basisvalues[0] = 1.00000000;
1418
basisvalues[1] = tmp0;
1419
for (unsigned int r = 0; r < 1; r++)
1420
{
1421
rr = (r + 1)*(r + 1 + 1)/2 + 1;
1422
ss = r*(r + 1)/2;
1423
basisvalues[rr] = basisvalues[ss]*(0.50000000 + r + Y*(1.50000000 + r));
1424
}// end loop over 'r'
1425
for (unsigned int r = 0; r < 2; r++)
1426
{
1427
for (unsigned int s = 0; s < 2 - r; s++)
1428
{
1429
rr = (r + s)*(r + s + 1)/2 + s;
1430
basisvalues[rr] *= std::sqrt((0.50000000 + r)*(1.00000000 + r + s));
1431
}// end loop over 's'
1432
}// end loop over 'r'
1433
1434
// Table(s) of coefficients.
1435
static const double coefficients0[3] = \
1436
{0.47140452, -0.28867513, -0.16666667};
1437
1438
// Tables of derivatives of the polynomial base (transpose).
1439
static const double dmats0[3][3] = \
1440
{{0.00000000, 0.00000000, 0.00000000},
1441
{4.89897949, 0.00000000, 0.00000000},
1442
{0.00000000, 0.00000000, 0.00000000}};
1443
1444
static const double dmats1[3][3] = \
1445
{{0.00000000, 0.00000000, 0.00000000},
1446
{2.44948974, 0.00000000, 0.00000000},
1447
{4.24264069, 0.00000000, 0.00000000}};
1448
1449
// Compute reference derivatives.
1450
// Declare pointer to array of derivatives on FIAT element.
1451
double *derivatives = new double[num_derivatives];
1452
for (unsigned int r = 0; r < num_derivatives; r++)
1453
{
1454
derivatives[r] = 0.00000000;
1455
}// end loop over 'r'
1456
1457
// Declare derivative matrix (of polynomial basis).
1458
double dmats[3][3] = \
1459
{{1.00000000, 0.00000000, 0.00000000},
1460
{0.00000000, 1.00000000, 0.00000000},
1461
{0.00000000, 0.00000000, 1.00000000}};
1462
1463
// Declare (auxiliary) derivative matrix (of polynomial basis).
1464
double dmats_old[3][3] = \
1465
{{1.00000000, 0.00000000, 0.00000000},
1466
{0.00000000, 1.00000000, 0.00000000},
1467
{0.00000000, 0.00000000, 1.00000000}};
1468
1469
// Loop possible derivatives.
1470
for (unsigned int r = 0; r < num_derivatives; r++)
1471
{
1472
// Resetting dmats values to compute next derivative.
1473
for (unsigned int t = 0; t < 3; t++)
1474
{
1475
for (unsigned int u = 0; u < 3; u++)
1476
{
1477
dmats[t][u] = 0.00000000;
1478
if (t == u)
1479
{
1480
dmats[t][u] = 1.00000000;
1481
}
1482
1483
}// end loop over 'u'
1484
}// end loop over 't'
1485
1486
// Looping derivative order to generate dmats.
1487
for (unsigned int s = 0; s < n; s++)
1488
{
1489
// Updating dmats_old with new values and resetting dmats.
1490
for (unsigned int t = 0; t < 3; t++)
1491
{
1492
for (unsigned int u = 0; u < 3; u++)
1493
{
1494
dmats_old[t][u] = dmats[t][u];
1495
dmats[t][u] = 0.00000000;
1496
}// end loop over 'u'
1497
}// end loop over 't'
1498
1499
// Update dmats using an inner product.
1500
if (combinations[r][s] == 0)
1501
{
1502
for (unsigned int t = 0; t < 3; t++)
1503
{
1504
for (unsigned int u = 0; u < 3; u++)
1505
{
1506
for (unsigned int tu = 0; tu < 3; tu++)
1507
{
1508
dmats[t][u] += dmats0[t][tu]*dmats_old[tu][u];
1509
}// end loop over 'tu'
1510
}// end loop over 'u'
1511
}// end loop over 't'
1512
}
1513
1514
if (combinations[r][s] == 1)
1515
{
1516
for (unsigned int t = 0; t < 3; t++)
1517
{
1518
for (unsigned int u = 0; u < 3; u++)
1519
{
1520
for (unsigned int tu = 0; tu < 3; tu++)
1521
{
1522
dmats[t][u] += dmats1[t][tu]*dmats_old[tu][u];
1523
}// end loop over 'tu'
1524
}// end loop over 'u'
1525
}// end loop over 't'
1526
}
1527
1528
}// end loop over 's'
1529
for (unsigned int s = 0; s < 3; s++)
1530
{
1531
for (unsigned int t = 0; t < 3; t++)
1532
{
1533
derivatives[r] += coefficients0[s]*dmats[s][t]*basisvalues[t];
1534
}// end loop over 't'
1535
}// end loop over 's'
1536
}// end loop over 'r'
1537
1538
// Transform derivatives back to physical element
1539
for (unsigned int r = 0; r < num_derivatives; r++)
1540
{
1541
for (unsigned int s = 0; s < num_derivatives; s++)
1542
{
1543
values[r] += transform[r][s]*derivatives[s];
1544
}// end loop over 's'
1545
}// end loop over 'r'
1546
1547
// Delete pointer to array of derivatives on FIAT element
1548
delete [] derivatives;
1549
1550
// Delete pointer to array of combinations of derivatives and transform
1551
for (unsigned int r = 0; r < num_derivatives; r++)
1552
{
1553
delete [] combinations[r];
1554
}// end loop over 'r'
1555
delete [] combinations;
1556
for (unsigned int r = 0; r < num_derivatives; r++)
1557
{
1558
delete [] transform[r];
1559
}// end loop over 'r'
1560
delete [] transform;
1561
break;
1562
}
1563
case 1:
1564
{
1565
1566
// Array of basisvalues.
1567
double basisvalues[3] = {0.00000000, 0.00000000, 0.00000000};
1568
1569
// Declare helper variables.
1570
unsigned int rr = 0;
1571
unsigned int ss = 0;
1572
double tmp0 = (1.00000000 + Y + 2.00000000*X)/2.00000000;
1573
1574
// Compute basisvalues.
1575
basisvalues[0] = 1.00000000;
1576
basisvalues[1] = tmp0;
1577
for (unsigned int r = 0; r < 1; r++)
1578
{
1579
rr = (r + 1)*(r + 1 + 1)/2 + 1;
1580
ss = r*(r + 1)/2;
1581
basisvalues[rr] = basisvalues[ss]*(0.50000000 + r + Y*(1.50000000 + r));
1582
}// end loop over 'r'
1583
for (unsigned int r = 0; r < 2; r++)
1584
{
1585
for (unsigned int s = 0; s < 2 - r; s++)
1586
{
1587
rr = (r + s)*(r + s + 1)/2 + s;
1588
basisvalues[rr] *= std::sqrt((0.50000000 + r)*(1.00000000 + r + s));
1589
}// end loop over 's'
1590
}// end loop over 'r'
1591
1592
// Table(s) of coefficients.
1593
static const double coefficients0[3] = \
1594
{0.47140452, 0.28867513, -0.16666667};
1595
1596
// Tables of derivatives of the polynomial base (transpose).
1597
static const double dmats0[3][3] = \
1598
{{0.00000000, 0.00000000, 0.00000000},
1599
{4.89897949, 0.00000000, 0.00000000},
1600
{0.00000000, 0.00000000, 0.00000000}};
1601
1602
static const double dmats1[3][3] = \
1603
{{0.00000000, 0.00000000, 0.00000000},
1604
{2.44948974, 0.00000000, 0.00000000},
1605
{4.24264069, 0.00000000, 0.00000000}};
1606
1607
// Compute reference derivatives.
1608
// Declare pointer to array of derivatives on FIAT element.
1609
double *derivatives = new double[num_derivatives];
1610
for (unsigned int r = 0; r < num_derivatives; r++)
1611
{
1612
derivatives[r] = 0.00000000;
1613
}// end loop over 'r'
1614
1615
// Declare derivative matrix (of polynomial basis).
1616
double dmats[3][3] = \
1617
{{1.00000000, 0.00000000, 0.00000000},
1618
{0.00000000, 1.00000000, 0.00000000},
1619
{0.00000000, 0.00000000, 1.00000000}};
1620
1621
// Declare (auxiliary) derivative matrix (of polynomial basis).
1622
double dmats_old[3][3] = \
1623
{{1.00000000, 0.00000000, 0.00000000},
1624
{0.00000000, 1.00000000, 0.00000000},
1625
{0.00000000, 0.00000000, 1.00000000}};
1626
1627
// Loop possible derivatives.
1628
for (unsigned int r = 0; r < num_derivatives; r++)
1629
{
1630
// Resetting dmats values to compute next derivative.
1631
for (unsigned int t = 0; t < 3; t++)
1632
{
1633
for (unsigned int u = 0; u < 3; u++)
1634
{
1635
dmats[t][u] = 0.00000000;
1636
if (t == u)
1637
{
1638
dmats[t][u] = 1.00000000;
1639
}
1640
1641
}// end loop over 'u'
1642
}// end loop over 't'
1643
1644
// Looping derivative order to generate dmats.
1645
for (unsigned int s = 0; s < n; s++)
1646
{
1647
// Updating dmats_old with new values and resetting dmats.
1648
for (unsigned int t = 0; t < 3; t++)
1649
{
1650
for (unsigned int u = 0; u < 3; u++)
1651
{
1652
dmats_old[t][u] = dmats[t][u];
1653
dmats[t][u] = 0.00000000;
1654
}// end loop over 'u'
1655
}// end loop over 't'
1656
1657
// Update dmats using an inner product.
1658
if (combinations[r][s] == 0)
1659
{
1660
for (unsigned int t = 0; t < 3; t++)
1661
{
1662
for (unsigned int u = 0; u < 3; u++)
1663
{
1664
for (unsigned int tu = 0; tu < 3; tu++)
1665
{
1666
dmats[t][u] += dmats0[t][tu]*dmats_old[tu][u];
1667
}// end loop over 'tu'
1668
}// end loop over 'u'
1669
}// end loop over 't'
1670
}
1671
1672
if (combinations[r][s] == 1)
1673
{
1674
for (unsigned int t = 0; t < 3; t++)
1675
{
1676
for (unsigned int u = 0; u < 3; u++)
1677
{
1678
for (unsigned int tu = 0; tu < 3; tu++)
1679
{
1680
dmats[t][u] += dmats1[t][tu]*dmats_old[tu][u];
1681
}// end loop over 'tu'
1682
}// end loop over 'u'
1683
}// end loop over 't'
1684
}
1685
1686
}// end loop over 's'
1687
for (unsigned int s = 0; s < 3; s++)
1688
{
1689
for (unsigned int t = 0; t < 3; t++)
1690
{
1691
derivatives[r] += coefficients0[s]*dmats[s][t]*basisvalues[t];
1692
}// end loop over 't'
1693
}// end loop over 's'
1694
}// end loop over 'r'
1695
1696
// Transform derivatives back to physical element
1697
for (unsigned int r = 0; r < num_derivatives; r++)
1698
{
1699
for (unsigned int s = 0; s < num_derivatives; s++)
1700
{
1701
values[r] += transform[r][s]*derivatives[s];
1702
}// end loop over 's'
1703
}// end loop over 'r'
1704
1705
// Delete pointer to array of derivatives on FIAT element
1706
delete [] derivatives;
1707
1708
// Delete pointer to array of combinations of derivatives and transform
1709
for (unsigned int r = 0; r < num_derivatives; r++)
1710
{
1711
delete [] combinations[r];
1712
}// end loop over 'r'
1713
delete [] combinations;
1714
for (unsigned int r = 0; r < num_derivatives; r++)
1715
{
1716
delete [] transform[r];
1717
}// end loop over 'r'
1718
delete [] transform;
1719
break;
1720
}
1721
case 2:
1722
{
1723
1724
// Array of basisvalues.
1725
double basisvalues[3] = {0.00000000, 0.00000000, 0.00000000};
1726
1727
// Declare helper variables.
1728
unsigned int rr = 0;
1729
unsigned int ss = 0;
1730
double tmp0 = (1.00000000 + Y + 2.00000000*X)/2.00000000;
1731
1732
// Compute basisvalues.
1733
basisvalues[0] = 1.00000000;
1734
basisvalues[1] = tmp0;
1735
for (unsigned int r = 0; r < 1; r++)
1736
{
1737
rr = (r + 1)*(r + 1 + 1)/2 + 1;
1738
ss = r*(r + 1)/2;
1739
basisvalues[rr] = basisvalues[ss]*(0.50000000 + r + Y*(1.50000000 + r));
1740
}// end loop over 'r'
1741
for (unsigned int r = 0; r < 2; r++)
1742
{
1743
for (unsigned int s = 0; s < 2 - r; s++)
1744
{
1745
rr = (r + s)*(r + s + 1)/2 + s;
1746
basisvalues[rr] *= std::sqrt((0.50000000 + r)*(1.00000000 + r + s));
1747
}// end loop over 's'
1748
}// end loop over 'r'
1749
1750
// Table(s) of coefficients.
1751
static const double coefficients0[3] = \
1752
{0.47140452, 0.00000000, 0.33333333};
1753
1754
// Tables of derivatives of the polynomial base (transpose).
1755
static const double dmats0[3][3] = \
1756
{{0.00000000, 0.00000000, 0.00000000},
1757
{4.89897949, 0.00000000, 0.00000000},
1758
{0.00000000, 0.00000000, 0.00000000}};
1759
1760
static const double dmats1[3][3] = \
1761
{{0.00000000, 0.00000000, 0.00000000},
1762
{2.44948974, 0.00000000, 0.00000000},
1763
{4.24264069, 0.00000000, 0.00000000}};
1764
1765
// Compute reference derivatives.
1766
// Declare pointer to array of derivatives on FIAT element.
1767
double *derivatives = new double[num_derivatives];
1768
for (unsigned int r = 0; r < num_derivatives; r++)
1769
{
1770
derivatives[r] = 0.00000000;
1771
}// end loop over 'r'
1772
1773
// Declare derivative matrix (of polynomial basis).
1774
double dmats[3][3] = \
1775
{{1.00000000, 0.00000000, 0.00000000},
1776
{0.00000000, 1.00000000, 0.00000000},
1777
{0.00000000, 0.00000000, 1.00000000}};
1778
1779
// Declare (auxiliary) derivative matrix (of polynomial basis).
1780
double dmats_old[3][3] = \
1781
{{1.00000000, 0.00000000, 0.00000000},
1782
{0.00000000, 1.00000000, 0.00000000},
1783
{0.00000000, 0.00000000, 1.00000000}};
1784
1785
// Loop possible derivatives.
1786
for (unsigned int r = 0; r < num_derivatives; r++)
1787
{
1788
// Resetting dmats values to compute next derivative.
1789
for (unsigned int t = 0; t < 3; t++)
1790
{
1791
for (unsigned int u = 0; u < 3; u++)
1792
{
1793
dmats[t][u] = 0.00000000;
1794
if (t == u)
1795
{
1796
dmats[t][u] = 1.00000000;
1797
}
1798
1799
}// end loop over 'u'
1800
}// end loop over 't'
1801
1802
// Looping derivative order to generate dmats.
1803
for (unsigned int s = 0; s < n; s++)
1804
{
1805
// Updating dmats_old with new values and resetting dmats.
1806
for (unsigned int t = 0; t < 3; t++)
1807
{
1808
for (unsigned int u = 0; u < 3; u++)
1809
{
1810
dmats_old[t][u] = dmats[t][u];
1811
dmats[t][u] = 0.00000000;
1812
}// end loop over 'u'
1813
}// end loop over 't'
1814
1815
// Update dmats using an inner product.
1816
if (combinations[r][s] == 0)
1817
{
1818
for (unsigned int t = 0; t < 3; t++)
1819
{
1820
for (unsigned int u = 0; u < 3; u++)
1821
{
1822
for (unsigned int tu = 0; tu < 3; tu++)
1823
{
1824
dmats[t][u] += dmats0[t][tu]*dmats_old[tu][u];
1825
}// end loop over 'tu'
1826
}// end loop over 'u'
1827
}// end loop over 't'
1828
}
1829
1830
if (combinations[r][s] == 1)
1831
{
1832
for (unsigned int t = 0; t < 3; t++)
1833
{
1834
for (unsigned int u = 0; u < 3; u++)
1835
{
1836
for (unsigned int tu = 0; tu < 3; tu++)
1837
{
1838
dmats[t][u] += dmats1[t][tu]*dmats_old[tu][u];
1839
}// end loop over 'tu'
1840
}// end loop over 'u'
1841
}// end loop over 't'
1842
}
1843
1844
}// end loop over 's'
1845
for (unsigned int s = 0; s < 3; s++)
1846
{
1847
for (unsigned int t = 0; t < 3; t++)
1848
{
1849
derivatives[r] += coefficients0[s]*dmats[s][t]*basisvalues[t];
1850
}// end loop over 't'
1851
}// end loop over 's'
1852
}// end loop over 'r'
1853
1854
// Transform derivatives back to physical element
1855
for (unsigned int r = 0; r < num_derivatives; r++)
1856
{
1857
for (unsigned int s = 0; s < num_derivatives; s++)
1858
{
1859
values[r] += transform[r][s]*derivatives[s];
1860
}// end loop over 's'
1861
}// end loop over 'r'
1862
1863
// Delete pointer to array of derivatives on FIAT element
1864
delete [] derivatives;
1865
1866
// Delete pointer to array of combinations of derivatives and transform
1867
for (unsigned int r = 0; r < num_derivatives; r++)
1868
{
1869
delete [] combinations[r];
1870
}// end loop over 'r'
1871
delete [] combinations;
1872
for (unsigned int r = 0; r < num_derivatives; r++)
1873
{
1874
delete [] transform[r];
1875
}// end loop over 'r'
1876
delete [] transform;
1877
break;
1878
}
1879
case 3:
1880
{
1881
1882
// Array of basisvalues.
1883
double basisvalues[3] = {0.00000000, 0.00000000, 0.00000000};
1884
1885
// Declare helper variables.
1886
unsigned int rr = 0;
1887
unsigned int ss = 0;
1888
double tmp0 = (1.00000000 + Y + 2.00000000*X)/2.00000000;
1889
1890
// Compute basisvalues.
1891
basisvalues[0] = 1.00000000;
1892
basisvalues[1] = tmp0;
1893
for (unsigned int r = 0; r < 1; r++)
1894
{
1895
rr = (r + 1)*(r + 1 + 1)/2 + 1;
1896
ss = r*(r + 1)/2;
1897
basisvalues[rr] = basisvalues[ss]*(0.50000000 + r + Y*(1.50000000 + r));
1898
}// end loop over 'r'
1899
for (unsigned int r = 0; r < 2; r++)
1900
{
1901
for (unsigned int s = 0; s < 2 - r; s++)
1902
{
1903
rr = (r + s)*(r + s + 1)/2 + s;
1904
basisvalues[rr] *= std::sqrt((0.50000000 + r)*(1.00000000 + r + s));
1905
}// end loop over 's'
1906
}// end loop over 'r'
1907
1908
// Table(s) of coefficients.
1909
static const double coefficients0[3] = \
1910
{0.47140452, -0.28867513, -0.16666667};
1911
1912
// Tables of derivatives of the polynomial base (transpose).
1913
static const double dmats0[3][3] = \
1914
{{0.00000000, 0.00000000, 0.00000000},
1915
{4.89897949, 0.00000000, 0.00000000},
1916
{0.00000000, 0.00000000, 0.00000000}};
1917
1918
static const double dmats1[3][3] = \
1919
{{0.00000000, 0.00000000, 0.00000000},
1920
{2.44948974, 0.00000000, 0.00000000},
1921
{4.24264069, 0.00000000, 0.00000000}};
1922
1923
// Compute reference derivatives.
1924
// Declare pointer to array of derivatives on FIAT element.
1925
double *derivatives = new double[num_derivatives];
1926
for (unsigned int r = 0; r < num_derivatives; r++)
1927
{
1928
derivatives[r] = 0.00000000;
1929
}// end loop over 'r'
1930
1931
// Declare derivative matrix (of polynomial basis).
1932
double dmats[3][3] = \
1933
{{1.00000000, 0.00000000, 0.00000000},
1934
{0.00000000, 1.00000000, 0.00000000},
1935
{0.00000000, 0.00000000, 1.00000000}};
1936
1937
// Declare (auxiliary) derivative matrix (of polynomial basis).
1938
double dmats_old[3][3] = \
1939
{{1.00000000, 0.00000000, 0.00000000},
1940
{0.00000000, 1.00000000, 0.00000000},
1941
{0.00000000, 0.00000000, 1.00000000}};
1942
1943
// Loop possible derivatives.
1944
for (unsigned int r = 0; r < num_derivatives; r++)
1945
{
1946
// Resetting dmats values to compute next derivative.
1947
for (unsigned int t = 0; t < 3; t++)
1948
{
1949
for (unsigned int u = 0; u < 3; u++)
1950
{
1951
dmats[t][u] = 0.00000000;
1952
if (t == u)
1953
{
1954
dmats[t][u] = 1.00000000;
1955
}
1956
1957
}// end loop over 'u'
1958
}// end loop over 't'
1959
1960
// Looping derivative order to generate dmats.
1961
for (unsigned int s = 0; s < n; s++)
1962
{
1963
// Updating dmats_old with new values and resetting dmats.
1964
for (unsigned int t = 0; t < 3; t++)
1965
{
1966
for (unsigned int u = 0; u < 3; u++)
1967
{
1968
dmats_old[t][u] = dmats[t][u];
1969
dmats[t][u] = 0.00000000;
1970
}// end loop over 'u'
1971
}// end loop over 't'
1972
1973
// Update dmats using an inner product.
1974
if (combinations[r][s] == 0)
1975
{
1976
for (unsigned int t = 0; t < 3; t++)
1977
{
1978
for (unsigned int u = 0; u < 3; u++)
1979
{
1980
for (unsigned int tu = 0; tu < 3; tu++)
1981
{
1982
dmats[t][u] += dmats0[t][tu]*dmats_old[tu][u];
1983
}// end loop over 'tu'
1984
}// end loop over 'u'
1985
}// end loop over 't'
1986
}
1987
1988
if (combinations[r][s] == 1)
1989
{
1990
for (unsigned int t = 0; t < 3; t++)
1991
{
1992
for (unsigned int u = 0; u < 3; u++)
1993
{
1994
for (unsigned int tu = 0; tu < 3; tu++)
1995
{
1996
dmats[t][u] += dmats1[t][tu]*dmats_old[tu][u];
1997
}// end loop over 'tu'
1998
}// end loop over 'u'
1999
}// end loop over 't'
2000
}
2001
2002
}// end loop over 's'
2003
for (unsigned int s = 0; s < 3; s++)
2004
{
2005
for (unsigned int t = 0; t < 3; t++)
2006
{
2007
derivatives[r] += coefficients0[s]*dmats[s][t]*basisvalues[t];
2008
}// end loop over 't'
2009
}// end loop over 's'
2010
}// end loop over 'r'
2011
2012
// Transform derivatives back to physical element
2013
for (unsigned int r = 0; r < num_derivatives; r++)
2014
{
2015
for (unsigned int s = 0; s < num_derivatives; s++)
2016
{
2017
values[num_derivatives + r] += transform[r][s]*derivatives[s];
2018
}// end loop over 's'
2019
}// end loop over 'r'
2020
2021
// Delete pointer to array of derivatives on FIAT element
2022
delete [] derivatives;
2023
2024
// Delete pointer to array of combinations of derivatives and transform
2025
for (unsigned int r = 0; r < num_derivatives; r++)
2026
{
2027
delete [] combinations[r];
2028
}// end loop over 'r'
2029
delete [] combinations;
2030
for (unsigned int r = 0; r < num_derivatives; r++)
2031
{
2032
delete [] transform[r];
2033
}// end loop over 'r'
2034
delete [] transform;
2035
break;
2036
}
2037
case 4:
2038
{
2039
2040
// Array of basisvalues.
2041
double basisvalues[3] = {0.00000000, 0.00000000, 0.00000000};
2042
2043
// Declare helper variables.
2044
unsigned int rr = 0;
2045
unsigned int ss = 0;
2046
double tmp0 = (1.00000000 + Y + 2.00000000*X)/2.00000000;
2047
2048
// Compute basisvalues.
2049
basisvalues[0] = 1.00000000;
2050
basisvalues[1] = tmp0;
2051
for (unsigned int r = 0; r < 1; r++)
2052
{
2053
rr = (r + 1)*(r + 1 + 1)/2 + 1;
2054
ss = r*(r + 1)/2;
2055
basisvalues[rr] = basisvalues[ss]*(0.50000000 + r + Y*(1.50000000 + r));
2056
}// end loop over 'r'
2057
for (unsigned int r = 0; r < 2; r++)
2058
{
2059
for (unsigned int s = 0; s < 2 - r; s++)
2060
{
2061
rr = (r + s)*(r + s + 1)/2 + s;
2062
basisvalues[rr] *= std::sqrt((0.50000000 + r)*(1.00000000 + r + s));
2063
}// end loop over 's'
2064
}// end loop over 'r'
2065
2066
// Table(s) of coefficients.
2067
static const double coefficients0[3] = \
2068
{0.47140452, 0.28867513, -0.16666667};
2069
2070
// Tables of derivatives of the polynomial base (transpose).
2071
static const double dmats0[3][3] = \
2072
{{0.00000000, 0.00000000, 0.00000000},
2073
{4.89897949, 0.00000000, 0.00000000},
2074
{0.00000000, 0.00000000, 0.00000000}};
2075
2076
static const double dmats1[3][3] = \
2077
{{0.00000000, 0.00000000, 0.00000000},
2078
{2.44948974, 0.00000000, 0.00000000},
2079
{4.24264069, 0.00000000, 0.00000000}};
2080
2081
// Compute reference derivatives.
2082
// Declare pointer to array of derivatives on FIAT element.
2083
double *derivatives = new double[num_derivatives];
2084
for (unsigned int r = 0; r < num_derivatives; r++)
2085
{
2086
derivatives[r] = 0.00000000;
2087
}// end loop over 'r'
2088
2089
// Declare derivative matrix (of polynomial basis).
2090
double dmats[3][3] = \
2091
{{1.00000000, 0.00000000, 0.00000000},
2092
{0.00000000, 1.00000000, 0.00000000},
2093
{0.00000000, 0.00000000, 1.00000000}};
2094
2095
// Declare (auxiliary) derivative matrix (of polynomial basis).
2096
double dmats_old[3][3] = \
2097
{{1.00000000, 0.00000000, 0.00000000},
2098
{0.00000000, 1.00000000, 0.00000000},
2099
{0.00000000, 0.00000000, 1.00000000}};
2100
2101
// Loop possible derivatives.
2102
for (unsigned int r = 0; r < num_derivatives; r++)
2103
{
2104
// Resetting dmats values to compute next derivative.
2105
for (unsigned int t = 0; t < 3; t++)
2106
{
2107
for (unsigned int u = 0; u < 3; u++)
2108
{
2109
dmats[t][u] = 0.00000000;
2110
if (t == u)
2111
{
2112
dmats[t][u] = 1.00000000;
2113
}
2114
2115
}// end loop over 'u'
2116
}// end loop over 't'
2117
2118
// Looping derivative order to generate dmats.
2119
for (unsigned int s = 0; s < n; s++)
2120
{
2121
// Updating dmats_old with new values and resetting dmats.
2122
for (unsigned int t = 0; t < 3; t++)
2123
{
2124
for (unsigned int u = 0; u < 3; u++)
2125
{
2126
dmats_old[t][u] = dmats[t][u];
2127
dmats[t][u] = 0.00000000;
2128
}// end loop over 'u'
2129
}// end loop over 't'
2130
2131
// Update dmats using an inner product.
2132
if (combinations[r][s] == 0)
2133
{
2134
for (unsigned int t = 0; t < 3; t++)
2135
{
2136
for (unsigned int u = 0; u < 3; u++)
2137
{
2138
for (unsigned int tu = 0; tu < 3; tu++)
2139
{
2140
dmats[t][u] += dmats0[t][tu]*dmats_old[tu][u];
2141
}// end loop over 'tu'
2142
}// end loop over 'u'
2143
}// end loop over 't'
2144
}
2145
2146
if (combinations[r][s] == 1)
2147
{
2148
for (unsigned int t = 0; t < 3; t++)
2149
{
2150
for (unsigned int u = 0; u < 3; u++)
2151
{
2152
for (unsigned int tu = 0; tu < 3; tu++)
2153
{
2154
dmats[t][u] += dmats1[t][tu]*dmats_old[tu][u];
2155
}// end loop over 'tu'
2156
}// end loop over 'u'
2157
}// end loop over 't'
2158
}
2159
2160
}// end loop over 's'
2161
for (unsigned int s = 0; s < 3; s++)
2162
{
2163
for (unsigned int t = 0; t < 3; t++)
2164
{
2165
derivatives[r] += coefficients0[s]*dmats[s][t]*basisvalues[t];
2166
}// end loop over 't'
2167
}// end loop over 's'
2168
}// end loop over 'r'
2169
2170
// Transform derivatives back to physical element
2171
for (unsigned int r = 0; r < num_derivatives; r++)
2172
{
2173
for (unsigned int s = 0; s < num_derivatives; s++)
2174
{
2175
values[num_derivatives + r] += transform[r][s]*derivatives[s];
2176
}// end loop over 's'
2177
}// end loop over 'r'
2178
2179
// Delete pointer to array of derivatives on FIAT element
2180
delete [] derivatives;
2181
2182
// Delete pointer to array of combinations of derivatives and transform
2183
for (unsigned int r = 0; r < num_derivatives; r++)
2184
{
2185
delete [] combinations[r];
2186
}// end loop over 'r'
2187
delete [] combinations;
2188
for (unsigned int r = 0; r < num_derivatives; r++)
2189
{
2190
delete [] transform[r];
2191
}// end loop over 'r'
2192
delete [] transform;
2193
break;
2194
}
2195
case 5:
2196
{
2197
2198
// Array of basisvalues.
2199
double basisvalues[3] = {0.00000000, 0.00000000, 0.00000000};
2200
2201
// Declare helper variables.
2202
unsigned int rr = 0;
2203
unsigned int ss = 0;
2204
double tmp0 = (1.00000000 + Y + 2.00000000*X)/2.00000000;
2205
2206
// Compute basisvalues.
2207
basisvalues[0] = 1.00000000;
2208
basisvalues[1] = tmp0;
2209
for (unsigned int r = 0; r < 1; r++)
2210
{
2211
rr = (r + 1)*(r + 1 + 1)/2 + 1;
2212
ss = r*(r + 1)/2;
2213
basisvalues[rr] = basisvalues[ss]*(0.50000000 + r + Y*(1.50000000 + r));
2214
}// end loop over 'r'
2215
for (unsigned int r = 0; r < 2; r++)
2216
{
2217
for (unsigned int s = 0; s < 2 - r; s++)
2218
{
2219
rr = (r + s)*(r + s + 1)/2 + s;
2220
basisvalues[rr] *= std::sqrt((0.50000000 + r)*(1.00000000 + r + s));
2221
}// end loop over 's'
2222
}// end loop over 'r'
2223
2224
// Table(s) of coefficients.
2225
static const double coefficients0[3] = \
2226
{0.47140452, 0.00000000, 0.33333333};
2227
2228
// Tables of derivatives of the polynomial base (transpose).
2229
static const double dmats0[3][3] = \
2230
{{0.00000000, 0.00000000, 0.00000000},
2231
{4.89897949, 0.00000000, 0.00000000},
2232
{0.00000000, 0.00000000, 0.00000000}};
2233
2234
static const double dmats1[3][3] = \
2235
{{0.00000000, 0.00000000, 0.00000000},
2236
{2.44948974, 0.00000000, 0.00000000},
2237
{4.24264069, 0.00000000, 0.00000000}};
2238
2239
// Compute reference derivatives.
2240
// Declare pointer to array of derivatives on FIAT element.
2241
double *derivatives = new double[num_derivatives];
2242
for (unsigned int r = 0; r < num_derivatives; r++)
2243
{
2244
derivatives[r] = 0.00000000;
2245
}// end loop over 'r'
2246
2247
// Declare derivative matrix (of polynomial basis).
2248
double dmats[3][3] = \
2249
{{1.00000000, 0.00000000, 0.00000000},
2250
{0.00000000, 1.00000000, 0.00000000},
2251
{0.00000000, 0.00000000, 1.00000000}};
2252
2253
// Declare (auxiliary) derivative matrix (of polynomial basis).
2254
double dmats_old[3][3] = \
2255
{{1.00000000, 0.00000000, 0.00000000},
2256
{0.00000000, 1.00000000, 0.00000000},
2257
{0.00000000, 0.00000000, 1.00000000}};
2258
2259
// Loop possible derivatives.
2260
for (unsigned int r = 0; r < num_derivatives; r++)
2261
{
2262
// Resetting dmats values to compute next derivative.
2263
for (unsigned int t = 0; t < 3; t++)
2264
{
2265
for (unsigned int u = 0; u < 3; u++)
2266
{
2267
dmats[t][u] = 0.00000000;
2268
if (t == u)
2269
{
2270
dmats[t][u] = 1.00000000;
2271
}
2272
2273
}// end loop over 'u'
2274
}// end loop over 't'
2275
2276
// Looping derivative order to generate dmats.
2277
for (unsigned int s = 0; s < n; s++)
2278
{
2279
// Updating dmats_old with new values and resetting dmats.
2280
for (unsigned int t = 0; t < 3; t++)
2281
{
2282
for (unsigned int u = 0; u < 3; u++)
2283
{
2284
dmats_old[t][u] = dmats[t][u];
2285
dmats[t][u] = 0.00000000;
2286
}// end loop over 'u'
2287
}// end loop over 't'
2288
2289
// Update dmats using an inner product.
2290
if (combinations[r][s] == 0)
2291
{
2292
for (unsigned int t = 0; t < 3; t++)
2293
{
2294
for (unsigned int u = 0; u < 3; u++)
2295
{
2296
for (unsigned int tu = 0; tu < 3; tu++)
2297
{
2298
dmats[t][u] += dmats0[t][tu]*dmats_old[tu][u];
2299
}// end loop over 'tu'
2300
}// end loop over 'u'
2301
}// end loop over 't'
2302
}
2303
2304
if (combinations[r][s] == 1)
2305
{
2306
for (unsigned int t = 0; t < 3; t++)
2307
{
2308
for (unsigned int u = 0; u < 3; u++)
2309
{
2310
for (unsigned int tu = 0; tu < 3; tu++)
2311
{
2312
dmats[t][u] += dmats1[t][tu]*dmats_old[tu][u];
2313
}// end loop over 'tu'
2314
}// end loop over 'u'
2315
}// end loop over 't'
2316
}
2317
2318
}// end loop over 's'
2319
for (unsigned int s = 0; s < 3; s++)
2320
{
2321
for (unsigned int t = 0; t < 3; t++)
2322
{
2323
derivatives[r] += coefficients0[s]*dmats[s][t]*basisvalues[t];
2324
}// end loop over 't'
2325
}// end loop over 's'
2326
}// end loop over 'r'
2327
2328
// Transform derivatives back to physical element
2329
for (unsigned int r = 0; r < num_derivatives; r++)
2330
{
2331
for (unsigned int s = 0; s < num_derivatives; s++)
2332
{
2333
values[num_derivatives + r] += transform[r][s]*derivatives[s];
2334
}// end loop over 's'
2335
}// end loop over 'r'
2336
2337
// Delete pointer to array of derivatives on FIAT element
2338
delete [] derivatives;
2339
2340
// Delete pointer to array of combinations of derivatives and transform
2341
for (unsigned int r = 0; r < num_derivatives; r++)
2342
{
2343
delete [] combinations[r];
2344
}// end loop over 'r'
2345
delete [] combinations;
2346
for (unsigned int r = 0; r < num_derivatives; r++)
2347
{
2348
delete [] transform[r];
2349
}// end loop over 'r'
2350
delete [] transform;
2351
break;
2352
}
1166
2353
}
1167
2354
1168
2355
}
Loggerhead is a web-based interface for Breezy
Version: 2.0.1