~ubuntu-branches/ubuntu/raring/ffc/raring

« back to all changes in this revision

Viewing changes to test/regression/references/r_quadrature/Conditional.h

  • Committer: Package Import Robot
  • Author(s): Johannes Ring
  • Date: 2011-10-26 17:52:20 UTC
  • mfrom: (1.1.10)
  • Revision ID: package-import@ubuntu.com-20111026175220-ope1dzqv4jn2b8pq
Tags: 1.0-beta2-1
* New upstream release. This release includes some performance
  improvements for evaluating basis functions. It also adds support
  for Bessel functions and error functions.
* debian/control: Bump version numbers for python-ufc, python-fiat,
  python-instant, and python-ufl.

Show diffs side-by-side

added added

removed removed

Lines of Context:
1
 
// This code conforms with the UFC specification version 2.0.0
2
 
// and was automatically generated by FFC version 0.9.10.
 
1
// This code conforms with the UFC specification version 2.0.2
 
2
// and was automatically generated by FFC version 1.0-beta+.
3
3
// 
4
4
// This code was generated with the following parameters:
5
5
// 
130
130
      double basisvalues[6] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
131
131
      
132
132
      // Declare helper variables.
133
 
      unsigned int rr = 0;
134
 
      unsigned int ss = 0;
135
 
      unsigned int tt = 0;
136
 
      double tmp5 = 0.0;
137
 
      double tmp6 = 0.0;
138
 
      double tmp7 = 0.0;
139
133
      double tmp0 = (1.0 + Y + 2.0*X)/2.0;
140
134
      double tmp1 = (1.0 - Y)/2.0;
141
135
      double tmp2 = tmp1*tmp1;
143
137
      // Compute basisvalues.
144
138
      basisvalues[0] = 1.0;
145
139
      basisvalues[1] = tmp0;
146
 
      for (unsigned int r = 1; r < 2; r++)
147
 
      {
148
 
        rr = (r + 1)*((r + 1) + 1)/2;
149
 
        ss = r*(r + 1)/2;
150
 
        tt = (r - 1)*((r - 1) + 1)/2;
151
 
        tmp5 = (1.0 + 2.0*r)/(1.0 + r);
152
 
        basisvalues[rr] = (basisvalues[ss]*tmp0*tmp5 - basisvalues[tt]*tmp2*r/(1.0 + r));
153
 
      }// end loop over 'r'
154
 
      for (unsigned int r = 0; r < 2; r++)
155
 
      {
156
 
        rr = (r + 1)*(r + 1 + 1)/2 + 1;
157
 
        ss = r*(r + 1)/2;
158
 
        basisvalues[rr] = basisvalues[ss]*(0.5 + r + Y*(1.5 + r));
159
 
      }// end loop over 'r'
160
 
      for (unsigned int r = 0; r < 1; r++)
161
 
      {
162
 
        for (unsigned int s = 1; s < 2 - r; s++)
163
 
        {
164
 
          rr = (r + s + 1)*(r + s + 1 + 1)/2 + s + 1;
165
 
          ss = (r + s)*(r + s + 1)/2 + s;
166
 
          tt = (r + s - 1)*(r + s - 1 + 1)/2 + s - 1;
167
 
          tmp5 = (2.0 + 2.0*r + 2.0*s)*(3.0 + 2.0*r + 2.0*s)/(2.0*(1.0 + s)*(2.0 + s + 2.0*r));
168
 
          tmp6 = (1.0 + 4.0*r*r + 4.0*r)*(2.0 + 2.0*r + 2.0*s)/(2.0*(1.0 + 2.0*r + 2.0*s)*(1.0 + s)*(2.0 + s + 2.0*r));
169
 
          tmp7 = (1.0 + s + 2.0*r)*(3.0 + 2.0*r + 2.0*s)*s/((1.0 + 2.0*r + 2.0*s)*(1.0 + s)*(2.0 + s + 2.0*r));
170
 
          basisvalues[rr] = (basisvalues[ss]*(tmp6 + Y*tmp5) - basisvalues[tt]*tmp7);
171
 
        }// end loop over 's'
172
 
      }// end loop over 'r'
173
 
      for (unsigned int r = 0; r < 3; r++)
174
 
      {
175
 
        for (unsigned int s = 0; s < 3 - r; s++)
176
 
        {
177
 
          rr = (r + s)*(r + s + 1)/2 + s;
178
 
          basisvalues[rr] *= std::sqrt((0.5 + r)*(1.0 + r + s));
179
 
        }// end loop over 's'
180
 
      }// end loop over 'r'
 
140
      basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
 
141
      basisvalues[2] = basisvalues[0]*(0.5 + 1.5*Y);
 
142
      basisvalues[4] = basisvalues[1]*(1.5 + 2.5*Y);
 
143
      basisvalues[5] = basisvalues[2]*(0.11111111 + Y*1.6666667) - basisvalues[0]*0.55555556;
 
144
      basisvalues[0] *= std::sqrt(0.5);
 
145
      basisvalues[2] *= std::sqrt(1.0);
 
146
      basisvalues[5] *= std::sqrt(1.5);
 
147
      basisvalues[1] *= std::sqrt(3.0);
 
148
      basisvalues[4] *= std::sqrt(4.5);
 
149
      basisvalues[3] *= std::sqrt(7.5);
181
150
      
182
151
      // Table(s) of coefficients.
183
152
      static const double coefficients0[6] = \
197
166
      double basisvalues[6] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
198
167
      
199
168
      // Declare helper variables.
200
 
      unsigned int rr = 0;
201
 
      unsigned int ss = 0;
202
 
      unsigned int tt = 0;
203
 
      double tmp5 = 0.0;
204
 
      double tmp6 = 0.0;
205
 
      double tmp7 = 0.0;
206
169
      double tmp0 = (1.0 + Y + 2.0*X)/2.0;
207
170
      double tmp1 = (1.0 - Y)/2.0;
208
171
      double tmp2 = tmp1*tmp1;
210
173
      // Compute basisvalues.
211
174
      basisvalues[0] = 1.0;
212
175
      basisvalues[1] = tmp0;
213
 
      for (unsigned int r = 1; r < 2; r++)
214
 
      {
215
 
        rr = (r + 1)*((r + 1) + 1)/2;
216
 
        ss = r*(r + 1)/2;
217
 
        tt = (r - 1)*((r - 1) + 1)/2;
218
 
        tmp5 = (1.0 + 2.0*r)/(1.0 + r);
219
 
        basisvalues[rr] = (basisvalues[ss]*tmp0*tmp5 - basisvalues[tt]*tmp2*r/(1.0 + r));
220
 
      }// end loop over 'r'
