51 #ifdef HAVE_STOKHOS_FORUQTK
58 template <
typename OrdinalType,
typename ValueType>
71 int order =
setup.basis.order();
76 int size =
setup.basis.size();
83 for (
int i=0; i<=
setup.p; i++)
84 n2[i] = 1.0/(2.0*i+1.0);
92 for (
int i=0; i<=
setup.p; i++) {
93 n1[i] =
setup.basis.norm_squared(i);
94 n2[i] = 1.0/(2.0*i+1.0);
107 for (
int i=0; i<=
setup.p; i++) {
109 for (
int j=0;
j<nqp;
j++)
110 n2[i] += w[
j]*v[
j][i]*v[
j][i];
124 for (
int i=0; i<=
setup.p; i++) {
126 for (
int k=0; k<nqp; k++)
127 mat(i,
j) += w[k]*v[k][i]*v[k][
j];
134 out <<
"\n Error, mat.normInf() < atol = " << mat.
normInf()
135 <<
" < " <<
setup.atol <<
": failed!\n";
136 out <<
"mat = " << mat << std::endl;
142 setup.basis.computeTripleProductTensor();
149 for (
int i=0; i<=
setup.p; i++) {
151 for (
int k=0; k<=
setup.p; k++) {
154 for (
int qp=0; qp<nqp; qp++)
155 c += w[qp]*v[qp][i]*v[qp][
j]*v[qp][k];
156 std::stringstream ss;
157 ss <<
"Cijk(" << i <<
"," <<
j <<
"," << k <<
")";
168 setup.basis.computeDerivDoubleProductTensor();
176 for (
int i=0; i<nqp; i++) {
179 setup.basis.evaluateBasesAndDerivatives(
x[i],
val[i], deriv[i]);
183 for (
int i=0; i<=
setup.p; i++) {
186 for (
int qp=0; qp<nqp; qp++)
187 b += w[qp]*deriv[qp][i]*
val[qp][
j];
188 std::stringstream ss;
189 ss <<
"Bij(" << i <<
"," <<
j <<
")";
199 setup.basis.computeDerivDoubleProductTensor();
206 deriv_coeffs[1][0] = 1.0;
210 deriv_coeffs[2][1] = 3.0;
212 for (
int k=3; k<=
setup.p; k++) {
214 deriv_coeffs[k][0] = 1.0/3.0*deriv_coeffs[k-1][1];
215 for (
int i=1; i<=k-3; i++)
216 deriv_coeffs[k][i] = i/(2.0*i - 1.0)*deriv_coeffs[k-1][i-1] +
217 (i+1.0)/(2.0*i + 3.0)*deriv_coeffs[k-1][i+1];
218 deriv_coeffs[k][k-2] = (k-2.0)/(2.0*k-5.0)*deriv_coeffs[k-1][k-3];
219 deriv_coeffs[k][k-1] = (k-1.0)/(2.0*k-3.0)*deriv_coeffs[k-1][k-2] + k;
223 for (
int i=0; i<=
setup.p; i++) {
225 double b = deriv_coeffs[i][
j]*
n[
j];
226 std::stringstream ss;
227 ss <<
"Bij(" << i <<
"," <<
j <<
")";
238 setup.basis.evaluateBases(
x, v1);
247 for (
int i=2; i<=
setup.p; i++)
248 v2[i] = (2.0*i-1.0)/i*
x*v2[i-1] - (i-1.0)/i*v2[i-2];
257 setup.basis.evaluateBasesAndDerivatives(
x, val1, deriv1);
269 for (
int i=2; i<=
setup.p; i++)
270 val2[i] = (2.0*i-1.0)/i*
x*val2[i-1] - (i-1.0)/i*val2[i-2];
275 for (
int i=2; i<=
setup.p; i++)
276 deriv2[i] = i*val2[i-1] +
x*deriv2[i-1];
289 for (
int i=0; i<=
setup.p; i++)
290 v1[i] =
setup.basis.evaluate(
x, i);
299 for (
int i=2; i<=
setup.p; i++)
300 v2[i] = (2.0*i-1.0)/i*
x*v2[i-1] - (i-1.0)/i*v2[i-2];
310 setup.basis.getRecurrenceCoefficients(a1, b1, c1, d1);
313 a2[0] = 0.0; b2[0] = 1.0; c2[0] = 1.0; d2[0] = 1.0;
314 for (
int i=1; i<=
setup.p; i++) {
317 c2[i] = (2.0*i+1.0)/(i+1);
336 setup.basis.getRecurrenceCoefficients(a1, b1, c1, d1);
344 for (
int i=0; i<=
setup.p; i++) {
346 for (
int j=0;
j<nqp;
j++)
347 a2[i] += w[
j]*
x[
j]*v[
j][i]*v[
j][i];
348 a2[i] = a2[i]*c1[i]/n1[i];
360 setup.basis.getRecurrenceCoefficients(a1, b1, c1, d1);
365 for (
int i=1; i<=
setup.p; i++) {
366 b2[i] = d1[i]*(c1[i]/c1[i-1])*(n1[i]/n1[i-1]);
372 #ifdef HAVE_STOKHOS_DAKOTA
374 int n = static_cast<int>(std::ceil((
setup.p+1)/2.0));
377 setup.basis.getQuadPoints(
setup.p, x1, w1, v1);
381 webbur::legendre_compute(n, &x2[0], &w2[0]);
383 for (
int i=0; i<
n; i++) {
385 v2[i].resize(
setup.p+1);
386 setup.basis.evaluateBases(x2[i], v2[i]);
393 for (
int i=0; i<
n; i++) {
394 std::stringstream ss1, ss2;
395 ss1 <<
"v1[" << i <<
"]";
396 ss2 <<
"v2[" << i <<
"]";
404 #ifdef HAVE_STOKHOS_FORUQTK
406 int n = static_cast<int>(std::ceil((
setup.p+1)/2.0));
409 setup.basis.getQuadPoints(
setup.p, x1, w1, v1);
415 double endpts[2] = {0.0, 0.0};
419 GAUSSQ_F77(&kind, &n, &alpha, &beta, &kpts, endpts, &b[0], &x2[0], &w2[0]);
421 for (
int i=0; i<
n; i++) {
423 v2[i].resize(
setup.p+1);
424 setup.basis.evaluateBases(x2[i], v2[i]);
431 for (
int i=0; i<
n; i++) {
432 std::stringstream ss1, ss2;
433 ss1 <<
"v1[" << i <<
"]";
434 ss2 <<
"v2[" << i <<
"]";