64
64
/// Return the topological dimension of the cell shape
65
virtual unsigned int topological_dimension() const
65
virtual std::size_t topological_dimension() const
70
70
/// Return the geometric dimension of the cell shape
71
virtual unsigned int geometric_dimension() const
71
virtual std::size_t geometric_dimension() const
76
76
/// Return the dimension of the finite element function space
77
virtual unsigned int space_dimension() const
77
virtual std::size_t space_dimension() const
82
82
/// Return the rank of the value space
83
virtual unsigned int value_rank() const
83
virtual std::size_t value_rank() const
88
88
/// Return the dimension of the value space for axis i
89
virtual unsigned int value_dimension(unsigned int i) const
89
virtual std::size_t value_dimension(std::size_t i) const
94
/// Evaluate basis function i at given point in cell
95
virtual void evaluate_basis(unsigned int i,
94
/// Evaluate basis function i at given point x in cell
95
virtual void evaluate_basis(std::size_t i,
97
const double* coordinates,
98
const ufc::cell& c) const
98
const double* vertex_coordinates,
99
int cell_orientation) const
100
// Extract vertex coordinates
102
// Compute Jacobian of affine map from reference cell
104
// Compute determinant of Jacobian
106
// Compute inverse of Jacobian
103
compute_jacobian_triangle_2d(J, vertex_coordinates);
105
// Compute Jacobian inverse and determinant
108
compute_jacobian_inverse_triangle_2d(K, detJ, J);
108
111
// Compute constants
110
113
// Get coordinates and map to the reference (FIAT) element
115
// Array of basisvalues.
118
// Array of basisvalues
116
119
double basisvalues[1] = {0.0};
118
// Declare helper variables.
121
// Declare helper variables
120
// Compute basisvalues.
123
// Compute basisvalues
121
124
basisvalues[0] = 1.0;
123
// Table(s) of coefficients.
126
// Table(s) of coefficients
124
127
static const double coefficients0[1] = \
128
131
for (unsigned int r = 0; r < 1; r++)
130
133
*values += coefficients0[r]*basisvalues[r];
131
134
}// end loop over 'r'
134
/// Evaluate all basis functions at given point in cell
137
/// Evaluate all basis functions at given point x in cell
135
138
virtual void evaluate_basis_all(double* values,
136
const double* coordinates,
137
const ufc::cell& c) const
140
const double* vertex_coordinates,
141
int cell_orientation) const
139
143
// Element is constant, calling evaluate_basis.
140
evaluate_basis(0, values, coordinates, c);
144
evaluate_basis(0, values, x, vertex_coordinates, cell_orientation);
143
/// Evaluate order n derivatives of basis function i at given point in cell
144
virtual void evaluate_basis_derivatives(unsigned int i,
147
/// Evaluate order n derivatives of basis function i at given point x in cell
148
virtual void evaluate_basis_derivatives(std::size_t i,
147
const double* coordinates,
148
const ufc::cell& c) const
152
const double* vertex_coordinates,
153
int cell_orientation) const
150
// Extract vertex coordinates
151
const double * const * x = c.coordinates;
153
// Compute Jacobian of affine map from reference cell
154
const double J_00 = x[1][0] - x[0][0];
155
const double J_01 = x[2][0] - x[0][0];
156
const double J_10 = x[1][1] - x[0][1];
157
const double J_11 = x[2][1] - x[0][1];
159
// Compute determinant of Jacobian
160
double detJ = J_00*J_11 - J_01*J_10;
162
// Compute inverse of Jacobian
163
const double K_00 = J_11 / detJ;
164
const double K_01 = -J_01 / detJ;
165
const double K_10 = -J_10 / detJ;
166
const double K_11 = J_00 / detJ;
157
compute_jacobian_triangle_2d(J, vertex_coordinates);
159
// Compute Jacobian inverse and determinant
162
compute_jacobian_inverse_triangle_2d(K, detJ, J);
168
165
// Compute constants
363
360
delete [] transform;
366
/// Evaluate order n derivatives of all basis functions at given point in cell
367
virtual void evaluate_basis_derivatives_all(unsigned int n,
363
/// Evaluate order n derivatives of all basis functions at given point x in cell
364
virtual void evaluate_basis_derivatives_all(std::size_t n,
369
const double* coordinates,
370
const ufc::cell& c) const
367
const double* vertex_coordinates,
368
int cell_orientation) const
372
370
// Element is constant, calling evaluate_basis_derivatives.
373
evaluate_basis_derivatives(0, n, values, coordinates, c);
371
evaluate_basis_derivatives(0, n, values, x, vertex_coordinates, cell_orientation);
376
374
/// Evaluate linear functional for dof i on the function f
377
virtual double evaluate_dof(unsigned int i,
375
virtual double evaluate_dof(std::size_t i,
378
376
const ufc::function& f,
377
const double* vertex_coordinates,
378
int cell_orientation,
379
379
const ufc::cell& c) const
381
// Declare variables for result of evaluation.
381
// Declare variables for result of evaluation
384
// Declare variable for physical coordinates.
384
// Declare variable for physical coordinates
386
const double * const * x = c.coordinates;
391
y[0] = 0.33333333*x[0][0] + 0.33333333*x[1][0] + 0.33333333*x[2][0];
392
y[1] = 0.33333333*x[0][1] + 0.33333333*x[1][1] + 0.33333333*x[2][1];
390
y[0] = 0.33333333*vertex_coordinates[0] + 0.33333333*vertex_coordinates[2] + 0.33333333*vertex_coordinates[4];
391
y[1] = 0.33333333*vertex_coordinates[1] + 0.33333333*vertex_coordinates[3] + 0.33333333*vertex_coordinates[5];
393
392
f.evaluate(vals, y, c);
1172
1176
/// Return the topological dimension of the cell shape
1173
virtual unsigned int topological_dimension() const
1177
virtual std::size_t topological_dimension() const
1178
1182
/// Return the geometric dimension of the cell shape
1179
virtual unsigned int geometric_dimension() const
1183
virtual std::size_t geometric_dimension() const
1184
1188
/// Return the dimension of the finite element function space
1185
virtual unsigned int space_dimension() const
1189
virtual std::size_t space_dimension() const
1190
1194
/// Return the rank of the value space
1191
virtual unsigned int value_rank() const
1195
virtual std::size_t value_rank() const
1196
1200
/// Return the dimension of the value space for axis i
1197
virtual unsigned int value_dimension(unsigned int i) const
1201
virtual std::size_t value_dimension(std::size_t i) const
1202
/// Evaluate basis function i at given point in cell
1203
virtual void evaluate_basis(unsigned int i,
1206
/// Evaluate basis function i at given point x in cell
1207
virtual void evaluate_basis(std::size_t i,
1204
1208
double* values,
1205
const double* coordinates,
1206
const ufc::cell& c) const
1210
const double* vertex_coordinates,
1211
int cell_orientation) const
1208
// Extract vertex coordinates
1209
const double * const * x = c.coordinates;
1211
// Compute Jacobian of affine map from reference cell
1212
const double J_00 = x[1][0] - x[0][0];
1213
const double J_01 = x[2][0] - x[0][0];
1214
const double J_10 = x[1][1] - x[0][1];
1215
const double J_11 = x[2][1] - x[0][1];
1217
// Compute determinant of Jacobian
1218
double detJ = J_00*J_11 - J_01*J_10;
1220
// Compute inverse of Jacobian
1215
compute_jacobian_triangle_2d(J, vertex_coordinates);
1217
// Compute Jacobian inverse and determinant
1220
compute_jacobian_inverse_triangle_2d(K, detJ, J);
1222
1223
// Compute constants
1223
const double C0 = x[1][0] + x[2][0];
1224
const double C1 = x[1][1] + x[2][1];
1224
const double C0 = vertex_coordinates[2] + vertex_coordinates[4];
1225
const double C1 = vertex_coordinates[3] + vertex_coordinates[5];
1226
1227
// Get coordinates and map to the reference (FIAT) element
1227
double X = (J_01*(C1 - 2.0*coordinates[1]) + J_11*(2.0*coordinates[0] - C0)) / detJ;
1228
double Y = (J_00*(2.0*coordinates[1] - C1) + J_10*(C0 - 2.0*coordinates[0])) / detJ;
1228
double X = (J[1]*(C1 - 2.0*x[1]) + J[3]*(2.0*x[0] - C0)) / detJ;
1229
double Y = (J[0]*(2.0*x[1] - C1) + J[2]*(C0 - 2.0*x[0])) / detJ;
1237
// Array of basisvalues.
1238
// Array of basisvalues
1238
1239
double basisvalues[6] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
1240
// Declare helper variables.
1241
// Declare helper variables
1241
1242
double tmp0 = (1.0 + Y + 2.0*X)/2.0;
1242
1243
double tmp1 = (1.0 - Y)/2.0;
1243
1244
double tmp2 = tmp1*tmp1;
1245
// Compute basisvalues.
1246
// Compute basisvalues
1246
1247
basisvalues[0] = 1.0;
1247
1248
basisvalues[1] = tmp0;
1248
1249
basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
1454
/// Evaluate all basis functions at given point in cell
1455
/// Evaluate all basis functions at given point x in cell
1455
1456
virtual void evaluate_basis_all(double* values,
1456
const double* coordinates,
1457
const ufc::cell& c) const
1458
const double* vertex_coordinates,
1459
int cell_orientation) const
1459
1461
// Helper variable to hold values of a single dof.
1460
1462
double dof_values = 0.0;
1462
// Loop dofs and call evaluate_basis.
1464
// Loop dofs and call evaluate_basis
1463
1465
for (unsigned int r = 0; r < 6; r++)
1465
evaluate_basis(r, &dof_values, coordinates, c);
1467
evaluate_basis(r, &dof_values, x, vertex_coordinates, cell_orientation);
1466
1468
values[r] = dof_values;
1467
1469
}// end loop over 'r'
1470
/// Evaluate order n derivatives of basis function i at given point in cell
1471
virtual void evaluate_basis_derivatives(unsigned int i,
1472
/// Evaluate order n derivatives of basis function i at given point x in cell
1473
virtual void evaluate_basis_derivatives(std::size_t i,
1473
1475
double* values,
1474
const double* coordinates,
1475
const ufc::cell& c) const
1477
const double* vertex_coordinates,
1478
int cell_orientation) const
1477
// Extract vertex coordinates
1478
const double * const * x = c.coordinates;
1480
// Compute Jacobian of affine map from reference cell
1481
const double J_00 = x[1][0] - x[0][0];
1482
const double J_01 = x[2][0] - x[0][0];
1483
const double J_10 = x[1][1] - x[0][1];
1484
const double J_11 = x[2][1] - x[0][1];
1486
// Compute determinant of Jacobian
1487
double detJ = J_00*J_11 - J_01*J_10;
1489
// Compute inverse of Jacobian
1490
const double K_00 = J_11 / detJ;
1491
const double K_01 = -J_01 / detJ;
1492
const double K_10 = -J_10 / detJ;
1493
const double K_11 = J_00 / detJ;
1482
compute_jacobian_triangle_2d(J, vertex_coordinates);
1484
// Compute Jacobian inverse and determinant
1487
compute_jacobian_inverse_triangle_2d(K, detJ, J);
1495
1490
// Compute constants
1496
const double C0 = x[1][0] + x[2][0];
1497
const double C1 = x[1][1] + x[2][1];
1491
const double C0 = vertex_coordinates[2] + vertex_coordinates[4];
1492
const double C1 = vertex_coordinates[3] + vertex_coordinates[5];
1499
1494
// Get coordinates and map to the reference (FIAT) element
1500
double X = (J_01*(C1 - 2.0*coordinates[1]) + J_11*(2.0*coordinates[0] - C0)) / detJ;
1501
double Y = (J_00*(2.0*coordinates[1] - C1) + J_10*(C0 - 2.0*coordinates[0])) / detJ;
1495
double X = (J[1]*(C1 - 2.0*x[1]) + J[3]*(2.0*x[0] - C0)) / detJ;
1496
double Y = (J[0]*(2.0*x[1] - C1) + J[2]*(C0 - 2.0*x[0])) / detJ;
1503
1498
// Compute number of derivatives.
1504
1499
unsigned int num_derivatives = 1;
2822
/// Evaluate basis function i at given point in cell
2823
virtual void evaluate_basis(unsigned int i,
2822
/// Evaluate basis function i at given point x in cell
2823
virtual void evaluate_basis(std::size_t i,
2824
2824
double* values,
2825
const double* coordinates,
2826
const ufc::cell& c) const
2826
const double* vertex_coordinates,
2827
int cell_orientation) const
2828
// Extract vertex coordinates
2829
const double * const * x = c.coordinates;
2831
// Compute Jacobian of affine map from reference cell
2832
const double J_00 = x[1][0] - x[0][0];
2833
const double J_01 = x[2][0] - x[0][0];
2834
const double J_10 = x[1][1] - x[0][1];
2835
const double J_11 = x[2][1] - x[0][1];
2837
// Compute determinant of Jacobian
2838
double detJ = J_00*J_11 - J_01*J_10;
2840
// Compute inverse of Jacobian
2831
compute_jacobian_triangle_2d(J, vertex_coordinates);
2833
// Compute Jacobian inverse and determinant
2836
compute_jacobian_inverse_triangle_2d(K, detJ, J);
2842
2839
// Compute constants
2843
const double C0 = x[1][0] + x[2][0];
2844
const double C1 = x[1][1] + x[2][1];
2840
const double C0 = vertex_coordinates[2] + vertex_coordinates[4];
2841
const double C1 = vertex_coordinates[3] + vertex_coordinates[5];
2846
2843
// Get coordinates and map to the reference (FIAT) element
2847
double X = (J_01*(C1 - 2.0*coordinates[1]) + J_11*(2.0*coordinates[0] - C0)) / detJ;
2848
double Y = (J_00*(2.0*coordinates[1] - C1) + J_10*(C0 - 2.0*coordinates[0])) / detJ;
2844
double X = (J[1]*(C1 - 2.0*x[1]) + J[3]*(2.0*x[0] - C0)) / detJ;
2845
double Y = (J[0]*(2.0*x[1] - C1) + J[2]*(C0 - 2.0*x[0])) / detJ;
2851
2848
values[0] = 0.0;
2852
2849
values[1] = 0.0;
2853
2850
values[2] = 0.0;
3136
3134
}// end loop over 'r'
3139
/// Evaluate order n derivatives of basis function i at given point in cell
3140
virtual void evaluate_basis_derivatives(unsigned int i,
3137
/// Evaluate order n derivatives of basis function i at given point x in cell
3138
virtual void evaluate_basis_derivatives(std::size_t i,
3142
3140
double* values,
3143
const double* coordinates,
3144
const ufc::cell& c) const
3142
const double* vertex_coordinates,
3143
int cell_orientation) const
3146
// Extract vertex coordinates
3147
const double * const * x = c.coordinates;
3149
// Compute Jacobian of affine map from reference cell
3150
const double J_00 = x[1][0] - x[0][0];
3151
const double J_01 = x[2][0] - x[0][0];
3152
const double J_10 = x[1][1] - x[0][1];
3153
const double J_11 = x[2][1] - x[0][1];
3155
// Compute determinant of Jacobian
3156
double detJ = J_00*J_11 - J_01*J_10;
3158
// Compute inverse of Jacobian
3159
const double K_00 = J_11 / detJ;
3160
const double K_01 = -J_01 / detJ;
3161
const double K_10 = -J_10 / detJ;
3162
const double K_11 = J_00 / detJ;
3147
compute_jacobian_triangle_2d(J, vertex_coordinates);
3149
// Compute Jacobian inverse and determinant
3152
compute_jacobian_inverse_triangle_2d(K, detJ, J);
3164
3155
// Compute constants
3165
const double C0 = x[1][0] + x[2][0];
3166
const double C1 = x[1][1] + x[2][1];
3156
const double C0 = vertex_coordinates[2] + vertex_coordinates[4];
3157
const double C1 = vertex_coordinates[3] + vertex_coordinates[5];
3168
3159
// Get coordinates and map to the reference (FIAT) element
3169
double X = (J_01*(C1 - 2.0*coordinates[1]) + J_11*(2.0*coordinates[0] - C0)) / detJ;
3170
double Y = (J_00*(2.0*coordinates[1] - C1) + J_10*(C0 - 2.0*coordinates[0])) / detJ;
3160
double X = (J[1]*(C1 - 2.0*x[1]) + J[3]*(2.0*x[0] - C0)) / detJ;
3161
double Y = (J[0]*(2.0*x[1] - C1) + J[2]*(C0 - 2.0*x[0])) / detJ;
3172
3163
// Compute number of derivatives.
3173
3164
unsigned int num_derivatives = 1;
4536
4528
/// Evaluate linear functional for dof i on the function f
4537
virtual double evaluate_dof(unsigned int i,
4529
virtual double evaluate_dof(std::size_t i,
4538
4530
const ufc::function& f,
4531
const double* vertex_coordinates,
4532
int cell_orientation,
4539
4533
const ufc::cell& c) const
4541
// Declare variables for result of evaluation.
4535
// Declare variables for result of evaluation
4542
4536
double vals[3];
4544
// Declare variable for physical coordinates.
4538
// Declare variable for physical coordinates
4546
const double * const * x = c.coordinates;
4551
y[0] = 0.33333333*x[0][0] + 0.33333333*x[1][0] + 0.33333333*x[2][0];
4552
y[1] = 0.33333333*x[0][1] + 0.33333333*x[1][1] + 0.33333333*x[2][1];
4544
y[0] = 0.33333333*vertex_coordinates[0] + 0.33333333*vertex_coordinates[2] + 0.33333333*vertex_coordinates[4];
4545
y[1] = 0.33333333*vertex_coordinates[1] + 0.33333333*vertex_coordinates[3] + 0.33333333*vertex_coordinates[5];
4553
4546
f.evaluate(vals, y, c);
4554
4547
return vals[0];
4559
y[0] = 0.33333333*x[0][0] + 0.33333333*x[1][0] + 0.33333333*x[2][0];
4560
y[1] = 0.33333333*x[0][1] + 0.33333333*x[1][1] + 0.33333333*x[2][1];
4552
y[0] = 0.33333333*vertex_coordinates[0] + 0.33333333*vertex_coordinates[2] + 0.33333333*vertex_coordinates[4];
4553
y[1] = 0.33333333*vertex_coordinates[1] + 0.33333333*vertex_coordinates[3] + 0.33333333*vertex_coordinates[5];
4561
4554
f.evaluate(vals, y, c);
4562
4555
return vals[1];
4560
y[0] = vertex_coordinates[0];
4561
y[1] = vertex_coordinates[1];
4569
4562
f.evaluate(vals, y, c);
4570
4563
return vals[2];
4568
y[0] = vertex_coordinates[2];
4569
y[1] = vertex_coordinates[3];
4577
4570
f.evaluate(vals, y, c);
4578
4571
return vals[2];
4576
y[0] = vertex_coordinates[4];
4577
y[1] = vertex_coordinates[5];
4585
4578
f.evaluate(vals, y, c);
4586
4579
return vals[2];
4591
y[0] = 0.5*x[1][0] + 0.5*x[2][0];
4592
y[1] = 0.5*x[1][1] + 0.5*x[2][1];
4584
y[0] = 0.5*vertex_coordinates[2] + 0.5*vertex_coordinates[4];
4585
y[1] = 0.5*vertex_coordinates[3] + 0.5*vertex_coordinates[5];
4593
4586
f.evaluate(vals, y, c);
4594
4587
return vals[2];
4599
y[0] = 0.5*x[0][0] + 0.5*x[2][0];
4600
y[1] = 0.5*x[0][1] + 0.5*x[2][1];
4592
y[0] = 0.5*vertex_coordinates[0] + 0.5*vertex_coordinates[4];
4593
y[1] = 0.5*vertex_coordinates[1] + 0.5*vertex_coordinates[5];
4601
4594
f.evaluate(vals, y, c);
4602
4595
return vals[2];
4607
y[0] = 0.5*x[0][0] + 0.5*x[1][0];
4608
y[1] = 0.5*x[0][1] + 0.5*x[1][1];
4600
y[0] = 0.5*vertex_coordinates[0] + 0.5*vertex_coordinates[2];
4601
y[1] = 0.5*vertex_coordinates[1] + 0.5*vertex_coordinates[3];
4609
4602
f.evaluate(vals, y, c);
4610
4603
return vals[2];
4618
4611
/// Evaluate linear functionals for all dofs on the function f
4619
4612
virtual void evaluate_dofs(double* values,
4620
4613
const ufc::function& f,
4614
const double* vertex_coordinates,
4615
int cell_orientation,
4621
4616
const ufc::cell& c) const
4623
// Declare variables for result of evaluation.
4618
// Declare variables for result of evaluation
4624
4619
double vals[3];
4626
// Declare variable for physical coordinates.
4621
// Declare variable for physical coordinates
4628
const double * const * x = c.coordinates;
4629
y[0] = 0.33333333*x[0][0] + 0.33333333*x[1][0] + 0.33333333*x[2][0];
4630
y[1] = 0.33333333*x[0][1] + 0.33333333*x[1][1] + 0.33333333*x[2][1];
4623
y[0] = 0.33333333*vertex_coordinates[0] + 0.33333333*vertex_coordinates[2] + 0.33333333*vertex_coordinates[4];
4624
y[1] = 0.33333333*vertex_coordinates[1] + 0.33333333*vertex_coordinates[3] + 0.33333333*vertex_coordinates[5];
4631
4625
f.evaluate(vals, y, c);
4632
4626
values[0] = vals[0];
4633
y[0] = 0.33333333*x[0][0] + 0.33333333*x[1][0] + 0.33333333*x[2][0];
4634
y[1] = 0.33333333*x[0][1] + 0.33333333*x[1][1] + 0.33333333*x[2][1];
4627
y[0] = 0.33333333*vertex_coordinates[0] + 0.33333333*vertex_coordinates[2] + 0.33333333*vertex_coordinates[4];
4628
y[1] = 0.33333333*vertex_coordinates[1] + 0.33333333*vertex_coordinates[3] + 0.33333333*vertex_coordinates[5];
4635
4629
f.evaluate(vals, y, c);
4636
4630
values[1] = vals[1];
4631
y[0] = vertex_coordinates[0];
4632
y[1] = vertex_coordinates[1];
4639
4633
f.evaluate(vals, y, c);
4640
4634
values[2] = vals[2];
4635
y[0] = vertex_coordinates[2];
4636
y[1] = vertex_coordinates[3];
4643
4637
f.evaluate(vals, y, c);
4644
4638
values[3] = vals[2];
4639
y[0] = vertex_coordinates[4];
4640
y[1] = vertex_coordinates[5];
4647
4641
f.evaluate(vals, y, c);
4648
4642
values[4] = vals[2];
4649
y[0] = 0.5*x[1][0] + 0.5*x[2][0];
4650
y[1] = 0.5*x[1][1] + 0.5*x[2][1];
4643
y[0] = 0.5*vertex_coordinates[2] + 0.5*vertex_coordinates[4];
4644
y[1] = 0.5*vertex_coordinates[3] + 0.5*vertex_coordinates[5];
4651
4645
f.evaluate(vals, y, c);
4652
4646
values[5] = vals[2];
4653
y[0] = 0.5*x[0][0] + 0.5*x[2][0];
4654
y[1] = 0.5*x[0][1] + 0.5*x[2][1];
4647
y[0] = 0.5*vertex_coordinates[0] + 0.5*vertex_coordinates[4];
4648
y[1] = 0.5*vertex_coordinates[1] + 0.5*vertex_coordinates[5];
4655
4649
f.evaluate(vals, y, c);
4656
4650
values[6] = vals[2];
4657
y[0] = 0.5*x[0][0] + 0.5*x[1][0];
4658
y[1] = 0.5*x[0][1] + 0.5*x[1][1];
4651
y[0] = 0.5*vertex_coordinates[0] + 0.5*vertex_coordinates[2];
4652
y[1] = 0.5*vertex_coordinates[1] + 0.5*vertex_coordinates[3];
4659
4653
f.evaluate(vals, y, c);
4660
4654
values[7] = vals[2];
4801
/// Evaluate basis function i at given point in cell
4802
virtual void evaluate_basis(unsigned int i,
4797
/// Evaluate basis function i at given point x in cell
4798
virtual void evaluate_basis(std::size_t i,
4803
4799
double* values,
4804
const double* coordinates,
4805
const ufc::cell& c) const
4801
const double* vertex_coordinates,
4802
int cell_orientation) const
4807
// Extract vertex coordinates
4808
const double * const * x = c.coordinates;
4810
// Compute Jacobian of affine map from reference cell
4811
const double J_00 = x[1][0] - x[0][0];
4812
const double J_01 = x[2][0] - x[0][0];
4813
const double J_10 = x[1][1] - x[0][1];
4814
const double J_11 = x[2][1] - x[0][1];
4816
// Compute determinant of Jacobian
4817
double detJ = J_00*J_11 - J_01*J_10;
4819
// Compute inverse of Jacobian
4806
compute_jacobian_triangle_2d(J, vertex_coordinates);
4808
// Compute Jacobian inverse and determinant
4811
compute_jacobian_inverse_triangle_2d(K, detJ, J);
4821
4815
// Compute constants
4822
const double C0 = x[1][0] + x[2][0];
4823
const double C1 = x[1][1] + x[2][1];
4816
const double C0 = vertex_coordinates[2] + vertex_coordinates[4];
4817
const double C1 = vertex_coordinates[3] + vertex_coordinates[5];
4825
4819
// Get coordinates and map to the reference (FIAT) element
4826
double X = (J_01*(C1 - 2.0*coordinates[1]) + J_11*(2.0*coordinates[0] - C0)) / detJ;
4827
double Y = (J_00*(2.0*coordinates[1] - C1) + J_10*(C0 - 2.0*coordinates[0])) / detJ;
4820
double X = (J[1]*(C1 - 2.0*x[1]) + J[3]*(2.0*x[0] - C0)) / detJ;
4821
double Y = (J[0]*(2.0*x[1] - C1) + J[2]*(C0 - 2.0*x[0])) / detJ;
4830
4824
values[0] = 0.0;
4831
4825
values[1] = 0.0;
4864
4858
basisvalues[7] *= std::sqrt(10.0);
4865
4859
basisvalues[6] *= std::sqrt(14.0);
4867
// Table(s) of coefficients.
4861
// Table(s) of coefficients
4868
4862
static const double coefficients0[10] = \
4869
4863
{0.18856181, 0.28867513, -0.16666667, 0.26082027, -0.20203051, 0.11664237, 0.1069045, -0.09035079, 0.069985421, -0.040406102};
4871
4865
static const double coefficients1[10] = \
4872
4866
{0.14142136, -0.038490018, 0.13333333, 0.069552071, -0.12121831, 0.089425816, 0.0, 0.06023386, -0.077761579, 0.053874802};
4874
// Compute value(s).
4875
4869
for (unsigned int r = 0; r < 10; r++)
4877
4871
values[0] += coefficients0[r]*basisvalues[r];
4878
4872
values[1] += coefficients1[r]*basisvalues[r];
4879
4873
}// end loop over 'r'
4881
// Using contravariant Piola transform to map values back to the physical element.
4875
// Using contravariant Piola transform to map values back to the physical element
4882
4876
const double tmp_ref0 = values[0];
4883
4877
const double tmp_ref1 = values[1];
4884
values[0] = (1.0/detJ)*(J_00*tmp_ref0 + J_01*tmp_ref1);
4885
values[1] = (1.0/detJ)*(J_10*tmp_ref0 + J_11*tmp_ref1);
4878
values[0] = (1.0/detJ)*(J[0]*tmp_ref0 + J[1]*tmp_ref1);
4879
values[1] = (1.0/detJ)*(J[2]*tmp_ref0 + J[3]*tmp_ref1);
4891
// Array of basisvalues.
4885
// Array of basisvalues
4892
4886
double basisvalues[10] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
4894
// Declare helper variables.
4888
// Declare helper variables
4895
4889
double tmp0 = (1.0 + Y + 2.0*X)/2.0;
4896
4890
double tmp1 = (1.0 - Y)/2.0;
4897
4891
double tmp2 = tmp1*tmp1;
4899
// Compute basisvalues.
4893
// Compute basisvalues
4900
4894
basisvalues[0] = 1.0;
4901
4895
basisvalues[1] = tmp0;
4902
4896
basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
4918
4912
basisvalues[7] *= std::sqrt(10.0);
4919
4913
basisvalues[6] *= std::sqrt(14.0);
4921
// Table(s) of coefficients.
4915
// Table(s) of coefficients
4922
4916
static const double coefficients0[10] = \
4923
4917
{-0.18856181, -0.25018512, 0.23333333, -0.10432811, 0.32324881, -0.25661321, 0.0, 0.12046772, -0.15552316, 0.1077496};
4925
4919
static const double coefficients1[10] = \
4926
4920
{-0.18856181, 0.076980036, -0.33333333, 0.0, 0.24243661, -0.34992711, 0.0, 0.0, 0.15552316, -0.16162441};
4928
// Compute value(s).
4929
4923
for (unsigned int r = 0; r < 10; r++)
4931
4925
values[0] += coefficients0[r]*basisvalues[r];
4932
4926
values[1] += coefficients1[r]*basisvalues[r];
4933
4927
}// end loop over 'r'
4935
// Using contravariant Piola transform to map values back to the physical element.
4929
// Using contravariant Piola transform to map values back to the physical element
4936
4930
const double tmp_ref0 = values[0];
4937
4931
const double tmp_ref1 = values[1];
4938
values[0] = (1.0/detJ)*(J_00*tmp_ref0 + J_01*tmp_ref1);
4939
values[1] = (1.0/detJ)*(J_10*tmp_ref0 + J_11*tmp_ref1);
4932
values[0] = (1.0/detJ)*(J[0]*tmp_ref0 + J[1]*tmp_ref1);
4933
values[1] = (1.0/detJ)*(J[2]*tmp_ref0 + J[3]*tmp_ref1);
4945
// Array of basisvalues.
4939
// Array of basisvalues
4946
4940
double basisvalues[10] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
4948
// Declare helper variables.
4942
// Declare helper variables
4949
4943
double tmp0 = (1.0 + Y + 2.0*X)/2.0;
4950
4944
double tmp1 = (1.0 - Y)/2.0;
4951
4945
double tmp2 = tmp1*tmp1;
4953
// Compute basisvalues.
4947
// Compute basisvalues
4954
4948
basisvalues[0] = 1.0;
4955
4949
basisvalues[1] = tmp0;
4956
4950
basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
4972
4966
basisvalues[7] *= std::sqrt(10.0);
4973
4967
basisvalues[6] *= std::sqrt(14.0);
4975
// Table(s) of coefficients.
4969
// Table(s) of coefficients
4976
4970
static const double coefficients0[10] = \
4977
4971
{0.14142136, 0.096225045, -0.1, 0.0, -0.067343503, 0.15163508, 0.0, 0.0, 0.077761579, -0.080812204};
4979
4973
static const double coefficients1[10] = \
4980
4974
{0.18856181, 0.0, 0.33333333, 0.0, 0.0, 0.34992711, 0.0, 0.0, 0.0, 0.16162441};
4982
// Compute value(s).
4983
4977
for (unsigned int r = 0; r < 10; r++)
4985
4979
values[0] += coefficients0[r]*basisvalues[r];
4986
4980
values[1] += coefficients1[r]*basisvalues[r];
4987
4981
}// end loop over 'r'
4989
// Using contravariant Piola transform to map values back to the physical element.
4983
// Using contravariant Piola transform to map values back to the physical element
4990
4984
const double tmp_ref0 = values[0];
4991
4985
const double tmp_ref1 = values[1];
4992
values[0] = (1.0/detJ)*(J_00*tmp_ref0 + J_01*tmp_ref1);
4993
values[1] = (1.0/detJ)*(J_10*tmp_ref0 + J_11*tmp_ref1);
4986
values[0] = (1.0/detJ)*(J[0]*tmp_ref0 + J[1]*tmp_ref1);
4987
values[1] = (1.0/detJ)*(J[2]*tmp_ref0 + J[3]*tmp_ref1);
4999
// Array of basisvalues.
4993
// Array of basisvalues
5000
4994
double basisvalues[10] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
5002
// Declare helper variables.
4996
// Declare helper variables
5003
4997
double tmp0 = (1.0 + Y + 2.0*X)/2.0;
5004
4998
double tmp1 = (1.0 - Y)/2.0;
5005
4999
double tmp2 = tmp1*tmp1;
5007
// Compute basisvalues.
5001
// Compute basisvalues
5008
5002
basisvalues[0] = 1.0;
5009
5003
basisvalues[1] = tmp0;
5010
5004
basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
5026
5020
basisvalues[7] *= std::sqrt(10.0);
5027
5021
basisvalues[6] *= std::sqrt(14.0);
5029
// Table(s) of coefficients.
5023
// Table(s) of coefficients
5030
5024
static const double coefficients0[10] = \
5031
5025
{0.32998316, -0.25018512, -0.033333333, 0.33037234, 0.32324881, 0.20606818, -0.1069045, -0.03011693, 0.0077761579, 0.013468701};
5033
5027
static const double coefficients1[10] = \
5034
5028
{-0.14142136, -0.038490018, -0.13333333, -0.069552071, -0.12121831, -0.089425816, 0.0, -0.06023386, -0.077761579, -0.053874802};
5036
// Compute value(s).
5037
5031
for (unsigned int r = 0; r < 10; r++)
5039
5033
values[0] += coefficients0[r]*basisvalues[r];
5040
5034
values[1] += coefficients1[r]*basisvalues[r];
5041
5035
}// end loop over 'r'
5043
// Using contravariant Piola transform to map values back to the physical element.
5037
// Using contravariant Piola transform to map values back to the physical element
5044
5038
const double tmp_ref0 = values[0];
5045
5039
const double tmp_ref1 = values[1];
5046
values[0] = (1.0/detJ)*(J_00*tmp_ref0 + J_01*tmp_ref1);
5047
values[1] = (1.0/detJ)*(J_10*tmp_ref0 + J_11*tmp_ref1);
5040
values[0] = (1.0/detJ)*(J[0]*tmp_ref0 + J[1]*tmp_ref1);
5041
values[1] = (1.0/detJ)*(J[2]*tmp_ref0 + J[3]*tmp_ref1);
5053
// Array of basisvalues.
5047
// Array of basisvalues
5054
5048
double basisvalues[10] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
5056
// Declare helper variables.
5050
// Declare helper variables
5057
5051
double tmp0 = (1.0 + Y + 2.0*X)/2.0;
5058
5052
double tmp1 = (1.0 - Y)/2.0;
5059
5053
double tmp2 = tmp1*tmp1;
5061
// Compute basisvalues.
5055
// Compute basisvalues
5062
5056
basisvalues[0] = 1.0;
5063
5057
basisvalues[1] = tmp0;
5064
5058
basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
5080
5074
basisvalues[7] *= std::sqrt(10.0);
5081
5075
basisvalues[6] *= std::sqrt(14.0);
5083
// Table(s) of coefficients.
5077
// Table(s) of coefficients
5084
5078
static const double coefficients0[10] = \
5085
5079
{-0.37712362, 0.17320508, -0.1, -0.10432811, -0.56568542, -0.60654032, 0.0, 0.12046772, 0.0, -0.053874802};
5087
5081
static const double coefficients1[10] = \
5088
5082
{0.18856181, 0.076980036, 0.33333333, 0.0, 0.24243661, 0.34992711, 0.0, 0.0, 0.15552316, 0.16162441};
5090
// Compute value(s).
5091
5085
for (unsigned int r = 0; r < 10; r++)
5093
5087
values[0] += coefficients0[r]*basisvalues[r];
5094
5088
values[1] += coefficients1[r]*basisvalues[r];
5095
5089
}// end loop over 'r'
5097
// Using contravariant Piola transform to map values back to the physical element.
5091
// Using contravariant Piola transform to map values back to the physical element
5098
5092
const double tmp_ref0 = values[0];
5099
5093
const double tmp_ref1 = values[1];
5100
values[0] = (1.0/detJ)*(J_00*tmp_ref0 + J_01*tmp_ref1);
5101
values[1] = (1.0/detJ)*(J_10*tmp_ref0 + J_11*tmp_ref1);
5094
values[0] = (1.0/detJ)*(J[0]*tmp_ref0 + J[1]*tmp_ref1);
5095
values[1] = (1.0/detJ)*(J[2]*tmp_ref0 + J[3]*tmp_ref1);
5107
// Array of basisvalues.
5101
// Array of basisvalues
5108
5102
double basisvalues[10] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
5110
// Declare helper variables.
5104
// Declare helper variables
5111
5105
double tmp0 = (1.0 + Y + 2.0*X)/2.0;
5112
5106
double tmp1 = (1.0 - Y)/2.0;
5113
5107
double tmp2 = tmp1*tmp1;
5115
// Compute basisvalues.
5109
// Compute basisvalues
5116
5110
basisvalues[0] = 1.0;
5117
5111
basisvalues[1] = tmp0;
5118
5112
basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
5134
5128
basisvalues[7] *= std::sqrt(10.0);
5135
5129
basisvalues[6] *= std::sqrt(14.0);
5137
// Table(s) of coefficients.
5131
// Table(s) of coefficients
5138
5132
static const double coefficients0[10] = \
5139
5133
{0.32998316, -0.096225045, 0.23333333, 0.0, 0.067343503, 0.50156219, 0.0, 0.0, -0.077761579, 0.080812204};
5141
5135
static const double coefficients1[10] = \
5142
5136
{-0.18856181, 0.0, -0.33333333, 0.0, 0.0, -0.34992711, 0.0, 0.0, 0.0, -0.16162441};
5144
// Compute value(s).
5145
5139
for (unsigned int r = 0; r < 10; r++)
5147
5141
values[0] += coefficients0[r]*basisvalues[r];
5148
5142
values[1] += coefficients1[r]*basisvalues[r];
5149
5143
}// end loop over 'r'
5151
// Using contravariant Piola transform to map values back to the physical element.
5145
// Using contravariant Piola transform to map values back to the physical element
5152
5146
const double tmp_ref0 = values[0];
5153
5147
const double tmp_ref1 = values[1];
5154
values[0] = (1.0/detJ)*(J_00*tmp_ref0 + J_01*tmp_ref1);
5155
values[1] = (1.0/detJ)*(J_10*tmp_ref0 + J_11*tmp_ref1);
5148
values[0] = (1.0/detJ)*(J[0]*tmp_ref0 + J[1]*tmp_ref1);
5149
values[1] = (1.0/detJ)*(J[2]*tmp_ref0 + J[3]*tmp_ref1);
5161
// Array of basisvalues.
5155
// Array of basisvalues
5162
5156
double basisvalues[10] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
5164
// Declare helper variables.
5158
// Declare helper variables
5165
5159
double tmp0 = (1.0 + Y + 2.0*X)/2.0;
5166
5160
double tmp1 = (1.0 - Y)/2.0;
5167
5161
double tmp2 = tmp1*tmp1;
5169
// Compute basisvalues.
5163
// Compute basisvalues
5170
5164
basisvalues[0] = 1.0;
5171
5165
basisvalues[1] = tmp0;
5172
5166
basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
5188
5182
basisvalues[7] *= std::sqrt(10.0);
5189
5183
basisvalues[6] *= std::sqrt(14.0);
5191
// Table(s) of coefficients.
5185
// Table(s) of coefficients
5192
5186
static const double coefficients0[10] = \
5193
5187
{0.14142136, 0.13471506, -0.033333333, 0.15649216, 0.053874802, 0.011664237, 0.1069045, 0.03011693, -0.0077761579, -0.013468701};
5195
5189
static const double coefficients1[10] = \
5196
5190
{-0.32998316, 0.15396007, 0.2, -0.41731242, -0.25590531, -0.12830661, 0.0, 0.06023386, 0.077761579, 0.053874802};
5198
// Compute value(s).
5199
5193
for (unsigned int r = 0; r < 10; r++)
5201
5195
values[0] += coefficients0[r]*basisvalues[r];
5202
5196
values[1] += coefficients1[r]*basisvalues[r];
5203
5197
}// end loop over 'r'
5205
// Using contravariant Piola transform to map values back to the physical element.
5199
// Using contravariant Piola transform to map values back to the physical element
5206
5200
const double tmp_ref0 = values[0];
5207
5201
const double tmp_ref1 = values[1];
5208
values[0] = (1.0/detJ)*(J_00*tmp_ref0 + J_01*tmp_ref1);
5209
values[1] = (1.0/detJ)*(J_10*tmp_ref0 + J_11*tmp_ref1);
5202
values[0] = (1.0/detJ)*(J[0]*tmp_ref0 + J[1]*tmp_ref1);
5203
values[1] = (1.0/detJ)*(J[2]*tmp_ref0 + J[3]*tmp_ref1);
5215
// Array of basisvalues.
5209
// Array of basisvalues
5216
5210
double basisvalues[10] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
5218
// Declare helper variables.
5212
// Declare helper variables
5219
5213
double tmp0 = (1.0 + Y + 2.0*X)/2.0;
5220
5214
double tmp1 = (1.0 - Y)/2.0;
5221
5215
double tmp2 = tmp1*tmp1;
5223
// Compute basisvalues.
5217
// Compute basisvalues
5224
5218
basisvalues[0] = 1.0;
5225
5219
basisvalues[1] = tmp0;
5226
5220
basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
5242
5236
basisvalues[7] *= std::sqrt(10.0);
5243
5237
basisvalues[6] *= std::sqrt(14.0);
5245
// Table(s) of coefficients.
5239
// Table(s) of coefficients
5246
5240
static const double coefficients0[10] = \
5247
5241
{-0.18856181, -0.32716515, 0.1, -0.41731242, 0.080812204, 0.023328474, -0.21380899, 0.06023386, 0.015552316, -0.026937401};
5249
5243
static const double coefficients1[10] = \
5250
5244
{0.37712362, 0.0, -0.2, 0.83462485, 0.0, -0.046656947, 0.0, -0.12046772, 0.0, 0.053874802};
5252
// Compute value(s).
5253
5247
for (unsigned int r = 0; r < 10; r++)
5255
5249
values[0] += coefficients0[r]*basisvalues[r];
5256
5250
values[1] += coefficients1[r]*basisvalues[r];
5257
5251
}// end loop over 'r'
5259
// Using contravariant Piola transform to map values back to the physical element.
5253
// Using contravariant Piola transform to map values back to the physical element
5260
5254
const double tmp_ref0 = values[0];
5261
5255
const double tmp_ref1 = values[1];
5262
values[0] = (1.0/detJ)*(J_00*tmp_ref0 + J_01*tmp_ref1);
5263
values[1] = (1.0/detJ)*(J_10*tmp_ref0 + J_11*tmp_ref1);
5256
values[0] = (1.0/detJ)*(J[0]*tmp_ref0 + J[1]*tmp_ref1);
5257
values[1] = (1.0/detJ)*(J[2]*tmp_ref0 + J[3]*tmp_ref1);
5269
// Array of basisvalues.
5263
// Array of basisvalues
5270
5264
double basisvalues[10] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
5272
// Declare helper variables.
5266
// Declare helper variables
5273
5267
double tmp0 = (1.0 + Y + 2.0*X)/2.0;
5274
5268
double tmp1 = (1.0 - Y)/2.0;
5275
5269
double tmp2 = tmp1*tmp1;
5277
// Compute basisvalues.
5271
// Compute basisvalues
5278
5272
basisvalues[0] = 1.0;
5279
5273
basisvalues[1] = tmp0;
5280
5274
basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
5296
5290
basisvalues[7] *= std::sqrt(10.0);
5297
5291
basisvalues[6] *= std::sqrt(14.0);
5299
// Table(s) of coefficients.
5293
// Table(s) of coefficients
5300
5294
static const double coefficients0[10] = \
5301
5295
{0.18856181, 0.28867513, -0.16666667, 0.26082027, -0.20203051, 0.11664237, 0.1069045, -0.09035079, 0.069985421, -0.040406102};
5303
5297
static const double coefficients1[10] = \
5304
5298
{-0.32998316, -0.15396007, 0.2, -0.41731242, 0.25590531, -0.12830661, 0.0, 0.06023386, -0.077761579, 0.053874802};
5306
// Compute value(s).
5307
5301
for (unsigned int r = 0; r < 10; r++)
5309
5303
values[0] += coefficients0[r]*basisvalues[r];
5310
5304
values[1] += coefficients1[r]*basisvalues[r];
5311
5305
}// end loop over 'r'
5313
// Using contravariant Piola transform to map values back to the physical element.
5307
// Using contravariant Piola transform to map values back to the physical element
5314
5308
const double tmp_ref0 = values[0];
5315
5309
const double tmp_ref1 = values[1];
5316
values[0] = (1.0/detJ)*(J_00*tmp_ref0 + J_01*tmp_ref1);
5317
values[1] = (1.0/detJ)*(J_10*tmp_ref0 + J_11*tmp_ref1);
5310
values[0] = (1.0/detJ)*(J[0]*tmp_ref0 + J[1]*tmp_ref1);
5311
values[1] = (1.0/detJ)*(J[2]*tmp_ref0 + J[3]*tmp_ref1);
5323
// Array of basisvalues.
5317
// Array of basisvalues
5324
5318
double basisvalues[10] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
5326
// Declare helper variables.
5320
// Declare helper variables
5327
5321
double tmp0 = (1.0 + Y + 2.0*X)/2.0;
5328
5322
double tmp1 = (1.0 - Y)/2.0;
5329
5323
double tmp2 = tmp1*tmp1;
5331
// Compute basisvalues.
5325
// Compute basisvalues
5332
5326
basisvalues[0] = 1.0;
5333
5327
basisvalues[1] = tmp0;
5334
5328
basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
5350
5344
basisvalues[7] *= std::sqrt(10.0);
5351
5345
basisvalues[6] *= std::sqrt(14.0);
5353
// Table(s) of coefficients.
5347
// Table(s) of coefficients
5354
5348
static const double coefficients0[10] = \
5355
5349
{0.32998316, -0.69282032, -0.53333333, -0.39992441, -0.026937401, 0.20606818, 0.32071349, -0.03011693, -0.023328474, 0.013468701};
5357
5351
static const double coefficients1[10] = \
5358
5352
{0.23570226, 0.038490018, 0.066666667, 0.20865621, 0.12121831, -0.081649658, 0.0, 0.18070158, 0.077761579, 0.0};
5360
// Compute value(s).
5361
5355
for (unsigned int r = 0; r < 10; r++)
5363
5357
values[0] += coefficients0[r]*basisvalues[r];
5364
5358
values[1] += coefficients1[r]*basisvalues[r];
5365
5359
}// end loop over 'r'
5367
// Using contravariant Piola transform to map values back to the physical element.
5361
// Using contravariant Piola transform to map values back to the physical element
5368
5362
const double tmp_ref0 = values[0];
5369
5363
const double tmp_ref1 = values[1];
5370
values[0] = (1.0/detJ)*(J_00*tmp_ref0 + J_01*tmp_ref1);
5371
values[1] = (1.0/detJ)*(J_10*tmp_ref0 + J_11*tmp_ref1);
5364
values[0] = (1.0/detJ)*(J[0]*tmp_ref0 + J[1]*tmp_ref1);
5365
values[1] = (1.0/detJ)*(J[2]*tmp_ref0 + J[3]*tmp_ref1);
5377
// Array of basisvalues.
5371
// Array of basisvalues
5378
5372
double basisvalues[10] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
5380
// Declare helper variables.
5374
// Declare helper variables
5381
5375
double tmp0 = (1.0 + Y + 2.0*X)/2.0;
5382
5376
double tmp1 = (1.0 - Y)/2.0;
5383
5377
double tmp2 = tmp1*tmp1;
5385
// Compute basisvalues.
5379
// Compute basisvalues
5386
5380
basisvalues[0] = 1.0;
5387
5381
basisvalues[1] = tmp0;
5388
5382
basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
5404
5398
basisvalues[7] *= std::sqrt(10.0);
5405
5399
basisvalues[6] *= std::sqrt(14.0);
5407
// Table(s) of coefficients.
5401
// Table(s) of coefficients
5408
5402
static const double coefficients0[10] = \
5409
5403
{0.56568542, 0.65433031, -0.46666667, -0.19126819, -0.094280904, 0.12441853, -0.32071349, 0.15058465, -0.054433105, 0.013468701};
5411
5405
static const double coefficients1[10] = \
5412
5406
{-0.23570226, 0.038490018, -0.066666667, -0.20865621, 0.12121831, 0.081649658, 0.0, -0.18070158, 0.077761579, 0.0};
5414
// Compute value(s).
5415
5409
for (unsigned int r = 0; r < 10; r++)
5417
5411
values[0] += coefficients0[r]*basisvalues[r];
5418
5412
values[1] += coefficients1[r]*basisvalues[r];
5419
5413
}// end loop over 'r'
5421
// Using contravariant Piola transform to map values back to the physical element.
5415
// Using contravariant Piola transform to map values back to the physical element
5422
5416
const double tmp_ref0 = values[0];
5423
5417
const double tmp_ref1 = values[1];
5424
values[0] = (1.0/detJ)*(J_00*tmp_ref0 + J_01*tmp_ref1);
5425
values[1] = (1.0/detJ)*(J_10*tmp_ref0 + J_11*tmp_ref1);
5418
values[0] = (1.0/detJ)*(J[0]*tmp_ref0 + J[1]*tmp_ref1);
5419
values[1] = (1.0/detJ)*(J[2]*tmp_ref0 + J[3]*tmp_ref1);
5431
// Array of basisvalues.
5425
// Array of basisvalues
5432
5426
double basisvalues[10] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
5434
// Declare helper variables.
5428
// Declare helper variables
5435
5429
double tmp0 = (1.0 + Y + 2.0*X)/2.0;
5436
5430
double tmp1 = (1.0 - Y)/2.0;
5437
5431
double tmp2 = tmp1*tmp1;
5439
// Compute basisvalues.
5433
// Compute basisvalues
5440
5434
basisvalues[0] = 1.0;
5441
5435
basisvalues[1] = tmp0;
5442
5436
basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
5458
5452
basisvalues[7] *= std::sqrt(10.0);
5459
5453
basisvalues[6] *= std::sqrt(14.0);
5461
// Table(s) of coefficients.
5455
// Table(s) of coefficients
5462
5456
static const double coefficients0[10] = \
5463
5457
{0.094280904, 0.076980036, 0.93333333, 0.20865621, 0.24243661, -0.443241, 0.0, -0.24093544, 0.15552316, -0.053874802};
5465
5459
static const double coefficients1[10] = \
5466
5460
{0.0, -0.15396007, 0.0, 0.0, -0.48487322, 0.0, 0.0, 0.0, -0.31104632, 0.0};
5468
// Compute value(s).
5469
5463
for (unsigned int r = 0; r < 10; r++)
5471
5465
values[0] += coefficients0[r]*basisvalues[r];
5472
5466
values[1] += coefficients1[r]*basisvalues[r];
5473
5467
}// end loop over 'r'
5475
// Using contravariant Piola transform to map values back to the physical element.
5469
// Using contravariant Piola transform to map values back to the physical element
5476
5470
const double tmp_ref0 = values[0];
5477
5471
const double tmp_ref1 = values[1];
5478
values[0] = (1.0/detJ)*(J_00*tmp_ref0 + J_01*tmp_ref1);
5479
values[1] = (1.0/detJ)*(J_10*tmp_ref0 + J_11*tmp_ref1);
5472
values[0] = (1.0/detJ)*(J[0]*tmp_ref0 + J[1]*tmp_ref1);
5473
values[1] = (1.0/detJ)*(J[2]*tmp_ref0 + J[3]*tmp_ref1);
5485
// Array of basisvalues.
5479
// Array of basisvalues
5486
5480
double basisvalues[10] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
5488
// Declare helper variables.
5482
// Declare helper variables
5489
5483
double tmp0 = (1.0 + Y + 2.0*X)/2.0;
5490
5484
double tmp1 = (1.0 - Y)/2.0;
5491
5485
double tmp2 = tmp1*tmp1;
5493
// Compute basisvalues.
5487
// Compute basisvalues
5494
5488
basisvalues[0] = 1.0;
5495
5489
basisvalues[1] = tmp0;
5496
5490
basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
5512
5506
basisvalues[7] *= std::sqrt(10.0);
5513
5507
basisvalues[6] *= std::sqrt(14.0);
5515
// Table(s) of coefficients.
5509
// Table(s) of coefficients
5516
5510
static const double coefficients0[10] = \
5517
5511
{0.23570226, 0.076980036, 0.0, 0.052164053, 0.24243661, 0.058321184, 0.1069045, 0.15058465, -0.0077761579, -0.067343503};
5519
5513
static const double coefficients1[10] = \
5520
5514
{0.32998316, -0.80829038, -0.33333333, 0.069552071, -0.39059232, -0.21384434, 0.0, 0.06023386, 0.23328474, 0.21549921};
5522
// Compute value(s).
5523
5517
for (unsigned int r = 0; r < 10; r++)
5525
5519
values[0] += coefficients0[r]*basisvalues[r];
5526
5520
values[1] += coefficients1[r]*basisvalues[r];
5527
5521
}// end loop over 'r'
5529
// Using contravariant Piola transform to map values back to the physical element.
5523
// Using contravariant Piola transform to map values back to the physical element
5530
5524
const double tmp_ref0 = values[0];
5531
5525
const double tmp_ref1 = values[1];
5532
values[0] = (1.0/detJ)*(J_00*tmp_ref0 + J_01*tmp_ref1);
5533
values[1] = (1.0/detJ)*(J_10*tmp_ref0 + J_11*tmp_ref1);
5526
values[0] = (1.0/detJ)*(J[0]*tmp_ref0 + J[1]*tmp_ref1);
5527
values[1] = (1.0/detJ)*(J[2]*tmp_ref0 + J[3]*tmp_ref1);
5539
// Array of basisvalues.
5533
// Array of basisvalues
5540
5534
double basisvalues[10] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
5542
// Declare helper variables.
5536
// Declare helper variables
5543
5537
double tmp0 = (1.0 + Y + 2.0*X)/2.0;
5544
5538
double tmp1 = (1.0 - Y)/2.0;
5545
5539
double tmp2 = tmp1*tmp1;
5547
// Compute basisvalues.
5541
// Compute basisvalues
5548
5542
basisvalues[0] = 1.0;
5549
5543
basisvalues[1] = tmp0;
5550
5544
basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
5566
5560
basisvalues[7] *= std::sqrt(10.0);
5567
5561
basisvalues[6] *= std::sqrt(14.0);
5569
// Table(s) of coefficients.
5563
// Table(s) of coefficients
5570
5564
static const double coefficients0[10] = \
5571
5565
{0.0, -0.076980036, -0.13333333, -0.31298432, -0.24243661, 0.27994168, -0.21380899, -0.06023386, 0.17107547, -0.13468701};
5573
5567
static const double coefficients1[10] = \
5574
5568
{0.094280904, 0.84678039, -0.4, -0.13910414, 0.51181062, -0.13219468, 0.0, -0.12046772, -0.15552316, 0.21549921};
5576
// Compute value(s).
5577
5571
for (unsigned int r = 0; r < 10; r++)
5579
5573
values[0] += coefficients0[r]*basisvalues[r];
5580
5574
values[1] += coefficients1[r]*basisvalues[r];
5581
5575
}// end loop over 'r'
5583
// Using contravariant Piola transform to map values back to the physical element.
5577
// Using contravariant Piola transform to map values back to the physical element
5584
5578
const double tmp_ref0 = values[0];
5585
5579
const double tmp_ref1 = values[1];
5586
values[0] = (1.0/detJ)*(J_00*tmp_ref0 + J_01*tmp_ref1);
5587
values[1] = (1.0/detJ)*(J_10*tmp_ref0 + J_11*tmp_ref1);
5580
values[0] = (1.0/detJ)*(J[0]*tmp_ref0 + J[1]*tmp_ref1);
5581
values[1] = (1.0/detJ)*(J[2]*tmp_ref0 + J[3]*tmp_ref1);
5593
// Array of basisvalues.
5587
// Array of basisvalues
5594
5588
double basisvalues[10] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
5596
// Declare helper variables.
5590
// Declare helper variables
5597
5591
double tmp0 = (1.0 + Y + 2.0*X)/2.0;
5598
5592
double tmp1 = (1.0 - Y)/2.0;
5599
5593
double tmp2 = tmp1*tmp1;
5601
// Compute basisvalues.
5595
// Compute basisvalues
5602
5596
basisvalues[0] = 1.0;
5603
5597
basisvalues[1] = tmp0;
5604
5598
basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
5620
5614
basisvalues[7] *= std::sqrt(10.0);
5621
5615
basisvalues[6] *= std::sqrt(14.0);
5623
// Table(s) of coefficients.
5617
// Table(s) of coefficients
5624
5618
static const double coefficients0[10] = \
5625
5619
{-0.23570226, -0.038490018, 0.066666667, 0.10432811, -0.12121831, -0.19829203, 0.0, -0.12046772, -0.077761579, 0.13468701};
5627
5621
static const double coefficients1[10] = \
5628
5622
{0.56568542, -0.076980036, 0.8, 0.0, -0.24243661, -0.046656947, 0.0, 0.0, -0.15552316, -0.32324881};
5630
// Compute value(s).
5631
5625
for (unsigned int r = 0; r < 10; r++)
5633
5627
values[0] += coefficients0[r]*basisvalues[r];
5634
5628
values[1] += coefficients1[r]*basisvalues[r];
5635
5629
}// end loop over 'r'
5637
// Using contravariant Piola transform to map values back to the physical element.
5631
// Using contravariant Piola transform to map values back to the physical element
5638
5632
const double tmp_ref0 = values[0];
5639
5633
const double tmp_ref1 = values[1];
5640
values[0] = (1.0/detJ)*(J_00*tmp_ref0 + J_01*tmp_ref1);
5641
values[1] = (1.0/detJ)*(J_10*tmp_ref0 + J_11*tmp_ref1);
5634
values[0] = (1.0/detJ)*(J[0]*tmp_ref0 + J[1]*tmp_ref1);
5635
values[1] = (1.0/detJ)*(J[2]*tmp_ref0 + J[3]*tmp_ref1);
5648
/// Evaluate all basis functions at given point in cell
5642
/// Evaluate all basis functions at given point x in cell
5649
5643
virtual void evaluate_basis_all(double* values,
5650
const double* coordinates,
5651
const ufc::cell& c) const
5645
const double* vertex_coordinates,
5646
int cell_orientation) const
5653
5648
// Helper variable to hold values of a single dof.
5654
5649
double dof_values[2] = {0.0, 0.0};
5656
// Loop dofs and call evaluate_basis.
5651
// Loop dofs and call evaluate_basis
5657
5652
for (unsigned int r = 0; r < 15; r++)
5659
evaluate_basis(r, dof_values, coordinates, c);
5654
evaluate_basis(r, dof_values, x, vertex_coordinates, cell_orientation);
5660
5655
for (unsigned int s = 0; s < 2; s++)
5662
5657
values[r*2 + s] = dof_values[s];
5664
5659
}// end loop over 'r'
5667
/// Evaluate order n derivatives of basis function i at given point in cell
5668
virtual void evaluate_basis_derivatives(unsigned int i,
5662
/// Evaluate order n derivatives of basis function i at given point x in cell
5663
virtual void evaluate_basis_derivatives(std::size_t i,
5670
5665
double* values,
5671
const double* coordinates,
5672
const ufc::cell& c) const
5667
const double* vertex_coordinates,
5668
int cell_orientation) const
5674
// Extract vertex coordinates
5675
const double * const * x = c.coordinates;
5677
// Compute Jacobian of affine map from reference cell
5678
const double J_00 = x[1][0] - x[0][0];
5679
const double J_01 = x[2][0] - x[0][0];
5680
const double J_10 = x[1][1] - x[0][1];
5681
const double J_11 = x[2][1] - x[0][1];
5683
// Compute determinant of Jacobian
5684
double detJ = J_00*J_11 - J_01*J_10;
5686
// Compute inverse of Jacobian
5687
const double K_00 = J_11 / detJ;
5688
const double K_01 = -J_01 / detJ;
5689
const double K_10 = -J_10 / detJ;
5690
const double K_11 = J_00 / detJ;
5672
compute_jacobian_triangle_2d(J, vertex_coordinates);
5674
// Compute Jacobian inverse and determinant
5677
compute_jacobian_inverse_triangle_2d(K, detJ, J);
5692
5681
// Compute constants
5693
const double C0 = x[1][0] + x[2][0];
5694
const double C1 = x[1][1] + x[2][1];
5682
const double C0 = vertex_coordinates[2] + vertex_coordinates[4];
5683
const double C1 = vertex_coordinates[3] + vertex_coordinates[5];
5696
5685
// Get coordinates and map to the reference (FIAT) element
5697
double X = (J_01*(C1 - 2.0*coordinates[1]) + J_11*(2.0*coordinates[0] - C0)) / detJ;
5698
double Y = (J_00*(2.0*coordinates[1] - C1) + J_10*(C0 - 2.0*coordinates[0])) / detJ;
5686
double X = (J[1]*(C1 - 2.0*x[1]) + J[3]*(2.0*x[0] - C0)) / detJ;
5687
double Y = (J[0]*(2.0*x[1] - C1) + J[2]*(C0 - 2.0*x[0])) / detJ;
5700
5689
// Compute number of derivatives.
5701
5690
unsigned int num_derivatives = 1;
8819
8959
/// Evaluate linear functional for dof i on the function f
8820
virtual double evaluate_dof(unsigned int i,
8960
virtual double evaluate_dof(std::size_t i,
8821
8961
const ufc::function& f,
8962
const double* vertex_coordinates,
8963
int cell_orientation,
8822
8964
const ufc::cell& c) const
8824
// Declare variables for result of evaluation.
8966
// Declare variables for result of evaluation
8825
8967
double vals[2];
8827
// Declare variable for physical coordinates.
8969
// Declare variable for physical coordinates
8830
// Extract vertex coordinates
8831
const double * const * x = c.coordinates;
8833
// Compute Jacobian of affine map from reference cell
8834
const double J_00 = x[1][0] - x[0][0];
8835
const double J_01 = x[2][0] - x[0][0];
8836
const double J_10 = x[1][1] - x[0][1];
8837
const double J_11 = x[2][1] - x[0][1];
8839
// Compute determinant of Jacobian
8840
double detJ = J_00*J_11 - J_01*J_10;
8842
// Compute inverse of Jacobian
8843
const double K_00 = J_11 / detJ;
8844
const double K_01 = -J_01 / detJ;
8845
const double K_10 = -J_10 / detJ;
8846
const double K_11 = J_00 / detJ;switch (i)
8975
compute_jacobian_triangle_2d(J, vertex_coordinates);
8978
// Compute Jacobian inverse and determinant
8981
compute_jacobian_inverse_triangle_2d(K, detJ, J);
8850
y[0] = 0.75*x[1][0] + 0.25*x[2][0];
8851
y[1] = 0.75*x[1][1] + 0.25*x[2][1];
8988
y[0] = 0.75*vertex_coordinates[2] + 0.25*vertex_coordinates[4];
8989
y[1] = 0.75*vertex_coordinates[3] + 0.25*vertex_coordinates[5];
8852
8990
f.evaluate(vals, y, c);
8853
result = (detJ*(K_00*vals[0] + K_01*vals[1])) + (detJ*(K_10*vals[0] + K_11*vals[1]));
8991
result = (detJ*(K[0]*vals[0] + K[1]*vals[1])) + (detJ*(K[2]*vals[0] + K[3]*vals[1]));
8859
y[0] = 0.5*x[1][0] + 0.5*x[2][0];
8860
y[1] = 0.5*x[1][1] + 0.5*x[2][1];
8997
y[0] = 0.5*vertex_coordinates[2] + 0.5*vertex_coordinates[4];
8998
y[1] = 0.5*vertex_coordinates[3] + 0.5*vertex_coordinates[5];
8861
8999
f.evaluate(vals, y, c);
8862
result = (detJ*(K_00*vals[0] + K_01*vals[1])) + (detJ*(K_10*vals[0] + K_11*vals[1]));
9000
result = (detJ*(K[0]*vals[0] + K[1]*vals[1])) + (detJ*(K[2]*vals[0] + K[3]*vals[1]));
8868
y[0] = 0.25*x[1][0] + 0.75*x[2][0];
8869
y[1] = 0.25*x[1][1] + 0.75*x[2][1];
9006
y[0] = 0.25*vertex_coordinates[2] + 0.75*vertex_coordinates[4];
9007
y[1] = 0.25*vertex_coordinates[3] + 0.75*vertex_coordinates[5];
8870
9008
f.evaluate(vals, y, c);
8871
result = (detJ*(K_00*vals[0] + K_01*vals[1])) + (detJ*(K_10*vals[0] + K_11*vals[1]));
9009
result = (detJ*(K[0]*vals[0] + K[1]*vals[1])) + (detJ*(K[2]*vals[0] + K[3]*vals[1]));
8877
y[0] = 0.75*x[0][0] + 0.25*x[2][0];
8878
y[1] = 0.75*x[0][1] + 0.25*x[2][1];
9015
y[0] = 0.75*vertex_coordinates[0] + 0.25*vertex_coordinates[4];
9016
y[1] = 0.75*vertex_coordinates[1] + 0.25*vertex_coordinates[5];
8879
9017
f.evaluate(vals, y, c);
8880
result = (detJ*(K_00*vals[0] + K_01*vals[1]));
9018
result = (detJ*(K[0]*vals[0] + K[1]*vals[1]));
8886
y[0] = 0.5*x[0][0] + 0.5*x[2][0];
8887
y[1] = 0.5*x[0][1] + 0.5*x[2][1];
9024
y[0] = 0.5*vertex_coordinates[0] + 0.5*vertex_coordinates[4];
9025
y[1] = 0.5*vertex_coordinates[1] + 0.5*vertex_coordinates[5];
8888
9026
f.evaluate(vals, y, c);
8889
result = (detJ*(K_00*vals[0] + K_01*vals[1]));
9027
result = (detJ*(K[0]*vals[0] + K[1]*vals[1]));
8895
y[0] = 0.25*x[0][0] + 0.75*x[2][0];
8896
y[1] = 0.25*x[0][1] + 0.75*x[2][1];
9033
y[0] = 0.25*vertex_coordinates[0] + 0.75*vertex_coordinates[4];
9034
y[1] = 0.25*vertex_coordinates[1] + 0.75*vertex_coordinates[5];
8897
9035
f.evaluate(vals, y, c);
8898
result = (detJ*(K_00*vals[0] + K_01*vals[1]));
9036
result = (detJ*(K[0]*vals[0] + K[1]*vals[1]));
8904
y[0] = 0.75*x[0][0] + 0.25*x[1][0];
8905
y[1] = 0.75*x[0][1] + 0.25*x[1][1];
9042
y[0] = 0.75*vertex_coordinates[0] + 0.25*vertex_coordinates[2];
9043
y[1] = 0.75*vertex_coordinates[1] + 0.25*vertex_coordinates[3];
8906
9044
f.evaluate(vals, y, c);
8907
result = (-1.0)*(detJ*(K_10*vals[0] + K_11*vals[1]));
9045
result = (-1.0)*(detJ*(K[2]*vals[0] + K[3]*vals[1]));
8913
y[0] = 0.5*x[0][0] + 0.5*x[1][0];
8914
y[1] = 0.5*x[0][1] + 0.5*x[1][1];
9051
y[0] = 0.5*vertex_coordinates[0] + 0.5*vertex_coordinates[2];
9052
y[1] = 0.5*vertex_coordinates[1] + 0.5*vertex_coordinates[3];
8915
9053
f.evaluate(vals, y, c);
8916
result = (-1.0)*(detJ*(K_10*vals[0] + K_11*vals[1]));
9054
result = (-1.0)*(detJ*(K[2]*vals[0] + K[3]*vals[1]));
8922
y[0] = 0.25*x[0][0] + 0.75*x[1][0];
8923
y[1] = 0.25*x[0][1] + 0.75*x[1][1];
9060
y[0] = 0.25*vertex_coordinates[0] + 0.75*vertex_coordinates[2];
9061
y[1] = 0.25*vertex_coordinates[1] + 0.75*vertex_coordinates[3];
8924
9062
f.evaluate(vals, y, c);
8925
result = (-1.0)*(detJ*(K_10*vals[0] + K_11*vals[1]));
9063
result = (-1.0)*(detJ*(K[2]*vals[0] + K[3]*vals[1]));
8931
y[0] = 0.5*x[0][0] + 0.25*x[1][0] + 0.25*x[2][0];
8932
y[1] = 0.5*x[0][1] + 0.25*x[1][1] + 0.25*x[2][1];
9069
y[0] = 0.5*vertex_coordinates[0] + 0.25*vertex_coordinates[2] + 0.25*vertex_coordinates[4];
9070
y[1] = 0.5*vertex_coordinates[1] + 0.25*vertex_coordinates[3] + 0.25*vertex_coordinates[5];
8933
9071
f.evaluate(vals, y, c);
8934
result = (detJ*(K_00*vals[0] + K_01*vals[1]));
9072
result = (detJ*(K[0]*vals[0] + K[1]*vals[1]));
8940
y[0] = 0.25*x[0][0] + 0.5*x[1][0] + 0.25*x[2][0];
8941
y[1] = 0.25*x[0][1] + 0.5*x[1][1] + 0.25*x[2][1];
9078
y[0] = 0.25*vertex_coordinates[0] + 0.5*vertex_coordinates[2] + 0.25*vertex_coordinates[4];
9079
y[1] = 0.25*vertex_coordinates[1] + 0.5*vertex_coordinates[3] + 0.25*vertex_coordinates[5];
8942
9080
f.evaluate(vals, y, c);
8943
result = (detJ*(K_00*vals[0] + K_01*vals[1]));
9081
result = (detJ*(K[0]*vals[0] + K[1]*vals[1]));
8949
y[0] = 0.25*x[0][0] + 0.25*x[1][0] + 0.5*x[2][0];
8950
y[1] = 0.25*x[0][1] + 0.25*x[1][1] + 0.5*x[2][1];
9087
y[0] = 0.25*vertex_coordinates[0] + 0.25*vertex_coordinates[2] + 0.5*vertex_coordinates[4];
9088
y[1] = 0.25*vertex_coordinates[1] + 0.25*vertex_coordinates[3] + 0.5*vertex_coordinates[5];
8951
9089
f.evaluate(vals, y, c);
8952
result = (detJ*(K_00*vals[0] + K_01*vals[1]));
9090
result = (detJ*(K[0]*vals[0] + K[1]*vals[1]));
8958
y[0] = 0.5*x[0][0] + 0.25*x[1][0] + 0.25*x[2][0];
8959
y[1] = 0.5*x[0][1] + 0.25*x[1][1] + 0.25*x[2][1];
9096
y[0] = 0.5*vertex_coordinates[0] + 0.25*vertex_coordinates[2] + 0.25*vertex_coordinates[4];
9097
y[1] = 0.5*vertex_coordinates[1] + 0.25*vertex_coordinates[3] + 0.25*vertex_coordinates[5];
8960
9098
f.evaluate(vals, y, c);
8961
result = (detJ*(K_10*vals[0] + K_11*vals[1]));
9099
result = (detJ*(K[2]*vals[0] + K[3]*vals[1]));
8967
y[0] = 0.25*x[0][0] + 0.5*x[1][0] + 0.25*x[2][0];
8968
y[1] = 0.25*x[0][1] + 0.5*x[1][1] + 0.25*x[2][1];
9105
y[0] = 0.25*vertex_coordinates[0] + 0.5*vertex_coordinates[2] + 0.25*vertex_coordinates[4];
9106
y[1] = 0.25*vertex_coordinates[1] + 0.5*vertex_coordinates[3] + 0.25*vertex_coordinates[5];
8969
9107
f.evaluate(vals, y, c);
8970
result = (detJ*(K_10*vals[0] + K_11*vals[1]));
9108
result = (detJ*(K[2]*vals[0] + K[3]*vals[1]));
8976
y[0] = 0.25*x[0][0] + 0.25*x[1][0] + 0.5*x[2][0];
8977
y[1] = 0.25*x[0][1] + 0.25*x[1][1] + 0.5*x[2][1];
9114
y[0] = 0.25*vertex_coordinates[0] + 0.25*vertex_coordinates[2] + 0.5*vertex_coordinates[4];
9115
y[1] = 0.25*vertex_coordinates[1] + 0.25*vertex_coordinates[3] + 0.5*vertex_coordinates[5];
8978
9116
f.evaluate(vals, y, c);
8979
result = (detJ*(K_10*vals[0] + K_11*vals[1]));
9117
result = (detJ*(K[2]*vals[0] + K[3]*vals[1]));
8988
9126
/// Evaluate linear functionals for all dofs on the function f
8989
9127
virtual void evaluate_dofs(double* values,
8990
9128
const ufc::function& f,
9129
const double* vertex_coordinates,
9130
int cell_orientation,
8991
9131
const ufc::cell& c) const
8993
// Declare variables for result of evaluation.
9133
// Declare variables for result of evaluation
8994
9134
double vals[2];
8996
// Declare variable for physical coordinates.
9136
// Declare variable for physical coordinates
8999
// Extract vertex coordinates
9000
const double * const * x = c.coordinates;
9002
// Compute Jacobian of affine map from reference cell
9003
const double J_00 = x[1][0] - x[0][0];
9004
const double J_01 = x[2][0] - x[0][0];
9005
const double J_10 = x[1][1] - x[0][1];
9006
const double J_11 = x[2][1] - x[0][1];
9008
// Compute determinant of Jacobian
9009
double detJ = J_00*J_11 - J_01*J_10;
9011
// Compute inverse of Jacobian
9012
const double K_00 = J_11 / detJ;
9013
const double K_01 = -J_01 / detJ;
9014
const double K_10 = -J_10 / detJ;
9015
const double K_11 = J_00 / detJ;y[0] = 0.75*x[1][0] + 0.25*x[2][0];
9016
y[1] = 0.75*x[1][1] + 0.25*x[2][1];
9142
compute_jacobian_triangle_2d(J, vertex_coordinates);
9145
// Compute Jacobian inverse and determinant
9148
compute_jacobian_inverse_triangle_2d(K, detJ, J);
9151
y[0] = 0.75*vertex_coordinates[2] + 0.25*vertex_coordinates[4];
9152
y[1] = 0.75*vertex_coordinates[3] + 0.25*vertex_coordinates[5];
9017
9153
f.evaluate(vals, y, c);
9018
result = (detJ*(K_00*vals[0] + K_01*vals[1])) + (detJ*(K_10*vals[0] + K_11*vals[1]));
9154
result = (detJ*(K[0]*vals[0] + K[1]*vals[1])) + (detJ*(K[2]*vals[0] + K[3]*vals[1]));
9019
9155
values[0] = result;
9020
y[0] = 0.5*x[1][0] + 0.5*x[2][0];
9021
y[1] = 0.5*x[1][1] + 0.5*x[2][1];
9156
y[0] = 0.5*vertex_coordinates[2] + 0.5*vertex_coordinates[4];
9157
y[1] = 0.5*vertex_coordinates[3] + 0.5*vertex_coordinates[5];
9022
9158
f.evaluate(vals, y, c);
9023
result = (detJ*(K_00*vals[0] + K_01*vals[1])) + (detJ*(K_10*vals[0] + K_11*vals[1]));
9159
result = (detJ*(K[0]*vals[0] + K[1]*vals[1])) + (detJ*(K[2]*vals[0] + K[3]*vals[1]));
9024
9160
values[1] = result;
9025
y[0] = 0.25*x[1][0] + 0.75*x[2][0];
9026
y[1] = 0.25*x[1][1] + 0.75*x[2][1];
9161
y[0] = 0.25*vertex_coordinates[2] + 0.75*vertex_coordinates[4];
9162
y[1] = 0.25*vertex_coordinates[3] + 0.75*vertex_coordinates[5];
9027
9163
f.evaluate(vals, y, c);
9028
result = (detJ*(K_00*vals[0] + K_01*vals[1])) + (detJ*(K_10*vals[0] + K_11*vals[1]));
9164
result = (detJ*(K[0]*vals[0] + K[1]*vals[1])) + (detJ*(K[2]*vals[0] + K[3]*vals[1]));
9029
9165
values[2] = result;
9030
y[0] = 0.75*x[0][0] + 0.25*x[2][0];
9031
y[1] = 0.75*x[0][1] + 0.25*x[2][1];
9166
y[0] = 0.75*vertex_coordinates[0] + 0.25*vertex_coordinates[4];
9167
y[1] = 0.75*vertex_coordinates[1] + 0.25*vertex_coordinates[5];
9032
9168
f.evaluate(vals, y, c);
9033
result = (detJ*(K_00*vals[0] + K_01*vals[1]));
9169
result = (detJ*(K[0]*vals[0] + K[1]*vals[1]));
9034
9170
values[3] = result;
9035
y[0] = 0.5*x[0][0] + 0.5*x[2][0];
9036
y[1] = 0.5*x[0][1] + 0.5*x[2][1];
9171
y[0] = 0.5*vertex_coordinates[0] + 0.5*vertex_coordinates[4];
9172
y[1] = 0.5*vertex_coordinates[1] + 0.5*vertex_coordinates[5];
9037
9173
f.evaluate(vals, y, c);
9038
result = (detJ*(K_00*vals[0] + K_01*vals[1]));
9174
result = (detJ*(K[0]*vals[0] + K[1]*vals[1]));
9039
9175
values[4] = result;
9040
y[0] = 0.25*x[0][0] + 0.75*x[2][0];
9041
y[1] = 0.25*x[0][1] + 0.75*x[2][1];
9176
y[0] = 0.25*vertex_coordinates[0] + 0.75*vertex_coordinates[4];
9177
y[1] = 0.25*vertex_coordinates[1] + 0.75*vertex_coordinates[5];
9042
9178
f.evaluate(vals, y, c);
9043
result = (detJ*(K_00*vals[0] + K_01*vals[1]));
9179
result = (detJ*(K[0]*vals[0] + K[1]*vals[1]));
9044
9180
values[5] = result;
9045
y[0] = 0.75*x[0][0] + 0.25*x[1][0];
9046
y[1] = 0.75*x[0][1] + 0.25*x[1][1];
9181
y[0] = 0.75*vertex_coordinates[0] + 0.25*vertex_coordinates[2];
9182
y[1] = 0.75*vertex_coordinates[1] + 0.25*vertex_coordinates[3];
9047
9183
f.evaluate(vals, y, c);
9048
result = (-1.0)*(detJ*(K_10*vals[0] + K_11*vals[1]));
9184
result = (-1.0)*(detJ*(K[2]*vals[0] + K[3]*vals[1]));
9049
9185
values[6] = result;
9050
y[0] = 0.5*x[0][0] + 0.5*x[1][0];
9051
y[1] = 0.5*x[0][1] + 0.5*x[1][1];
9186
y[0] = 0.5*vertex_coordinates[0] + 0.5*vertex_coordinates[2];
9187
y[1] = 0.5*vertex_coordinates[1] + 0.5*vertex_coordinates[3];
9052
9188
f.evaluate(vals, y, c);
9053
result = (-1.0)*(detJ*(K_10*vals[0] + K_11*vals[1]));
9189
result = (-1.0)*(detJ*(K[2]*vals[0] + K[3]*vals[1]));
9054
9190
values[7] = result;
9055
y[0] = 0.25*x[0][0] + 0.75*x[1][0];
9056
y[1] = 0.25*x[0][1] + 0.75*x[1][1];
9191
y[0] = 0.25*vertex_coordinates[0] + 0.75*vertex_coordinates[2];
9192
y[1] = 0.25*vertex_coordinates[1] + 0.75*vertex_coordinates[3];
9057
9193
f.evaluate(vals, y, c);
9058
result = (-1.0)*(detJ*(K_10*vals[0] + K_11*vals[1]));
9194
result = (-1.0)*(detJ*(K[2]*vals[0] + K[3]*vals[1]));
9059
9195
values[8] = result;
9060
y[0] = 0.5*x[0][0] + 0.25*x[1][0] + 0.25*x[2][0];
9061
y[1] = 0.5*x[0][1] + 0.25*x[1][1] + 0.25*x[2][1];
9196
y[0] = 0.5*vertex_coordinates[0] + 0.25*vertex_coordinates[2] + 0.25*vertex_coordinates[4];
9197
y[1] = 0.5*vertex_coordinates[1] + 0.25*vertex_coordinates[3] + 0.25*vertex_coordinates[5];
9062
9198
f.evaluate(vals, y, c);
9063
result = (detJ*(K_00*vals[0] + K_01*vals[1]));
9199
result = (detJ*(K[0]*vals[0] + K[1]*vals[1]));
9064
9200
values[9] = result;
9065
y[0] = 0.25*x[0][0] + 0.5*x[1][0] + 0.25*x[2][0];
9066
y[1] = 0.25*x[0][1] + 0.5*x[1][1] + 0.25*x[2][1];
9201
y[0] = 0.25*vertex_coordinates[0] + 0.5*vertex_coordinates[2] + 0.25*vertex_coordinates[4];
9202
y[1] = 0.25*vertex_coordinates[1] + 0.5*vertex_coordinates[3] + 0.25*vertex_coordinates[5];
9067
9203
f.evaluate(vals, y, c);
9068
result = (detJ*(K_00*vals[0] + K_01*vals[1]));
9204
result = (detJ*(K[0]*vals[0] + K[1]*vals[1]));
9069
9205
values[10] = result;
9070
y[0] = 0.25*x[0][0] + 0.25*x[1][0] + 0.5*x[2][0];
9071
y[1] = 0.25*x[0][1] + 0.25*x[1][1] + 0.5*x[2][1];
9206
y[0] = 0.25*vertex_coordinates[0] + 0.25*vertex_coordinates[2] + 0.5*vertex_coordinates[4];
9207
y[1] = 0.25*vertex_coordinates[1] + 0.25*vertex_coordinates[3] + 0.5*vertex_coordinates[5];
9072
9208
f.evaluate(vals, y, c);
9073
result = (detJ*(K_00*vals[0] + K_01*vals[1]));
9209
result = (detJ*(K[0]*vals[0] + K[1]*vals[1]));
9074
9210
values[11] = result;
9075
y[0] = 0.5*x[0][0] + 0.25*x[1][0] + 0.25*x[2][0];
9076
y[1] = 0.5*x[0][1] + 0.25*x[1][1] + 0.25*x[2][1];
9211
y[0] = 0.5*vertex_coordinates[0] + 0.25*vertex_coordinates[2] + 0.25*vertex_coordinates[4];
9212
y[1] = 0.5*vertex_coordinates[1] + 0.25*vertex_coordinates[3] + 0.25*vertex_coordinates[5];
9077
9213
f.evaluate(vals, y, c);
9078
result = (detJ*(K_10*vals[0] + K_11*vals[1]));
9214
result = (detJ*(K[2]*vals[0] + K[3]*vals[1]));
9079
9215
values[12] = result;
9080
y[0] = 0.25*x[0][0] + 0.5*x[1][0] + 0.25*x[2][0];
9081
y[1] = 0.25*x[0][1] + 0.5*x[1][1] + 0.25*x[2][1];
9216
y[0] = 0.25*vertex_coordinates[0] + 0.5*vertex_coordinates[2] + 0.25*vertex_coordinates[4];
9217
y[1] = 0.25*vertex_coordinates[1] + 0.5*vertex_coordinates[3] + 0.25*vertex_coordinates[5];
9082
9218
f.evaluate(vals, y, c);
9083
result = (detJ*(K_10*vals[0] + K_11*vals[1]));
9219
result = (detJ*(K[2]*vals[0] + K[3]*vals[1]));
9084
9220
values[13] = result;
9085
y[0] = 0.25*x[0][0] + 0.25*x[1][0] + 0.5*x[2][0];
9086
y[1] = 0.25*x[0][1] + 0.25*x[1][1] + 0.5*x[2][1];
9221
y[0] = 0.25*vertex_coordinates[0] + 0.25*vertex_coordinates[2] + 0.5*vertex_coordinates[4];
9222
y[1] = 0.25*vertex_coordinates[1] + 0.25*vertex_coordinates[3] + 0.5*vertex_coordinates[5];
9087
9223
f.evaluate(vals, y, c);
9088
result = (detJ*(K_10*vals[0] + K_11*vals[1]));
9224
result = (detJ*(K[2]*vals[0] + K[3]*vals[1]));
9089
9225
values[14] = result;
9092
9228
/// Interpolate vertex values from dof values
9093
9229
virtual void interpolate_vertex_values(double* vertex_values,
9094
9230
const double* dof_values,
9231
const double* vertex_coordinates,
9232
int cell_orientation,
9095
9233
const ufc::cell& c) const
9097
// Extract vertex coordinates
9098
const double * const * x = c.coordinates;
9100
// Compute Jacobian of affine map from reference cell
9101
const double J_00 = x[1][0] - x[0][0];
9102
const double J_01 = x[2][0] - x[0][0];
9103
const double J_10 = x[1][1] - x[0][1];
9104
const double J_11 = x[2][1] - x[0][1];
9106
// Compute determinant of Jacobian
9107
double detJ = J_00*J_11 - J_01*J_10;
9109
// Compute inverse of Jacobian
9237
compute_jacobian_triangle_2d(J, vertex_coordinates);
9240
// Compute Jacobian inverse and determinant
9243
compute_jacobian_inverse_triangle_2d(K, detJ, J);
9110
9247
// Evaluate function and change variables
9111
vertex_values[0] = dof_values[3]*(1.0/detJ)*J_00*3.0 + dof_values[4]*((1.0/detJ)*(J_00*(-3.0))) + dof_values[5]*(1.0/detJ)*J_00 + dof_values[6]*((1.0/detJ)*(J_01*(-3.0))) + dof_values[7]*(1.0/detJ)*J_01*3.0 + dof_values[8]*((1.0/detJ)*(J_01*(-1.0)));
9112
vertex_values[2] = dof_values[0]*(1.0/detJ)*J_00*3.0 + dof_values[1]*((1.0/detJ)*(J_00*(-3.0))) + dof_values[2]*(1.0/detJ)*J_00 + dof_values[6]*((1.0/detJ)*(J_00 + J_01*(-1.0))) + dof_values[7]*((1.0/detJ)*(J_00*(-3.0) + J_01*3.0)) + dof_values[8]*((1.0/detJ)*(J_00*3.0 + J_01*(-3.0)));
9113
vertex_values[4] = dof_values[0]*(1.0/detJ)*J_01 + dof_values[1]*((1.0/detJ)*(J_01*(-3.0))) + dof_values[2]*(1.0/detJ)*J_01*3.0 + dof_values[3]*((1.0/detJ)*(J_00 + J_01*(-1.0))) + dof_values[4]*((1.0/detJ)*(J_00*(-3.0) + J_01*3.0)) + dof_values[5]*((1.0/detJ)*(J_00*3.0 + J_01*(-3.0)));
9114
vertex_values[1] = dof_values[3]*(1.0/detJ)*J_10*3.0 + dof_values[4]*((1.0/detJ)*(J_10*(-3.0))) + dof_values[5]*(1.0/detJ)*J_10 + dof_values[6]*((1.0/detJ)*(J_11*(-3.0))) + dof_values[7]*(1.0/detJ)*J_11*3.0 + dof_values[8]*((1.0/detJ)*(J_11*(-1.0)));
9115
vertex_values[3] = dof_values[0]*(1.0/detJ)*J_10*3.0 + dof_values[1]*((1.0/detJ)*(J_10*(-3.0))) + dof_values[2]*(1.0/detJ)*J_10 + dof_values[6]*((1.0/detJ)*(J_10 + J_11*(-1.0))) + dof_values[7]*((1.0/detJ)*(J_10*(-3.0) + J_11*3.0)) + dof_values[8]*((1.0/detJ)*(J_10*3.0 + J_11*(-3.0)));
9116
vertex_values[5] = dof_values[0]*(1.0/detJ)*J_11 + dof_values[1]*((1.0/detJ)*(J_11*(-3.0))) + dof_values[2]*(1.0/detJ)*J_11*3.0 + dof_values[3]*((1.0/detJ)*(J_10 + J_11*(-1.0))) + dof_values[4]*((1.0/detJ)*(J_10*(-3.0) + J_11*3.0)) + dof_values[5]*((1.0/detJ)*(J_10*3.0 + J_11*(-3.0)));
9248
vertex_values[0] = dof_values[3]*(1.0/detJ)*J[0]*3.0 + dof_values[4]*((1.0/detJ)*(J[0]*(-3.0))) + dof_values[5]*(1.0/detJ)*J[0] + dof_values[6]*((1.0/detJ)*(J[1]*(-3.0))) + dof_values[7]*(1.0/detJ)*J[1]*3.0 + dof_values[8]*((1.0/detJ)*(J[1]*(-1.0)));
9249
vertex_values[2] = dof_values[0]*(1.0/detJ)*J[0]*3.0 + dof_values[1]*((1.0/detJ)*(J[0]*(-3.0))) + dof_values[2]*(1.0/detJ)*J[0] + dof_values[6]*((1.0/detJ)*(J[0] + J[1]*(-1.0))) + dof_values[7]*((1.0/detJ)*(J[0]*(-3.0) + J[1]*3.0)) + dof_values[8]*((1.0/detJ)*(J[0]*3.0 + J[1]*(-3.0)));
9250
vertex_values[4] = dof_values[0]*(1.0/detJ)*J[1] + dof_values[1]*((1.0/detJ)*(J[1]*(-3.0))) + dof_values[2]*(1.0/detJ)*J[1]*3.0 + dof_values[3]*((1.0/detJ)*(J[0] + J[1]*(-1.0))) + dof_values[4]*((1.0/detJ)*(J[0]*(-3.0) + J[1]*3.0)) + dof_values[5]*((1.0/detJ)*(J[0]*3.0 + J[1]*(-3.0)));
9251
vertex_values[1] = dof_values[3]*(1.0/detJ)*J[2]*3.0 + dof_values[4]*((1.0/detJ)*(J[2]*(-3.0))) + dof_values[5]*(1.0/detJ)*J[2] + dof_values[6]*((1.0/detJ)*(J[3]*(-3.0))) + dof_values[7]*(1.0/detJ)*J[3]*3.0 + dof_values[8]*((1.0/detJ)*(J[3]*(-1.0)));
9252
vertex_values[3] = dof_values[0]*(1.0/detJ)*J[2]*3.0 + dof_values[1]*((1.0/detJ)*(J[2]*(-3.0))) + dof_values[2]*(1.0/detJ)*J[2] + dof_values[6]*((1.0/detJ)*(J[2] + J[3]*(-1.0))) + dof_values[7]*((1.0/detJ)*(J[2]*(-3.0) + J[3]*3.0)) + dof_values[8]*((1.0/detJ)*(J[2]*3.0 + J[3]*(-3.0)));
9253
vertex_values[5] = dof_values[0]*(1.0/detJ)*J[3] + dof_values[1]*((1.0/detJ)*(J[3]*(-3.0))) + dof_values[2]*(1.0/detJ)*J[3]*3.0 + dof_values[3]*((1.0/detJ)*(J[2] + J[3]*(-1.0))) + dof_values[4]*((1.0/detJ)*(J[2]*(-3.0) + J[3]*3.0)) + dof_values[5]*((1.0/detJ)*(J[2]*3.0 + J[3]*(-3.0)));
9119
9256
/// Map coordinate xhat from reference cell to coordinate x in cell
9224
/// Evaluate basis function i at given point in cell
9225
virtual void evaluate_basis(unsigned int i,
9361
/// Evaluate basis function i at given point x in cell
9362
virtual void evaluate_basis(std::size_t i,
9226
9363
double* values,
9227
const double* coordinates,
9228
const ufc::cell& c) const
9365
const double* vertex_coordinates,
9366
int cell_orientation) const
9230
// Extract vertex coordinates
9231
const double * const * x = c.coordinates;
9233
// Compute Jacobian of affine map from reference cell
9234
const double J_00 = x[1][0] - x[0][0];
9235
const double J_01 = x[2][0] - x[0][0];
9236
const double J_10 = x[1][1] - x[0][1];
9237
const double J_11 = x[2][1] - x[0][1];
9239
// Compute determinant of Jacobian
9240
double detJ = J_00*J_11 - J_01*J_10;
9242
// Compute inverse of Jacobian
9370
compute_jacobian_triangle_2d(J, vertex_coordinates);
9372
// Compute Jacobian inverse and determinant
9375
compute_jacobian_inverse_triangle_2d(K, detJ, J);
9244
9379
// Compute constants
9245
const double C0 = x[1][0] + x[2][0];
9246
const double C1 = x[1][1] + x[2][1];
9380
const double C0 = vertex_coordinates[2] + vertex_coordinates[4];
9381
const double C1 = vertex_coordinates[3] + vertex_coordinates[5];
9248
9383
// Get coordinates and map to the reference (FIAT) element
9249
double X = (J_01*(C1 - 2.0*coordinates[1]) + J_11*(2.0*coordinates[0] - C0)) / detJ;
9250
double Y = (J_00*(2.0*coordinates[1] - C1) + J_10*(C0 - 2.0*coordinates[0])) / detJ;
9384
double X = (J[1]*(C1 - 2.0*x[1]) + J[3]*(2.0*x[0] - C0)) / detJ;
9385
double Y = (J[0]*(2.0*x[1] - C1) + J[2]*(C0 - 2.0*x[0])) / detJ;
9253
9388
values[0] = 0.0;
9254
9389
values[1] = 0.0;
9255
9390
values[2] = 0.0;
9550
9685
basisvalues[7] *= std::sqrt(10.0);
9551
9686
basisvalues[6] *= std::sqrt(14.0);
9553
// Table(s) of coefficients.
9688
// Table(s) of coefficients
9554
9689
static const double coefficients0[10] = \
9555
9690
{0.18856181, 0.28867513, -0.16666667, 0.26082027, -0.20203051, 0.11664237, 0.1069045, -0.09035079, 0.069985421, -0.040406102};
9557
9692
static const double coefficients1[10] = \
9558
9693
{0.14142136, -0.038490018, 0.13333333, 0.069552071, -0.12121831, 0.089425816, 0.0, 0.06023386, -0.077761579, 0.053874802};
9560
// Compute value(s).
9561
9696
for (unsigned int r = 0; r < 10; r++)
9563
9698
values[3] += coefficients0[r]*basisvalues[r];
9564
9699
values[4] += coefficients1[r]*basisvalues[r];
9565
9700
}// end loop over 'r'
9567
// Using contravariant Piola transform to map values back to the physical element.
9702
// Using contravariant Piola transform to map values back to the physical element
9568
9703
const double tmp_ref0 = values[3];
9569
9704
const double tmp_ref1 = values[4];
9570
values[3] = (1.0/detJ)*(J_00*tmp_ref0 + J_01*tmp_ref1);
9571
values[4] = (1.0/detJ)*(J_10*tmp_ref0 + J_11*tmp_ref1);
9705
values[3] = (1.0/detJ)*(J[0]*tmp_ref0 + J[1]*tmp_ref1);
9706
values[4] = (1.0/detJ)*(J[2]*tmp_ref0 + J[3]*tmp_ref1);
9577
// Array of basisvalues.
9712
// Array of basisvalues
9578
9713
double basisvalues[10] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
9580
// Declare helper variables.
9715
// Declare helper variables
9581
9716
double tmp0 = (1.0 + Y + 2.0*X)/2.0;
9582
9717
double tmp1 = (1.0 - Y)/2.0;
9583
9718
double tmp2 = tmp1*tmp1;
9585
// Compute basisvalues.
9720
// Compute basisvalues
9586
9721
basisvalues[0] = 1.0;
9587
9722
basisvalues[1] = tmp0;
9588
9723
basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
9604
9739
basisvalues[7] *= std::sqrt(10.0);
9605
9740
basisvalues[6] *= std::sqrt(14.0);
9607
// Table(s) of coefficients.
9742
// Table(s) of coefficients
9608
9743
static const double coefficients0[10] = \
9609
9744
{-0.18856181, -0.25018512, 0.23333333, -0.10432811, 0.32324881, -0.25661321, 0.0, 0.12046772, -0.15552316, 0.1077496};
9611
9746
static const double coefficients1[10] = \
9612
9747
{-0.18856181, 0.076980036, -0.33333333, 0.0, 0.24243661, -0.34992711, 0.0, 0.0, 0.15552316, -0.16162441};
9614
// Compute value(s).
9615
9750
for (unsigned int r = 0; r < 10; r++)
9617
9752
values[3] += coefficients0[r]*basisvalues[r];
9618
9753
values[4] += coefficients1[r]*basisvalues[r];
9619
9754
}// end loop over 'r'
9621
// Using contravariant Piola transform to map values back to the physical element.
9756
// Using contravariant Piola transform to map values back to the physical element
9622
9757
const double tmp_ref0 = values[3];
9623
9758
const double tmp_ref1 = values[4];
9624
values[3] = (1.0/detJ)*(J_00*tmp_ref0 + J_01*tmp_ref1);
9625
values[4] = (1.0/detJ)*(J_10*tmp_ref0 + J_11*tmp_ref1);
9759
values[3] = (1.0/detJ)*(J[0]*tmp_ref0 + J[1]*tmp_ref1);
9760
values[4] = (1.0/detJ)*(J[2]*tmp_ref0 + J[3]*tmp_ref1);
9631
// Array of basisvalues.
9766
// Array of basisvalues
9632
9767
double basisvalues[10] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
9634
// Declare helper variables.
9769
// Declare helper variables
9635
9770
double tmp0 = (1.0 + Y + 2.0*X)/2.0;
9636
9771
double tmp1 = (1.0 - Y)/2.0;
9637
9772
double tmp2 = tmp1*tmp1;
9639
// Compute basisvalues.
9774
// Compute basisvalues
9640
9775
basisvalues[0] = 1.0;
9641
9776
basisvalues[1] = tmp0;
9642
9777
basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
9658
9793
basisvalues[7] *= std::sqrt(10.0);
9659
9794
basisvalues[6] *= std::sqrt(14.0);
9661
// Table(s) of coefficients.
9796
// Table(s) of coefficients
9662
9797
static const double coefficients0[10] = \
9663
9798
{0.14142136, 0.096225045, -0.1, 0.0, -0.067343503, 0.15163508, 0.0, 0.0, 0.077761579, -0.080812204};
9665
9800
static const double coefficients1[10] = \
9666
9801
{0.18856181, 0.0, 0.33333333, 0.0, 0.0, 0.34992711, 0.0, 0.0, 0.0, 0.16162441};
9668
// Compute value(s).
9669
9804
for (unsigned int r = 0; r < 10; r++)
9671
9806
values[3] += coefficients0[r]*basisvalues[r];
9672
9807
values[4] += coefficients1[r]*basisvalues[r];
9673
9808
}// end loop over 'r'
9675
// Using contravariant Piola transform to map values back to the physical element.
9810
// Using contravariant Piola transform to map values back to the physical element
9676
9811
const double tmp_ref0 = values[3];
9677
9812
const double tmp_ref1 = values[4];
9678
values[3] = (1.0/detJ)*(J_00*tmp_ref0 + J_01*tmp_ref1);
9679
values[4] = (1.0/detJ)*(J_10*tmp_ref0 + J_11*tmp_ref1);
9813
values[3] = (1.0/detJ)*(J[0]*tmp_ref0 + J[1]*tmp_ref1);
9814
values[4] = (1.0/detJ)*(J[2]*tmp_ref0 + J[3]*tmp_ref1);
9685
// Array of basisvalues.
9820
// Array of basisvalues
9686
9821
double basisvalues[10] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
9688
// Declare helper variables.
9823
// Declare helper variables
9689
9824
double tmp0 = (1.0 + Y + 2.0*X)/2.0;
9690
9825
double tmp1 = (1.0 - Y)/2.0;
9691
9826
double tmp2 = tmp1*tmp1;
9693
// Compute basisvalues.
9828
// Compute basisvalues
9694
9829
basisvalues[0] = 1.0;
9695
9830
basisvalues[1] = tmp0;
9696
9831
basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
9712
9847
basisvalues[7] *= std::sqrt(10.0);
9713
9848
basisvalues[6] *= std::sqrt(14.0);
9715
// Table(s) of coefficients.
9850
// Table(s) of coefficients
9716
9851
static const double coefficients0[10] = \
9717
9852
{0.32998316, -0.25018512, -0.033333333, 0.33037234, 0.32324881, 0.20606818, -0.1069045, -0.03011693, 0.0077761579, 0.013468701};
9719
9854
static const double coefficients1[10] = \
9720
9855
{-0.14142136, -0.038490018, -0.13333333, -0.069552071, -0.12121831, -0.089425816, 0.0, -0.06023386, -0.077761579, -0.053874802};
9722
// Compute value(s).
9723
9858
for (unsigned int r = 0; r < 10; r++)
9725
9860
values[3] += coefficients0[r]*basisvalues[r];
9726
9861
values[4] += coefficients1[r]*basisvalues[r];
9727
9862
}// end loop over 'r'
9729
// Using contravariant Piola transform to map values back to the physical element.
9864
// Using contravariant Piola transform to map values back to the physical element
9730
9865
const double tmp_ref0 = values[3];
9731
9866
const double tmp_ref1 = values[4];
9732
values[3] = (1.0/detJ)*(J_00*tmp_ref0 + J_01*tmp_ref1);
9733
values[4] = (1.0/detJ)*(J_10*tmp_ref0 + J_11*tmp_ref1);
9867
values[3] = (1.0/detJ)*(J[0]*tmp_ref0 + J[1]*tmp_ref1);
9868
values[4] = (1.0/detJ)*(J[2]*tmp_ref0 + J[3]*tmp_ref1);
9739
// Array of basisvalues.
9874
// Array of basisvalues
9740
9875
double basisvalues[10] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
9742
// Declare helper variables.
9877
// Declare helper variables
9743
9878
double tmp0 = (1.0 + Y + 2.0*X)/2.0;
9744
9879
double tmp1 = (1.0 - Y)/2.0;
9745
9880
double tmp2 = tmp1*tmp1;
9747
// Compute basisvalues.
9882
// Compute basisvalues
9748
9883
basisvalues[0] = 1.0;
9749
9884
basisvalues[1] = tmp0;
9750
9885
basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
9766
9901
basisvalues[7] *= std::sqrt(10.0);
9767
9902
basisvalues[6] *= std::sqrt(14.0);
9769
// Table(s) of coefficients.
9904
// Table(s) of coefficients
9770
9905
static const double coefficients0[10] = \
9771
9906
{-0.37712362, 0.17320508, -0.1, -0.10432811, -0.56568542, -0.60654032, 0.0, 0.12046772, 0.0, -0.053874802};
9773
9908
static const double coefficients1[10] = \
9774
9909
{0.18856181, 0.076980036, 0.33333333, 0.0, 0.24243661, 0.34992711, 0.0, 0.0, 0.15552316, 0.16162441};
9776
// Compute value(s).
9777
9912
for (unsigned int r = 0; r < 10; r++)
9779
9914
values[3] += coefficients0[r]*basisvalues[r];
9780
9915
values[4] += coefficients1[r]*basisvalues[r];
9781
9916
}// end loop over 'r'
9783
// Using contravariant Piola transform to map values back to the physical element.
9918
// Using contravariant Piola transform to map values back to the physical element
9784
9919
const double tmp_ref0 = values[3];
9785
9920
const double tmp_ref1 = values[4];
9786
values[3] = (1.0/detJ)*(J_00*tmp_ref0 + J_01*tmp_ref1);
9787
values[4] = (1.0/detJ)*(J_10*tmp_ref0 + J_11*tmp_ref1);
9921
values[3] = (1.0/detJ)*(J[0]*tmp_ref0 + J[1]*tmp_ref1);
9922
values[4] = (1.0/detJ)*(J[2]*tmp_ref0 + J[3]*tmp_ref1);
9793
// Array of basisvalues.
9928
// Array of basisvalues
9794
9929
double basisvalues[10] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
9796
// Declare helper variables.
9931
// Declare helper variables
9797
9932
double tmp0 = (1.0 + Y + 2.0*X)/2.0;
9798
9933
double tmp1 = (1.0 - Y)/2.0;
9799
9934
double tmp2 = tmp1*tmp1;
9801
// Compute basisvalues.
9936
// Compute basisvalues
9802
9937
basisvalues[0] = 1.0;
9803
9938
basisvalues[1] = tmp0;
9804
9939
basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
9820
9955
basisvalues[7] *= std::sqrt(10.0);
9821
9956
basisvalues[6] *= std::sqrt(14.0);
9823
// Table(s) of coefficients.
9958
// Table(s) of coefficients
9824
9959
static const double coefficients0[10] = \
9825
9960
{0.32998316, -0.096225045, 0.23333333, 0.0, 0.067343503, 0.50156219, 0.0, 0.0, -0.077761579, 0.080812204};
9827
9962
static const double coefficients1[10] = \
9828
9963
{-0.18856181, 0.0, -0.33333333, 0.0, 0.0, -0.34992711, 0.0, 0.0, 0.0, -0.16162441};
9830
// Compute value(s).
9831
9966
for (unsigned int r = 0; r < 10; r++)
9833
9968
values[3] += coefficients0[r]*basisvalues[r];
9834
9969
values[4] += coefficients1[r]*basisvalues[r];
9835
9970
}// end loop over 'r'
9837
// Using contravariant Piola transform to map values back to the physical element.
9972
// Using contravariant Piola transform to map values back to the physical element
9838
9973
const double tmp_ref0 = values[3];
9839
9974
const double tmp_ref1 = values[4];
9840
values[3] = (1.0/detJ)*(J_00*tmp_ref0 + J_01*tmp_ref1);
9841
values[4] = (1.0/detJ)*(J_10*tmp_ref0 + J_11*tmp_ref1);
9975
values[3] = (1.0/detJ)*(J[0]*tmp_ref0 + J[1]*tmp_ref1);
9976
values[4] = (1.0/detJ)*(J[2]*tmp_ref0 + J[3]*tmp_ref1);
9847
// Array of basisvalues.
9982
// Array of basisvalues
9848
9983
double basisvalues[10] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
9850
// Declare helper variables.
9985
// Declare helper variables
9851
9986
double tmp0 = (1.0 + Y + 2.0*X)/2.0;
9852
9987
double tmp1 = (1.0 - Y)/2.0;
9853
9988
double tmp2 = tmp1*tmp1;
9855
// Compute basisvalues.
9990
// Compute basisvalues
9856
9991
basisvalues[0] = 1.0;
9857
9992
basisvalues[1] = tmp0;
9858
9993
basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
9874
10009
basisvalues[7] *= std::sqrt(10.0);
9875
10010
basisvalues[6] *= std::sqrt(14.0);
9877
// Table(s) of coefficients.
10012
// Table(s) of coefficients
9878
10013
static const double coefficients0[10] = \
9879
10014
{0.14142136, 0.13471506, -0.033333333, 0.15649216, 0.053874802, 0.011664237, 0.1069045, 0.03011693, -0.0077761579, -0.013468701};
9881
10016
static const double coefficients1[10] = \
9882
10017
{-0.32998316, 0.15396007, 0.2, -0.41731242, -0.25590531, -0.12830661, 0.0, 0.06023386, 0.077761579, 0.053874802};
9884
// Compute value(s).
10019
// Compute value(s)
9885
10020
for (unsigned int r = 0; r < 10; r++)
9887
10022
values[3] += coefficients0[r]*basisvalues[r];
9888
10023
values[4] += coefficients1[r]*basisvalues[r];
9889
10024
}// end loop over 'r'
9891
// Using contravariant Piola transform to map values back to the physical element.
10026
// Using contravariant Piola transform to map values back to the physical element
9892
10027
const double tmp_ref0 = values[3];
9893
10028
const double tmp_ref1 = values[4];
9894
values[3] = (1.0/detJ)*(J_00*tmp_ref0 + J_01*tmp_ref1);
9895
values[4] = (1.0/detJ)*(J_10*tmp_ref0 + J_11*tmp_ref1);
10029
values[3] = (1.0/detJ)*(J[0]*tmp_ref0 + J[1]*tmp_ref1);
10030
values[4] = (1.0/detJ)*(J[2]*tmp_ref0 + J[3]*tmp_ref1);
9901
// Array of basisvalues.
10036
// Array of basisvalues
9902
10037
double basisvalues[10] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
9904
// Declare helper variables.
10039
// Declare helper variables
9905
10040
double tmp0 = (1.0 + Y + 2.0*X)/2.0;
9906
10041
double tmp1 = (1.0 - Y)/2.0;
9907
10042
double tmp2 = tmp1*tmp1;
9909
// Compute basisvalues.
10044
// Compute basisvalues
9910
10045
basisvalues[0] = 1.0;
9911
10046
basisvalues[1] = tmp0;
9912
10047
basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
9928
10063
basisvalues[7] *= std::sqrt(10.0);
9929
10064
basisvalues[6] *= std::sqrt(14.0);
9931
// Table(s) of coefficients.
10066
// Table(s) of coefficients
9932
10067
static const double coefficients0[10] = \
9933
10068
{-0.18856181, -0.32716515, 0.1, -0.41731242, 0.080812204, 0.023328474, -0.21380899, 0.06023386, 0.015552316, -0.026937401};
9935
10070
static const double coefficients1[10] = \
9936
10071
{0.37712362, 0.0, -0.2, 0.83462485, 0.0, -0.046656947, 0.0, -0.12046772, 0.0, 0.053874802};
9938
// Compute value(s).
10073
// Compute value(s)
9939
10074
for (unsigned int r = 0; r < 10; r++)
9941
10076
values[3] += coefficients0[r]*basisvalues[r];
9942
10077
values[4] += coefficients1[r]*basisvalues[r];
9943
10078
}// end loop over 'r'
9945
// Using contravariant Piola transform to map values back to the physical element.
10080
// Using contravariant Piola transform to map values back to the physical element
9946
10081
const double tmp_ref0 = values[3];
9947
10082
const double tmp_ref1 = values[4];
9948
values[3] = (1.0/detJ)*(J_00*tmp_ref0 + J_01*tmp_ref1);
9949
values[4] = (1.0/detJ)*(J_10*tmp_ref0 + J_11*tmp_ref1);
10083
values[3] = (1.0/detJ)*(J[0]*tmp_ref0 + J[1]*tmp_ref1);
10084
values[4] = (1.0/detJ)*(J[2]*tmp_ref0 + J[3]*tmp_ref1);
9955
// Array of basisvalues.
10090
// Array of basisvalues
9956
10091
double basisvalues[10] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
9958
// Declare helper variables.
10093
// Declare helper variables
9959
10094
double tmp0 = (1.0 + Y + 2.0*X)/2.0;
9960
10095
double tmp1 = (1.0 - Y)/2.0;
9961
10096
double tmp2 = tmp1*tmp1;
9963
// Compute basisvalues.
10098
// Compute basisvalues
9964
10099
basisvalues[0] = 1.0;
9965
10100
basisvalues[1] = tmp0;
9966
10101
basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
9982
10117
basisvalues[7] *= std::sqrt(10.0);
9983
10118
basisvalues[6] *= std::sqrt(14.0);
9985
// Table(s) of coefficients.
10120
// Table(s) of coefficients
9986
10121
static const double coefficients0[10] = \
9987
10122
{0.18856181, 0.28867513, -0.16666667, 0.26082027, -0.20203051, 0.11664237, 0.1069045, -0.09035079, 0.069985421, -0.040406102};
9989
10124
static const double coefficients1[10] = \
9990
10125
{-0.32998316, -0.15396007, 0.2, -0.41731242, 0.25590531, -0.12830661, 0.0, 0.06023386, -0.077761579, 0.053874802};
9992
// Compute value(s).
10127
// Compute value(s)
9993
10128
for (unsigned int r = 0; r < 10; r++)
9995
10130
values[3] += coefficients0[r]*basisvalues[r];
9996
10131
values[4] += coefficients1[r]*basisvalues[r];
9997
10132
}// end loop over 'r'
9999
// Using contravariant Piola transform to map values back to the physical element.
10134
// Using contravariant Piola transform to map values back to the physical element
10000
10135
const double tmp_ref0 = values[3];
10001
10136
const double tmp_ref1 = values[4];
10002
values[3] = (1.0/detJ)*(J_00*tmp_ref0 + J_01*tmp_ref1);
10003
values[4] = (1.0/detJ)*(J_10*tmp_ref0 + J_11*tmp_ref1);
10137
values[3] = (1.0/detJ)*(J[0]*tmp_ref0 + J[1]*tmp_ref1);
10138
values[4] = (1.0/detJ)*(J[2]*tmp_ref0 + J[3]*tmp_ref1);
10009
// Array of basisvalues.
10144
// Array of basisvalues
10010
10145
double basisvalues[10] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
10012
// Declare helper variables.
10147
// Declare helper variables
10013
10148
double tmp0 = (1.0 + Y + 2.0*X)/2.0;
10014
10149
double tmp1 = (1.0 - Y)/2.0;
10015
10150
double tmp2 = tmp1*tmp1;
10017
// Compute basisvalues.
10152
// Compute basisvalues
10018
10153
basisvalues[0] = 1.0;
10019
10154
basisvalues[1] = tmp0;
10020
10155
basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
10036
10171
basisvalues[7] *= std::sqrt(10.0);
10037
10172
basisvalues[6] *= std::sqrt(14.0);
10039
// Table(s) of coefficients.
10174
// Table(s) of coefficients
10040
10175
static const double coefficients0[10] = \
10041
10176
{0.32998316, -0.69282032, -0.53333333, -0.39992441, -0.026937401, 0.20606818, 0.32071349, -0.03011693, -0.023328474, 0.013468701};
10043
10178
static const double coefficients1[10] = \
10044
10179
{0.23570226, 0.038490018, 0.066666667, 0.20865621, 0.12121831, -0.081649658, 0.0, 0.18070158, 0.077761579, 0.0};
10046
// Compute value(s).
10181
// Compute value(s)
10047
10182
for (unsigned int r = 0; r < 10; r++)
10049
10184
values[3] += coefficients0[r]*basisvalues[r];
10050
10185
values[4] += coefficients1[r]*basisvalues[r];
10051
10186
}// end loop over 'r'
10053
// Using contravariant Piola transform to map values back to the physical element.
10188
// Using contravariant Piola transform to map values back to the physical element
10054
10189
const double tmp_ref0 = values[3];
10055
10190
const double tmp_ref1 = values[4];
10056
values[3] = (1.0/detJ)*(J_00*tmp_ref0 + J_01*tmp_ref1);
10057
values[4] = (1.0/detJ)*(J_10*tmp_ref0 + J_11*tmp_ref1);
10191
values[3] = (1.0/detJ)*(J[0]*tmp_ref0 + J[1]*tmp_ref1);
10192
values[4] = (1.0/detJ)*(J[2]*tmp_ref0 + J[3]*tmp_ref1);
10063
// Array of basisvalues.
10198
// Array of basisvalues
10064
10199
double basisvalues[10] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
10066
// Declare helper variables.
10201
// Declare helper variables
10067
10202
double tmp0 = (1.0 + Y + 2.0*X)/2.0;
10068
10203
double tmp1 = (1.0 - Y)/2.0;
10069
10204
double tmp2 = tmp1*tmp1;
10071
// Compute basisvalues.
10206
// Compute basisvalues
10072
10207
basisvalues[0] = 1.0;
10073
10208
basisvalues[1] = tmp0;
10074
10209
basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
10090
10225
basisvalues[7] *= std::sqrt(10.0);
10091
10226
basisvalues[6] *= std::sqrt(14.0);
10093
// Table(s) of coefficients.
10228
// Table(s) of coefficients
10094
10229
static const double coefficients0[10] = \
10095
10230
{0.56568542, 0.65433031, -0.46666667, -0.19126819, -0.094280904, 0.12441853, -0.32071349, 0.15058465, -0.054433105, 0.013468701};
10097
10232
static const double coefficients1[10] = \
10098
10233
{-0.23570226, 0.038490018, -0.066666667, -0.20865621, 0.12121831, 0.081649658, 0.0, -0.18070158, 0.077761579, 0.0};
10100
// Compute value(s).
10235
// Compute value(s)
10101
10236
for (unsigned int r = 0; r < 10; r++)
10103
10238
values[3] += coefficients0[r]*basisvalues[r];
10104
10239
values[4] += coefficients1[r]*basisvalues[r];
10105
10240
}// end loop over 'r'
10107
// Using contravariant Piola transform to map values back to the physical element.
10242
// Using contravariant Piola transform to map values back to the physical element
10108
10243
const double tmp_ref0 = values[3];
10109
10244
const double tmp_ref1 = values[4];
10110
values[3] = (1.0/detJ)*(J_00*tmp_ref0 + J_01*tmp_ref1);
10111
values[4] = (1.0/detJ)*(J_10*tmp_ref0 + J_11*tmp_ref1);
10245
values[3] = (1.0/detJ)*(J[0]*tmp_ref0 + J[1]*tmp_ref1);
10246
values[4] = (1.0/detJ)*(J[2]*tmp_ref0 + J[3]*tmp_ref1);
10117
// Array of basisvalues.
10252
// Array of basisvalues
10118
10253
double basisvalues[10] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
10120
// Declare helper variables.
10255
// Declare helper variables
10121
10256
double tmp0 = (1.0 + Y + 2.0*X)/2.0;
10122
10257
double tmp1 = (1.0 - Y)/2.0;
10123
10258
double tmp2 = tmp1*tmp1;
10125
// Compute basisvalues.
10260
// Compute basisvalues
10126
10261
basisvalues[0] = 1.0;
10127
10262
basisvalues[1] = tmp0;
10128
10263
basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
10144
10279
basisvalues[7] *= std::sqrt(10.0);
10145
10280
basisvalues[6] *= std::sqrt(14.0);
10147
// Table(s) of coefficients.
10282
// Table(s) of coefficients
10148
10283
static const double coefficients0[10] = \
10149
10284
{0.094280904, 0.076980036, 0.93333333, 0.20865621, 0.24243661, -0.443241, 0.0, -0.24093544, 0.15552316, -0.053874802};
10151
10286
static const double coefficients1[10] = \
10152
10287
{0.0, -0.15396007, 0.0, 0.0, -0.48487322, 0.0, 0.0, 0.0, -0.31104632, 0.0};
10154
// Compute value(s).
10289
// Compute value(s)
10155
10290
for (unsigned int r = 0; r < 10; r++)
10157
10292
values[3] += coefficients0[r]*basisvalues[r];
10158
10293
values[4] += coefficients1[r]*basisvalues[r];
10159
10294
}// end loop over 'r'
10161
// Using contravariant Piola transform to map values back to the physical element.
10296
// Using contravariant Piola transform to map values back to the physical element
10162
10297
const double tmp_ref0 = values[3];
10163
10298
const double tmp_ref1 = values[4];
10164
values[3] = (1.0/detJ)*(J_00*tmp_ref0 + J_01*tmp_ref1);
10165
values[4] = (1.0/detJ)*(J_10*tmp_ref0 + J_11*tmp_ref1);
10299
values[3] = (1.0/detJ)*(J[0]*tmp_ref0 + J[1]*tmp_ref1);
10300
values[4] = (1.0/detJ)*(J[2]*tmp_ref0 + J[3]*tmp_ref1);
10171
// Array of basisvalues.
10306
// Array of basisvalues
10172
10307
double basisvalues[10] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
10174
// Declare helper variables.
10309
// Declare helper variables
10175
10310
double tmp0 = (1.0 + Y + 2.0*X)/2.0;
10176
10311
double tmp1 = (1.0 - Y)/2.0;
10177
10312
double tmp2 = tmp1*tmp1;
10179
// Compute basisvalues.
10314
// Compute basisvalues
10180
10315
basisvalues[0] = 1.0;
10181
10316
basisvalues[1] = tmp0;
10182
10317
basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
10198
10333
basisvalues[7] *= std::sqrt(10.0);
10199
10334
basisvalues[6] *= std::sqrt(14.0);
10201
// Table(s) of coefficients.
10336
// Table(s) of coefficients
10202
10337
static const double coefficients0[10] = \
10203
10338
{0.23570226, 0.076980036, 0.0, 0.052164053, 0.24243661, 0.058321184, 0.1069045, 0.15058465, -0.0077761579, -0.067343503};
10205
10340
static const double coefficients1[10] = \
10206
10341
{0.32998316, -0.80829038, -0.33333333, 0.069552071, -0.39059232, -0.21384434, 0.0, 0.06023386, 0.23328474, 0.21549921};
10208
// Compute value(s).
10343
// Compute value(s)
10209
10344
for (unsigned int r = 0; r < 10; r++)
10211
10346
values[3] += coefficients0[r]*basisvalues[r];
10212
10347
values[4] += coefficients1[r]*basisvalues[r];
10213
10348
}// end loop over 'r'
10215
// Using contravariant Piola transform to map values back to the physical element.
10350
// Using contravariant Piola transform to map values back to the physical element
10216
10351
const double tmp_ref0 = values[3];
10217
10352
const double tmp_ref1 = values[4];
10218
values[3] = (1.0/detJ)*(J_00*tmp_ref0 + J_01*tmp_ref1);
10219
values[4] = (1.0/detJ)*(J_10*tmp_ref0 + J_11*tmp_ref1);
10353
values[3] = (1.0/detJ)*(J[0]*tmp_ref0 + J[1]*tmp_ref1);
10354
values[4] = (1.0/detJ)*(J[2]*tmp_ref0 + J[3]*tmp_ref1);
10225
// Array of basisvalues.
10360
// Array of basisvalues
10226
10361
double basisvalues[10] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
10228
// Declare helper variables.
10363
// Declare helper variables
10229
10364
double tmp0 = (1.0 + Y + 2.0*X)/2.0;
10230
10365
double tmp1 = (1.0 - Y)/2.0;
10231
10366
double tmp2 = tmp1*tmp1;
10233
// Compute basisvalues.
10368
// Compute basisvalues
10234
10369
basisvalues[0] = 1.0;
10235
10370
basisvalues[1] = tmp0;
10236
10371
basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
10252
10387
basisvalues[7] *= std::sqrt(10.0);
10253
10388
basisvalues[6] *= std::sqrt(14.0);
10255
// Table(s) of coefficients.
10390
// Table(s) of coefficients
10256
10391
static const double coefficients0[10] = \
10257
10392
{0.0, -0.076980036, -0.13333333, -0.31298432, -0.24243661, 0.27994168, -0.21380899, -0.06023386, 0.17107547, -0.13468701};
10259
10394
static const double coefficients1[10] = \
10260
10395
{0.094280904, 0.84678039, -0.4, -0.13910414, 0.51181062, -0.13219468, 0.0, -0.12046772, -0.15552316, 0.21549921};
10262
// Compute value(s).
10397
// Compute value(s)
10263
10398
for (unsigned int r = 0; r < 10; r++)
10265
10400
values[3] += coefficients0[r]*basisvalues[r];
10266
10401
values[4] += coefficients1[r]*basisvalues[r];
10267
10402
}// end loop over 'r'
10269
// Using contravariant Piola transform to map values back to the physical element.
10404
// Using contravariant Piola transform to map values back to the physical element
10270
10405
const double tmp_ref0 = values[3];
10271
10406
const double tmp_ref1 = values[4];
10272
values[3] = (1.0/detJ)*(J_00*tmp_ref0 + J_01*tmp_ref1);
10273
values[4] = (1.0/detJ)*(J_10*tmp_ref0 + J_11*tmp_ref1);
10407
values[3] = (1.0/detJ)*(J[0]*tmp_ref0 + J[1]*tmp_ref1);
10408
values[4] = (1.0/detJ)*(J[2]*tmp_ref0 + J[3]*tmp_ref1);
10279
// Array of basisvalues.
10414
// Array of basisvalues
10280
10415
double basisvalues[10] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
10282
// Declare helper variables.
10417
// Declare helper variables
10283
10418
double tmp0 = (1.0 + Y + 2.0*X)/2.0;
10284
10419
double tmp1 = (1.0 - Y)/2.0;
10285
10420
double tmp2 = tmp1*tmp1;
10287
// Compute basisvalues.
10422
// Compute basisvalues
10288
10423
basisvalues[0] = 1.0;
10289
10424
basisvalues[1] = tmp0;
10290
10425
basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
10306
10441
basisvalues[7] *= std::sqrt(10.0);
10307
10442
basisvalues[6] *= std::sqrt(14.0);
10309
// Table(s) of coefficients.
10444
// Table(s) of coefficients
10310
10445
static const double coefficients0[10] = \
10311
10446
{-0.23570226, -0.038490018, 0.066666667, 0.10432811, -0.12121831, -0.19829203, 0.0, -0.12046772, -0.077761579, 0.13468701};
10313
10448
static const double coefficients1[10] = \
10314
10449
{0.56568542, -0.076980036, 0.8, 0.0, -0.24243661, -0.046656947, 0.0, 0.0, -0.15552316, -0.32324881};
10316
// Compute value(s).
10451
// Compute value(s)
10317
10452
for (unsigned int r = 0; r < 10; r++)
10319
10454
values[3] += coefficients0[r]*basisvalues[r];
10320
10455
values[4] += coefficients1[r]*basisvalues[r];
10321
10456
}// end loop over 'r'
10323
// Using contravariant Piola transform to map values back to the physical element.
10458
// Using contravariant Piola transform to map values back to the physical element
10324
10459
const double tmp_ref0 = values[3];
10325
10460
const double tmp_ref1 = values[4];
10326
values[3] = (1.0/detJ)*(J_00*tmp_ref0 + J_01*tmp_ref1);
10327
values[4] = (1.0/detJ)*(J_10*tmp_ref0 + J_11*tmp_ref1);
10461
values[3] = (1.0/detJ)*(J[0]*tmp_ref0 + J[1]*tmp_ref1);
10462
values[4] = (1.0/detJ)*(J[2]*tmp_ref0 + J[3]*tmp_ref1);
10334
/// Evaluate all basis functions at given point in cell
10469
/// Evaluate all basis functions at given point x in cell
10335
10470
virtual void evaluate_basis_all(double* values,
10336
const double* coordinates,
10337
const ufc::cell& c) const
10472
const double* vertex_coordinates,
10473
int cell_orientation) const
10339
10475
// Helper variable to hold values of a single dof.
10340
10476
double dof_values[5] = {0.0, 0.0, 0.0, 0.0, 0.0};
10342
// Loop dofs and call evaluate_basis.
10478
// Loop dofs and call evaluate_basis
10343
10479
for (unsigned int r = 0; r < 23; r++)
10345
evaluate_basis(r, dof_values, coordinates, c);
10481
evaluate_basis(r, dof_values, x, vertex_coordinates, cell_orientation);
10346
10482
for (unsigned int s = 0; s < 5; s++)
10348
10484
values[r*5 + s] = dof_values[s];
10350
10486
}// end loop over 'r'
10353
/// Evaluate order n derivatives of basis function i at given point in cell
10354
virtual void evaluate_basis_derivatives(unsigned int i,
10489
/// Evaluate order n derivatives of basis function i at given point x in cell
10490
virtual void evaluate_basis_derivatives(std::size_t i,
10356
10492
double* values,
10357
const double* coordinates,
10358
const ufc::cell& c) const
10494
const double* vertex_coordinates,
10495
int cell_orientation) const
10360
// Extract vertex coordinates
10361
const double * const * x = c.coordinates;
10363
// Compute Jacobian of affine map from reference cell
10364
const double J_00 = x[1][0] - x[0][0];
10365
const double J_01 = x[2][0] - x[0][0];
10366
const double J_10 = x[1][1] - x[0][1];
10367
const double J_11 = x[2][1] - x[0][1];
10369
// Compute determinant of Jacobian
10370
double detJ = J_00*J_11 - J_01*J_10;
10372
// Compute inverse of Jacobian
10373
const double K_00 = J_11 / detJ;
10374
const double K_01 = -J_01 / detJ;
10375
const double K_10 = -J_10 / detJ;
10376
const double K_11 = J_00 / detJ;
10497
// Compute Jacobian
10499
compute_jacobian_triangle_2d(J, vertex_coordinates);
10501
// Compute Jacobian inverse and determinant
10504
compute_jacobian_inverse_triangle_2d(K, detJ, J);
10378
10508
// Compute constants
10379
const double C0 = x[1][0] + x[2][0];
10380
const double C1 = x[1][1] + x[2][1];
10509
const double C0 = vertex_coordinates[2] + vertex_coordinates[4];
10510
const double C1 = vertex_coordinates[3] + vertex_coordinates[5];
10382
10512
// Get coordinates and map to the reference (FIAT) element
10383
double X = (J_01*(C1 - 2.0*coordinates[1]) + J_11*(2.0*coordinates[0] - C0)) / detJ;
10384
double Y = (J_00*(2.0*coordinates[1] - C1) + J_10*(C0 - 2.0*coordinates[0])) / detJ;
10513
double X = (J[1]*(C1 - 2.0*x[1]) + J[3]*(2.0*x[0] - C0)) / detJ;
10514
double Y = (J[0]*(2.0*x[1] - C1) + J[2]*(C0 - 2.0*x[0])) / detJ;
10386
10516
// Compute number of derivatives.
10387
10517
unsigned int num_derivatives = 1;
14765
15046
/// Evaluate linear functional for dof i on the function f
14766
virtual double evaluate_dof(unsigned int i,
15047
virtual double evaluate_dof(std::size_t i,
14767
15048
const ufc::function& f,
15049
const double* vertex_coordinates,
15050
int cell_orientation,
14768
15051
const ufc::cell& c) const
14770
// Declare variables for result of evaluation.
15053
// Declare variables for result of evaluation
14771
15054
double vals[5];
14773
// Declare variable for physical coordinates.
15056
// Declare variable for physical coordinates
14775
15059
double result;
14776
// Extract vertex coordinates
14777
const double * const * x = c.coordinates;
14779
// Compute Jacobian of affine map from reference cell
14780
const double J_00 = x[1][0] - x[0][0];
14781
const double J_01 = x[2][0] - x[0][0];
14782
const double J_10 = x[1][1] - x[0][1];
14783
const double J_11 = x[2][1] - x[0][1];
14785
// Compute determinant of Jacobian
14786
double detJ = J_00*J_11 - J_01*J_10;
14788
// Compute inverse of Jacobian
14789
const double K_00 = J_11 / detJ;
14790
const double K_01 = -J_01 / detJ;
14791
const double K_10 = -J_10 / detJ;
14792
const double K_11 = J_00 / detJ;switch (i)
15060
// Compute Jacobian
15062
compute_jacobian_triangle_2d(J, vertex_coordinates);
15065
// Compute Jacobian inverse and determinant
15068
compute_jacobian_inverse_triangle_2d(K, detJ, J);
14796
y[0] = 0.33333333*x[0][0] + 0.33333333*x[1][0] + 0.33333333*x[2][0];
14797
y[1] = 0.33333333*x[0][1] + 0.33333333*x[1][1] + 0.33333333*x[2][1];
15075
y[0] = 0.33333333*vertex_coordinates[0] + 0.33333333*vertex_coordinates[2] + 0.33333333*vertex_coordinates[4];
15076
y[1] = 0.33333333*vertex_coordinates[1] + 0.33333333*vertex_coordinates[3] + 0.33333333*vertex_coordinates[5];
14798
15077
f.evaluate(vals, y, c);
14799
15078
return vals[0];
14804
y[0] = 0.33333333*x[0][0] + 0.33333333*x[1][0] + 0.33333333*x[2][0];
14805
y[1] = 0.33333333*x[0][1] + 0.33333333*x[1][1] + 0.33333333*x[2][1];
15083
y[0] = 0.33333333*vertex_coordinates[0] + 0.33333333*vertex_coordinates[2] + 0.33333333*vertex_coordinates[4];
15084
y[1] = 0.33333333*vertex_coordinates[1] + 0.33333333*vertex_coordinates[3] + 0.33333333*vertex_coordinates[5];
14806
15085
f.evaluate(vals, y, c);
14807
15086
return vals[1];
15091
y[0] = vertex_coordinates[0];
15092
y[1] = vertex_coordinates[1];
14814
15093
f.evaluate(vals, y, c);
14815
15094
return vals[2];
15099
y[0] = vertex_coordinates[2];
15100
y[1] = vertex_coordinates[3];
14822
15101
f.evaluate(vals, y, c);
14823
15102
return vals[2];
15107
y[0] = vertex_coordinates[4];
15108
y[1] = vertex_coordinates[5];
14830
15109
f.evaluate(vals, y, c);
14831
15110
return vals[2];
14836
y[0] = 0.5*x[1][0] + 0.5*x[2][0];
14837
y[1] = 0.5*x[1][1] + 0.5*x[2][1];
15115
y[0] = 0.5*vertex_coordinates[2] + 0.5*vertex_coordinates[4];
15116
y[1] = 0.5*vertex_coordinates[3] + 0.5*vertex_coordinates[5];
14838
15117
f.evaluate(vals, y, c);
14839
15118
return vals[2];
14844
y[0] = 0.5*x[0][0] + 0.5*x[2][0];
14845
y[1] = 0.5*x[0][1] + 0.5*x[2][1];
15123
y[0] = 0.5*vertex_coordinates[0] + 0.5*vertex_coordinates[4];
15124
y[1] = 0.5*vertex_coordinates[1] + 0.5*vertex_coordinates[5];
14846
15125
f.evaluate(vals, y, c);
14847
15126
return vals[2];
14852
y[0] = 0.5*x[0][0] + 0.5*x[1][0];
14853
y[1] = 0.5*x[0][1] + 0.5*x[1][1];
15131
y[0] = 0.5*vertex_coordinates[0] + 0.5*vertex_coordinates[2];
15132
y[1] = 0.5*vertex_coordinates[1] + 0.5*vertex_coordinates[3];
14854
15133
f.evaluate(vals, y, c);
14855
15134
return vals[2];
14860
y[0] = 0.75*x[1][0] + 0.25*x[2][0];
14861
y[1] = 0.75*x[1][1] + 0.25*x[2][1];
15139
y[0] = 0.75*vertex_coordinates[2] + 0.25*vertex_coordinates[4];
15140
y[1] = 0.75*vertex_coordinates[3] + 0.25*vertex_coordinates[5];
14862
15141
f.evaluate(vals, y, c);
14863
result = (detJ*(K_00*vals[3] + K_01*vals[4])) + (detJ*(K_10*vals[3] + K_11*vals[4]));
15142
result = (detJ*(K[0]*vals[3] + K[1]*vals[4])) + (detJ*(K[2]*vals[3] + K[3]*vals[4]));
14864
15143
return result;
14869
y[0] = 0.5*x[1][0] + 0.5*x[2][0];
14870
y[1] = 0.5*x[1][1] + 0.5*x[2][1];
15148
y[0] = 0.5*vertex_coordinates[2] + 0.5*vertex_coordinates[4];
15149
y[1] = 0.5*vertex_coordinates[3] + 0.5*vertex_coordinates[5];
14871
15150
f.evaluate(vals, y, c);
14872
result = (detJ*(K_00*vals[3] + K_01*vals[4])) + (detJ*(K_10*vals[3] + K_11*vals[4]));
15151
result = (detJ*(K[0]*vals[3] + K[1]*vals[4])) + (detJ*(K[2]*vals[3] + K[3]*vals[4]));
14873
15152
return result;
14878
y[0] = 0.25*x[1][0] + 0.75*x[2][0];
14879
y[1] = 0.25*x[1][1] + 0.75*x[2][1];
15157
y[0] = 0.25*vertex_coordinates[2] + 0.75*vertex_coordinates[4];
15158
y[1] = 0.25*vertex_coordinates[3] + 0.75*vertex_coordinates[5];
14880
15159
f.evaluate(vals, y, c);
14881
result = (detJ*(K_00*vals[3] + K_01*vals[4])) + (detJ*(K_10*vals[3] + K_11*vals[4]));
15160
result = (detJ*(K[0]*vals[3] + K[1]*vals[4])) + (detJ*(K[2]*vals[3] + K[3]*vals[4]));
14882
15161
return result;
14887
y[0] = 0.75*x[0][0] + 0.25*x[2][0];
14888
y[1] = 0.75*x[0][1] + 0.25*x[2][1];
15166
y[0] = 0.75*vertex_coordinates[0] + 0.25*vertex_coordinates[4];
15167
y[1] = 0.75*vertex_coordinates[1] + 0.25*vertex_coordinates[5];
14889
15168
f.evaluate(vals, y, c);
14890
result = (detJ*(K_00*vals[3] + K_01*vals[4]));
15169
result = (detJ*(K[0]*vals[3] + K[1]*vals[4]));
14891
15170
return result;
14896
y[0] = 0.5*x[0][0] + 0.5*x[2][0];
14897
y[1] = 0.5*x[0][1] + 0.5*x[2][1];
15175
y[0] = 0.5*vertex_coordinates[0] + 0.5*vertex_coordinates[4];
15176
y[1] = 0.5*vertex_coordinates[1] + 0.5*vertex_coordinates[5];
14898
15177
f.evaluate(vals, y, c);
14899
result = (detJ*(K_00*vals[3] + K_01*vals[4]));
15178
result = (detJ*(K[0]*vals[3] + K[1]*vals[4]));
14900
15179
return result;
14905
y[0] = 0.25*x[0][0] + 0.75*x[2][0];
14906
y[1] = 0.25*x[0][1] + 0.75*x[2][1];
15184
y[0] = 0.25*vertex_coordinates[0] + 0.75*vertex_coordinates[4];
15185
y[1] = 0.25*vertex_coordinates[1] + 0.75*vertex_coordinates[5];
14907
15186
f.evaluate(vals, y, c);
14908
result = (detJ*(K_00*vals[3] + K_01*vals[4]));
15187
result = (detJ*(K[0]*vals[3] + K[1]*vals[4]));
14909
15188
return result;
14914
y[0] = 0.75*x[0][0] + 0.25*x[1][0];
14915
y[1] = 0.75*x[0][1] + 0.25*x[1][1];
15193
y[0] = 0.75*vertex_coordinates[0] + 0.25*vertex_coordinates[2];
15194
y[1] = 0.75*vertex_coordinates[1] + 0.25*vertex_coordinates[3];
14916
15195
f.evaluate(vals, y, c);
14917
result = (-1.0)*(detJ*(K_10*vals[3] + K_11*vals[4]));
15196
result = (-1.0)*(detJ*(K[2]*vals[3] + K[3]*vals[4]));
14918
15197
return result;
14923
y[0] = 0.5*x[0][0] + 0.5*x[1][0];
14924
y[1] = 0.5*x[0][1] + 0.5*x[1][1];
15202
y[0] = 0.5*vertex_coordinates[0] + 0.5*vertex_coordinates[2];
15203
y[1] = 0.5*vertex_coordinates[1] + 0.5*vertex_coordinates[3];
14925
15204
f.evaluate(vals, y, c);
14926
result = (-1.0)*(detJ*(K_10*vals[3] + K_11*vals[4]));
15205
result = (-1.0)*(detJ*(K[2]*vals[3] + K[3]*vals[4]));
14927
15206
return result;
14932
y[0] = 0.25*x[0][0] + 0.75*x[1][0];
14933
y[1] = 0.25*x[0][1] + 0.75*x[1][1];
15211
y[0] = 0.25*vertex_coordinates[0] + 0.75*vertex_coordinates[2];
15212
y[1] = 0.25*vertex_coordinates[1] + 0.75*vertex_coordinates[3];
14934
15213
f.evaluate(vals, y, c);
14935
result = (-1.0)*(detJ*(K_10*vals[3] + K_11*vals[4]));
15214
result = (-1.0)*(detJ*(K[2]*vals[3] + K[3]*vals[4]));
14936
15215
return result;
14941
y[0] = 0.5*x[0][0] + 0.25*x[1][0] + 0.25*x[2][0];
14942
y[1] = 0.5*x[0][1] + 0.25*x[1][1] + 0.25*x[2][1];
15220
y[0] = 0.5*vertex_coordinates[0] + 0.25*vertex_coordinates[2] + 0.25*vertex_coordinates[4];
15221
y[1] = 0.5*vertex_coordinates[1] + 0.25*vertex_coordinates[3] + 0.25*vertex_coordinates[5];
14943
15222
f.evaluate(vals, y, c);
14944
result = (detJ*(K_00*vals[3] + K_01*vals[4]));
15223
result = (detJ*(K[0]*vals[3] + K[1]*vals[4]));
14945
15224
return result;
14950
y[0] = 0.25*x[0][0] + 0.5*x[1][0] + 0.25*x[2][0];
14951
y[1] = 0.25*x[0][1] + 0.5*x[1][1] + 0.25*x[2][1];
15229
y[0] = 0.25*vertex_coordinates[0] + 0.5*vertex_coordinates[2] + 0.25*vertex_coordinates[4];
15230
y[1] = 0.25*vertex_coordinates[1] + 0.5*vertex_coordinates[3] + 0.25*vertex_coordinates[5];
14952
15231
f.evaluate(vals, y, c);
14953
result = (detJ*(K_00*vals[3] + K_01*vals[4]));
15232
result = (detJ*(K[0]*vals[3] + K[1]*vals[4]));
14954
15233
return result;
14959
y[0] = 0.25*x[0][0] + 0.25*x[1][0] + 0.5*x[2][0];
14960
y[1] = 0.25*x[0][1] + 0.25*x[1][1] + 0.5*x[2][1];
15238
y[0] = 0.25*vertex_coordinates[0] + 0.25*vertex_coordinates[2] + 0.5*vertex_coordinates[4];
15239
y[1] = 0.25*vertex_coordinates[1] + 0.25*vertex_coordinates[3] + 0.5*vertex_coordinates[5];
14961
15240
f.evaluate(vals, y, c);
14962
result = (detJ*(K_00*vals[3] + K_01*vals[4]));
15241
result = (detJ*(K[0]*vals[3] + K[1]*vals[4]));
14963
15242
return result;
14968
y[0] = 0.5*x[0][0] + 0.25*x[1][0] + 0.25*x[2][0];
14969
y[1] = 0.5*x[0][1] + 0.25*x[1][1] + 0.25*x[2][1];
15247
y[0] = 0.5*vertex_coordinates[0] + 0.25*vertex_coordinates[2] + 0.25*vertex_coordinates[4];
15248
y[1] = 0.5*vertex_coordinates[1] + 0.25*vertex_coordinates[3] + 0.25*vertex_coordinates[5];
14970
15249
f.evaluate(vals, y, c);
14971
result = (detJ*(K_10*vals[3] + K_11*vals[4]));
15250
result = (detJ*(K[2]*vals[3] + K[3]*vals[4]));
14972
15251
return result;
14977
y[0] = 0.25*x[0][0] + 0.5*x[1][0] + 0.25*x[2][0];
14978
y[1] = 0.25*x[0][1] + 0.5*x[1][1] + 0.25*x[2][1];
15256
y[0] = 0.25*vertex_coordinates[0] + 0.5*vertex_coordinates[2] + 0.25*vertex_coordinates[4];
15257
y[1] = 0.25*vertex_coordinates[1] + 0.5*vertex_coordinates[3] + 0.25*vertex_coordinates[5];
14979
15258
f.evaluate(vals, y, c);
14980
result = (detJ*(K_10*vals[3] + K_11*vals[4]));
15259
result = (detJ*(K[2]*vals[3] + K[3]*vals[4]));
14981
15260
return result;
14986
y[0] = 0.25*x[0][0] + 0.25*x[1][0] + 0.5*x[2][0];
14987
y[1] = 0.25*x[0][1] + 0.25*x[1][1] + 0.5*x[2][1];
15265
y[0] = 0.25*vertex_coordinates[0] + 0.25*vertex_coordinates[2] + 0.5*vertex_coordinates[4];
15266
y[1] = 0.25*vertex_coordinates[1] + 0.25*vertex_coordinates[3] + 0.5*vertex_coordinates[5];
14988
15267
f.evaluate(vals, y, c);
14989
result = (detJ*(K_10*vals[3] + K_11*vals[4]));
15268
result = (detJ*(K[2]*vals[3] + K[3]*vals[4]));
14990
15269
return result;
14998
15277
/// Evaluate linear functionals for all dofs on the function f
14999
15278
virtual void evaluate_dofs(double* values,
15000
15279
const ufc::function& f,
15280
const double* vertex_coordinates,
15281
int cell_orientation,
15001
15282
const ufc::cell& c) const
15003
// Declare variables for result of evaluation.
15284
// Declare variables for result of evaluation
15004
15285
double vals[5];
15006
// Declare variable for physical coordinates.
15287
// Declare variable for physical coordinates
15008
15290
double result;
15009
// Extract vertex coordinates
15010
const double * const * x = c.coordinates;
15012
// Compute Jacobian of affine map from reference cell
15013
const double J_00 = x[1][0] - x[0][0];
15014
const double J_01 = x[2][0] - x[0][0];
15015
const double J_10 = x[1][1] - x[0][1];
15016
const double J_11 = x[2][1] - x[0][1];
15018
// Compute determinant of Jacobian
15019
double detJ = J_00*J_11 - J_01*J_10;
15021
// Compute inverse of Jacobian
15022
const double K_00 = J_11 / detJ;
15023
const double K_01 = -J_01 / detJ;
15024
const double K_10 = -J_10 / detJ;
15025
const double K_11 = J_00 / detJ;y[0] = 0.33333333*x[0][0] + 0.33333333*x[1][0] + 0.33333333*x[2][0];
15026
y[1] = 0.33333333*x[0][1] + 0.33333333*x[1][1] + 0.33333333*x[2][1];
15291
// Compute Jacobian
15293
compute_jacobian_triangle_2d(J, vertex_coordinates);
15296
// Compute Jacobian inverse and determinant
15299
compute_jacobian_inverse_triangle_2d(K, detJ, J);
15302
y[0] = 0.33333333*vertex_coordinates[0] + 0.33333333*vertex_coordinates[2] + 0.33333333*vertex_coordinates[4];
15303
y[1] = 0.33333333*vertex_coordinates[1] + 0.33333333*vertex_coordinates[3] + 0.33333333*vertex_coordinates[5];
15027
15304
f.evaluate(vals, y, c);
15028
15305
values[0] = vals[0];
15029
y[0] = 0.33333333*x[0][0] + 0.33333333*x[1][0] + 0.33333333*x[2][0];
15030
y[1] = 0.33333333*x[0][1] + 0.33333333*x[1][1] + 0.33333333*x[2][1];
15306
y[0] = 0.33333333*vertex_coordinates[0] + 0.33333333*vertex_coordinates[2] + 0.33333333*vertex_coordinates[4];
15307
y[1] = 0.33333333*vertex_coordinates[1] + 0.33333333*vertex_coordinates[3] + 0.33333333*vertex_coordinates[5];
15031
15308
f.evaluate(vals, y, c);
15032
15309
values[1] = vals[1];
15310
y[0] = vertex_coordinates[0];
15311
y[1] = vertex_coordinates[1];
15035
15312
f.evaluate(vals, y, c);
15036
15313
values[2] = vals[2];
15314
y[0] = vertex_coordinates[2];
15315
y[1] = vertex_coordinates[3];
15039
15316
f.evaluate(vals, y, c);
15040
15317
values[3] = vals[2];
15318
y[0] = vertex_coordinates[4];
15319
y[1] = vertex_coordinates[5];
15043
15320
f.evaluate(vals, y, c);
15044
15321
values[4] = vals[2];
15045
y[0] = 0.5*x[1][0] + 0.5*x[2][0];
15046
y[1] = 0.5*x[1][1] + 0.5*x[2][1];
15322
y[0] = 0.5*vertex_coordinates[2] + 0.5*vertex_coordinates[4];
15323
y[1] = 0.5*vertex_coordinates[3] + 0.5*vertex_coordinates[5];
15047
15324
f.evaluate(vals, y, c);
15048
15325
values[5] = vals[2];
15049
y[0] = 0.5*x[0][0] + 0.5*x[2][0];
15050
y[1] = 0.5*x[0][1] + 0.5*x[2][1];
15326
y[0] = 0.5*vertex_coordinates[0] + 0.5*vertex_coordinates[4];
15327
y[1] = 0.5*vertex_coordinates[1] + 0.5*vertex_coordinates[5];
15051
15328
f.evaluate(vals, y, c);
15052
15329
values[6] = vals[2];
15053
y[0] = 0.5*x[0][0] + 0.5*x[1][0];
15054
y[1] = 0.5*x[0][1] + 0.5*x[1][1];
15330
y[0] = 0.5*vertex_coordinates[0] + 0.5*vertex_coordinates[2];
15331
y[1] = 0.5*vertex_coordinates[1] + 0.5*vertex_coordinates[3];
15055
15332
f.evaluate(vals, y, c);
15056
15333
values[7] = vals[2];
15057
y[0] = 0.75*x[1][0] + 0.25*x[2][0];
15058
y[1] = 0.75*x[1][1] + 0.25*x[2][1];
15334
y[0] = 0.75*vertex_coordinates[2] + 0.25*vertex_coordinates[4];
15335
y[1] = 0.75*vertex_coordinates[3] + 0.25*vertex_coordinates[5];
15059
15336
f.evaluate(vals, y, c);
15060
result = (detJ*(K_00*vals[3] + K_01*vals[4])) + (detJ*(K_10*vals[3] + K_11*vals[4]));
15337
result = (detJ*(K[0]*vals[3] + K[1]*vals[4])) + (detJ*(K[2]*vals[3] + K[3]*vals[4]));
15061
15338
values[8] = result;
15062
y[0] = 0.5*x[1][0] + 0.5*x[2][0];
15063
y[1] = 0.5*x[1][1] + 0.5*x[2][1];
15339
y[0] = 0.5*vertex_coordinates[2] + 0.5*vertex_coordinates[4];
15340
y[1] = 0.5*vertex_coordinates[3] + 0.5*vertex_coordinates[5];
15064
15341
f.evaluate(vals, y, c);
15065
result = (detJ*(K_00*vals[3] + K_01*vals[4])) + (detJ*(K_10*vals[3] + K_11*vals[4]));
15342
result = (detJ*(K[0]*vals[3] + K[1]*vals[4])) + (detJ*(K[2]*vals[3] + K[3]*vals[4]));
15066
15343
values[9] = result;
15067
y[0] = 0.25*x[1][0] + 0.75*x[2][0];
15068
y[1] = 0.25*x[1][1] + 0.75*x[2][1];
15344
y[0] = 0.25*vertex_coordinates[2] + 0.75*vertex_coordinates[4];
15345
y[1] = 0.25*vertex_coordinates[3] + 0.75*vertex_coordinates[5];
15069
15346
f.evaluate(vals, y, c);
15070
result = (detJ*(K_00*vals[3] + K_01*vals[4])) + (detJ*(K_10*vals[3] + K_11*vals[4]));
15347
result = (detJ*(K[0]*vals[3] + K[1]*vals[4])) + (detJ*(K[2]*vals[3] + K[3]*vals[4]));
15071
15348
values[10] = result;
15072
y[0] = 0.75*x[0][0] + 0.25*x[2][0];
15073
y[1] = 0.75*x[0][1] + 0.25*x[2][1];
15349
y[0] = 0.75*vertex_coordinates[0] + 0.25*vertex_coordinates[4];
15350
y[1] = 0.75*vertex_coordinates[1] + 0.25*vertex_coordinates[5];
15074
15351
f.evaluate(vals, y, c);
15075
result = (detJ*(K_00*vals[3] + K_01*vals[4]));
15352
result = (detJ*(K[0]*vals[3] + K[1]*vals[4]));
15076
15353
values[11] = result;
15077
y[0] = 0.5*x[0][0] + 0.5*x[2][0];
15078
y[1] = 0.5*x[0][1] + 0.5*x[2][1];
15354
y[0] = 0.5*vertex_coordinates[0] + 0.5*vertex_coordinates[4];
15355
y[1] = 0.5*vertex_coordinates[1] + 0.5*vertex_coordinates[5];
15079
15356
f.evaluate(vals, y, c);
15080
result = (detJ*(K_00*vals[3] + K_01*vals[4]));
15357
result = (detJ*(K[0]*vals[3] + K[1]*vals[4]));
15081
15358
values[12] = result;
15082
y[0] = 0.25*x[0][0] + 0.75*x[2][0];
15083
y[1] = 0.25*x[0][1] + 0.75*x[2][1];
15359
y[0] = 0.25*vertex_coordinates[0] + 0.75*vertex_coordinates[4];
15360
y[1] = 0.25*vertex_coordinates[1] + 0.75*vertex_coordinates[5];
15084
15361
f.evaluate(vals, y, c);
15085
result = (detJ*(K_00*vals[3] + K_01*vals[4]));
15362
result = (detJ*(K[0]*vals[3] + K[1]*vals[4]));
15086
15363
values[13] = result;
15087
y[0] = 0.75*x[0][0] + 0.25*x[1][0];
15088
y[1] = 0.75*x[0][1] + 0.25*x[1][1];
15364
y[0] = 0.75*vertex_coordinates[0] + 0.25*vertex_coordinates[2];
15365
y[1] = 0.75*vertex_coordinates[1] + 0.25*vertex_coordinates[3];
15089
15366
f.evaluate(vals, y, c);
15090
result = (-1.0)*(detJ*(K_10*vals[3] + K_11*vals[4]));
15367
result = (-1.0)*(detJ*(K[2]*vals[3] + K[3]*vals[4]));
15091
15368
values[14] = result;
15092
y[0] = 0.5*x[0][0] + 0.5*x[1][0];
15093
y[1] = 0.5*x[0][1] + 0.5*x[1][1];
15369
y[0] = 0.5*vertex_coordinates[0] + 0.5*vertex_coordinates[2];
15370
y[1] = 0.5*vertex_coordinates[1] + 0.5*vertex_coordinates[3];
15094
15371
f.evaluate(vals, y, c);
15095
result = (-1.0)*(detJ*(K_10*vals[3] + K_11*vals[4]));
15372
result = (-1.0)*(detJ*(K[2]*vals[3] + K[3]*vals[4]));
15096
15373
values[15] = result;
15097
y[0] = 0.25*x[0][0] + 0.75*x[1][0];
15098
y[1] = 0.25*x[0][1] + 0.75*x[1][1];
15374
y[0] = 0.25*vertex_coordinates[0] + 0.75*vertex_coordinates[2];
15375
y[1] = 0.25*vertex_coordinates[1] + 0.75*vertex_coordinates[3];
15099
15376
f.evaluate(vals, y, c);
15100
result = (-1.0)*(detJ*(K_10*vals[3] + K_11*vals[4]));
15377
result = (-1.0)*(detJ*(K[2]*vals[3] + K[3]*vals[4]));
15101
15378
values[16] = result;
15102
y[0] = 0.5*x[0][0] + 0.25*x[1][0] + 0.25*x[2][0];
15103
y[1] = 0.5*x[0][1] + 0.25*x[1][1] + 0.25*x[2][1];
15379
y[0] = 0.5*vertex_coordinates[0] + 0.25*vertex_coordinates[2] + 0.25*vertex_coordinates[4];
15380
y[1] = 0.5*vertex_coordinates[1] + 0.25*vertex_coordinates[3] + 0.25*vertex_coordinates[5];
15104
15381
f.evaluate(vals, y, c);
15105
result = (detJ*(K_00*vals[3] + K_01*vals[4]));
15382
result = (detJ*(K[0]*vals[3] + K[1]*vals[4]));
15106
15383
values[17] = result;
15107
y[0] = 0.25*x[0][0] + 0.5*x[1][0] + 0.25*x[2][0];
15108
y[1] = 0.25*x[0][1] + 0.5*x[1][1] + 0.25*x[2][1];
15384
y[0] = 0.25*vertex_coordinates[0] + 0.5*vertex_coordinates[2] + 0.25*vertex_coordinates[4];
15385
y[1] = 0.25*vertex_coordinates[1] + 0.5*vertex_coordinates[3] + 0.25*vertex_coordinates[5];
15109
15386
f.evaluate(vals, y, c);
15110
result = (detJ*(K_00*vals[3] + K_01*vals[4]));
15387
result = (detJ*(K[0]*vals[3] + K[1]*vals[4]));
15111
15388
values[18] = result;
15112
y[0] = 0.25*x[0][0] + 0.25*x[1][0] + 0.5*x[2][0];
15113
y[1] = 0.25*x[0][1] + 0.25*x[1][1] + 0.5*x[2][1];
15389
y[0] = 0.25*vertex_coordinates[0] + 0.25*vertex_coordinates[2] + 0.5*vertex_coordinates[4];
15390
y[1] = 0.25*vertex_coordinates[1] + 0.25*vertex_coordinates[3] + 0.5*vertex_coordinates[5];
15114
15391
f.evaluate(vals, y, c);
15115
result = (detJ*(K_00*vals[3] + K_01*vals[4]));
15392
result = (detJ*(K[0]*vals[3] + K[1]*vals[4]));
15116
15393
values[19] = result;
15117
y[0] = 0.5*x[0][0] + 0.25*x[1][0] + 0.25*x[2][0];
15118
y[1] = 0.5*x[0][1] + 0.25*x[1][1] + 0.25*x[2][1];
15394
y[0] = 0.5*vertex_coordinates[0] + 0.25*vertex_coordinates[2] + 0.25*vertex_coordinates[4];
15395
y[1] = 0.5*vertex_coordinates[1] + 0.25*vertex_coordinates[3] + 0.25*vertex_coordinates[5];
15119
15396
f.evaluate(vals, y, c);
15120
result = (detJ*(K_10*vals[3] + K_11*vals[4]));
15397
result = (detJ*(K[2]*vals[3] + K[3]*vals[4]));
15121
15398
values[20] = result;
15122
y[0] = 0.25*x[0][0] + 0.5*x[1][0] + 0.25*x[2][0];
15123
y[1] = 0.25*x[0][1] + 0.5*x[1][1] + 0.25*x[2][1];
15399
y[0] = 0.25*vertex_coordinates[0] + 0.5*vertex_coordinates[2] + 0.25*vertex_coordinates[4];
15400
y[1] = 0.25*vertex_coordinates[1] + 0.5*vertex_coordinates[3] + 0.25*vertex_coordinates[5];
15124
15401
f.evaluate(vals, y, c);
15125
result = (detJ*(K_10*vals[3] + K_11*vals[4]));
15402
result = (detJ*(K[2]*vals[3] + K[3]*vals[4]));
15126
15403
values[21] = result;
15127
y[0] = 0.25*x[0][0] + 0.25*x[1][0] + 0.5*x[2][0];
15128
y[1] = 0.25*x[0][1] + 0.25*x[1][1] + 0.5*x[2][1];
15404
y[0] = 0.25*vertex_coordinates[0] + 0.25*vertex_coordinates[2] + 0.5*vertex_coordinates[4];
15405
y[1] = 0.25*vertex_coordinates[1] + 0.25*vertex_coordinates[3] + 0.5*vertex_coordinates[5];
15129
15406
f.evaluate(vals, y, c);
15130
result = (detJ*(K_10*vals[3] + K_11*vals[4]));
15407
result = (detJ*(K[2]*vals[3] + K[3]*vals[4]));
15131
15408
values[22] = result;
15134
15411
/// Interpolate vertex values from dof values
15135
15412
virtual void interpolate_vertex_values(double* vertex_values,
15136
15413
const double* dof_values,
15414
const double* vertex_coordinates,
15415
int cell_orientation,
15137
15416
const ufc::cell& c) const
15139
// Extract vertex coordinates
15140
const double * const * x = c.coordinates;
15142
// Compute Jacobian of affine map from reference cell
15143
const double J_00 = x[1][0] - x[0][0];
15144
const double J_01 = x[2][0] - x[0][0];
15145
const double J_10 = x[1][1] - x[0][1];
15146
const double J_11 = x[2][1] - x[0][1];
15148
// Compute determinant of Jacobian
15149
double detJ = J_00*J_11 - J_01*J_10;
15151
// Compute inverse of Jacobian
15418
// Compute Jacobian
15420
compute_jacobian_triangle_2d(J, vertex_coordinates);
15423
// Compute Jacobian inverse and determinant
15426
compute_jacobian_inverse_triangle_2d(K, detJ, J);
15152
15430
// Evaluate function and change variables
15153
15431
vertex_values[0] = dof_values[0];
15154
15432
vertex_values[5] = dof_values[0];
15273
15548
return false;
15276
/// Initialize dofmap for mesh (return true iff init_cell() is needed)
15277
virtual bool init_mesh(const ufc::mesh& m)
15279
_global_dimension = m.num_entities[2];
15283
/// Initialize dofmap for given cell
15284
virtual void init_cell(const ufc::mesh& m,
15285
const ufc::cell& c)
15290
/// Finish initialization of dofmap for cells
15291
virtual void init_cell_finalize()
15296
15551
/// Return the topological dimension of the associated cell shape
15297
virtual unsigned int topological_dimension() const
15552
virtual std::size_t topological_dimension() const
15302
15557
/// Return the geometric dimension of the associated cell shape
15303
virtual unsigned int geometric_dimension() const
15558
virtual std::size_t geometric_dimension() const
15308
15563
/// Return the dimension of the global finite element function space
15309
virtual unsigned int global_dimension() const
15564
virtual std::size_t global_dimension(const std::vector<std::size_t>&
15565
num_global_entities) const
15311
return _global_dimension;
15567
return num_global_entities[2];
15314
15570
/// Return the dimension of the local finite element function space for a cell
15315
virtual unsigned int local_dimension(const ufc::cell& c) const
15571
virtual std::size_t local_dimension(const ufc::cell& c) const
15320
15576
/// Return the maximum dimension of the local finite element function space
15321
virtual unsigned int max_local_dimension() const
15577
virtual std::size_t max_local_dimension() const
15326
15582
/// Return the number of dofs on each cell facet
15327
virtual unsigned int num_facet_dofs() const
15583
virtual std::size_t num_facet_dofs() const
15332
15588
/// Return the number of dofs associated with each cell entity of dimension d
15333
virtual unsigned int num_entity_dofs(unsigned int d) const
15589
virtual std::size_t num_entity_dofs(std::size_t d) const
15505
15756
return false;
15508
/// Initialize dofmap for mesh (return true iff init_cell() is needed)
15509
virtual bool init_mesh(const ufc::mesh& m)
15511
_global_dimension = 2*m.num_entities[2];
15515
/// Initialize dofmap for given cell
15516
virtual void init_cell(const ufc::mesh& m,
15517
const ufc::cell& c)
15522
/// Finish initialization of dofmap for cells
15523
virtual void init_cell_finalize()
15528
15759
/// Return the topological dimension of the associated cell shape
15529
virtual unsigned int topological_dimension() const
15760
virtual std::size_t topological_dimension() const
15534
15765
/// Return the geometric dimension of the associated cell shape
15535
virtual unsigned int geometric_dimension() const
15766
virtual std::size_t geometric_dimension() const
15540
15771
/// Return the dimension of the global finite element function space
15541
virtual unsigned int global_dimension() const
15772
virtual std::size_t global_dimension(const std::vector<std::size_t>&
15773
num_global_entities) const
15543
return _global_dimension;
15775
return 2*num_global_entities[2];
15546
15778
/// Return the dimension of the local finite element function space for a cell
15547
virtual unsigned int local_dimension(const ufc::cell& c) const
15779
virtual std::size_t local_dimension(const ufc::cell& c) const
15552
15784
/// Return the maximum dimension of the local finite element function space
15553
virtual unsigned int max_local_dimension() const
15785
virtual std::size_t max_local_dimension() const
15558
15790
/// Return the number of dofs on each cell facet
15559
virtual unsigned int num_facet_dofs() const
15791
virtual std::size_t num_facet_dofs() const
15564
15796
/// Return the number of dofs associated with each cell entity of dimension d
15565
virtual unsigned int num_entity_dofs(unsigned int d) const
15797
virtual std::size_t num_entity_dofs(std::size_t d) const
15758
15985
return false;
15761
/// Initialize dofmap for mesh (return true iff init_cell() is needed)
15762
virtual bool init_mesh(const ufc::mesh& m)
15764
_global_dimension = m.num_entities[0] + m.num_entities[1];
15768
/// Initialize dofmap for given cell
15769
virtual void init_cell(const ufc::mesh& m,
15770
const ufc::cell& c)
15775
/// Finish initialization of dofmap for cells
15776
virtual void init_cell_finalize()
15781
15988
/// Return the topological dimension of the associated cell shape
15782
virtual unsigned int topological_dimension() const
15989
virtual std::size_t topological_dimension() const
15787
15994
/// Return the geometric dimension of the associated cell shape
15788
virtual unsigned int geometric_dimension() const
15995
virtual std::size_t geometric_dimension() const
15793
16000
/// Return the dimension of the global finite element function space
15794
virtual unsigned int global_dimension() const
16001
virtual std::size_t global_dimension(const std::vector<std::size_t>&
16002
num_global_entities) const
15796
return _global_dimension;
16004
return num_global_entities[0] + num_global_entities[1];
15799
16007
/// Return the dimension of the local finite element function space for a cell
15800
virtual unsigned int local_dimension(const ufc::cell& c) const
16008
virtual std::size_t local_dimension(const ufc::cell& c) const
15805
16013
/// Return the maximum dimension of the local finite element function space
15806
virtual unsigned int max_local_dimension() const
16014
virtual std::size_t max_local_dimension() const
15811
16019
/// Return the number of dofs on each cell facet
15812
virtual unsigned int num_facet_dofs() const
16020
virtual std::size_t num_facet_dofs() const
15817
16025
/// Return the number of dofs associated with each cell entity of dimension d
15818
virtual unsigned int num_entity_dofs(unsigned int d) const
16026
virtual std::size_t num_entity_dofs(std::size_t d) const
15965
16173
/// Tabulate the coordinates of all dofs on a cell
15966
virtual void tabulate_coordinates(double** coordinates,
15967
const ufc::cell& c) const
16174
virtual void tabulate_coordinates(double** dof_coordinates,
16175
const double* vertex_coordinates) const
15969
const double * const * x = c.coordinates;
15971
coordinates[0][0] = x[0][0];
15972
coordinates[0][1] = x[0][1];
15973
coordinates[1][0] = x[1][0];
15974
coordinates[1][1] = x[1][1];
15975
coordinates[2][0] = x[2][0];
15976
coordinates[2][1] = x[2][1];
15977
coordinates[3][0] = 0.5*x[1][0] + 0.5*x[2][0];
15978
coordinates[3][1] = 0.5*x[1][1] + 0.5*x[2][1];
15979
coordinates[4][0] = 0.5*x[0][0] + 0.5*x[2][0];
15980
coordinates[4][1] = 0.5*x[0][1] + 0.5*x[2][1];
15981
coordinates[5][0] = 0.5*x[0][0] + 0.5*x[1][0];
15982
coordinates[5][1] = 0.5*x[0][1] + 0.5*x[1][1];
16177
dof_coordinates[0][0] = vertex_coordinates[0];
16178
dof_coordinates[0][1] = vertex_coordinates[1];
16179
dof_coordinates[1][0] = vertex_coordinates[2];
16180
dof_coordinates[1][1] = vertex_coordinates[3];
16181
dof_coordinates[2][0] = vertex_coordinates[4];
16182
dof_coordinates[2][1] = vertex_coordinates[5];
16183
dof_coordinates[3][0] = 0.5*vertex_coordinates[2] + 0.5*vertex_coordinates[4];
16184
dof_coordinates[3][1] = 0.5*vertex_coordinates[3] + 0.5*vertex_coordinates[5];
16185
dof_coordinates[4][0] = 0.5*vertex_coordinates[0] + 0.5*vertex_coordinates[4];
16186
dof_coordinates[4][1] = 0.5*vertex_coordinates[1] + 0.5*vertex_coordinates[5];
16187
dof_coordinates[5][0] = 0.5*vertex_coordinates[0] + 0.5*vertex_coordinates[2];
16188
dof_coordinates[5][1] = 0.5*vertex_coordinates[1] + 0.5*vertex_coordinates[3];
15985
16191
/// Return the number of sub dofmaps (for a mixed element)
15986
virtual unsigned int num_sub_dofmaps() const
16192
virtual std::size_t num_sub_dofmaps() const
15991
16197
/// Create a new dofmap for sub dofmap i (for a mixed element)
15992
virtual ufc::dofmap* create_sub_dofmap(unsigned int i) const
16198
virtual ufc::dofmap* create_sub_dofmap(std::size_t i) const
16055
16258
return false;
16058
/// Initialize dofmap for mesh (return true iff init_cell() is needed)
16059
virtual bool init_mesh(const ufc::mesh& m)
16061
_global_dimension = m.num_entities[0] + m.num_entities[1] + 2*m.num_entities[2];
16065
/// Initialize dofmap for given cell
16066
virtual void init_cell(const ufc::mesh& m,
16067
const ufc::cell& c)
16072
/// Finish initialization of dofmap for cells
16073
virtual void init_cell_finalize()
16078
16261
/// Return the topological dimension of the associated cell shape
16079
virtual unsigned int topological_dimension() const
16262
virtual std::size_t topological_dimension() const
16084
16267
/// Return the geometric dimension of the associated cell shape
16085
virtual unsigned int geometric_dimension() const
16268
virtual std::size_t geometric_dimension() const
16090
16273
/// Return the dimension of the global finite element function space
16091
virtual unsigned int global_dimension() const
16274
virtual std::size_t global_dimension(const std::vector<std::size_t>&
16275
num_global_entities) const
16093
return _global_dimension;
16277
return num_global_entities[0] + num_global_entities[1] + 2*num_global_entities[2];
16096
16280
/// Return the dimension of the local finite element function space for a cell
16097
virtual unsigned int local_dimension(const ufc::cell& c) const
16281
virtual std::size_t local_dimension(const ufc::cell& c) const
16102
16286
/// Return the maximum dimension of the local finite element function space
16103
virtual unsigned int max_local_dimension() const
16287
virtual std::size_t max_local_dimension() const
16108
16292
/// Return the number of dofs on each cell facet
16109
virtual unsigned int num_facet_dofs() const
16293
virtual std::size_t num_facet_dofs() const
16114
16298
/// Return the number of dofs associated with each cell entity of dimension d
16115
virtual unsigned int num_entity_dofs(unsigned int d) const
16299
virtual std::size_t num_entity_dofs(std::size_t d) const
16272
16456
/// Tabulate the coordinates of all dofs on a cell
16273
virtual void tabulate_coordinates(double** coordinates,
16274
const ufc::cell& c) const
16457
virtual void tabulate_coordinates(double** dof_coordinates,
16458
const double* vertex_coordinates) const
16276
const double * const * x = c.coordinates;
16278
coordinates[0][0] = 0.33333333*x[0][0] + 0.33333333*x[1][0] + 0.33333333*x[2][0];
16279
coordinates[0][1] = 0.33333333*x[0][1] + 0.33333333*x[1][1] + 0.33333333*x[2][1];
16280
coordinates[1][0] = 0.33333333*x[0][0] + 0.33333333*x[1][0] + 0.33333333*x[2][0];
16281
coordinates[1][1] = 0.33333333*x[0][1] + 0.33333333*x[1][1] + 0.33333333*x[2][1];
16282
coordinates[2][0] = x[0][0];
16283
coordinates[2][1] = x[0][1];
16284
coordinates[3][0] = x[1][0];
16285
coordinates[3][1] = x[1][1];
16286
coordinates[4][0] = x[2][0];
16287
coordinates[4][1] = x[2][1];
16288
coordinates[5][0] = 0.5*x[1][0] + 0.5*x[2][0];
16289
coordinates[5][1] = 0.5*x[1][1] + 0.5*x[2][1];
16290
coordinates[6][0] = 0.5*x[0][0] + 0.5*x[2][0];
16291
coordinates[6][1] = 0.5*x[0][1] + 0.5*x[2][1];
16292
coordinates[7][0] = 0.5*x[0][0] + 0.5*x[1][0];
16293
coordinates[7][1] = 0.5*x[0][1] + 0.5*x[1][1];
16460
dof_coordinates[0][0] = 0.33333333*vertex_coordinates[0] + 0.33333333*vertex_coordinates[2] + 0.33333333*vertex_coordinates[4];
16461
dof_coordinates[0][1] = 0.33333333*vertex_coordinates[1] + 0.33333333*vertex_coordinates[3] + 0.33333333*vertex_coordinates[5];
16462
dof_coordinates[1][0] = 0.33333333*vertex_coordinates[0] + 0.33333333*vertex_coordinates[2] + 0.33333333*vertex_coordinates[4];
16463
dof_coordinates[1][1] = 0.33333333*vertex_coordinates[1] + 0.33333333*vertex_coordinates[3] + 0.33333333*vertex_coordinates[5];
16464
dof_coordinates[2][0] = vertex_coordinates[0];
16465
dof_coordinates[2][1] = vertex_coordinates[1];
16466
dof_coordinates[3][0] = vertex_coordinates[2];
16467
dof_coordinates[3][1] = vertex_coordinates[3];
16468
dof_coordinates[4][0] = vertex_coordinates[4];
16469
dof_coordinates[4][1] = vertex_coordinates[5];
16470
dof_coordinates[5][0] = 0.5*vertex_coordinates[2] + 0.5*vertex_coordinates[4];
16471
dof_coordinates[5][1] = 0.5*vertex_coordinates[3] + 0.5*vertex_coordinates[5];
16472
dof_coordinates[6][0] = 0.5*vertex_coordinates[0] + 0.5*vertex_coordinates[4];
16473
dof_coordinates[6][1] = 0.5*vertex_coordinates[1] + 0.5*vertex_coordinates[5];
16474
dof_coordinates[7][0] = 0.5*vertex_coordinates[0] + 0.5*vertex_coordinates[2];
16475
dof_coordinates[7][1] = 0.5*vertex_coordinates[1] + 0.5*vertex_coordinates[3];
16296
16478
/// Return the number of sub dofmaps (for a mixed element)
16297
virtual unsigned int num_sub_dofmaps() const
16479
virtual std::size_t num_sub_dofmaps() const
16302
16484
/// Create a new dofmap for sub dofmap i (for a mixed element)
16303
virtual ufc::dofmap* create_sub_dofmap(unsigned int i) const
16485
virtual ufc::dofmap* create_sub_dofmap(std::size_t i) const
16380
16559
return false;
16383
/// Initialize dofmap for mesh (return true iff init_cell() is needed)
16384
virtual bool init_mesh(const ufc::mesh& m)
16386
_global_dimension = 3*m.num_entities[1] + 6*m.num_entities[2];
16390
/// Initialize dofmap for given cell
16391
virtual void init_cell(const ufc::mesh& m,
16392
const ufc::cell& c)
16397
/// Finish initialization of dofmap for cells
16398
virtual void init_cell_finalize()
16403
16562
/// Return the topological dimension of the associated cell shape
16404
virtual unsigned int topological_dimension() const
16563
virtual std::size_t topological_dimension() const
16409
16568
/// Return the geometric dimension of the associated cell shape
16410
virtual unsigned int geometric_dimension() const
16569
virtual std::size_t geometric_dimension() const
16415
16574
/// Return the dimension of the global finite element function space
16416
virtual unsigned int global_dimension() const
16575
virtual std::size_t global_dimension(const std::vector<std::size_t>&
16576
num_global_entities) const
16418
return _global_dimension;
16578
return 3*num_global_entities[1] + 6*num_global_entities[2];
16421
16581
/// Return the dimension of the local finite element function space for a cell
16422
virtual unsigned int local_dimension(const ufc::cell& c) const
16582
virtual std::size_t local_dimension(const ufc::cell& c) const
16427
16587
/// Return the maximum dimension of the local finite element function space
16428
virtual unsigned int max_local_dimension() const
16588
virtual std::size_t max_local_dimension() const
16433
16593
/// Return the number of dofs on each cell facet
16434
virtual unsigned int num_facet_dofs() const
16594
virtual std::size_t num_facet_dofs() const
16439
16599
/// Return the number of dofs associated with each cell entity of dimension d
16440
virtual unsigned int num_entity_dofs(unsigned int d) const
16600
virtual std::size_t num_entity_dofs(std::size_t d) const
16589
16749
/// Tabulate the coordinates of all dofs on a cell
16590
virtual void tabulate_coordinates(double** coordinates,
16591
const ufc::cell& c) const
16750
virtual void tabulate_coordinates(double** dof_coordinates,
16751
const double* vertex_coordinates) const
16593
const double * const * x = c.coordinates;
16595
coordinates[0][0] = 0.75*x[1][0] + 0.25*x[2][0];
16596
coordinates[0][1] = 0.75*x[1][1] + 0.25*x[2][1];
16597
coordinates[1][0] = 0.5*x[1][0] + 0.5*x[2][0];
16598
coordinates[1][1] = 0.5*x[1][1] + 0.5*x[2][1];
16599
coordinates[2][0] = 0.25*x[1][0] + 0.75*x[2][0];
16600
coordinates[2][1] = 0.25*x[1][1] + 0.75*x[2][1];
16601
coordinates[3][0] = 0.75*x[0][0] + 0.25*x[2][0];
16602
coordinates[3][1] = 0.75*x[0][1] + 0.25*x[2][1];
16603
coordinates[4][0] = 0.5*x[0][0] + 0.5*x[2][0];
16604
coordinates[4][1] = 0.5*x[0][1] + 0.5*x[2][1];
16605
coordinates[5][0] = 0.25*x[0][0] + 0.75*x[2][0];
16606
coordinates[5][1] = 0.25*x[0][1] + 0.75*x[2][1];
16607
coordinates[6][0] = 0.75*x[0][0] + 0.25*x[1][0];
16608
coordinates[6][1] = 0.75*x[0][1] + 0.25*x[1][1];
16609
coordinates[7][0] = 0.5*x[0][0] + 0.5*x[1][0];
16610
coordinates[7][1] = 0.5*x[0][1] + 0.5*x[1][1];
16611
coordinates[8][0] = 0.25*x[0][0] + 0.75*x[1][0];
16612
coordinates[8][1] = 0.25*x[0][1] + 0.75*x[1][1];
16613
coordinates[9][0] = 0.5*x[0][0] + 0.25*x[1][0] + 0.25*x[2][0];
16614
coordinates[9][1] = 0.5*x[0][1] + 0.25*x[1][1] + 0.25*x[2][1];
16615
coordinates[10][0] = 0.25*x[0][0] + 0.5*x[1][0] + 0.25*x[2][0];
16616
coordinates[10][1] = 0.25*x[0][1] + 0.5*x[1][1] + 0.25*x[2][1];
16617
coordinates[11][0] = 0.25*x[0][0] + 0.25*x[1][0] + 0.5*x[2][0];
16618
coordinates[11][1] = 0.25*x[0][1] + 0.25*x[1][1] + 0.5*x[2][1];
16619
coordinates[12][0] = 0.5*x[0][0] + 0.25*x[1][0] + 0.25*x[2][0];
16620
coordinates[12][1] = 0.5*x[0][1] + 0.25*x[1][1] + 0.25*x[2][1];
16621
coordinates[13][0] = 0.25*x[0][0] + 0.5*x[1][0] + 0.25*x[2][0];
16622
coordinates[13][1] = 0.25*x[0][1] + 0.5*x[1][1] + 0.25*x[2][1];
16623
coordinates[14][0] = 0.25*x[0][0] + 0.25*x[1][0] + 0.5*x[2][0];
16624
coordinates[14][1] = 0.25*x[0][1] + 0.25*x[1][1] + 0.5*x[2][1];
16753
dof_coordinates[0][0] = 0.75*vertex_coordinates[2] + 0.25*vertex_coordinates[4];
16754
dof_coordinates[0][1] = 0.75*vertex_coordinates[3] + 0.25*vertex_coordinates[5];
16755
dof_coordinates[1][0] = 0.5*vertex_coordinates[2] + 0.5*vertex_coordinates[4];
16756
dof_coordinates[1][1] = 0.5*vertex_coordinates[3] + 0.5*vertex_coordinates[5];
16757
dof_coordinates[2][0] = 0.25*vertex_coordinates[2] + 0.75*vertex_coordinates[4];
16758
dof_coordinates[2][1] = 0.25*vertex_coordinates[3] + 0.75*vertex_coordinates[5];
16759
dof_coordinates[3][0] = 0.75*vertex_coordinates[0] + 0.25*vertex_coordinates[4];
16760
dof_coordinates[3][1] = 0.75*vertex_coordinates[1] + 0.25*vertex_coordinates[5];
16761
dof_coordinates[4][0] = 0.5*vertex_coordinates[0] + 0.5*vertex_coordinates[4];
16762
dof_coordinates[4][1] = 0.5*vertex_coordinates[1] + 0.5*vertex_coordinates[5];
16763
dof_coordinates[5][0] = 0.25*vertex_coordinates[0] + 0.75*vertex_coordinates[4];
16764
dof_coordinates[5][1] = 0.25*vertex_coordinates[1] + 0.75*vertex_coordinates[5];
16765
dof_coordinates[6][0] = 0.75*vertex_coordinates[0] + 0.25*vertex_coordinates[2];
16766
dof_coordinates[6][1] = 0.75*vertex_coordinates[1] + 0.25*vertex_coordinates[3];
16767
dof_coordinates[7][0] = 0.5*vertex_coordinates[0] + 0.5*vertex_coordinates[2];
16768
dof_coordinates[7][1] = 0.5*vertex_coordinates[1] + 0.5*vertex_coordinates[3];
16769
dof_coordinates[8][0] = 0.25*vertex_coordinates[0] + 0.75*vertex_coordinates[2];
16770
dof_coordinates[8][1] = 0.25*vertex_coordinates[1] + 0.75*vertex_coordinates[3];
16771
dof_coordinates[9][0] = 0.5*vertex_coordinates[0] + 0.25*vertex_coordinates[2] + 0.25*vertex_coordinates[4];
16772
dof_coordinates[9][1] = 0.5*vertex_coordinates[1] + 0.25*vertex_coordinates[3] + 0.25*vertex_coordinates[5];
16773
dof_coordinates[10][0] = 0.25*vertex_coordinates[0] + 0.5*vertex_coordinates[2] + 0.25*vertex_coordinates[4];
16774
dof_coordinates[10][1] = 0.25*vertex_coordinates[1] + 0.5*vertex_coordinates[3] + 0.25*vertex_coordinates[5];
16775
dof_coordinates[11][0] = 0.25*vertex_coordinates[0] + 0.25*vertex_coordinates[2] + 0.5*vertex_coordinates[4];
16776
dof_coordinates[11][1] = 0.25*vertex_coordinates[1] + 0.25*vertex_coordinates[3] + 0.5*vertex_coordinates[5];
16777
dof_coordinates[12][0] = 0.5*vertex_coordinates[0] + 0.25*vertex_coordinates[2] + 0.25*vertex_coordinates[4];
16778
dof_coordinates[12][1] = 0.5*vertex_coordinates[1] + 0.25*vertex_coordinates[3] + 0.25*vertex_coordinates[5];
16779
dof_coordinates[13][0] = 0.25*vertex_coordinates[0] + 0.5*vertex_coordinates[2] + 0.25*vertex_coordinates[4];
16780
dof_coordinates[13][1] = 0.25*vertex_coordinates[1] + 0.5*vertex_coordinates[3] + 0.25*vertex_coordinates[5];
16781
dof_coordinates[14][0] = 0.25*vertex_coordinates[0] + 0.25*vertex_coordinates[2] + 0.5*vertex_coordinates[4];
16782
dof_coordinates[14][1] = 0.25*vertex_coordinates[1] + 0.25*vertex_coordinates[3] + 0.5*vertex_coordinates[5];
16627
16785
/// Return the number of sub dofmaps (for a mixed element)
16628
virtual unsigned int num_sub_dofmaps() const
16786
virtual std::size_t num_sub_dofmaps() const
16633
16791
/// Create a new dofmap for sub dofmap i (for a mixed element)
16634
virtual ufc::dofmap* create_sub_dofmap(unsigned int i) const
16792
virtual ufc::dofmap* create_sub_dofmap(std::size_t i) const
16697
16852
return false;
16700
/// Initialize dofmap for mesh (return true iff init_cell() is needed)
16701
virtual bool init_mesh(const ufc::mesh& m)
16703
_global_dimension = m.num_entities[0] + 4*m.num_entities[1] + 8*m.num_entities[2];
16707
/// Initialize dofmap for given cell
16708
virtual void init_cell(const ufc::mesh& m,
16709
const ufc::cell& c)
16714
/// Finish initialization of dofmap for cells
16715
virtual void init_cell_finalize()
16720
16855
/// Return the topological dimension of the associated cell shape
16721
virtual unsigned int topological_dimension() const
16856
virtual std::size_t topological_dimension() const
16726
16861
/// Return the geometric dimension of the associated cell shape
16727
virtual unsigned int geometric_dimension() const
16862
virtual std::size_t geometric_dimension() const
16732
16867
/// Return the dimension of the global finite element function space
16733
virtual unsigned int global_dimension() const
16868
virtual std::size_t global_dimension(const std::vector<std::size_t>&
16869
num_global_entities) const
16735
return _global_dimension;
16871
return num_global_entities[0] + 4*num_global_entities[1] + 8*num_global_entities[2];
16738
16874
/// Return the dimension of the local finite element function space for a cell
16739
virtual unsigned int local_dimension(const ufc::cell& c) const
16875
virtual std::size_t local_dimension(const ufc::cell& c) const
16744
16880
/// Return the maximum dimension of the local finite element function space
16745
virtual unsigned int max_local_dimension() const
16881
virtual std::size_t max_local_dimension() const
16750
16886
/// Return the number of dofs on each cell facet
16751
virtual unsigned int num_facet_dofs() const
16887
virtual std::size_t num_facet_dofs() const
16756
16892
/// Return the number of dofs associated with each cell entity of dimension d
16757
virtual unsigned int num_entity_dofs(unsigned int d) const
16893
virtual std::size_t num_entity_dofs(std::size_t d) const
16955
17091
/// Tabulate the coordinates of all dofs on a cell
16956
virtual void tabulate_coordinates(double** coordinates,
16957
const ufc::cell& c) const
17092
virtual void tabulate_coordinates(double** dof_coordinates,
17093
const double* vertex_coordinates) const
16959
const double * const * x = c.coordinates;
16961
coordinates[0][0] = 0.33333333*x[0][0] + 0.33333333*x[1][0] + 0.33333333*x[2][0];
16962
coordinates[0][1] = 0.33333333*x[0][1] + 0.33333333*x[1][1] + 0.33333333*x[2][1];
16963
coordinates[1][0] = 0.33333333*x[0][0] + 0.33333333*x[1][0] + 0.33333333*x[2][0];
16964
coordinates[1][1] = 0.33333333*x[0][1] + 0.33333333*x[1][1] + 0.33333333*x[2][1];
16965
coordinates[2][0] = x[0][0];
16966
coordinates[2][1] = x[0][1];
16967
coordinates[3][0] = x[1][0];
16968
coordinates[3][1] = x[1][1];
16969
coordinates[4][0] = x[2][0];
16970
coordinates[4][1] = x[2][1];
16971
coordinates[5][0] = 0.5*x[1][0] + 0.5*x[2][0];
16972
coordinates[5][1] = 0.5*x[1][1] + 0.5*x[2][1];
16973
coordinates[6][0] = 0.5*x[0][0] + 0.5*x[2][0];
16974
coordinates[6][1] = 0.5*x[0][1] + 0.5*x[2][1];
16975
coordinates[7][0] = 0.5*x[0][0] + 0.5*x[1][0];
16976
coordinates[7][1] = 0.5*x[0][1] + 0.5*x[1][1];
16977
coordinates[8][0] = 0.75*x[1][0] + 0.25*x[2][0];
16978
coordinates[8][1] = 0.75*x[1][1] + 0.25*x[2][1];
16979
coordinates[9][0] = 0.5*x[1][0] + 0.5*x[2][0];
16980
coordinates[9][1] = 0.5*x[1][1] + 0.5*x[2][1];
16981
coordinates[10][0] = 0.25*x[1][0] + 0.75*x[2][0];
16982
coordinates[10][1] = 0.25*x[1][1] + 0.75*x[2][1];
16983
coordinates[11][0] = 0.75*x[0][0] + 0.25*x[2][0];
16984
coordinates[11][1] = 0.75*x[0][1] + 0.25*x[2][1];
16985
coordinates[12][0] = 0.5*x[0][0] + 0.5*x[2][0];
16986
coordinates[12][1] = 0.5*x[0][1] + 0.5*x[2][1];
16987
coordinates[13][0] = 0.25*x[0][0] + 0.75*x[2][0];
16988
coordinates[13][1] = 0.25*x[0][1] + 0.75*x[2][1];
16989
coordinates[14][0] = 0.75*x[0][0] + 0.25*x[1][0];
16990
coordinates[14][1] = 0.75*x[0][1] + 0.25*x[1][1];
16991
coordinates[15][0] = 0.5*x[0][0] + 0.5*x[1][0];
16992
coordinates[15][1] = 0.5*x[0][1] + 0.5*x[1][1];
16993
coordinates[16][0] = 0.25*x[0][0] + 0.75*x[1][0];
16994
coordinates[16][1] = 0.25*x[0][1] + 0.75*x[1][1];
16995
coordinates[17][0] = 0.5*x[0][0] + 0.25*x[1][0] + 0.25*x[2][0];
16996
coordinates[17][1] = 0.5*x[0][1] + 0.25*x[1][1] + 0.25*x[2][1];
16997
coordinates[18][0] = 0.25*x[0][0] + 0.5*x[1][0] + 0.25*x[2][0];
16998
coordinates[18][1] = 0.25*x[0][1] + 0.5*x[1][1] + 0.25*x[2][1];
16999
coordinates[19][0] = 0.25*x[0][0] + 0.25*x[1][0] + 0.5*x[2][0];
17000
coordinates[19][1] = 0.25*x[0][1] + 0.25*x[1][1] + 0.5*x[2][1];
17001
coordinates[20][0] = 0.5*x[0][0] + 0.25*x[1][0] + 0.25*x[2][0];
17002
coordinates[20][1] = 0.5*x[0][1] + 0.25*x[1][1] + 0.25*x[2][1];
17003
coordinates[21][0] = 0.25*x[0][0] + 0.5*x[1][0] + 0.25*x[2][0];
17004
coordinates[21][1] = 0.25*x[0][1] + 0.5*x[1][1] + 0.25*x[2][1];
17005
coordinates[22][0] = 0.25*x[0][0] + 0.25*x[1][0] + 0.5*x[2][0];
17006
coordinates[22][1] = 0.25*x[0][1] + 0.25*x[1][1] + 0.5*x[2][1];
17095
dof_coordinates[0][0] = 0.33333333*vertex_coordinates[0] + 0.33333333*vertex_coordinates[2] + 0.33333333*vertex_coordinates[4];
17096
dof_coordinates[0][1] = 0.33333333*vertex_coordinates[1] + 0.33333333*vertex_coordinates[3] + 0.33333333*vertex_coordinates[5];
17097
dof_coordinates[1][0] = 0.33333333*vertex_coordinates[0] + 0.33333333*vertex_coordinates[2] + 0.33333333*vertex_coordinates[4];
17098
dof_coordinates[1][1] = 0.33333333*vertex_coordinates[1] + 0.33333333*vertex_coordinates[3] + 0.33333333*vertex_coordinates[5];
17099
dof_coordinates[2][0] = vertex_coordinates[0];
17100
dof_coordinates[2][1] = vertex_coordinates[1];
17101
dof_coordinates[3][0] = vertex_coordinates[2];
17102
dof_coordinates[3][1] = vertex_coordinates[3];
17103
dof_coordinates[4][0] = vertex_coordinates[4];
17104
dof_coordinates[4][1] = vertex_coordinates[5];
17105
dof_coordinates[5][0] = 0.5*vertex_coordinates[2] + 0.5*vertex_coordinates[4];
17106
dof_coordinates[5][1] = 0.5*vertex_coordinates[3] + 0.5*vertex_coordinates[5];
17107
dof_coordinates[6][0] = 0.5*vertex_coordinates[0] + 0.5*vertex_coordinates[4];
17108
dof_coordinates[6][1] = 0.5*vertex_coordinates[1] + 0.5*vertex_coordinates[5];
17109
dof_coordinates[7][0] = 0.5*vertex_coordinates[0] + 0.5*vertex_coordinates[2];
17110
dof_coordinates[7][1] = 0.5*vertex_coordinates[1] + 0.5*vertex_coordinates[3];
17111
dof_coordinates[8][0] = 0.75*vertex_coordinates[2] + 0.25*vertex_coordinates[4];
17112
dof_coordinates[8][1] = 0.75*vertex_coordinates[3] + 0.25*vertex_coordinates[5];
17113
dof_coordinates[9][0] = 0.5*vertex_coordinates[2] + 0.5*vertex_coordinates[4];
17114
dof_coordinates[9][1] = 0.5*vertex_coordinates[3] + 0.5*vertex_coordinates[5];
17115
dof_coordinates[10][0] = 0.25*vertex_coordinates[2] + 0.75*vertex_coordinates[4];
17116
dof_coordinates[10][1] = 0.25*vertex_coordinates[3] + 0.75*vertex_coordinates[5];
17117
dof_coordinates[11][0] = 0.75*vertex_coordinates[0] + 0.25*vertex_coordinates[4];
17118
dof_coordinates[11][1] = 0.75*vertex_coordinates[1] + 0.25*vertex_coordinates[5];
17119
dof_coordinates[12][0] = 0.5*vertex_coordinates[0] + 0.5*vertex_coordinates[4];
17120
dof_coordinates[12][1] = 0.5*vertex_coordinates[1] + 0.5*vertex_coordinates[5];
17121
dof_coordinates[13][0] = 0.25*vertex_coordinates[0] + 0.75*vertex_coordinates[4];
17122
dof_coordinates[13][1] = 0.25*vertex_coordinates[1] + 0.75*vertex_coordinates[5];
17123
dof_coordinates[14][0] = 0.75*vertex_coordinates[0] + 0.25*vertex_coordinates[2];
17124
dof_coordinates[14][1] = 0.75*vertex_coordinates[1] + 0.25*vertex_coordinates[3];
17125
dof_coordinates[15][0] = 0.5*vertex_coordinates[0] + 0.5*vertex_coordinates[2];
17126
dof_coordinates[15][1] = 0.5*vertex_coordinates[1] + 0.5*vertex_coordinates[3];
17127
dof_coordinates[16][0] = 0.25*vertex_coordinates[0] + 0.75*vertex_coordinates[2];
17128
dof_coordinates[16][1] = 0.25*vertex_coordinates[1] + 0.75*vertex_coordinates[3];
17129
dof_coordinates[17][0] = 0.5*vertex_coordinates[0] + 0.25*vertex_coordinates[2] + 0.25*vertex_coordinates[4];
17130
dof_coordinates[17][1] = 0.5*vertex_coordinates[1] + 0.25*vertex_coordinates[3] + 0.25*vertex_coordinates[5];
17131
dof_coordinates[18][0] = 0.25*vertex_coordinates[0] + 0.5*vertex_coordinates[2] + 0.25*vertex_coordinates[4];
17132
dof_coordinates[18][1] = 0.25*vertex_coordinates[1] + 0.5*vertex_coordinates[3] + 0.25*vertex_coordinates[5];
17133
dof_coordinates[19][0] = 0.25*vertex_coordinates[0] + 0.25*vertex_coordinates[2] + 0.5*vertex_coordinates[4];
17134
dof_coordinates[19][1] = 0.25*vertex_coordinates[1] + 0.25*vertex_coordinates[3] + 0.5*vertex_coordinates[5];
17135
dof_coordinates[20][0] = 0.5*vertex_coordinates[0] + 0.25*vertex_coordinates[2] + 0.25*vertex_coordinates[4];
17136
dof_coordinates[20][1] = 0.5*vertex_coordinates[1] + 0.25*vertex_coordinates[3] + 0.25*vertex_coordinates[5];
17137
dof_coordinates[21][0] = 0.25*vertex_coordinates[0] + 0.5*vertex_coordinates[2] + 0.25*vertex_coordinates[4];
17138
dof_coordinates[21][1] = 0.25*vertex_coordinates[1] + 0.5*vertex_coordinates[3] + 0.25*vertex_coordinates[5];
17139
dof_coordinates[22][0] = 0.25*vertex_coordinates[0] + 0.25*vertex_coordinates[2] + 0.5*vertex_coordinates[4];
17140
dof_coordinates[22][1] = 0.25*vertex_coordinates[1] + 0.25*vertex_coordinates[3] + 0.5*vertex_coordinates[5];
17009
17143
/// Return the number of sub dofmaps (for a mixed element)
17010
virtual unsigned int num_sub_dofmaps() const
17144
virtual std::size_t num_sub_dofmaps() const
17015
17149
/// Create a new dofmap for sub dofmap i (for a mixed element)
17016
virtual ufc::dofmap* create_sub_dofmap(unsigned int i) const
17150
virtual ufc::dofmap* create_sub_dofmap(std::size_t i) const
17173
17291
/// Return a string identifying the form
17174
17292
virtual const char* signature() const
17176
return "Form([Integral(Indexed(Argument(MixedElement(*[MixedElement(*[VectorElement('Discontinuous Lagrange', Cell('triangle', Space(2)), 0, 2, None), FiniteElement('Lagrange', Cell('triangle', Space(2)), 2, None)], **{'value_shape': (3,) }), FiniteElement('Raviart-Thomas', Cell('triangle', Space(2)), 3, None)], **{'value_shape': (5,) }), 0), MultiIndex((FixedIndex(0),), {})), Measure('cell', 0, None))])";
17294
return "4dbd6ffce91104cca6203350dd5899d91933e7a0ff03b5e3d1e43c46be855704ea0938c37d2dc1542a93ea929629248d53f538a8c8fd87383de2bb660568b10b";
17179
17297
/// Return the rank of the global tensor (r)
17180
virtual unsigned int rank() const
17298
virtual std::size_t rank() const
17185
17303
/// Return the number of coefficients (n)
17186
virtual unsigned int num_coefficients() const
17304
virtual std::size_t num_coefficients() const
17191
17309
/// Return the number of cell domains
17192
virtual unsigned int num_cell_domains() const
17310
virtual std::size_t num_cell_domains() const
17197
17315
/// Return the number of exterior facet domains
17198
virtual unsigned int num_exterior_facet_domains() const
17316
virtual std::size_t num_exterior_facet_domains() const
17203
17321
/// Return the number of interior facet domains
17204
virtual unsigned int num_interior_facet_domains() const
17322
virtual std::size_t num_interior_facet_domains() const
17327
/// Return the number of point domains
17328
virtual std::size_t num_point_domains() const
17333
/// Return whether the form has any cell integrals
17334
virtual bool has_cell_integrals() const
17339
/// Return whether the form has any exterior facet integrals
17340
virtual bool has_exterior_facet_integrals() const
17345
/// Return whether the form has any interior facet integrals
17346
virtual bool has_interior_facet_integrals() const
17351
/// Return whether the form has any point integrals
17352
virtual bool has_point_integrals() const
17209
17357
/// Create a new finite element for argument function i
17210
virtual ufc::finite_element* create_finite_element(unsigned int i) const
17358
virtual ufc::finite_element* create_finite_element(std::size_t i) const