221
 
      for (unsigned int r = 0; r < 2; r++)
222
 
      {
223
 
        rr = (r + 1)*(r + 1 + 1)/2 + 1;
224
 
        ss = r*(r + 1)/2;
225
 
        basisvalues[rr] = basisvalues[ss]*(0.5 + r + Y*(1.5 + r));
226
 
      }// end loop over 'r'
227
 
      for (unsigned int r = 0; r < 1; r++)
228
 
      {
229
 
        for (unsigned int s = 1; s < 2 - r; s++)
230
 
        {
231
 
          rr = (r + s + 1)*(r + s + 1 + 1)/2 + s + 1;
232
 
          ss = (r + s)*(r + s + 1)/2 + s;
233
 
          tt = (r + s - 1)*(r + s - 1 + 1)/2 + s - 1;
234
 
          tmp5 = (2.0 + 2.0*r + 2.0*s)*(3.0 + 2.0*r + 2.0*s)/(2.0*(1.0 + s)*(2.0 + s + 2.0*r));
235
 
          tmp6 = (1.0 + 4.0*r*r + 4.0*r)*(2.0 + 2.0*r + 2.0*s)/(2.0*(1.0 + 2.0*r + 2.0*s)*(1.0 + s)*(2.0 + s + 2.0*r));
236
 
          tmp7 = (1.0 + s + 2.0*r)*(3.0 + 2.0*r + 2.0*s)*s/((1.0 + 2.0*r + 2.0*s)*(1.0 + s)*(2.0 + s + 2.0*r));
237
 
          basisvalues[rr] = (basisvalues[ss]*(tmp6 + Y*tmp5) - basisvalues[tt]*tmp7);
238
 
        }// end loop over 's'
239
 
      }// end loop over 'r'
240
 
      for (unsigned int r = 0; r < 3; r++)
241
 
      {
242
 
        for (unsigned int s = 0; s < 3 - r; s++)
243
 
        {
244
 
          rr = (r + s)*(r + s + 1)/2 + s;
245
 
          basisvalues[rr] *= std::sqrt((0.5 + r)*(1.0 + r + s));
246
 
        }// end loop over 's'
247
 
      }// end loop over 'r'
 
176
      basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
 
177
      basisvalues[2] = basisvalues[0]*(0.5 + 1.5*Y);
 
178
      basisvalues[4] = basisvalues[1]*(1.5 + 2.5*Y);
 
179
      basisvalues[5] = basisvalues[2]*(0.11111111 + Y*1.6666667) - basisvalues[0]*0.55555556;
 
180
      basisvalues[0] *= std::sqrt(0.5);
 
181
      basisvalues[2] *= std::sqrt(1.0);
 
182
      basisvalues[5] *= std::sqrt(1.5);
 
183
      basisvalues[1] *= std::sqrt(3.0);
 
184
      basisvalues[4] *= std::sqrt(4.5);
 
185
      basisvalues[3] *= std::sqrt(7.5);
248
186
      
249
187
      // Table(s) of coefficients.
250
188
      static const double coefficients0[6] = \
264
202
      double basisvalues[6] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
265
203
      
266
204
      // Declare helper variables.
267
 
      unsigned int rr = 0;
268
 
      unsigned int ss = 0;
269
 
      unsigned int tt = 0;
270
 
      double tmp5 = 0.0;
271
 
      double tmp6 = 0.0;
272
 
      double tmp7 = 0.0;
273
205
      double tmp0 = (1.0 + Y + 2.0*X)/2.0;
274
206
      double tmp1 = (1.0 - Y)/2.0;
275
207
      double tmp2 = tmp1*tmp1;
277
209
      // Compute basisvalues.
278
210
      basisvalues[0] = 1.0;
279
211
      basisvalues[1] = tmp0;
280
 
      for (unsigned int r = 1; r < 2; r++)
281
 
      {
282
 
        rr = (r + 1)*((r + 1) + 1)/2;
283
 
        ss = r*(r + 1)/2;
284
 
        tt = (r - 1)*((r - 1) + 1)/2;
285
 
        tmp5 = (1.0 + 2.0*r)/(1.0 + r);
286
 
        basisvalues[rr] = (basisvalues[ss]*tmp0*tmp5 - basisvalues[tt]*tmp2*r/(1.0 + r));
287
 
      }// end loop over 'r'
288
 
      for (unsigned int r = 0; r < 2; r++)
289
 
      {
290
 
        rr = (r + 1)*(r + 1 + 1)/2 + 1;
291
 
        ss = r*(r + 1)/2;
292
 
        basisvalues[rr] = basisvalues[ss]*(0.5 + r + Y*(1.5 + r));
293
 
      }// end loop over 'r'
294
 
      for (unsigned int r = 0; r < 1; r++)
295
 
      {
296
 
        for (unsigned int s = 1; s < 2 - r; s++)
297
 
        {
298
 
          rr = (r + s + 1)*(r + s + 1 + 1)/2 + s + 1;
299
 
          ss = (r + s)*(r + s + 1)/2 + s;
300
 
          tt = (r + s - 1)*(r + s - 1 + 1)/2 + s - 1;
301
 
          tmp5 = (2.0 + 2.0*r + 2.0*s)*(3.0 + 2.0*r + 2.0*s)/(2.0*(1.0 + s)*(2.0 + s + 2.0*r));
302
 
          tmp6 = (1.0 + 4.0*r*r + 4.0*r)*(2.0 + 2.0*r + 2.0*s)/(2.0*(1.0 + 2.0*r + 2.0*s)*(1.0 + s)*(2.0 + s + 2.0*r));
303
 
          tmp7 = (1.0 + s + 2.0*r)*(3.0 + 2.0*r + 2.0*s)*s/((1.0 + 2.0*r + 2.0*s)*(1.0 + s)*(2.0 + s + 2.0*r));
304
 
          basisvalues[rr] = (basisvalues[ss]*(tmp6 + Y*tmp5) - basisvalues[tt]*tmp7);
305
 
        }// end loop over 's'
306
 
      }// end loop over 'r'
307
 
      for (unsigned int r = 0; r < 3; r++)
308
 
      {
309
 
        for (unsigned int s = 0; s < 3 - r; s++)
310
 
        {
311
 
          rr = (r + s)*(r + s + 1)/2 + s;
312
 
          basisvalues[rr] *= std::sqrt((0.5 + r)*(1.0 + r + s));
313
 
        }// end loop over 's'
314
 
      }// end loop over 'r'
 
212
      basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
 
213
      basisvalues[2] = basisvalues[0]*(0.5 + 1.5*Y);
 
214
      basisvalues[4] = basisvalues[1]*(1.5 + 2.5*Y);
 
