Stokhos Package Browser (Single Doxygen Collection) Version of the Day
Loading...
Searching...
No Matches
Stokhos_ForUQTKOrthogPolyExpansionImp.hpp
Go to the documentation of this file.
1// @HEADER
2// ***********************************************************************
3//
4// Stokhos Package
5// Copyright (2009) Sandia Corporation
6//
7// Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8// license for use of this work by or on behalf of the U.S. Government.
9//
10// Redistribution and use in source and binary forms, with or without
11// modification, are permitted provided that the following conditions are
12// met:
13//
14// 1. Redistributions of source code must retain the above copyright
15// notice, this list of conditions and the following disclaimer.
16//
17// 2. Redistributions in binary form must reproduce the above copyright
18// notice, this list of conditions and the following disclaimer in the
19// documentation and/or other materials provided with the distribution.
20//
21// 3. Neither the name of the Corporation nor the names of the
22// contributors may be used to endorse or promote products derived from
23// this software without specific prior written permission.
24//
25// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36//
37// Questions? Contact Eric T. Phipps (etphipp@sandia.gov).
38//
39// ***********************************************************************
40// @HEADER
41
42#include "Teuchos_Assert.hpp"
44#include "Teuchos_ConfigDefs.hpp"
45
46#define UQ_PREP_F77 F77_FUNC_(uq_prep,UQ_PREP)
47#define UQ_PROD2_F77 F77_FUNC_(uq_prod2,UQ_PROD2)
48#define UQ_DIV_F77 F77_FUNC_(uq_div,UQ_DIV)
49#define UQ_EXP_F77 F77_FUNC_(uq_exp,UQ_EXP)
50#define UQ_LOG_F77 F77_FUNC_(uq_log,UQ_LOG)
51#define UQ_SQRT_F77 F77_FUNC_(uq_sqrt,UQ_SQRT)
52#define UQ_EXP_INT_F77 F77_FUNC_(uq_exp_int,UQ_EXP)
53#define UQ_LOG_INT_F77 F77_FUNC_(uq_log_int,UQ_LOG)
54
55extern "C" {
56 void UQ_PREP_F77(int*, int*, int*);
57 void UQ_PROD2_F77(const double*, const double*, double*, int*);
58 void UQ_DIV_F77(const double*, const double*, double*, int*);
59 void UQ_EXP_F77(const double*, double*, int*, int*, double*, int*);
60 void UQ_LOG_F77(const double*, double*, int*, int*, double*, int*);
61 void UQ_SQRT_F77(const double*, double*, int*, int*);
62 void UQ_EXP_INT_F77(const double*, double*, int*);
63 void UQ_LOG_INT_F77(const double*, double*, int*);
64}
65
66template <typename ordinal_type, typename value_type>
67Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
68ForUQTKOrthogPolyExpansion(
69 const Teuchos::RCP<const Stokhos::OrthogPolyBasis<ordinal_type, value_type> >& basis_,
70 const Teuchos::RCP<const Stokhos::Sparse3Tensor<ordinal_type, value_type> >& Cijk_,
71 EXPANSION_METHOD method_,
72 value_type rtol_) :
73 OrthogPolyExpansionBase<ordinal_type, value_type, node_type>(basis_, Cijk_),
74 rtol(rtol_),
75 method(method_)
76{
77 order = this->basis->order();
78 dim = this->basis->dimension();
79 int nup;
80 UQ_PREP_F77(&order, &dim, &nup);
81 sz = nup+1;
82}
83
84template <typename ordinal_type, typename value_type>
85void
86Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
88 const value_type& val)
89{
90 OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::timesEqual(c,val);
91}
92
93template <typename ordinal_type, typename value_type>
94void
95Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
97 const value_type& val)
98{
99 OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::divideEqual(c,val);
100}
101
102template <typename ordinal_type, typename value_type>
103void
104Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
105timesEqual(
108{
109 ordinal_type p = c.size();
110 ordinal_type xp = x.size();
111 ordinal_type pc;
112 if (p > 1 && xp > 1)
113 pc = sz;
114 else
115 pc = p*xp;
116 TEUCHOS_TEST_FOR_EXCEPTION(sz < pc, std::logic_error,
117 "Stokhos::ForUQTKOrthogPolyExpansion::timesEqual()" <<
118 ": Expansion size (" << sz <<
119 ") is too small for computation.");
120 if (c.size() != pc)
121 c.resize(pc);
122
123 value_type* cc = c.coeff();
124 const value_type* xc = x.coeff();
125
126 if (p > 1 && xp > 1) {
127 TEUCHOS_TEST_FOR_EXCEPTION(pc != xp, std::logic_error,
128 "Stokhos::ForUQTKOrthogPolyExpansion::timesEqual()"
129 << ": Arguments have incompatible sizes: "
130 << "x.size() = " << xp << ", c.size() = " << pc << ".");
131
132 // Copy c coefficients into temporary array
134
135 int nup = pc-1;
136 UQ_PROD2_F77(cc, xc, tc, &nup);
138 }
139 else if (p > 1) {
140 for (ordinal_type i=0; i<p; i++)
141 cc[i] *= xc[0];
142 }
143 else if (xp > 1) {
144 for (ordinal_type i=1; i<xp; i++)
145 cc[i] = cc[0]*xc[i];
146 cc[0] *= xc[0];
147 }
148 else {
149 cc[0] *= xc[0];
150 }
151}
152
153template <typename ordinal_type, typename value_type>
154void
155Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
157 const OrthogPolyApprox<ordinal_type, value_type, node_type>& x)
158{
159 ordinal_type p = c.size();
160 ordinal_type xp = x.size();
161 ordinal_type pc;
162 if (xp > 1)
163 pc = sz;
164 else
165 pc = p;
166 TEUCHOS_TEST_FOR_EXCEPTION(sz < pc, std::logic_error,
167 "Stokhos::ForUQTKOrthogPolyExpansion::divideEqual()" <<
168 ": Expansion size (" << sz <<
169 ") is too small for computation.");
170 if (c.size() != pc)
171 c.resize(pc);
172
173 value_type* cc = c.coeff();
174 const value_type* xc = x.coeff();
175
176 if (xp > 1) {
177 TEUCHOS_TEST_FOR_EXCEPTION(pc != xp, std::logic_error,
178 "Stokhos::ForUQTKOrthogPolyExpansion::divideEqual()"
179 << ": Arguments have incompatible sizes: "
180 << "x.size() = " << xp << ", c.size() = " << pc << ".");
181
182 // Copy c coefficients into temporary array
184
185 int nup = pc-1;
186 UQ_DIV_F77(tc, xc, cc, &nup);
187
188 }
189 else {
190 for (ordinal_type i=0; i<pc; i++)
191 cc[i] /= xc[0];
192 }
193}
194
195template <typename ordinal_type, typename value_type>
196void
197Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
201{
202 ordinal_type pa = a.size();
203 ordinal_type pb = b.size();
204 ordinal_type pc;
205 if (pa > 1 && pb > 1)
206 pc = sz;
207 else
208 pc = pa*pb;
209 TEUCHOS_TEST_FOR_EXCEPTION(sz < pc, std::logic_error,
210 "Stokhos::ForUQTKOrthogPolyExpansion::times()" <<
211 ": Expansion size (" << sz <<
212 ") is too small for computation.");
213 if (c.size() != pc)
214 c.resize(pc);
215
216 const value_type* ca = a.coeff();
217 const value_type* cb = b.coeff();
218 value_type* cc = c.coeff();
219
220 if (pa > 1 && pb > 1) {
221 TEUCHOS_TEST_FOR_EXCEPTION(pa != pc || pb != pc, std::logic_error,
222 "Stokhos::ForUQTKOrthogPolyExpansion::times()"
223 << ": Arguments have incompatible sizes: "
224 << "a.size() = " << pa << ", b.size() = " << pb
225 << ", required size = " << pc << ".");
226
227 int nup = pc-1;
228 UQ_PROD2_F77(ca, cb, cc, &nup);
229 }
230 else if (pa > 1) {
231 for (ordinal_type i=0; i<pc; i++)
232 cc[i] = ca[i]*cb[0];
233 }
234 else if (pb > 1) {
235 for (ordinal_type i=0; i<pc; i++)
236 cc[i] = ca[0]*cb[i];
237 }
238 else {
239 cc[0] = ca[0]*cb[0];
240 }
241}
242
243template <typename ordinal_type, typename value_type>
244void
245Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
247 const value_type& a,
249{
250 OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::times(c,a,b);
251}
252
253template <typename ordinal_type, typename value_type>
254void
255Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
258 const value_type& b)
259{
260 OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::times(c,a,b);
261}
262
263template <typename ordinal_type, typename value_type>
264void
265Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
269{
270 ordinal_type pa = a.size();
271 ordinal_type pb = b.size();
272 ordinal_type pc;
273 if (pb > 1)
274 pc = sz;
275 else
276 pc = pa;
277 TEUCHOS_TEST_FOR_EXCEPTION(sz < pc, std::logic_error,
278 "Stokhos::ForUQTKOrthogPolyExpansion::divide()" <<
279 ": Expansion size (" << sz <<
280 ") is too small for computation.");
281 if (c.size() != pc)
282 c.resize(pc);
283
284 const value_type* ca = a.coeff();
285 const value_type* cb = b.coeff();
286 value_type* cc = c.coeff();
287
288 if (pb > 1) {
289 TEUCHOS_TEST_FOR_EXCEPTION(pa != pc || pb != pc, std::logic_error,
290 "Stokhos::ForUQTKOrthogPolyExpansion::divide()"
291 << ": Arguments have incompatible sizes: "
292 << "a.size() = " << pa << ", b.size() = " << pb
293 << ", required size = " << pc << ".");
294
295 int nup = pc-1;
296 UQ_DIV_F77(ca, cb, cc, &nup);
297 }
298 else {
299 for (ordinal_type i=0; i<pa; i++)
300 cc[i] = ca[i]/cb[0];
301 }
302}
303
304template <typename ordinal_type, typename value_type>
305void
306Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
308 const value_type& a,
310{
311 ordinal_type pb = b.size();
312 ordinal_type pc;
313 if (pb > 1)
314 pc = sz;
315 else
316 pc = 1;
317 if (c.size() != pc)
318 c.resize(pc);
319
320 const value_type* cb = b.coeff();
321 value_type* cc = c.coeff();
322
323 if (pb > 1) {
324 TEUCHOS_TEST_FOR_EXCEPTION(pb != pc, std::logic_error,
325 "Stokhos::ForUQTKOrthogPolyExpansion::divide()"
326 << ": Arguments have incompatible sizes: "
327 << "b.size() = " << pb
328 << ", required size = " << pc << ".");
329
331 ca[0] = a;
332 int nup = pc-1;
333 UQ_DIV_F77(ca, cb, cc, &nup);
334 }
335 else
336 cc[0] = a / cb[0];
337}
338
339template <typename ordinal_type, typename value_type>
340void
341Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
344 const value_type& b)
345{
346 OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::divide(c,a,b);
347}
348
349template <typename ordinal_type, typename value_type>
350void
351Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
354{
355 ordinal_type pa = a.size();
356 ordinal_type pc;
357 if (pa > 1)
358 pc = sz;
359 else
360 pc = 1;
361 if (c.size() != pc)
362 c.resize(pc);
363
364 const value_type* ca = a.coeff();
365 value_type* cc = c.coeff();
366
367 if (pa > 1) {
368 TEUCHOS_TEST_FOR_EXCEPTION(pa != pc, std::logic_error,
369 "Stokhos::ForUQTKOrthogPolyExpansion::exp()"
370 << ": Arguments have incompatible sizes: "
371 << "a.size() = " << pa << ", c.size() = " << pc
372 << ".");
373 int nup = pc-1;
374 if (method == TAYLOR) {
375 int nrm = 1;
376 UQ_EXP_F77(ca, cc, &dim, &nup, &rtol, &nrm);
377 }
378 else
379 UQ_EXP_INT_F77(ca, cc, &nup);
380 }
381 else
382 cc[0] = std::exp(ca[0]);
383}
384
385template <typename ordinal_type, typename value_type>
386void
387Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
390{
391 ordinal_type pa = a.size();
392 ordinal_type pc;
393 if (pa > 1)
394 pc = sz;
395 else
396 pc = 1;
397 if (c.size() != pc)
398 c.resize(pc);
399
400 const value_type* ca = a.coeff();
401 value_type* cc = c.coeff();
402
403 if (pa > 1) {
404 TEUCHOS_TEST_FOR_EXCEPTION(pa != pc, std::logic_error,
405 "Stokhos::ForUQTKOrthogPolyExpansion::log()"
406 << ": Arguments have incompatible sizes: "
407 << "a.size() = " << pa << ", c.size() = " << pc
408 << ".");
409 int nup = pc-1;
410 if (method == TAYLOR) {
411 int nrm = 1;
412 UQ_LOG_F77(ca, cc, &dim, &nup, &rtol, &nrm);
413 }
414 else
415 UQ_LOG_INT_F77(ca, cc, &nup);
416 }
417 else
418 cc[0] = std::log(ca[0]);
419}
420
421template <typename ordinal_type, typename value_type>
422void
423Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
426{
427 if (a.size() > 1) {
428 log(c,a);
429 divide(c,c,std::log(10.0));
430 }
431 else {
432 if (c.size() != 1)
433 c.resize(1);
434 c[0] = std::log10(a[0]);
435 }
436}
437
438template <typename ordinal_type, typename value_type>
439void
440Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
443{
444 ordinal_type pa = a.size();
445 ordinal_type pc;
446 if (pa > 1)
447 pc = sz;
448 else
449 pc = 1;
450 if (c.size() != pc)
451 c.resize(pc);
452
453 const value_type* ca = a.coeff();
454 value_type* cc = c.coeff();
455
456 if (pa > 1) {
457 TEUCHOS_TEST_FOR_EXCEPTION(pa != pc, std::logic_error,
458 "Stokhos::ForUQTKOrthogPolyExpansion::sqrt()"
459 << ": Arguments have incompatible sizes: "
460 << "a.size() = " << pa << ", c.size() = " << pc
461 << ".");
462 int iguess = 0;
463 int nup = pc-1;
464 UQ_SQRT_F77(ca, cc, &nup, &iguess);
465 }
466 else
467 cc[0] = std::sqrt(ca[0]);
468}
469
470template <typename ordinal_type, typename value_type>
471void
472Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
475{
476 if (a.size() > 1) {
477 log(c,a);
478 timesEqual(c,value_type(1.0/3.0));
479 exp(c,c);
480 }
481 else {
482 if (c.size() != 1)
483 c.resize(1);
484 c[0] = std::cbrt(a[0]);
485 }
486}
487
488template <typename ordinal_type, typename value_type>
489void
490Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
494{
495 if (a.size() > 1 || b.size() > 1) {
496 log(c,a);
497 timesEqual(c,b);
498 exp(c,c);
499 }
500 else {
501 if (c.size() != 1)
502 c.resize(1);
503 c[0] = std::pow(a[0], b[0]);
504 }
505}
506
507template <typename ordinal_type, typename value_type>
508void
509Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
511 const value_type& a,
513{
514 if (b.size() > 1) {
515 times(c,std::log(a),b);
516 exp(c,c);
517 }
518 else {
519 if (c.size() != 1)
520 c.resize(1);
521 c[0] = std::pow(a, b[0]);
522 }
523}
524
525template <typename ordinal_type, typename value_type>
526void
527Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
530 const value_type& b)
531{
532 if (a.size() > 1) {
533 log(c,a);
534 timesEqual(c,b);
535 exp(c,c);
536 }
537 else {
538 if (c.size() != 1)
539 c.resize(1);
540 c[0] = std::pow(a[0], b);
541 }
542}
543
544template <typename ordinal_type, typename value_type>
545void
546Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
549{
550 if (a.size() == 1) {
551 if (s.size() != 1)
552 s.resize(1);
553 s[0] = std::sin(a[0]);
554 }
555 else
556 TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
557 "Stokhos::ForUQTKOrthogPolyExpansion::sin()"
558 << ": Method not implemented!");
559}
560
561template <typename ordinal_type, typename value_type>
562void
563Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
566{
567 if (a.size() == 1) {
568 if (c.size() != 1)
569 c.resize(1);
570 c[0] = std::cos(a[0]);
571 }
572 else
573 TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
574 "Stokhos::ForUQTKOrthogPolyExpansion::cos()"
575 << ": Method not implemented!");
576}
577
578template <typename ordinal_type, typename value_type>
579void
580Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
583{
584 if (a.size() == 1) {
585 if (t.size() != 1)
586 t.resize(1);
587 t[0] = std::tan(a[0]);
588 }
589 else
590 TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
591 "Stokhos::ForUQTKOrthogPolyExpansion::tan()"
592 << ": Method not implemented!");
593}
594
595template <typename ordinal_type, typename value_type>
596void
597Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
600{
601 if (a.size() > 1) {
602 // sinh(x) = (exp(x) - exp(-x))/2.0
604 timesEqual(t, -1.0);
605 exp(s, a);
606 exp(t, t);
607 this->minusEqual(s, t);
608 divideEqual(s, 2.0);
609 }
610 else {
611 if (s.size() != 1)
612 s.resize(1);
613 s[0] = std::sinh(a[0]);
614 }
615}
616
617template <typename ordinal_type, typename value_type>
618void
619Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
622{
623 if (a.size() > 1) {
624 // cosh(x) = (exp(x) + exp(-x))/2.0
626 timesEqual(t, -1.0);
627 exp(c, a);
628 exp(t, t);
629 this->plusEqual(c, t);
630 divideEqual(c, 2.0);
631 }
632 else {
633 if (c.size() != 1)
634 c.resize(1);
635 c[0] = std::cosh(a[0]);
636 }
637}
638
639template <typename ordinal_type, typename value_type>
640void
641Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
644{
645 if (a.size() > 1) {
646 // tanh(x) = (exp(x) - exp(-x))/(exp(x) + exp(-x))
649 timesEqual(s, -1.0);
650 exp(s, s);
651 exp(c, a);
652 this->minus(t, c, s);
653 this->plusEqual(c, s);
654 divideEqual(t, c);
655 }
656 else {
657 if (t.size() != 1)
658 t.resize(1);
659 t[0] = std::tanh(a[0]);
660 }
661}
662
663template <typename ordinal_type, typename value_type>
664void
665Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
668{
669 if (a.size() == 1) {
670 if (c.size() != 1)
671 c.resize(1);
672 c[0] = std::acos(a[0]);
673 }
674 else
675 TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
676 "Stokhos::ForUQTKOrthogPolyExpansion::acos()"
677 << ": Method not implemented!");
678}
679
680template <typename ordinal_type, typename value_type>
681void
682Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
685{
686 if (a.size() == 1) {
687 if (c.size() != 1)
688 c.resize(1);
689 c[0] = std::asin(a[0]);
690 }
691 else
692 TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
693 "Stokhos::ForUQTKOrthogPolyExpansion::asin()"
694 << ": Method not implemented!");
695}
696
697template <typename ordinal_type, typename value_type>
698void
699Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
702{
703 if (a.size() == 1) {
704 if (c.size() != 1)
705 c.resize(1);
706 c[0] = std::atan(a[0]);
707 }
708 else
709 TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
710 "Stokhos::ForUQTKOrthogPolyExpansion::atan()"
711 << ": Method not implemented!");
712}
713
714template <typename ordinal_type, typename value_type>
715void
716Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
720{
721 if (a.size() == 1 && b.size() == 1) {
722 if (c.size() != 1)
723 c.resize(1);
724 c[0] = std::atan2(a[0], b[0]);
725 }
726 else
727 TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
728 "Stokhos::ForUQTKOrthogPolyExpansion::atan2()"
729 << ": Method not implemented!");
730}
731
732template <typename ordinal_type, typename value_type>
733void
734Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
736 const value_type& a,
738{
739 if (b.size() == 1) {
740 if (c.size() != 1)
741 c.resize(1);
742 c[0] = std::atan2(a, b[0]);
743 }
744 else
745 TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
746 "Stokhos::ForUQTKOrthogPolyExpansion::atan2()"
747 << ": Method not implemented!");
748}
749
750template <typename ordinal_type, typename value_type>
751void
752Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
755 const value_type& b)
756{
757 if (a.size() == 1) {
758 if (c.size() != 1)
759 c.resize(1);
760 c[0] = std::atan2(a[0], b);
761 }
762 else
763 TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
764 "Stokhos::ForUQTKOrthogPolyExpansion::atan2()"
765 << ": Method not implemented!");
766}
767
768template <typename ordinal_type, typename value_type>
769void
770Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
773{
774 if (a.size() == 1) {
775 if (c.size() != 1)
776 c.resize(1);
777 c[0] = std::log(a[0]+std::sqrt(a[0]*a[0]-value_type(1.0)));
778 }
779 else
780 TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
781 "Stokhos::ForUQTKOrthogPolyExpansion::acosh()"
782 << ": Method not implemented!");
783}
784
785template <typename ordinal_type, typename value_type>
786void
787Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
790{
791 if (a.size() == 1) {
792 if (c.size() != 1)
793 c.resize(1);
794 c[0] = std::log(a[0]+std::sqrt(a[0]*a[0]+value_type(1.0)));
795 }
796 else
797 TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
798 "Stokhos::ForUQTKOrthogPolyExpansion::asinh()"
799 << ": Method not implemented!");
800}
801
802template <typename ordinal_type, typename value_type>
803void
804Stokhos::ForUQTKOrthogPolyExpansion<ordinal_type, value_type>::
807{
808 if (a.size() == 1) {
809 if (c.size() != 1)
810 c.resize(1);
811 c[0] = 0.5*std::log((value_type(1.0)+a[0])/(value_type(1.0)-a[0]));
812 }
813 else
814 TEUCHOS_TEST_FOR_EXCEPTION(true, std::logic_error,
815 "Stokhos::ForUQTKOrthogPolyExpansion::atanh()"
816 << ": Method not implemented!");
817}
expr val()
Kokkos::Serial node_type
Class to store coefficients of a projection onto an orthogonal polynomial basis.
void resize(ordinal_type sz)
Resize coefficient array (coefficients are preserved)
pointer coeff()
Return coefficient array.
ordinal_type size() const
Return size.
Abstract base class for multivariate orthogonal polynomials.
Data structure storing a sparse 3-tensor C(i,j,k) in a a compressed format.
KOKKOS_INLINE_FUNCTION PCE< Storage > log(const PCE< Storage > &a)
KOKKOS_INLINE_FUNCTION PCE< Storage > exp(const PCE< Storage > &a)
const IndexType const IndexType const IndexType const IndexType const ValueType const ValueType * x
Definition: csr_vector.h:265
static void copy(const T *src, T *dest, int sz)
Copy array from src to dest of length sz.
static T * get_and_fill(int sz)
Get memory for new array of length sz and fill with zeros.