Stokhos Package Browser (Single Doxygen Collection)  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Stokhos_ForUQTKOrthogPolyExpansionImp.hpp
Go to the documentation of this file.
1 // @HEADER
2 // *****************************************************************************
3 // Stokhos Package
4 //
5 // Copyright 2009 NTESS and the Stokhos contributors.
6 // SPDX-License-Identifier: BSD-3-Clause
7 // *****************************************************************************
8 // @HEADER
9 
10 #include "Teuchos_Assert.hpp"
11 #include "Stokhos_DynamicArrayTraits.hpp"
12 #include "Teuchos_ConfigDefs.hpp"
13 
14 #define UQ_PREP_F77 F77_FUNC_(uq_prep,UQ_PREP)
15 #define UQ_PROD2_F77 F77_FUNC_(uq_prod2,UQ_PROD2)
16 #define UQ_DIV_F77 F77_FUNC_(uq_div,UQ_DIV)
17 #define UQ_EXP_F77 F77_FUNC_(uq_exp,UQ_EXP)
18 #define UQ_LOG_F77 F77_FUNC_(uq_log,UQ_LOG)
19 #define UQ_SQRT_F77 F77_FUNC_(uq_sqrt,UQ_SQRT)
20 #define UQ_EXP_INT_F77 F77_FUNC_(uq_exp_int,UQ_EXP)
21 #define UQ_LOG_INT_F77 F77_FUNC_(uq_log_int,UQ_LOG)
22 
23 extern "C" {
24  void UQ_PREP_F77(int*, int*, int*);
25  void UQ_PROD2_F77(const double*, const double*, double*, int*);
26  void UQ_DIV_F77(const double*, const double*, double*, int*);
27  void UQ_EXP_F77(const double*, double*, int*, int*, double*, int*);
28  void UQ_LOG_F77(const double*, double*, int*, int*, double*, int*);
29  void UQ_SQRT_F77(const double*, double*, int*, int*);
30  void UQ_EXP_INT_F77(const double*, double*, int*);
31  void UQ_LOG_INT_F77(const double*, double*, int*);
32 }
33 
34 template <typename ordinal_type, typename value_type>
35 Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
36 ForUQTKOrthogPolyExpansion(
37  const Teuchos::RCP<const Stokhos::OrthogPolyBasis<ordinal_type, value_type> >& basis_,
38  const Teuchos::RCP<const Stokhos::Sparse3Tensor<ordinal_type, value_type> >& Cijk_,
39  EXPANSION_METHOD method_,
40  value_type rtol_) :
41  OrthogPolyExpansionBase<ordinal_type, value_type, node_type>(basis_, Cijk_),
42  rtol(rtol_),
43  method(method_)
44 {
45  order = this->basis->order();
46  dim = this->basis->dimension();
47  int nup;
48  UQ_PREP_F77(&order, &dim, &nup);
49  sz = nup+1;
50 }
51 
52 template <typename ordinal_type, typename value_type>
53 void
54 Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
55 timesEqual(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
56  const value_type& val)
57 {
58  OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::timesEqual(c,val);
59 }
60 
61 template <typename ordinal_type, typename value_type>
62 void
63 Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
64 divideEqual(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
65  const value_type& val)
66 {
67  OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::divideEqual(c,val);
68 }
69 
70 template <typename ordinal_type, typename value_type>
71 void
72 Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
73 timesEqual(
74  Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
75  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& x)
76 {
77  ordinal_type p = c.size();
78  ordinal_type xp = x.size();
79  ordinal_type pc;
80  if (p > 1 && xp > 1)
81  pc = sz;
82  else
83  pc = p*xp;
84  TEUCHOS_TEST_FOR_EXCEPTION(sz < pc, std::logic_error,
85  "Stokhos::ForUQTKOrthogPolyExpansion::timesEqual()" <<
86  ": Expansion size (" << sz <<
87  ") is too small for computation.");
88  if (c.size() != pc)
89  c.resize(pc);
90 
91  value_type* cc = c.coeff();
92  const value_type* xc = x.coeff();
93 
94  if (p > 1 && xp > 1) {
95  TEUCHOS_TEST_FOR_EXCEPTION(pc != xp, std::logic_error,
96  "Stokhos::ForUQTKOrthogPolyExpansion::timesEqual()"
97  << ": Arguments have incompatible sizes: "
98  << "x.size() = " << xp << ", c.size() = " << pc << ".");
99 
100  // Copy c coefficients into temporary array
101  value_type* tc = Stokhos::ds_array<value_type>::get_and_fill(cc,pc);
102 
103  int nup = pc-1;
104  UQ_PROD2_F77(cc, xc, tc, &nup);
105  Stokhos::ds_array<value_type>::copy(tc, cc, p);
106  }
107  else if (p > 1) {
108  for (ordinal_type i=0; i<p; i++)
109  cc[i] *= xc[0];
110  }
111  else if (xp > 1) {
112  for (ordinal_type i=1; i<xp; i++)
113  cc[i] = cc[0]*xc[i];
114  cc[0] *= xc[0];
115  }
116  else {
117  cc[0] *= xc[0];
118  }
119 }
120 
121 template <typename ordinal_type, typename value_type>
122 void
123 Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
124 divideEqual(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
125  const OrthogPolyApprox<ordinal_type, value_type, node_type>& x)
126 {
127  ordinal_type p = c.size();
128  ordinal_type xp = x.size();
129  ordinal_type pc;
130  if (xp > 1)
131  pc = sz;
132  else
133  pc = p;
134  TEUCHOS_TEST_FOR_EXCEPTION(sz < pc, std::logic_error,
135  "Stokhos::ForUQTKOrthogPolyExpansion::divideEqual()" <<
136  ": Expansion size (" << sz <<
137  ") is too small for computation.");
138  if (c.size() != pc)
139  c.resize(pc);
140 
141  value_type* cc = c.coeff();
142  const value_type* xc = x.coeff();
143 
144  if (xp > 1) {
145  TEUCHOS_TEST_FOR_EXCEPTION(pc != xp, std::logic_error,
146  "Stokhos::ForUQTKOrthogPolyExpansion::divideEqual()"
147  << ": Arguments have incompatible sizes: "
148  << "x.size() = " << xp << ", c.size() = " << pc << ".");
149 
150  // Copy c coefficients into temporary array
151  value_type* tc = Stokhos::ds_array<value_type>::get_and_fill(cc,pc);
152 
153  int nup = pc-1;
154  UQ_DIV_F77(tc, xc, cc, &nup);
155 
156  }
157  else {
158  for (ordinal_type i=0; i<pc; i++)
159  cc[i] /= xc[0];
160  }
161 }
162 
163 template <typename ordinal_type, typename value_type>
164 void
165 Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
166 times(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
167  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a,
168  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& b)
169 {
170  ordinal_type pa = a.size();
171  ordinal_type pb = b.size();
172  ordinal_type pc;
173  if (pa > 1 && pb > 1)
174  pc = sz;
175  else
176  pc = pa*pb;
177  TEUCHOS_TEST_FOR_EXCEPTION(sz < pc, std::logic_error,
178  "Stokhos::ForUQTKOrthogPolyExpansion::times()" <<
179  ": Expansion size (" << sz <<
180  ") is too small for computation.");
181  if (c.size() != pc)
182  c.resize(pc);
183 
184  const value_type* ca = a.coeff();
185  const value_type* cb = b.coeff();
186  value_type* cc = c.coeff();
187 
188  if (pa > 1 && pb > 1) {
189  TEUCHOS_TEST_FOR_EXCEPTION(pa != pc || pb != pc, std::logic_error,
190  "Stokhos::ForUQTKOrthogPolyExpansion::times()"
191  << ": Arguments have incompatible sizes: "
192  << "a.size() = " << pa << ", b.size() = " << pb
193  << ", required size = " << pc << ".");
194 
195  int nup = pc-1;
196  UQ_PROD2_F77(ca, cb, cc, &nup);
197  }
198  else if (pa > 1) {
199  for (ordinal_type i=0; i<pc; i++)
200  cc[i] = ca[i]*cb[0];
201  }
202  else if (pb > 1) {
203  for (ordinal_type i=0; i<pc; i++)
204  cc[i] = ca[0]*cb[i];
205  }
206  else {
207  cc[0] = ca[0]*cb[0];
208  }
209 }
210 
211 template <typename ordinal_type, typename value_type>
212 void
213 Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
214 times(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
215  const value_type& a,
216  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& b)
217 {
218  OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::times(c,a,b);
219 }
220 
221 template <typename ordinal_type, typename value_type>
222 void
223 Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
224 times(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
225  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a,
226  const value_type& b)
227 {
228  OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::times(c,a,b);
229 }
230 
231 template <typename ordinal_type, typename value_type>
232 void
233 Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
234 divide(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
235  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a,
236  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& b)
237 {
238  ordinal_type pa = a.size();
239  ordinal_type pb = b.size();
240  ordinal_type pc;
241  if (pb > 1)
242  pc = sz;
243  else
244  pc = pa;
245  TEUCHOS_TEST_FOR_EXCEPTION(sz < pc, std::logic_error,
246  "Stokhos::ForUQTKOrthogPolyExpansion::divide()" <<
247  ": Expansion size (" << sz <<
248  ") is too small for computation.");
249  if (c.size() != pc)
250  c.resize(pc);
251 
252  const value_type* ca = a.coeff();
253  const value_type* cb = b.coeff();
254  value_type* cc = c.coeff();
255 
256  if (pb > 1) {
257  TEUCHOS_TEST_FOR_EXCEPTION(pa != pc || pb != pc, std::logic_error,
258  "Stokhos::ForUQTKOrthogPolyExpansion::divide()"
259  << ": Arguments have incompatible sizes: "
260  << "a.size() = " << pa << ", b.size() = " << pb
261  << ", required size = " << pc << ".");
262 
263  int nup = pc-1;
264  UQ_DIV_F77(ca, cb, cc, &nup);
265  }
266  else {
267  for (ordinal_type i=0; i<pa; i++)
268  cc[i] = ca[i]/cb[0];
269  }
270 }
271 
272 template <typename ordinal_type, typename value_type>
273 void
274 Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
275 divide(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
276  const value_type& a,
277  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& b)
278 {
279  ordinal_type pb = b.size();
280  ordinal_type pc;
281  if (pb > 1)
282  pc = sz;
283  else
284  pc = 1;
285  if (c.size() != pc)
286  c.resize(pc);
287 
288  const value_type* cb = b.coeff();
289  value_type* cc = c.coeff();
290 
291  if (pb > 1) {
292  TEUCHOS_TEST_FOR_EXCEPTION(pb != pc, std::logic_error,
293  "Stokhos::ForUQTKOrthogPolyExpansion::divide()"
294  << ": Arguments have incompatible sizes: "
295  << "b.size() = " << pb
296  << ", required size = " << pc << ".");
297 
298  value_type* ca = Stokhos::ds_array<value_type>::get_and_fill(pc);
299  ca[0] = a;
300  int nup = pc-1;
301  UQ_DIV_F77(ca, cb, cc, &nup);
302  }
303  else
304  cc[0] = a / cb[0];
305 }
306 
307 template <typename ordinal_type, typename value_type>
308 void
309 Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
310 divide(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
311  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a,
312  const value_type& b)
313 {
314  OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::divide(c,a,b);
315 }
316 
317 template <typename ordinal_type, typename value_type>
318 void
319 Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
320 exp(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
321  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a)
322 {
323  ordinal_type pa = a.size();
324  ordinal_type pc;
325  if (pa > 1)
326  pc = sz;
327  else
328  pc = 1;
329  if (c.size() != pc)
330  c.resize(pc);
331 
332  const value_type* ca = a.coeff();
333  value_type* cc = c.coeff();
334 
335  if (pa > 1) {
336  TEUCHOS_TEST_FOR_EXCEPTION(pa != pc, std::logic_error,
337  "Stokhos::ForUQTKOrthogPolyExpansion::exp()"
338  << ": Arguments have incompatible sizes: "
339  << "a.size() = " << pa << ", c.size() = " << pc
340  << ".");
341  int nup = pc-1;
342  if (method == TAYLOR) {
343  int nrm = 1;
344  UQ_EXP_F77(ca, cc, &dim, &nup, &rtol, &nrm);
345  }
346  else
347  UQ_EXP_INT_F77(ca, cc, &nup);
348  }
349  else
350  cc[0] = std::exp(ca[0]);
351 }
352 
353 template <typename ordinal_type, typename value_type>
354 void
355 Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
356 log(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
357  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a)
358 {
359  ordinal_type pa = a.size();
360  ordinal_type pc;
361  if (pa > 1)
362  pc = sz;
363  else
364  pc = 1;
365  if (c.size() != pc)
366  c.resize(pc);
367 
368  const value_type* ca = a.coeff();
369  value_type* cc = c.coeff();
370 
371  if (pa > 1) {
372  TEUCHOS_TEST_FOR_EXCEPTION(pa != pc, std::logic_error,
373  "Stokhos::ForUQTKOrthogPolyExpansion::log()"
374  << ": Arguments have incompatible sizes: "
375  << "a.size() = " << pa << ", c.size() = " << pc
376  << ".");
377  int nup = pc-1;
378  if (method == TAYLOR) {
379  int nrm = 1;
380  UQ_LOG_F77(ca, cc, &dim, &nup, &rtol, &nrm);
381  }
382  else
383  UQ_LOG_INT_F77(ca, cc, &nup);
384  }
385  else
386  cc[0] = std::log(ca[0]);
387 }
388 
389 template <typename ordinal_type, typename value_type>
390 void
391 Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
392 log10(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
393  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a)
394 {
395  if (a.size() > 1) {
396  log(c,a);
397  divide(c,c,std::log(10.0));
398  }
399  else {
400  if (c.size() != 1)
401  c.resize(1);
402  c[0] = std::log10(a[0]);
403  }
404 }
405 
406 template <typename ordinal_type, typename value_type>
407 void
408 Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
409 sqrt(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
410  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a)
411 {
412  ordinal_type pa = a.size();
413  ordinal_type pc;
414  if (pa > 1)
415  pc = sz;
416  else
417  pc = 1;
418  if (c.size() != pc)
419  c.resize(pc);
420 
421  const value_type* ca = a.coeff();
422  value_type* cc = c.coeff();
423 
424  if (pa > 1) {
425  TEUCHOS_TEST_FOR_EXCEPTION(pa != pc, std::logic_error,
426  "Stokhos::ForUQTKOrthogPolyExpansion::sqrt()"
427  << ": Arguments have incompatible sizes: "
428  << "a.size() = " << pa << ", c.size() = " << pc
429  << ".");
430  int iguess = 0;
431  int nup = pc-1;
432  UQ_SQRT_F77(ca, cc, &nup, &iguess);
433  }
434  else
435  cc[0] = std::sqrt(ca[0]);
436 }
437 
438 template <typename ordinal_type, typename value_type>
439 void
440 Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
441 cbrt(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
442  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a)
443 {
444  if (a.size() > 1) {
445  log(c,a);
446  timesEqual(c,value_type(1.0/3.0));
447  exp(c,c);
448  }
449  else {
450  if (c.size() != 1)
451  c.resize(1);
452  c[0] = std::cbrt(a[0]);
453  }
454 }
455 
456 template <typename ordinal_type, typename value_type>
457 void
458 Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
459 pow(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
460  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a,
461  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& b)
462 {
463  if (a.size() > 1 || b.size() > 1) {
464  log(c,a);
465  timesEqual(c,b);
466  exp(c,c);
467  }
468  else {
469  if (c.size() != 1)
470  c.resize(1);
471  c[0] = std::pow(a[0], b[0]);
472  }
473 }
474 
475 template <typename ordinal_type, typename value_type>
476 void
477 Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
478 pow(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
479  const value_type& a,
480  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& b)
481 {
482  if (b.size() > 1) {
483  times(c,std::log(a),b);
484  exp(c,c);
485  }
486  else {
487  if (c.size() != 1)
488  c.resize(1);
489  c[0] = std::pow(a, b[0]);
490  }
491 }
492 
493 template <typename ordinal_type, typename value_type>
494 void
495 Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
496 pow(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
497  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a,
498  const value_type& b)
499 {
500  if (a.size() > 1) {
501  log(c,a);
502  timesEqual(c,b);
503  exp(c,c);
504  }
505  else {
506  if (c.size() != 1)
507  c.resize(1);
508  c[0] = std::pow(a[0], b);
509  }
510 }
511 
512 template <typename ordinal_type, typename value_type>
513 void
514 Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
515 sin(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& s,
516  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a)
517 {
518  if (a.size() == 1) {
519  if (s.size() != 1)
520  s.resize(1);
521  s[0] = std::sin(a[0]);
522  }
523  else
524  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
525  "Stokhos::ForUQTKOrthogPolyExpansion::sin()"
526  << ": Method not implemented!");
527 }
528 
529 template <typename ordinal_type, typename value_type>
530 void
531 Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
532 cos(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
533  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a)
534 {
535  if (a.size() == 1) {
536  if (c.size() != 1)
537  c.resize(1);
538  c[0] = std::cos(a[0]);
539  }
540  else
541  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
542  "Stokhos::ForUQTKOrthogPolyExpansion::cos()"
543  << ": Method not implemented!");
544 }
545 
546 template <typename ordinal_type, typename value_type>
547 void
548 Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
549 tan(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& t,
550  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a)
551 {
552  if (a.size() == 1) {
553  if (t.size() != 1)
554  t.resize(1);
555  t[0] = std::tan(a[0]);
556  }
557  else
558  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
559  "Stokhos::ForUQTKOrthogPolyExpansion::tan()"
560  << ": Method not implemented!");
561 }
562 
563 template <typename ordinal_type, typename value_type>
564 void
565 Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
566 sinh(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& s,
567  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a)
568 {
569  if (a.size() > 1) {
570  // sinh(x) = (exp(x) - exp(-x))/2.0
571  Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type> t(a);
572  timesEqual(t, -1.0);
573  exp(s, a);
574  exp(t, t);
575  this->minusEqual(s, t);
576  divideEqual(s, 2.0);
577  }
578  else {
579  if (s.size() != 1)
580  s.resize(1);
581  s[0] = std::sinh(a[0]);
582  }
583 }
584 
585 template <typename ordinal_type, typename value_type>
586 void
587 Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
588 cosh(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
589  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a)
590 {
591  if (a.size() > 1) {
592  // cosh(x) = (exp(x) + exp(-x))/2.0
593  Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type> t(a);
594  timesEqual(t, -1.0);
595  exp(c, a);
596  exp(t, t);
597  this->plusEqual(c, t);
598  divideEqual(c, 2.0);
599  }
600  else {
601  if (c.size() != 1)
602  c.resize(1);
603  c[0] = std::cosh(a[0]);
604  }
605 }
606 
607 template <typename ordinal_type, typename value_type>
608 void
609 Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
610 tanh(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& t,
611  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a)
612 {
613  if (a.size() > 1) {
614  // tanh(x) = (exp(x) - exp(-x))/(exp(x) + exp(-x))
615  Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type> s(a);
616  Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type> c(a);
617  timesEqual(s, -1.0);
618  exp(s, s);
619  exp(c, a);
620  this->minus(t, c, s);
621  this->plusEqual(c, s);
622  divideEqual(t, c);
623  }
624  else {
625  if (t.size() != 1)
626  t.resize(1);
627  t[0] = std::tanh(a[0]);
628  }
629 }
630 
631 template <typename ordinal_type, typename value_type>
632 void
633 Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
634 acos(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
635  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a)
636 {
637  if (a.size() == 1) {
638  if (c.size() != 1)
639  c.resize(1);
640  c[0] = std::acos(a[0]);
641  }
642  else
643  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
644  "Stokhos::ForUQTKOrthogPolyExpansion::acos()"
645  << ": Method not implemented!");
646 }
647 
648 template <typename ordinal_type, typename value_type>
649 void
650 Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
651 asin(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
652  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a)
653 {
654  if (a.size() == 1) {
655  if (c.size() != 1)
656  c.resize(1);
657  c[0] = std::asin(a[0]);
658  }
659  else
660  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
661  "Stokhos::ForUQTKOrthogPolyExpansion::asin()"
662  << ": Method not implemented!");
663 }
664 
665 template <typename ordinal_type, typename value_type>
666 void
667 Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
668 atan(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
669  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a)
670 {
671  if (a.size() == 1) {
672  if (c.size() != 1)
673  c.resize(1);
674  c[0] = std::atan(a[0]);
675  }
676  else
677  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
678  "Stokhos::ForUQTKOrthogPolyExpansion::atan()"
679  << ": Method not implemented!");
680 }
681 
682 template <typename ordinal_type, typename value_type>
683 void
684 Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
685 atan2(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
686  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a,
687  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& b)
688 {
689  if (a.size() == 1 && b.size() == 1) {
690  if (c.size() != 1)
691  c.resize(1);
692  c[0] = std::atan2(a[0], b[0]);
693  }
694  else
695  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
696  "Stokhos::ForUQTKOrthogPolyExpansion::atan2()"
697  << ": Method not implemented!");
698 }
699 
700 template <typename ordinal_type, typename value_type>
701 void
702 Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
703 atan2(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
704  const value_type& a,
705  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& b)
706 {
707  if (b.size() == 1) {
708  if (c.size() != 1)
709  c.resize(1);
710  c[0] = std::atan2(a, b[0]);
711  }
712  else
713  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
714  "Stokhos::ForUQTKOrthogPolyExpansion::atan2()"
715  << ": Method not implemented!");
716 }
717 
718 template <typename ordinal_type, typename value_type>
719 void
720 Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
721 atan2(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
722  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a,
723  const value_type& b)
724 {
725  if (a.size() == 1) {
726  if (c.size() != 1)
727  c.resize(1);
728  c[0] = std::atan2(a[0], b);
729  }
730  else
731  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
732  "Stokhos::ForUQTKOrthogPolyExpansion::atan2()"
733  << ": Method not implemented!");
734 }
735 
736 template <typename ordinal_type, typename value_type>
737 void
738 Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
739 acosh(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
740  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a)
741 {
742  if (a.size() == 1) {
743  if (c.size() != 1)
744  c.resize(1);
745  c[0] = std::log(a[0]+std::sqrt(a[0]*a[0]-value_type(1.0)));
746  }
747  else
748  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
749  "Stokhos::ForUQTKOrthogPolyExpansion::acosh()"
750  << ": Method not implemented!");
751 }
752 
753 template <typename ordinal_type, typename value_type>
754 void
755 Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
756 asinh(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
757  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a)
758 {
759  if (a.size() == 1) {
760  if (c.size() != 1)
761  c.resize(1);
762  c[0] = std::log(a[0]+std::sqrt(a[0]*a[0]+value_type(1.0)));
763  }
764  else
765  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
766  "Stokhos::ForUQTKOrthogPolyExpansion::asinh()"
767  << ": Method not implemented!");
768 }
769 
770 template <typename ordinal_type, typename value_type>
771 void
772 Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
773 atanh(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
774  const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a)
775 {
776  if (a.size() == 1) {
777  if (c.size() != 1)
778  c.resize(1);
779  c[0] = 0.5*std::log((value_type(1.0)+a[0])/(value_type(1.0)-a[0]));
780  }
781  else
782  TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
783  "Stokhos::ForUQTKOrthogPolyExpansion::atanh()"
784  << ": Method not implemented!");
785 }
Sacado::UQ::sqrt
KOKKOS_INLINE_FUNCTION PCE< Storage > sqrt(const PCE< Storage > &a)
Definition: Sacado_UQ_PCE_Imp.hpp:847
UQ_LOG_INT_F77
#define UQ_LOG_INT_F77
Definition: Stokhos_ForUQTKOrthogPolyExpansionImp.hpp:21
Sacado::UQ::tan
KOKKOS_INLINE_FUNCTION PCE< Storage > tan(const PCE< Storage > &a)
Definition: Sacado_UQ_PCE_Imp.hpp:968
Sacado::UQ::sinh
KOKKOS_INLINE_FUNCTION PCE< Storage > sinh(const PCE< Storage > &a)
Definition: Sacado_UQ_PCE_Imp.hpp:985
Stokhos::OrthogPolyApprox::resize
void resize(ordinal_type sz)
Resize coefficient array (coefficients are preserved)
Definition: Stokhos_OrthogPolyApproxImp.hpp:118
atanh
atanh(expr.val())
Stokhos::Sparse3Tensor
Data structure storing a sparse 3-tensor C(i,j,k) in a a compressed format.
Definition: Stokhos_Sparse3Tensor.hpp:24
Sacado::UQ::pow
KOKKOS_INLINE_FUNCTION PCE< Storage > pow(const PCE< Storage > &a, const PCE< Storage > &b)
Definition: Sacado_UQ_PCE_Imp.hpp:881
TEUCHOS_TEST_FOR_EXCEPTION
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
UQ_EXP_INT_F77
#define UQ_EXP_INT_F77
Definition: Stokhos_ForUQTKOrthogPolyExpansionImp.hpp:20
Teuchos_ConfigDefs.hpp
Stokhos::OrthogPolyApprox::coeff
pointer coeff()
Return coefficient array.
Definition: Stokhos_OrthogPolyApproxImp.hpp:134
Sacado::UQ::tanh
KOKKOS_INLINE_FUNCTION PCE< Storage > tanh(const PCE< Storage > &a)
Definition: Sacado_UQ_PCE_Imp.hpp:1019
Sacado::UQ::cbrt
KOKKOS_INLINE_FUNCTION PCE< Storage > cbrt(const PCE< Storage > &a)
Definition: Sacado_UQ_PCE_Imp.hpp:864
Stokhos::ds_array::get_and_fill
static T * get_and_fill(int sz)
Get memory for new array of length sz and fill with zeros.
Definition: Stokhos_DynamicArrayTraits.hpp:48
Sacado::UQ::acos
KOKKOS_INLINE_FUNCTION PCE< Storage > acos(const PCE< Storage > &a)
Definition: Sacado_UQ_PCE_Imp.hpp:1036
atan2
atan2(expr1.val(), expr2.val())
node_type
Kokkos::Serial node_type
Definition: Stokhos_SacadoTraitsUnitTest.cpp:29
asinh
asinh(expr.val())
Stokhos::OrthogPolyBasis
Abstract base class for multivariate orthogonal polynomials.
Definition: Stokhos_OrthogPolyBasis.hpp:40
Sacado::UQ::cosh
KOKKOS_INLINE_FUNCTION PCE< Storage > cosh(const PCE< Storage > &a)
Definition: Sacado_UQ_PCE_Imp.hpp:1002
UQ_PREP_F77
#define UQ_PREP_F77
Definition: Stokhos_ForUQTKOrthogPolyExpansionImp.hpp:14
Sacado::UQ::atan
KOKKOS_INLINE_FUNCTION PCE< Storage > atan(const PCE< Storage > &a)
Definition: Sacado_UQ_PCE_Imp.hpp:1070
Sacado::UQ::exp
KOKKOS_INLINE_FUNCTION PCE< Storage > exp(const PCE< Storage > &a)
Definition: Sacado_UQ_PCE_Imp.hpp:796
Stokhos::OrthogPolyApprox
Class to store coefficients of a projection onto an orthogonal polynomial basis.
Definition: Stokhos_OrthogPolyApprox.hpp:29
val
expr val()
acosh
acosh(expr.val())
Stokhos::ds_array::copy
static void copy(const T *src, T *dest, int sz)
Copy array from src to dest of length sz.
Definition: Stokhos_DynamicArrayTraits.hpp:69
Stokhos_DynamicArrayTraits.hpp
UQ_PROD2_F77
#define UQ_PROD2_F77
Definition: Stokhos_ForUQTKOrthogPolyExpansionImp.hpp:15
Sacado::UQ::sin
KOKKOS_INLINE_FUNCTION PCE< Storage > sin(const PCE< Storage > &a)
Definition: Sacado_UQ_PCE_Imp.hpp:934
Stokhos::OrthogPolyApprox::size
ordinal_type size() const
Return size.
Definition: Stokhos_OrthogPolyApproxImp.hpp:126
UQ_SQRT_F77
#define UQ_SQRT_F77
Definition: Stokhos_ForUQTKOrthogPolyExpansionImp.hpp:19
Sacado::UQ::log
KOKKOS_INLINE_FUNCTION PCE< Storage > log(const PCE< Storage > &a)
Definition: Sacado_UQ_PCE_Imp.hpp:813
Teuchos_Assert.hpp
Sacado::UQ::log10
KOKKOS_INLINE_FUNCTION PCE< Storage > log10(const PCE< Storage > &a)
Definition: Sacado_UQ_PCE_Imp.hpp:830
UQ_LOG_F77
#define UQ_LOG_F77
Definition: Stokhos_ForUQTKOrthogPolyExpansionImp.hpp:18
Sacado::UQ::asin
KOKKOS_INLINE_FUNCTION PCE< Storage > asin(const PCE< Storage > &a)
Definition: Sacado_UQ_PCE_Imp.hpp:1053
Sacado::UQ::cos
KOKKOS_INLINE_FUNCTION PCE< Storage > cos(const PCE< Storage > &a)
Definition: Sacado_UQ_PCE_Imp.hpp:951
UQ_DIV_F77
#define UQ_DIV_F77
Definition: Stokhos_ForUQTKOrthogPolyExpansionImp.hpp:16
UQ_EXP_F77
#define UQ_EXP_F77
Definition: Stokhos_ForUQTKOrthogPolyExpansionImp.hpp:17
Generated on Fri Nov 22 2024 09:18:12 for Stokhos Package Browser (Single Doxygen Collection) by   doxygen 1.8.5