215
      basisvalues[5] = basisvalues[2]*(0.11111111 + Y*1.6666667) - basisvalues[0]*0.55555556;
 
216
      basisvalues[0] *= std::sqrt(0.5);
 
217
      basisvalues[2] *= std::sqrt(1.0);
 
218
      basisvalues[5] *= std::sqrt(1.5);
 
219
      basisvalues[1] *= std::sqrt(3.0);
 
220
      basisvalues[4] *= std::sqrt(4.5);
 
221
      basisvalues[3] *= std::sqrt(7.5);
315
222
      
316
223
      // Table(s) of coefficients.
317
224
      static const double coefficients0[6] = \
331
238
      double basisvalues[6] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
332
239
      
333
240
      // Declare helper variables.
334
 
      unsigned int rr = 0;
335
 
      unsigned int ss = 0;
336
 
      unsigned int tt = 0;
337
 
      double tmp5 = 0.0;
338
 
      double tmp6 = 0.0;
339
 
      double tmp7 = 0.0;
340
241
      double tmp0 = (1.0 + Y + 2.0*X)/2.0;
341
242
      double tmp1 = (1.0 - Y)/2.0;
342
243
      double tmp2 = tmp1*tmp1;
344
245
      // Compute basisvalues.
345
246
      basisvalues[0] = 1.0;
346
247
      basisvalues[1] = tmp0;
347
 
      for (unsigned int r = 1; r < 2; r++)
348
 
      {
349
 
        rr = (r + 1)*((r + 1) + 1)/2;
350
 
        ss = r*(r + 1)/2;
351
 
        tt = (r - 1)*((r - 1) + 1)/2;
352
 
        tmp5 = (1.0 + 2.0*r)/(1.0 + r);
353
 
        basisvalues[rr] = (basisvalues[ss]*tmp0*tmp5 - basisvalues[tt]*tmp2*r/(1.0 + r));
354
 
      }// end loop over 'r'
355
 
      for (unsigned int r = 0; r < 2; r++)
356
 
      {
357
 
        rr = (r + 1)*(r + 1 + 1)/2 + 1;
358
 
        ss = r*(r + 1)/2;
359
 
        basisvalues[rr] = basisvalues[ss]*(0.5 + r + Y*(1.5 + r));
360
 
      }// end loop over 'r'
361
 
      for (unsigned int r = 0; r < 1; r++)
362
 
      {
363
 
        for (unsigned int s = 1; s < 2 - r; s++)
364
 
        {
365
 
          rr = (r + s + 1)*(r + s + 1 + 1)/2 + s + 1;
366
 
          ss = (r + s)*(r + s + 1)/2 + s;
367
 
          tt = (r + s - 1)*(r + s - 1 + 1)/2 + s - 1;
368
 
          tmp5 = (2.0 + 2.0*r + 2.0*s)*(3.0 + 2.0*r + 2.0*s)/(2.0*(1.0 + s)*(2.0 + s + 2.0*r));
369
 
          tmp6 = (1.0 + 4.0*r*r + 4.0*r)*(2.0 + 2.0*r + 2.0*s)/(2.0*(1.0 + 2.0*r + 2.0*s)*(1.0 + s)*(2.0 + s + 2.0*r));
370
 
          tmp7 = (1.0 + s + 2.0*r)*(3.0 + 2.0*r + 2.0*s)*s/((1.0 + 2.0*r + 2.0*s)*(1.0 + s)*(2.0 + s + 2.0*r));
371
 
          basisvalues[rr] = (basisvalues[ss]*(tmp6 + Y*tmp5) - basisvalues[tt]*tmp7);
372
 
        }// end loop over 's'
373
 
      }// end loop over 'r'
374
 
      for (unsigned int r = 0; r < 3; r++)
375
 
      {
376
 
        for (unsigned int s = 0; s < 3 - r; s++)
377
 
        {
378
 
          rr = (r + s)*(r + s + 1)/2 + s;
379
 
          basisvalues[rr] *= std::sqrt((0.5 + r)*(1.0 + r + s));
380
 
        }// end loop over 's'
381
 
      }// end loop over 'r'
 
248
      basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
 
249
      basisvalues[2] = basisvalues[0]*(0.5 + 1.5*Y);
 
250
      basisvalues[4] = basisvalues[1]*(1.5 + 2.5*Y);
 
251
      basisvalues[5] = basisvalues[2]*(0.11111111 + Y*1.6666667) - basisvalues[0]*0.55555556;
 
252
      basisvalues[0] *= std::sqrt(0.5);
 
253
      basisvalues[2] *= std::sqrt(1.0);
 
254
      basisvalues[5] *= std::sqrt(1.5);
 
255
      basisvalues[1] *= std::sqrt(3.0);
 
256
      basisvalues[4] *= std::sqrt(4.5);
 
257
      basisvalues[3] *= std::sqrt(7.5);
382
258
      
383
259
      // Table(s) of coefficients.
384
260
      static const double coefficients0[6] = \
398
274
      double basisvalues[6] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
399
275
      
400
276
      // Declare helper variables.
401
 
      unsigned int rr = 0;
402
 
      unsigned int ss = 0;
403
 
      unsigned int tt = 0;
404
 
      double tmp5 = 0.0;
405
 
      double tmp6 = 0.0;
406
 
      double tmp7 = 0.0;
407
277
      double tmp0 = (1.0 + Y + 2.0*X)/2.0;
408
278
      double tmp1 = (1.0 - Y)/2.0;
409
279
      double tmp2 = tmp1*tmp1;
411
281
      // Compute basisvalues.
412
282
      basisvalues[0] = 1.0;
413
283
      basisvalues[1] = tmp0;
414
 
      for (unsigned int r = 1; r < 2; r++)
415
 
      {
416
 
        rr = (r + 1)*((r + 1) + 1)/2;
417
 
        ss = r*(r + 1)/2;
418
 
        tt = (r - 1)*((r - 1) + 1)/2;
419
 
        tmp5 = (1.0 + 2.0*r)/(1.0 + r);
420
 
        basisvalues[rr] = (basisvalues[ss]*tmp0*tmp5 - basisvalues[tt]*tmp2*r/(1.0 + r));
421
 
      }// end loop over 'r'
422
 
      for (unsigned int r = 0; r < 2; r++)
423
 
      {
424
 
        rr = (r + 1)*(r + 1 + 1)/2 + 1;
425
 
        ss = r*(r + 1)/2;
426
 
        basisvalues[rr] = basisvalues[ss]*(0.5 + r + Y*(1.5 + r));
427
 
      }// end loop over 'r'
428
 
      for (unsigned int r = 0; r < 1; r++)
