Stokhos Package Browser (Single Doxygen Collection)  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Stokhos_ConstantOrthogPolyExpansionImp.hpp
Go to the documentation of this file.
1 // $Id$
2 // $Source$
3 // @HEADER
4 // ***********************************************************************
5 //
6 // Stokhos Package
7 // Copyright (2009) Sandia Corporation
8 //
9 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
10 // license for use of this work by or on behalf of the U.S. Government.
11 //
12 // Redistribution and use in source and binary forms, with or without
13 // modification, are permitted provided that the following conditions are
14 // met:
15 //
16 // 1. Redistributions of source code must retain the above copyright
17 // notice, this list of conditions and the following disclaimer.
18 //
19 // 2. Redistributions in binary form must reproduce the above copyright
20 // notice, this list of conditions and the following disclaimer in the
21 // documentation and/or other materials provided with the distribution.
22 //
23 // 3. Neither the name of the Corporation nor the names of the
24 // contributors may be used to endorse or promote products derived from
25 // this software without specific prior written permission.
26 //
27 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
28 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
29 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
30 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
31 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
32 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
33 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
34 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
35 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
36 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
37 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
38 //
39 // Questions? Contact Eric T. Phipps (etphipp@sandia.gov).
40 //
41 // ***********************************************************************
42 // @HEADER
43 
44 template <typename ordinal_type, typename value_type>
47 {
48 }
49 
50 template <typename ordinal_type, typename value_type>
51 void
56 {
57  if (c.size() < 1)
58  c.resize(1);
59  c[0] = -a[0];
60 }
61 
62 template <typename ordinal_type, typename value_type>
63 void
66  const value_type& val)
67 {
68  if (c.size() < 1)
69  c.resize(1);
70  c[0] += val;
71 }
72 
73 template <typename ordinal_type, typename value_type>
74 void
77  const value_type& val)
78 {
79  if (c.size() < 1)
80  c.resize(1);
81  c[0] -= val;
82 }
83 
84 template <typename ordinal_type, typename value_type>
85 void
88  const value_type& val)
89 {
90  if (c.size() < 1)
91  c.resize(1);
92  c[0] *= val;
93 }
94 
95 template <typename ordinal_type, typename value_type>
96 void
99  const value_type& val)
100 {
101  if (c.size() < 1)
102  c.resize(1);
103  c[0] /= val;
104 }
105 
106 template <typename ordinal_type, typename value_type>
107 void
112 {
113  if (c.size() < 1)
114  c.resize(1);
115  c[0] += x[0];
116 }
117 
118 template <typename ordinal_type, typename value_type>
119 void
124 {
125  if (c.size() < 1)
126  c.resize(1);
127  c[0] -= x[0];
128 }
129 
130 template <typename ordinal_type, typename value_type>
131 void
136 {
137  if (c.size() < 1)
138  c.resize(1);
139  c[0] *= x[0];
140 }
141 
142 template <typename ordinal_type, typename value_type>
143 void
147 {
148  if (c.size() < 1)
149  c.resize(1);
150  c[0] /= x[0];
151 }
152 
153 template <typename ordinal_type, typename value_type>
154 void
159 {
160  if (c.size() < 1)
161  c.resize(1);
162  c[0] = a[0] + b[0];
163 }
164 
165 template <typename ordinal_type, typename value_type>
166 void
169  const value_type& a,
171 {
172  if (c.size() < 1)
173  c.resize(1);
174  c[0] = a + b[0];
175 }
176 
177 template <typename ordinal_type, typename value_type>
178 void
182  const value_type& b)
183 {
184  if (c.size() < 1)
185  c.resize(1);
186  c[0] = a[0] + b;
187 }
188 
189 template <typename ordinal_type, typename value_type>
190 void
195 {
196  if (c.size() < 1)
197  c.resize(1);
198  c[0] = a[0] - b[0];
199 }
200 
201 template <typename ordinal_type, typename value_type>
202 void
205  const value_type& a,
207 {
208  if (c.size() < 1)
209  c.resize(1);
210  c[0] = a - b[0];
211 }
212 
213 template <typename ordinal_type, typename value_type>
214 void
218  const value_type& b)
219 {
220  if (c.size() < 1)
221  c.resize(1);
222  c[0] = a[0] - b;
223 }
224 
225 template <typename ordinal_type, typename value_type>
226 void
231 {
232  if (c.size() < 1)
233  c.resize(1);
234  c[0] = a[0] * b[0];
235 }
236 
237 template <typename ordinal_type, typename value_type>
238 void
241  const value_type& a,
243 {
244  if (c.size() < 1)
245  c.resize(1);
246  c[0] = a * b[0];
247 }
248 
249 template <typename ordinal_type, typename value_type>
250 void
254  const value_type& b)
255 {
256  if (c.size() < 1)
257  c.resize(1);
258  c[0] = a[0] * b;
259 }
260 
261 template <typename ordinal_type, typename value_type>
262 void
267 {
268  if (c.size() < 1)
269  c.resize(1);
270  c[0] = a[0] / b[0];
271 }
272 
273 template <typename ordinal_type, typename value_type>
274 void
277  const value_type& a,
279 {
280  if (c.size() < 1)
281  c.resize(1);
282  c[0] = a / b[0];
283 }
284 
285 template <typename ordinal_type, typename value_type>
286 void
290  const value_type& b)
291 {
292  if (c.size() < 1)
293  c.resize(1);
294  c[0] = a[0] / b;
295 }
296 
297 template <typename ordinal_type, typename value_type>
298 void
302 {
303  if (c.size() < 1)
304  c.resize(1);
305  c[0] = std::exp(a[0]);
306 }
307 
308 template <typename ordinal_type, typename value_type>
309 void
313 {
314  if (c.size() < 1)
315  c.resize(1);
316  c[0] = std::log(a[0]);
317 }
318 
319 template <typename ordinal_type, typename value_type>
320 void
324 {
325  if (c.size() < 1)
326  c.resize(1);
327  c[0] = std::log10(a[0]);
328 }
329 
330 template <typename ordinal_type, typename value_type>
331 void
335 {
336  if (c.size() < 1)
337  c.resize(1);
338  c[0] = std::sqrt(a[0]);
339 }
340 
341 template <typename ordinal_type, typename value_type>
342 void
346 {
347  if (c.size() < 1)
348  c.resize(1);
349  c[0] = std::cbrt(a[0]);
350 }
351 
352 template <typename ordinal_type, typename value_type>
353 void
358 {
359  if (c.size() < 1)
360  c.resize(1);
361  c[0] = std::pow(a[0], b[0]);
362 }
363 
364 template <typename ordinal_type, typename value_type>
365 void
368  const value_type& a,
370 {
371  if (c.size() < 1)
372  c.resize(1);
373  c[0] = std::pow(a, b[0]);
374 }
375 
376 template <typename ordinal_type, typename value_type>
377 void
381  const value_type& b)
382 {
383  if (c.size() < 1)
384  c.resize(1);
385  c[0] = std::pow(a[0], b);
386 }
387 
388 template <typename ordinal_type, typename value_type>
389 void
393 {
394  if (s.size() < 1)
395  s.resize(1);
396  s[0] = std::sin(a[0]);
397 }
398 
399 template <typename ordinal_type, typename value_type>
400 void
404 {
405  if (c.size() < 1)
406  c.resize(1);
407  c[0] = std::cos(a[0]);
408 }
409 
410 template <typename ordinal_type, typename value_type>
411 void
415 {
416  if (t.size() < 1)
417  t.resize(1);
418  t[0] = std::tan(a[0]);
419 }
420 
421 template <typename ordinal_type, typename value_type>
422 void
426 {
427  if (s.size() < 1)
428  s.resize(1);
429  s[0] = std::sinh(a[0]);
430 }
431 
432 template <typename ordinal_type, typename value_type>
433 void
437 {
438  if (c.size() < 1)
439  c.resize(1);
440  c[0] = std::cosh(a[0]);
441 }
442 
443 template <typename ordinal_type, typename value_type>
444 void
448 {
449  if (t.size() < 1)
450  t.resize(1);
451  t[0] = std::tanh(a[0]);
452 }
453 
454 template <typename ordinal_type, typename value_type>
455 void
459 {
460  if (c.size() < 1)
461  c.resize(1);
462  c[0] = std::acos(a[0]);
463 }
464 
465 template <typename ordinal_type, typename value_type>
466 void
470 {
471  if (c.size() < 1)
472  c.resize(1);
473  c[0] = std::asin(a[0]);
474 }
475 
476 template <typename ordinal_type, typename value_type>
477 void
481 {
482  if (c.size() < 1)
483  c.resize(1);
484  c[0] = std::atan(a[0]);
485 }
486 
487 template <typename ordinal_type, typename value_type>
488 void
493 {
494  if (c.size() < 1)
495  c.resize(1);
496  c[0] = std::atan2(a[0], b[0]);
497 }
498 
499 template <typename ordinal_type, typename value_type>
500 void
503  const value_type& a,
505 {
506  if (c.size() < 1)
507  c.resize(1);
508  c[0] = std::atan2(a, b[0]);
509 }
510 
511 template <typename ordinal_type, typename value_type>
512 void
516  const value_type& b)
517 {
518  c[0] = std::atan2(a[0], b);
519 }
520 
521 template <typename ordinal_type, typename value_type>
522 void
526 {
527  if (c.size() < 1)
528  c.resize(1);
529  c[0] = std::log(a[0]+std::sqrt(a[0]*a[0]-value_type(1.0)));
530 }
531 
532 template <typename ordinal_type, typename value_type>
533 void
537 {
538  if (c.size() < 1)
539  c.resize(1);
540  c[0] = std::log(a[0]+std::sqrt(a[0]*a[0]+value_type(1.0)));
541 }
542 
543 template <typename ordinal_type, typename value_type>
544 void
548 {
549  if (c.size() < 1)
550  c.resize(1);
551  c[0] = 0.5*std::log((value_type(1.0)+a[0])/(value_type(1.0)-a[0]));
552 }
553 
554 template <typename ordinal_type, typename value_type>
555 void
559 {
560  if (c.size() < 1)
561  c.resize(1);
562  c[0] = std::fabs(a[0]);
563 }
564 
565 template <typename ordinal_type, typename value_type>
566 void
570 {
571  if (c.size() < 1)
572  c.resize(1);
573  c[0] = std::fabs(a[0]);
574 }
575 
576 template <typename ordinal_type, typename value_type>
577 void
582 {
583  if (c.size() < 1)
584  c.resize(1);
585  c[0] = std::max(a[0], b[0]);
586 }
587 
588 template <typename ordinal_type, typename value_type>
589 void
592  const value_type& a,
594 {
595  if (c.size() < 1)
596  c.resize(1);
597  c[0] = std::max(a, b[0]);
598 }
599 
600 template <typename ordinal_type, typename value_type>
601 void
605  const value_type& b)
606 {
607  if (c.size() < 1)
608  c.resize(1);
609  c[0] = std::max(a[0], b);
610 }
611 
612 template <typename ordinal_type, typename value_type>
613 void
618 {
619  if (c.size() < 1)
620  c.resize(1);
621  c[0] = std::min(a[0], b[0]);
622 }
623 
624 template <typename ordinal_type, typename value_type>
625 void
628  const value_type& a,
630 {
631  if (c.size() < 1)
632  c.resize(1);
633  c[0] = std::min(a, b[0]);
634 }
635 
636 template <typename ordinal_type, typename value_type>
637 void
641  const value_type& b)
642 {
643  if (c.size() < 1)
644  c.resize(1);
645  c[0] = std::min(a[0], b);
646 }
KOKKOS_INLINE_FUNCTION PCE< Storage > sqrt(const PCE< Storage > &a)
KOKKOS_INLINE_FUNCTION PCE< Storage > fabs(const PCE< Storage > &a)
KOKKOS_INLINE_FUNCTION PCE< Storage > tan(const PCE< Storage > &a)
void cbrt(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
void timesEqual(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const value_type &x)
KOKKOS_INLINE_FUNCTION PCE< Storage > sinh(const PCE< Storage > &a)
void resize(ordinal_type sz)
Resize coefficient array (coefficients are preserved)
KOKKOS_INLINE_FUNCTION PCE< Storage > pow(const PCE< Storage > &a, const PCE< Storage > &b)
void divide(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a, const OrthogPolyApprox< ordinal_type, value_type, node_type > &b)
void acosh(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
void cosh(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
void cos(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
void sin(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
void pow(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a, const OrthogPolyApprox< ordinal_type, value_type, node_type > &b)
KOKKOS_INLINE_FUNCTION PCE< Storage > tanh(const PCE< Storage > &a)
KOKKOS_INLINE_FUNCTION PCE< Storage > cbrt(const PCE< Storage > &a)
void sqrt(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
void unaryMinus(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
void plus(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a, const OrthogPolyApprox< ordinal_type, value_type, node_type > &b)
KOKKOS_INLINE_FUNCTION PCE< Storage > acos(const PCE< Storage > &a)
atan2(expr1.val(), expr2.val())
KOKKOS_INLINE_FUNCTION PCE< Storage > min(const typename PCE< Storage >::value_type &a, const PCE< Storage > &b)
void sinh(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
void minusEqual(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const value_type &x)
void asin(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
void max(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a, const OrthogPolyApprox< ordinal_type, value_type, node_type > &b)
KOKKOS_INLINE_FUNCTION PCE< Storage > max(const typename PCE< Storage >::value_type &a, const PCE< Storage > &b)
KOKKOS_INLINE_FUNCTION PCE< Storage > cosh(const PCE< Storage > &a)
void log10(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
void tanh(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
void atan(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
void divideEqual(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const value_type &x)
void plusEqual(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const value_type &x)
void fabs(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
void asinh(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
KOKKOS_INLINE_FUNCTION PCE< Storage > atan(const PCE< Storage > &a)
KOKKOS_INLINE_FUNCTION PCE< Storage > exp(const PCE< Storage > &a)
void min(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a, const OrthogPolyApprox< ordinal_type, value_type, node_type > &b)
void times(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a, const OrthogPolyApprox< ordinal_type, value_type, node_type > &b)
Class to store coefficients of a projection onto an orthogonal polynomial basis.
void tan(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
expr val()
void exp(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
void abs(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
KOKKOS_INLINE_FUNCTION PCE< Storage > sin(const PCE< Storage > &a)
ordinal_type size() const
Return size.
void atanh(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
void atan2(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a, const OrthogPolyApprox< ordinal_type, value_type, node_type > &b)
void acos(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
KOKKOS_INLINE_FUNCTION PCE< Storage > log(const PCE< Storage > &a)
KOKKOS_INLINE_FUNCTION PCE< Storage > log10(const PCE< Storage > &a)
KOKKOS_INLINE_FUNCTION PCE< Storage > asin(const PCE< Storage > &a)
void minus(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a, const OrthogPolyApprox< ordinal_type, value_type, node_type > &b)
KOKKOS_INLINE_FUNCTION PCE< Storage > cos(const PCE< Storage > &a)
void log(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)