429
 
      {
430
 
        for (unsigned int s = 1; s < 2 - r; s++)
431
 
        {
432
 
          rr = (r + s + 1)*(r + s + 1 + 1)/2 + s + 1;
433
 
          ss = (r + s)*(r + s + 1)/2 + s;
434
 
          tt = (r + s - 1)*(r + s - 1 + 1)/2 + s - 1;
435
 
          tmp5 = (2.0 + 2.0*r + 2.0*s)*(3.0 + 2.0*r + 2.0*s)/(2.0*(1.0 + s)*(2.0 + s + 2.0*r));
436
 
          tmp6 = (1.0 + 4.0*r*r + 4.0*r)*(2.0 + 2.0*r + 2.0*s)/(2.0*(1.0 + 2.0*r + 2.0*s)*(1.0 + s)*(2.0 + s + 2.0*r));
437
 
          tmp7 = (1.0 + s + 2.0*r)*(3.0 + 2.0*r + 2.0*s)*s/((1.0 + 2.0*r + 2.0*s)*(1.0 + s)*(2.0 + s + 2.0*r));
438
 
          basisvalues[rr] = (basisvalues[ss]*(tmp6 + Y*tmp5) - basisvalues[tt]*tmp7);
439
 
        }// end loop over 's'
440
 
      }// end loop over 'r'
441
 
      for (unsigned int r = 0; r < 3; r++)
442
 
      {
443
 
        for (unsigned int s = 0; s < 3 - r; s++)
444
 
        {
445
 
          rr = (r + s)*(r + s + 1)/2 + s;
446
 
          basisvalues[rr] *= std::sqrt((0.5 + r)*(1.0 + r + s));
447
 
        }// end loop over 's'
448
 
      }// end loop over 'r'
 
284
      basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
 
285
      basisvalues[2] = basisvalues[0]*(0.5 + 1.5*Y);
 
286
      basisvalues[4] = basisvalues[1]*(1.5 + 2.5*Y);
 
287
      basisvalues[5] = basisvalues[2]*(0.11111111 + Y*1.6666667) - basisvalues[0]*0.55555556;
 
288
      basisvalues[0] *= std::sqrt(0.5);
 
289
      basisvalues[2] *= std::sqrt(1.0);
 
290
      basisvalues[5] *= std::sqrt(1.5);
 
291
      basisvalues[1] *= std::sqrt(3.0);
 
292
      basisvalues[4] *= std::sqrt(4.5);
 
293
      basisvalues[3] *= std::sqrt(7.5);
449
294
      
450
295
      // Table(s) of coefficients.
451
296
      static const double coefficients0[6] = \
465
310
      double basisvalues[6] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
466
311
      
467
312
      // Declare helper variables.
468
 
      unsigned int rr = 0;
469
 
      unsigned int ss = 0;
470
 
      unsigned int tt = 0;
471
 
      double tmp5 = 0.0;
472
 
      double tmp6 = 0.0;
473
 
      double tmp7 = 0.0;
474
313
      double tmp0 = (1.0 + Y + 2.0*X)/2.0;
475
314
      double tmp1 = (1.0 - Y)/2.0;
476
315
      double tmp2 = tmp1*tmp1;
478
317
      // Compute basisvalues.
479
318
      basisvalues[0] = 1.0;
480
319
      basisvalues[1] = tmp0;
481
 
      for (unsigned int r = 1; r < 2; r++)
482
 
      {
483
 
        rr = (r + 1)*((r + 1) + 1)/2;
484
 
        ss = r*(r + 1)/2;
485
 
        tt = (r - 1)*((r - 1) + 1)/2;
486
 
        tmp5 = (1.0 + 2.0*r)/(1.0 + r);
487
 
        basisvalues[rr] = (basisvalues[ss]*tmp0*tmp5 - basisvalues[tt]*tmp2*r/(1.0 + r));
488
 
      }// end loop over 'r'
489
 
      for (unsigned int r = 0; r < 2; r++)
490
 
      {
491
 
        rr = (r + 1)*(r + 1 + 1)/2 + 1;
492
 
        ss = r*(r + 1)/2;
493
 
        basisvalues[rr] = basisvalues[ss]*(0.5 + r + Y*(1.5 + r));
494
 
      }// end loop over 'r'
495
 
      for (unsigned int r = 0; r < 1; r++)
496
 
      {
497
 
        for (unsigned int s = 1; s < 2 - r; s++)
498
 
        {
499
 
          rr = (r + s + 1)*(r + s + 1 + 1)/2 + s + 1;
500
 
          ss = (r + s)*(r + s + 1)/2 + s;
501
 
          tt = (r + s - 1)*(r + s - 1 + 1)/2 + s - 1;
502
 
          tmp5 = (2.0 + 2.0*r + 2.0*s)*(3.0 + 2.0*r + 2.0*s)/(2.0*(1.0 + s)*(2.0 + s + 2.0*r));
503
 
          tmp6 = (1.0 + 4.0*r*r + 4.0*r)*(2.0 + 2.0*r + 2.0*s)/(2.0*(1.0 + 2.0*r + 2.0*s)*(1.0 + s)*(2.0 + s + 2.0*r));
504
 
          tmp7 = (1.0 + s + 2.0*r)*(3.0 + 2.0*r + 2.0*s)*s/((1.0 + 2.0*r + 2.0*s)*(1.0 + s)*(2.0 + s + 2.0*r));
505
 
          basisvalues[rr] = (basisvalues[ss]*(tmp6 + Y*tmp5) - basisvalues[tt]*tmp7);
506
 
        }// end loop over 's'
507
 
      }// end loop over 'r'
508
 
      for (unsigned int r = 0; r < 3; r++)
509
 
      {
510
 
        for (unsigned int s = 0; s < 3 - r; s++)
511
 
        {
512
 
          rr = (r + s)*(r + s + 1)/2 + s;
513
 
          basisvalues[rr] *= std::sqrt((0.5 + r)*(1.0 + r + s));
514
 
        }// end loop over 's'
515
 
      }// end loop over 'r'
 
320
      basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
 
321
      basisvalues[2] = basisvalues[0]*(0.5 + 1.5*Y);
 
322
      basisvalues[4] = basisvalues[1]*(1.5 + 2.5*Y);
 
323
      basisvalues[5] = basisvalues[2]*(0.11111111 + Y*1.6666667) - basisvalues[0]*0.55555556;
 
324
      basisvalues[0] *= std::sqrt(0.5);
 
325
      basisvalues[2] *= std::sqrt(1.0);
 
326
      basisvalues[5] *= std::sqrt(1.5);
 
327
      basisvalues[1] *= std::sqrt(3.0);
 
328
      basisvalues[4] *= std::sqrt(4.5);
 
329
      basisvalues[3] *= std::sqrt(7.5);
516
330
      
517
331
      // Table(s) of coefficients.
518
332
      static const double coefficients0[6] = \
651
465
      double basisvalues[6] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
652
466
      
653
467
      // Declare helper variables.
654
 
      unsigned int rr = 0;
655
 
      unsigned int ss = 0;
656
 
      unsigned int tt = 0;
657
 
      double tmp5 = 0.0;
658
 
      double tmp6 = 0.0;
659
 
      double tmp7 = 0.0;
660
468
      double tmp0 = (1.0 + Y + 2.0*X)/2.0;
661
469
      double tmp1 = (1.0 - Y)/2.0;
662
470
      double tmp2 = tmp1*tmp1;
664
472
      // Compute basisvalues.
665
473
      basisvalues[0] = 1.0;
666
474
      basisvalues[1] = tmp0;
667
 
      for (unsigned int r = 1; r < 2; r++)
668
 
      {
669
 
        rr = (r + 1)*((r + 1) + 1)/2;
670
 
        ss = r*(r + 1)/2;
671
 
        tt = (r - 1)*((r - 1) + 1)/2;
672
 
        tmp5 = (1.0 + 2.0*r)/(1.0 + r);
673
 
        basisvalues[rr] = (basisvalues[ss]*tmp0*tmp5 - basisvalues[tt]*tmp2*r/(1.0 + r));
674
 
      }// end loop over 'r'
675
 
      for (unsigned int r = 0; r < 2; r++)
676
 
      {
677
 
        rr = (r + 1)*(r + 1 + 1)/2 + 1;
678
 
        ss = r*(r + 1)/2;
679
 
        basisvalues[rr] = basisvalues[ss]*(0.5 + r + Y*(1.5 + r));
680
 
      }// end loop over 'r'
681
 
      for (unsigned int r = 0; r < 1; r++)
682
 
      {
683
 
        for (unsigned int s = 1; s < 2 - r; s++)
684
 
        {
685
 
          rr = (r + s + 1)*(r + s + 1 + 1)/2 + s + 1;
686
 
          ss = (r + s)*(r + s + 1)/2 + s;
687
 
          tt = (r + s - 1)*(r + s - 1 + 1)/2 + s - 1;
688
 
          tmp5 = (2.0 + 2.0*r + 2.0*s)*(3.0 + 2.0*r + 2.0*s)/(2.0*(1.0 + s)*(2.0 + s + 2.0*r));
689
 
          tmp6 = (1.0 + 4.0*r*r + 4.0*r)*(2.0 + 2.0*r + 2.0*s)/(2.0*(1.0 + 2.0*r + 2.0*s)*(1.0 + s)*(2.0 + s + 2.0*r));
690
 
          tmp7 = (1.0 + s + 2.0*r)*(3.0 + 2.0*r + 2.0*s)*s/((1.0 + 2.0*r + 2.0*s)*(1.0 + s)*(2.0 + s + 2.0*r));
691
 
          basisvalues[rr] = (basisvalues[ss]*(tmp6 + Y*tmp5) - basisvalues[tt]*tmp7);
692
 
        }// end loop over 's'
693
 
      }// end loop over 'r'
694
 
      for (unsigned int r = 0; r < 3; r++)
695
 
      {
696
 
        for (unsigned int s = 0; s < 3 - r; s++)
697
 
        {
698
 
          rr = (r + s)*(r + s + 1)/2 + s;
699
 
          basisvalues[rr] *= std::sqrt((0.5 + r)*(1.0 + r + s));
700
 
        }// end loop over 's'
701
 
      }// end loop over 'r'
 
475
      basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
 
476
      basisvalues[2] = basisvalues[0]*(0.5 + 1.5*Y);
 
477
      basisvalues[4] = basisvalues[1]*(1.5 + 2.5*Y);
 
478
      basisvalues[5] = basisvalues[2]*(0.11111111 + Y*1.6666667) - basisvalues[0]*0.55555556;
 
479
      basisvalues[0] *= std::sqrt(0.5);
 
480
      basisvalues[2] *= std::sqrt(1.0);
 
481
      basisvalues[5] *= std::sqrt(1.5);
 
482
      basisvalues[1] *= std::sqrt(3.0);
 
483
      basisvalues[4] *= std::sqrt(4.5);
 
484
      basisvalues[3] *= std::sqrt(7.5);
702
485
      
703
486
      // Table(s) of coefficients.
704
487
      static const double coefficients0[6] = \
848
631
      double basisvalues[6] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
849
632
      
850
633
      // Declare helper variables.
851
 
      unsigned int rr = 0;
852
 
      unsigned int ss = 0;
853
 
      unsigned int tt = 0;
854
 
      double tmp5 = 0.0;
855
 
      double tmp6 = 0.0;
856
 
      double tmp7 = 0.0;
857
634
      double tmp0 = (1.0 + Y + 2.0*X)/2.0;
858
635
      double tmp1 = (1.0 - Y)/2.0;
859
636
      double tmp2 = tmp1*tmp1;
861
638
      // Compute basisvalues.
862
639
      basisvalues[0] = 1.0;
863
640
      basisvalues[1] = tmp0;
864
 
      for (unsigned int r = 1; r < 2; r++)
865
 
      {
866
 
        rr = (r + 1)*((r + 1) + 1)/2;
867
 
        ss = r*(r + 1)/2;
868
 
        tt = (r - 1)*((r - 1) + 1)/2;
869
 
        tmp5 = (1.0 + 2.0*r)/(1.0 + r);
870
 
        basisvalues[rr] = (basisvalues[ss]*tmp0*tmp5 - basisvalues[tt]*tmp2*r/(1.0 + r));
871
 
      }// end loop over 'r'
872
 
      for (unsigned int r = 0; r < 2; r++)
873
 
      {
874
 
        rr = (r + 1)*(r + 1 + 1)/2 + 1;
875
 
        ss = r*(r + 1)/2;
876
 
        basisvalues[rr] = basisvalues[ss]*(0.5 + r + Y*(1.5 + r));
877
 
      }// end loop over 'r'
878
 
      for (unsigned int r = 0; r < 1; r++)
879
 
      {
880
 
        for (unsigned int s = 1; s < 2 - r; s++)
881
 
        {
882
 
          rr = (r + s + 1)*(r + s + 1 + 1)/2 + s + 1;
883
 
          ss = (r + s)*(r + s + 1)/2 + s;
884
 
          tt = (r + s - 1)*(r + s - 1 + 1)/2 + s - 1;
885
 
          tmp5 = (2.0 + 2.0*r + 2.0*s)*(3.0 + 2.0*r + 2.0*s)/(2.0*(1.0 + s)*(2.0 + s + 2.0*r));
886
 
          tmp6 = (1.0 + 4.0*r*r + 4.0*r)*(2.0 + 2.0*r + 2.0*s)/(2.0*(1.0 + 2.0*r + 2.0*s)*(1.0 + s)*(2.0 + s + 2.0*r));
887
 
          tmp7 = (1.0 + s + 2.0*r)*(3.0 + 2.0*r + 2.0*s)*s/((1.0 + 2.0*r + 2.0*s)*(1.0 + s)*(2.0 + s + 2.0*r));
888
 
          basisvalues[rr] = (basisvalues[ss]*(tmp6 + Y*tmp5) - basisvalues[tt]*tmp7);
889
 
        }// end loop over 's'
890
 
      }// end loop over 'r'
891
 
      for (unsigned int r = 0; r < 3; r++)
892
 
      {
893
 
        for (unsigned int s = 0; s < 3 - r; s++)
894
 
        {
895
 
          rr = (r + s)*(r + s + 1)/2 + s;
896
 
          basisvalues[rr] *= std::sqrt((0.5 + r)*(1.0 + r + s));
897
 
        }// end loop over 's'
898
 
      }// end loop over 'r'
 
641
      basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
 
642
      basisvalues[2] = basisvalues[0]*(0.5 + 1.5*Y);
 
643
      basisvalues[4] = basisvalues[1]*(1.5 + 2.5*Y);
 
644
      basisvalues[5] = basisvalues[2]*(0.11111111 + Y*1.6666667) - basisvalues[0]*0.55555556;
 
645
      basisvalues[0] *= std::sqrt(0.5);
 
646
      basisvalues[2] *= std::sqrt(1.0);
 
647
      basisvalues[5] *= std::sqrt(1.5);
 
648
      basisvalues[1] *= std::sqrt(3.0);
 
649
      basisvalues[4] *= std::sqrt(4.5);
 
650
      basisvalues[3] *= std::sqrt(7.5);
899
651
      
900
652
      // Table(s) of coefficients.
901
653
      static const double coefficients0[6] = \
1045
797
      double basisvalues[6] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
1046
798
      
1047
799
      // Declare helper variables.
1048
 
      unsigned int rr = 0;
1049
 
      unsigned int ss = 0;
1050
 
      unsigned int tt = 0;
1051
 
      double tmp5 = 0.0;
1052
 
      double tmp6 = 0.0;
1053
 
      double tmp7 = 0.0;
1054
800
      double tmp0 = (1.0 + Y + 2.0*X)/2.0;
1055
801
      double tmp1 = (1.0 - Y)/2.0;
1056
802
      double tmp2 = tmp1*tmp1;
1058
804
      // Compute basisvalues.
1059
805
      basisvalues[0] = 1.0;
1060
806
      basisvalues[1] = tmp0;
1061
 
      for (unsigned int r = 1; r < 2; r++)
1062
 
      {
1063
 
        rr = (r + 1)*((r + 1) + 1)/2;
1064
 
        ss = r*(r + 1)/2;
1065
 
        tt = (r - 1)*((r - 1) + 1)/2;
1066
 
        tmp5 = (1.0 + 2.0*r)/(1.0 + r);
1067
 
        basisvalues[rr] = (basisvalues[ss]*tmp0*tmp5 - basisvalues[tt]*tmp2*r/(1.0 + r));
1068
 
      }// end loop over 'r'
1069
 
      for (unsigned int r = 0; r < 2; r++)
1070
 
      {
1071
 
        rr = (r + 1)*(r + 1 + 1)/2 + 1;
1072
 
        ss = r*(r + 1)/2;
1073
 
        basisvalues[rr] = basisvalues[ss]*(0.5 + r + Y*(1.5 + r));
1074
 
      }// end loop over 'r'
1075
 
      for (unsigned int r = 0; r < 1; r++)
1076
 
      {
1077
 
        for (unsigned int s = 1; s < 2 - r; s++)
1078
 
        {
1079
 
          rr = (r + s + 1)*(r + s + 1 + 1)/2 + s + 1;
1080
 
          ss = (r + s)*(r + s + 1)/2 + s;
1081
 
          tt = (r + s - 1)*(r + s - 1 + 1)/2 + s - 1;
1082
 
          tmp5 = (2.0 + 2.0*r + 2.0*s)*(3.0 + 2.0*r + 2.0*s)/(2.0*(1.0 + s)*(2.0 + s + 2.0*r));
1083
 
          tmp6 = (1.0 + 4.0*r*r + 4.0*r)*(2.0 + 2.0*r + 2.0*s)/(2.0*(1.0 + 2.0*r + 2.0*s)*(1.0 + s)*(2.0 + s + 2.0*r));
1084
 
          tmp7 = (1.0 + s + 2.0*r)*(3.0 + 2.0*r + 2.0*s)*s/((1.0 + 2.0*r + 2.0*s)*(1.0 + s)*(2.0 + s + 2.0*r));
1085
 
          basisvalues[rr] = (basisvalues[ss]*(tmp6 + Y*tmp5) - basisvalues[tt]*tmp7);
1086
 
        }// end loop over 's'
1087
 
      }// end loop over 'r'
1088
 
      for (unsigned int r = 0; r < 3; r++)
1089
 
      {
1090
 
        for (unsigned int s = 0; s < 3 - r; s++)
1091
 
        {
1092
 
          rr = (r + s)*(r + s + 1)/2 + s;
1093
 
          basisvalues[rr] *= std::sqrt((0.5 + r)*(1.0 + r + s));
1094
 
        }// end loop over 's'
1095
 
      }// end loop over 'r'
 
807
      basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
 
808
      basisvalues[2] = basisvalues[0]*(0.5 + 1.5*Y);
 
809
      basisvalues[4] = basisvalues[1]*(1.5 + 2.5*Y);
 
810
      basisvalues[5] = basisvalues[2]*(0.11111111 + Y*1.6666667) - basisvalues[0]*0.55555556;
 
811
      basisvalues[0] *= std::sqrt(0.5);
 
812
      basisvalues[2] *= std::sqrt(1.0);
 
813
      basisvalues[5] *= std::sqrt(1.5);
 
814
      basisvalues[1] *= std::sqrt(3.0);
 
815
      basisvalues[4] *= std::sqrt(4.5);
 
816
      basisvalues[3] *= std::sqrt(7.5);
1096
817
      
1097
818
      // Table(s) of coefficients.
1098
819
      static const double coefficients0[6] = \
1242
963
      double basisvalues[6] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
1243
964
      
1244
965
      // Declare helper variables.
1245
 
      unsigned int rr = 0;
1246
 
      unsigned int ss = 0;
1247
 
      unsigned int tt = 0;
1248
 
      double tmp5 = 0.0;
1249
 
      double tmp6 = 0.0;
1250
 
      double tmp7 = 0.0;
1251
966
      double tmp0 = (1.0 + Y + 2.0*X)/2.0;
1252
967
      double tmp1 = (1.0 - Y)/2.0;
1253
968
      double tmp2 = tmp1*tmp1;
1255
970
      // Compute basisvalues.
1256
971
      basisvalues[0] = 1.0;
1257
972
      basisvalues[1] = tmp0;
1258
 
      for (unsigned int r = 1; r < 2; r++)
1259
 
      {
1260
 
        rr = (r + 1)*((r + 1) + 1)/2;
1261
 
        ss = r*(r + 1)/2;
1262
 
        tt = (r - 1)*((r - 1) + 1)/2;
1263
 
        tmp5 = (1.0 + 2.0*r)/(1.0 + r);
1264
 
        basisvalues[rr] = (basisvalues[ss]*tmp0*tmp5 - basisvalues[tt]*tmp2*r/(1.0 + r));
1265
 
      }// end loop over 'r'
1266
 
      for (unsigned int r = 0; r < 2; r++)
1267
 
      {
1268
 
        rr = (r + 1)*(r + 1 + 1)/2 + 1;
1269
 
        ss = r*(r + 1)/2;
1270
 
        basisvalues[rr] = basisvalues[ss]*(0.5 + r + Y*(1.5 + r));
1271
 
      }// end loop over 'r'
1272
 
      for (unsigned int r = 0; r < 1; r++)
1273
 
      {
1274
 
        for (unsigned int s = 1; s < 2 - r; s++)
1275
 
        {
1276
 
          rr = (r + s + 1)*(r + s + 1 + 1)/2 + s + 1;
1277
 
          ss = (r + s)*(r + s + 1)/2 + s;
1278
 
          tt = (r + s - 1)*(r + s - 1 + 1)/2 + s - 1;
1279
 
          tmp5 = (2.0 + 2.0*r + 2.0*s)*(3.0 + 2.0*r + 2.0*s)/(2.0*(1.0 + s)*(2.0 + s + 2.0*r));
1280
 
          tmp6 = (1.0 + 4.0*r*r + 4.0*r)*(2.0 + 2.0*r + 2.0*s)/(2.0*(1.0 + 2.0*r + 2.0*s)*(1.0 + s)*(2.0 + s + 2.0*r));
1281
 
          tmp7 = (1.0 + s + 2.0*r)*(3.0 + 2.0*r + 2.0*s)*s/((1.0 + 2.0*r + 2.0*s)*(1.0 + s)*(2.0 + s + 2.0*r));
1282
 
          basisvalues[rr] = (basisvalues[ss]*(tmp6 + Y*tmp5) - basisvalues[tt]*tmp7);
1283
 
        }// end loop over 's'
1284
 
      }// end loop over 'r'
1285
 
      for (unsigned int r = 0; r < 3; r++)
1286
 
      {
1287
 
        for (unsigned int s = 0; s < 3 - r; s++)
1288
 
        {
1289
 
          rr = (r + s)*(r + s + 1)/2 + s;
1290
 
          basisvalues[rr] *= std::sqrt((0.5 + r)*(1.0 + r + s));
1291
 
        }// end loop over 's'
1292
 
      }// end loop over 'r'
 
973
      basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
 
974
      basisvalues[2] = basisvalues[0]*(0.5 + 1.5*Y);
 
975
      basisvalues[4] = basisvalues[1]*(1.5 + 2.5*Y);
 
976
      basisvalues[5] = basisvalues[2]*(0.11111111 + Y*1.6666667) - basisvalues[0]*0.55555556;
 
977
      basisvalues[0] *= std::sqrt(0.5);
 
978
      basisvalues[2] *= std::sqrt(1.0);
 
979
      basisvalues[5] *= std::sqrt(1.5);
 
980
      basisvalues[1] *= std::sqrt(3.0);
 
981
      basisvalues[4] *= std::sqrt(4.5);
 
982
      basisvalues[3] *= std::sqrt(7.5);
1293
983
      
1294
984
      // Table(s) of coefficients.
1295
985
      static const double coefficients0[6] = \
1439
1129
      double basisvalues[6] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
1440
1130
      
1441
1131
      // Declare helper variables.
1442
 
      unsigned int rr = 0;
1443
 
      unsigned int ss = 0;
1444
 
      unsigned int tt = 0;
1445
 
      double tmp5 = 0.0;
1446
 
      double tmp6 = 0.0;
1447
 
      double tmp7 = 0.0;
1448
1132
      double tmp0 = (1.0 + Y + 2.0*X)/2.0;
1449
1133
      double tmp1 = (1.0 - Y)/2.0;
1450
1134
      double tmp2 = tmp1*tmp1;
1452
1136
      // Compute basisvalues.
1453
1137
      basisvalues[0] = 1.0;
1454
1138
      basisvalues[1] = tmp0;
1455
 
      for (unsigned int r = 1; r < 2; r++)
1456
 
      {
1457
 
        rr = (r + 1)*((r + 1) + 1)/2;
1458
 
        ss = r*(r + 1)/2;
1459
 
        tt = (r - 1)*((r - 1) + 1)/2;
1460
 
        tmp5 = (1.0 + 2.0*r)/(1.0 + r);
1461
 
        basisvalues[rr] = (basisvalues[ss]*tmp0*tmp5 - basisvalues[tt]*tmp2*r/(1.0 + r));
1462
 
      }// end loop over 'r'
1463
 
      for (unsigned int r = 0; r < 2; r++)
1464
 
      {
1465
 
        rr = (r + 1)*(r + 1 + 1)/2 + 1;
1466
 
        ss = r*(r + 1)/2;
1467
 
        basisvalues[rr] = basisvalues[ss]*(0.5 + r + Y*(1.5 + r));
1468
 
      }// end loop over 'r'
1469
 
      for (unsigned int r = 0; r < 1; r++)
1470
 
      {
1471
 
        for (unsigned int s = 1; s < 2 - r; s++)
1472
 
        {
1473
 
          rr = (r + s + 1)*(r + s + 1 + 1)/2 + s + 1;
1474
 
          ss = (r + s)*(r + s + 1)/2 + s;
1475
 
          tt = (r + s - 1)*(r + s - 1 + 1)/2 + s - 1;
1476
 
          tmp5 = (2.0 + 2.0*r + 2.0*s)*(3.0 + 2.0*r + 2.0*s)/(2.0*(1.0 + s)*(2.0 + s + 2.0*r));
1477
 
          tmp6 = (1.0 + 4.0*r*r + 4.0*r)*(2.0 + 2.0*r + 2.0*s)/(2.0*(1.0 + 2.0*r + 2.0*s)*(1.0 + s)*(2.0 + s + 2.0*r));
1478
 
          tmp7 = (1.0 + s + 2.0*r)*(3.0 + 2.0*r + 2.0*s)*s/((1.0 + 2.0*r + 2.0*s)*(1.0 + s)*(2.0 + s + 2.0*r));
1479
 
          basisvalues[rr] = (basisvalues[ss]*(tmp6 + Y*tmp5) - basisvalues[tt]*tmp7);
1480
 
        }// end loop over 's'
1481
 
      }// end loop over 'r'
1482
 
      for (unsigned int r = 0; r < 3; r++)
1483
 
      {
1484
 
        for (unsigned int s = 0; s < 3 - r; s++)
1485
 
        {
1486
 
          rr = (r + s)*(r + s + 1)/2 + s;
1487
 
          basisvalues[rr] *= std::sqrt((0.5 + r)*(1.0 + r + s));
1488
 
        }// end loop over 's'
1489
 
      }// end loop over 'r'
 
1139
      basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
 
1140
      basisvalues[2] = basisvalues[0]*(0.5 + 1.5*Y);
 
1141
      basisvalues[4] = basisvalues[1]*(1.5 + 2.5*Y);
 
1142
      basisvalues[5] = basisvalues[2]*(0.11111111 + Y*1.6666667) - basisvalues[0]*0.55555556;
 
1143
      basisvalues[0] *= std::sqrt(0.5);
 
1144
      basisvalues[2] *= std::sqrt(1.0);
 
1145
      basisvalues[5] *= std::sqrt(1.5);
 
1146
      basisvalues[1] *= std::sqrt(3.0);
 
1147
      basisvalues[4] *= std::sqrt(4.5);
 
1148
      basisvalues[3] *= std::sqrt(7.5);
1490
1149
      
1491
1150
      // Table(s) of coefficients.
1492
1151
      static const double coefficients0[6] = \
1636
1295
      double basisvalues[6] = {0.0, 0.0, 0.0, 0.0, 0.0, 0.0};
1637
1296
      
1638
1297
      // Declare helper variables.
1639
 
      unsigned int rr = 0;
1640
 
      unsigned int ss = 0;
1641
 
      unsigned int tt = 0;
1642
 
      double tmp5 = 0.0;
1643
 
      double tmp6 = 0.0;
1644
 
      double tmp7 = 0.0;
1645
1298
      double tmp0 = (1.0 + Y + 2.0*X)/2.0;
1646
1299
      double tmp1 = (1.0 - Y)/2.0;
1647
1300
      double tmp2 = tmp1*tmp1;
1649
1302
      // Compute basisvalues.
1650
1303
      basisvalues[0] = 1.0;
1651
1304
      basisvalues[1] = tmp0;
1652
 
      for (unsigned int r = 1; r < 2; r++)
1653
 
      {
1654
 
        rr = (r + 1)*((r + 1) + 1)/2;
1655
 
        ss = r*(r + 1)/2;
1656
 
        tt = (r - 1)*((r - 1) + 1)/2;
1657
 
        tmp5 = (1.0 + 2.0*r)/(1.0 + r);
1658
 
        basisvalues[rr] = (basisvalues[ss]*tmp0*tmp5 - basisvalues[tt]*tmp2*r/(1.0 + r));
1659
 
      }// end loop over 'r'
1660
 
      for (unsigned int r = 0; r < 2; r++)
1661
 
      {
1662
 
        rr = (r + 1)*(r + 1 + 1)/2 + 1;
1663
 
        ss = r*(r + 1)/2;
1664
 
        basisvalues[rr] = basisvalues[ss]*(0.5 + r + Y*(1.5 + r));
1665
 
      }// end loop over 'r'
1666
 
      for (unsigned int r = 0; r < 1; r++)
1667
 
      {
1668
 
        for (unsigned int s = 1; s < 2 - r; s++)
1669
 
        {
1670
 
          rr = (r + s + 1)*(r + s + 1 + 1)/2 + s + 1;
1671
 
          ss = (r + s)*(r + s + 1)/2 + s;
1672
 
          tt = (r + s - 1)*(r + s - 1 + 1)/2 + s - 1;
1673
 
          tmp5 = (2.0 + 2.0*r + 2.0*s)*(3.0 + 2.0*r + 2.0*s)/(2.0*(1.0 + s)*(2.0 + s + 2.0*r));
1674
 
          tmp6 = (1.0 + 4.0*r*r + 4.0*r)*(2.0 + 2.0*r + 2.0*s)/(2.0*(1.0 + 2.0*r + 2.0*s)*(1.0 + s)*(2.0 + s + 2.0*r));
1675
 
          tmp7 = (1.0 + s + 2.0*r)*(3.0 + 2.0*r + 2.0*s)*s/((1.0 + 2.0*r + 2.0*s)*(1.0 + s)*(2.0 + s + 2.0*r));
1676
 
          basisvalues[rr] = (basisvalues[ss]*(tmp6 + Y*tmp5) - basisvalues[tt]*tmp7);
1677
 
        }// end loop over 's'
1678
 
      }// end loop over 'r'
1679
 
      for (unsigned int r = 0; r < 3; r++)
1680
 
      {
1681
 
        for (unsigned int s = 0; s < 3 - r; s++)
1682
 
        {
1683
 
          rr = (r + s)*(r + s + 1)/2 + s;
1684
 
          basisvalues[rr] *= std::sqrt((0.5 + r)*(1.0 + r + s));
1685
 
        }// end loop over 's'
1686
 
      }// end loop over 'r'
 
1305
      basisvalues[3] = basisvalues[1]*1.5*tmp0 - basisvalues[0]*0.5*tmp2;
 
1306
      basisvalues[2] = basisvalues[0]*(0.5 + 1.5*Y);
 
1307
      basisvalues[4] = basisvalues[1]*(1.5 + 2.5*Y);
 
1308
      basisvalues[5] = basisvalues[2]*(0.11111111 + Y*1.6666667) - basisvalues[0]*0.55555556;
 
1309
      basisvalues[0] *= std::sqrt(0.5);
 
1310
      basisvalues[2] *= std::sqrt(1.0);
 
1311
      basisvalues[5] *= std::sqrt(1.5);
 
1312
      basisvalues[1] *= std::sqrt(3.0);
 
1313
      basisvalues[4] *= std::sqrt(4.5);
 
1314
      basisvalues[3] *= std::sqrt(7.5);
1687
1315
      
1688
1316
      // Table(s) of coefficients.
1689
1317
      static const double coefficients0[6] = \