Stokhos Package Browser (Single Doxygen Collection) Version of the Day
Loading...
Searching...
No Matches
Stokhos_OrthogPolyExpansionBaseImp.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 "Stokhos_DenseDirectDivisionExpansionStrategy.hpp"
43#include "Stokhos_SPDDenseDirectDivisionExpansionStrategy.hpp"
44#include "Stokhos_MeanBasedDivisionExpansionStrategy.hpp"
45#include "Stokhos_GMRESDivisionExpansionStrategy.hpp"
46#include "Stokhos_CGDivisionExpansionStrategy.hpp"
47#include "Teuchos_Assert.hpp"
48#include "Teuchos_TimeMonitor.hpp"
49#include "Stokhos_DynamicArrayTraits.hpp"
50
51template <typename ordinal_type, typename value_type, typename node_type>
52Stokhos::OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::
53OrthogPolyExpansionBase(
54 const Teuchos::RCP<const Stokhos::OrthogPolyBasis<ordinal_type, value_type> >& basis_,
55 const Teuchos::RCP<const Stokhos::Sparse3Tensor<ordinal_type, value_type> >& Cijk_,
56 const Teuchos::RCP<Teuchos::ParameterList>& params_) :
57 basis(basis_),
58 Cijk(Cijk_),
59 params(params_),
60 sz(basis->size())
61{
62 if (params == Teuchos::null)
63 params = Teuchos::rcp(new Teuchos::ParameterList);
64
65 // Create division strategy
66 std::string name = params->get("Division Strategy","Dense Direct");
67 double tol = params->get("Division Tolerance", 1e-6);
68 int prec_iter = params->get("prec_iter", 1);
69 int max_it = params->get("max_it_div", 50);
70 std::string prec = params->get("Prec Strategy","None");
71 int PrecNum;
72 if (prec == "None")
73 PrecNum=0;
74 else if (prec == "Diag")
75 PrecNum=1;
76 else if (prec == "Jacobi")
77 PrecNum=2;
78 else if (prec == "GS")
79 PrecNum=3;
80 else if (prec == "Schur")
81 PrecNum=4;
82 else
83 PrecNum=-1;
84 std::string linopt = params->get("Prec option", "whole");
85 int linear = 0;
86 if (linopt == "linear")
87 linear = 1;
88
89 std::string schuropt = params->get("Schur option", "diag");
90 int diag = 0;
91 if (schuropt == "full")
92 diag = 1;
93
94 int equil = params->get("Equilibrate", 0);
95
96
97 if (name == "Dense Direct")
98 division_strategy =
99 Teuchos::rcp(new DenseDirectDivisionExpansionStrategy<ordinal_type,value_type,node_type>(this->basis, this->Cijk));
100 else if (name == "SPD Dense Direct")
101 division_strategy =
102 Teuchos::rcp(new SPDDenseDirectDivisionExpansionStrategy<ordinal_type,value_type,node_type>(this->basis, this->Cijk));
103 else if (name == "Mean-Based")
104 division_strategy =
105 Teuchos::rcp(new MeanBasedDivisionExpansionStrategy<ordinal_type,value_type,node_type>());
106 else if (name == "GMRES")
107 division_strategy =
108 Teuchos::rcp(new GMRESDivisionExpansionStrategy<ordinal_type,value_type,node_type>(this->basis, this->Cijk, prec_iter, tol, PrecNum, max_it, linear, diag, equil));
109 else if (name == "CG")
110 division_strategy =
111 Teuchos::rcp(new CGDivisionExpansionStrategy<ordinal_type,value_type,node_type>(this->basis, this->Cijk, prec_iter, tol, PrecNum, max_it, linear, diag, equil));
112
113// TEUCHOS_TEST_FOR_EXCEPTION(
114// true, std::logic_error,
115// "Invalid division strategy name" << name);
116}
117
118template <typename ordinal_type, typename value_type, typename node_type>
119void
120Stokhos::OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::
121unaryMinus(
122 Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
123 const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a)
124{
125#ifdef STOKHOS_TEUCHOS_TIME_MONITOR
126 TEUCHOS_FUNC_TIME_MONITOR("Stokhos::OrthogPolyExpansionBase::unaryMinus(OPA)");
127#endif
128 ordinal_type pc = a.size();
129 if (c.size() != pc)
130 c.resize(pc);
131
132 value_type* cc = c.coeff();
133 const value_type* ca = a.coeff();
134
135 for (ordinal_type i=0; i<pc; i++)
136 {
137 cc[i] = -ca[i];
138 }
139}
140template <typename ordinal_type, typename value_type, typename node_type>
141void
142Stokhos::OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::
143plusEqual(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
144 const value_type& val)
145{
146 c[0] += val;
147}
148
149template <typename ordinal_type, typename value_type, typename node_type>
150void
151Stokhos::OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::
152minusEqual(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
153 const value_type& val)
154{
155 c[0] -= val;
156}
157
158template <typename ordinal_type, typename value_type, typename node_type>
159void
160Stokhos::OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::
161timesEqual(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
162 const value_type& val)
163{
164#ifdef STOKHOS_TEUCHOS_TIME_MONITOR
165 TEUCHOS_FUNC_TIME_MONITOR("Stokhos::OrthogPolyExpansionBase::timesEqual(const)");
166#endif
167 ordinal_type pc = c.size();
168 value_type* cc = c.coeff();
169 for (ordinal_type i=0; i<pc; i++)
170 cc[i] *= val;
171}
172
173template <typename ordinal_type, typename value_type, typename node_type>
174void
175Stokhos::OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::
176divideEqual(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
177 const value_type& val)
178{
179#ifdef STOKHOS_TEUCHOS_TIME_MONITOR
180 TEUCHOS_FUNC_TIME_MONITOR("Stokhos::OrthogPolyExpansionBase::divideEqual(const)");
181#endif
182 ordinal_type pc = c.size();
183 value_type* cc=c.coeff();
184 for (ordinal_type i=0; i<pc; i++){
185 cc[i] /= val;
186}
187}
188template <typename ordinal_type, typename value_type, typename node_type>
189void
190Stokhos::OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::
191plusEqual(
192 Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
193 const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& x)
194{
195#ifdef STOKHOS_TEUCHOS_TIME_MONITOR
196 TEUCHOS_FUNC_TIME_MONITOR("Stokhos::OrthogPolyExpansionBase::plusEqual(OPA)");
197#endif
198 ordinal_type xp = x.size();
199 if (c.size() < xp)
200 c.resize(xp);
201
202 value_type* cc = c.coeff();
203 const value_type* xc = x.coeff();
204 for (ordinal_type i=0; i<xp; i++)
205 cc[i] += xc[i];
206}
207
208template <typename ordinal_type, typename value_type, typename node_type>
209void
210Stokhos::OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::
211minusEqual(
212 Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
213 const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& x)
214{
215#ifdef STOKHOS_TEUCHOS_TIME_MONITOR
216 TEUCHOS_FUNC_TIME_MONITOR("Stokhos::OrthogPolyExpansionBase::minusEqual(OPA)");
217#endif
218 ordinal_type xp = x.size();
219 if (c.size() < xp)
220 c.resize(xp);
221
222 value_type* cc = c.coeff();
223 const value_type* xc = x.coeff();
224 for (ordinal_type i=0; i<xp; i++)
225 cc[i] -= xc[i];
226}
227
228template <typename ordinal_type, typename value_type, typename node_type>
229void
230Stokhos::OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::
231timesEqual(
232 Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
233 const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& x)
234{
235#ifdef STOKHOS_TEUCHOS_TIME_MONITOR
236 TEUCHOS_FUNC_TIME_MONITOR("Stokhos::OrthogPolyExpansionBase::timesEqual(OPA)");
237#endif
238 ordinal_type p = c.size();
239 ordinal_type xp = x.size();
240 ordinal_type pc;
241 if (p > 1 && xp > 1)
242 pc = sz;
243 else
244 pc = p*xp;
245 TEUCHOS_TEST_FOR_EXCEPTION(sz < pc, std::logic_error,
246 "Stokhos::OrthogPolyExpansionBase::timesEqual()" <<
247 ": Expansion size (" << sz <<
248 ") is too small for computation.");
249 if (c.size() != pc)
250 c.resize(pc);
251
252 value_type* cc = c.coeff();
253 const value_type* xc = x.coeff();
254
255 if (p > 1 && xp > 1) {
256 // Copy c coefficients into temporary array
257 value_type* tc = Stokhos::ds_array<value_type>::get_and_fill(cc,p);
258
259 typename Cijk_type::i_iterator i_begin = Cijk->i_begin();
260 typename Cijk_type::i_iterator i_end = Cijk->i_end();
261 if (pc < Cijk->num_i())
262 i_end = Cijk->find_i(pc);
263 ordinal_type k_lim = p;
264 ordinal_type j_lim = xp;
265 const value_type* kc = tc;
266 const value_type* jc = xc;
267 if (xp < p) {
268 k_lim = xp;
269 j_lim = p;
270 kc = xc;
271 jc = tc;
272 }
273
274 value_type tmp, cijk;
275 ordinal_type i,j,k;
276 for (typename Cijk_type::i_iterator i_it=i_begin; i_it!=i_end; ++i_it) {
277 i = index(i_it);
278 tmp = value_type(0.0);
279 for (typename Cijk_type::ik_iterator k_it = Cijk->k_begin(i_it);
280 k_it != Cijk->k_end(i_it); ++k_it) {
281 k = index(k_it);
282 if (k < k_lim) {
283 for (typename Cijk_type::ikj_iterator j_it = Cijk->j_begin(k_it);
284 j_it != Cijk->j_end(k_it); ++j_it) {
285 j = index(j_it);
286 cijk = value(j_it);
287 if (j < j_lim)
288 tmp += cijk*kc[k]*jc[j];
289 }
290 }
291 }
292 cc[i] = tmp / basis->norm_squared(i);
293 }
294
295 Stokhos::ds_array<value_type>::destroy_and_release(tc, p);
296 }
297 else if (p > 1) {
298 for (ordinal_type i=0; i<p; i++)
299 cc[i] *= xc[0];
300 }
301 else if (xp > 1) {
302 for (ordinal_type i=1; i<xp; i++)
303 cc[i] = cc[0]*xc[i];
304 cc[0] *= xc[0];
305 }
306 else {
307 cc[0] *= xc[0];
308 }
309}
310
311template <typename ordinal_type, typename value_type, typename node_type>
312void
320
321template <typename ordinal_type, typename value_type, typename node_type>
322void
323Stokhos::OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::
324plus(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
325 const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a,
326 const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& b)
327{
328#ifdef STOKHOS_TEUCHOS_TIME_MONITOR
329 TEUCHOS_FUNC_TIME_MONITOR("Stokhos::OrthogPolyExpansionBase::plus(OPA,OPA)");
330#endif
331 ordinal_type pa = a.size();
332 ordinal_type pb = b.size();
333 ordinal_type pc = pa > pb ? pa : pb;
334 if (c.size() != pc)
335 c.resize(pc);
336
337 const value_type* ca = a.coeff();
338 const value_type* cb = b.coeff();
339 value_type* cc = c.coeff();
340
341 if (pa > pb) {
342 for (ordinal_type i=0; i<pb; i++)
343 cc[i] = ca[i] + cb[i];
344 for (ordinal_type i=pb; i<pc; i++)
345 cc[i] = ca[i];
346 }
347 else {
348 for (ordinal_type i=0; i<pa; i++)
349 cc[i] = ca[i] + cb[i];
350 for (ordinal_type i=pa; i<pc; i++)
351 cc[i] = cb[i];
352 }
353}
354
355template <typename ordinal_type, typename value_type, typename node_type>
356void
357Stokhos::OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::
358plus(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
359 const value_type& a,
360 const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& b)
361{
362#ifdef STOKHOS_TEUCHOS_TIME_MONITOR
363 TEUCHOS_FUNC_TIME_MONITOR("Stokhos::OrthogPolyExpansionBase::plus(const,OPA)");
364#endif
365 ordinal_type pc = b.size();
366 if (c.size() != pc)
367 c.resize(pc);
368
369 const value_type* cb = b.coeff();
370 value_type* cc = c.coeff();
371
372 cc[0] = a + cb[0];
373 for (ordinal_type i=1; i<pc; i++)
374 cc[i] = cb[i];
375}
376
377template <typename ordinal_type, typename value_type, typename node_type>
378void
379Stokhos::OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::
380plus(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
381 const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a,
382 const value_type& b)
383{
384#ifdef STOKHOS_TEUCHOS_TIME_MONITOR
385 TEUCHOS_FUNC_TIME_MONITOR("Stokhos::OrthogPolyExpansionBase::plus(OPA,const)");
386#endif
387 ordinal_type pc = a.size();
388 if (c.size() != pc)
389 c.resize(pc);
390
391 const value_type* ca = a.coeff();
392 value_type* cc = c.coeff();
393
394 cc[0] = ca[0] + b;
395 for (ordinal_type i=1; i<pc; i++)
396 cc[i] = ca[i];
397}
398
399template <typename ordinal_type, typename value_type, typename node_type>
400void
401Stokhos::OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::
402minus(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
403 const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a,
404 const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& b)
405{
406#ifdef STOKHOS_TEUCHOS_TIME_MONITOR
407 TEUCHOS_FUNC_TIME_MONITOR("Stokhos::OrthogPolyExpansionBase::minus(OPA,OPA)");
408#endif
409 ordinal_type pa = a.size();
410 ordinal_type pb = b.size();
411 ordinal_type pc = pa > pb ? pa : pb;
412 if (c.size() != pc)
413 c.resize(pc);
414
415 const value_type* ca = a.coeff();
416 const value_type* cb = b.coeff();
417 value_type* cc = c.coeff();
418
419 if (pa > pb) {
420 for (ordinal_type i=0; i<pb; i++)
421 cc[i] = ca[i] - cb[i];
422 for (ordinal_type i=pb; i<pc; i++)
423 cc[i] = ca[i];
424 }
425 else {
426 for (ordinal_type i=0; i<pa; i++)
427 cc[i] = ca[i] - cb[i];
428 for (ordinal_type i=pa; i<pc; i++)
429 cc[i] = -cb[i];
430 }
431}
432
433template <typename ordinal_type, typename value_type, typename node_type>
434void
435Stokhos::OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::
436minus(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
437 const value_type& a,
438 const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& b)
439{
440#ifdef STOKHOS_TEUCHOS_TIME_MONITOR
441 TEUCHOS_FUNC_TIME_MONITOR("Stokhos::OrthogPolyExpansionBase::minus(const,OPA)");
442#endif
443 ordinal_type pc = b.size();
444 if (c.size() != pc)
445 c.resize(pc);
446
447 const value_type* cb = b.coeff();
448 value_type* cc = c.coeff();
449
450 cc[0] = a - cb[0];
451 for (ordinal_type i=1; i<pc; i++)
452 cc[i] = -cb[i];
453}
454
455template <typename ordinal_type, typename value_type, typename node_type>
456void
457Stokhos::OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::
458minus(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
459 const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a,
460 const value_type& b)
461{
462#ifdef STOKHOS_TEUCHOS_TIME_MONITOR
463 TEUCHOS_FUNC_TIME_MONITOR("Stokhos::OrthogPolyExpansionBase::minus(OPA,const)");
464#endif
465 ordinal_type pc = a.size();
466 if (c.size() != pc)
467 c.resize(pc);
468
469 const value_type* ca = a.coeff();
470 value_type* cc = c.coeff();
471
472 cc[0] = ca[0] - b;
473 for (ordinal_type i=1; i<pc; i++)
474 cc[i] = ca[i];
475}
476
477template <typename ordinal_type, typename value_type, typename node_type>
478void
479Stokhos::OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::
480times(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
481 const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a,
482 const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& b)
483{
484#ifdef STOKHOS_TEUCHOS_TIME_MONITOR
485 TEUCHOS_FUNC_TIME_MONITOR("Stokhos::OrthogPolyExpansionBase::times(OPA,OPA)");
486#endif
487 ordinal_type pa = a.size();
488 ordinal_type pb = b.size();
489 ordinal_type pc;
490 if (pa > 1 && pb > 1)
491 pc = sz;
492 else
493 pc = pa*pb;
494 TEUCHOS_TEST_FOR_EXCEPTION(sz < pc, std::logic_error,
495 "Stokhos::OrthogPolyExpansionBase::times()" <<
496 ": Expansion size (" << sz <<
497 ") is too small for computation.");
498 if (c.size() != pc)
499 c.resize(pc);
500
501 const value_type* ca = a.coeff();
502 const value_type* cb = b.coeff();
503 value_type* cc = c.coeff();
504
505 if (pa > 1 && pb > 1) {
506 typename Cijk_type::i_iterator i_begin = Cijk->i_begin();
507 typename Cijk_type::i_iterator i_end = Cijk->i_end();
508 if (pc < Cijk->num_i())
509 i_end = Cijk->find_i(pc);
510 ordinal_type k_lim = pa;
511 ordinal_type j_lim = pb;
512 const value_type* kc = ca;
513 const value_type* jc = cb;
514 if (pb < pa) {
515 k_lim = pb;
516 j_lim = pa;
517 kc = cb;
518 jc = ca;
519 }
520 value_type tmp, cijk;
521 ordinal_type i,j,k;
522 for (typename Cijk_type::i_iterator i_it=i_begin; i_it!=i_end; ++i_it) {
523 i = index(i_it);
524 tmp = value_type(0.0);
525 for (typename Cijk_type::ik_iterator k_it = Cijk->k_begin(i_it);
526 k_it != Cijk->k_end(i_it); ++k_it) {
527 k = index(k_it);
528 if (k < k_lim) {
529 for (typename Cijk_type::ikj_iterator j_it = Cijk->j_begin(k_it);
530 j_it != Cijk->j_end(k_it); ++j_it) {
531 j = index(j_it);
532 cijk = value(j_it);
533 if (j < j_lim)
534 tmp += cijk*kc[k]*jc[j];
535 }
536 }
537 }
538 cc[i] = tmp / basis->norm_squared(i);
539 }
540 }
541 else if (pa > 1) {
542 for (ordinal_type i=0; i<pc; i++)
543 cc[i] = ca[i]*cb[0];
544 }
545 else if (pb > 1) {
546 for (ordinal_type i=0; i<pc; i++)
547 cc[i] = ca[0]*cb[i];
548 }
549 else {
550 cc[0] = ca[0]*cb[0];
551 }
552}
553
554template <typename ordinal_type, typename value_type, typename node_type>
555void
556Stokhos::OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::
557times(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
558 const value_type& a,
559 const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& b)
560{
561#ifdef STOKHOS_TEUCHOS_TIME_MONITOR
562 TEUCHOS_FUNC_TIME_MONITOR("Stokhos::OrthogPolyExpansionBase::times(const,OPA)");
563#endif
564 ordinal_type pc = b.size();
565 if (c.size() != pc)
566 c.resize(pc);
567
568 const value_type* cb = b.coeff();
569 value_type* cc = c.coeff();
570
571 for (ordinal_type i=0; i<pc; i++)
572 cc[i] = a*cb[i];
573}
574
575template <typename ordinal_type, typename value_type, typename node_type>
576void
577Stokhos::OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::
578times(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
579 const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a,
580 const value_type& b)
581{
582#ifdef STOKHOS_TEUCHOS_TIME_MONITOR
583 TEUCHOS_FUNC_TIME_MONITOR("Stokhos::OrthogPolyExpansionBase::times(OPA,const)");
584#endif
585 ordinal_type pc = a.size();
586 if (c.size() != pc)
587 c.resize(pc);
588
589 const value_type* ca = a.coeff();
590 value_type* cc = c.coeff();
591
592 for (ordinal_type i=0; i<pc; i++)
593 cc[i] = ca[i]*b;
594}
595
596template <typename ordinal_type, typename value_type, typename node_type>
597void
605
606template <typename ordinal_type, typename value_type, typename node_type>
607void
616
617template <typename ordinal_type, typename value_type, typename node_type>
618void
619Stokhos::OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::
620divide(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
621 const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a,
622 const value_type& b)
623{
624#ifdef STOKHOS_TEUCHOS_TIME_MONITOR
625 TEUCHOS_FUNC_TIME_MONITOR("Stokhos::OrthogPolyExpansionBase::divide(OPA,const)");
626#endif
627 ordinal_type pc = a.size();
628 if (c.size() != pc)
629 c.resize(pc);
630
631 const value_type* ca = a.coeff();
632 value_type* cc = c.coeff();
633
634 for (ordinal_type i=0; i<pc; i++)
635 cc[i] = ca[i]/b;
636}
637
638template <typename ordinal_type, typename value_type, typename node_type>
639void
640Stokhos::OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::
641fabs(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
642 const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a)
643{
644#ifdef STOKHOS_TEUCHOS_TIME_MONITOR
645 TEUCHOS_FUNC_TIME_MONITOR("Stokhos::OrthogPolyExpansionBase::fabs(OPA)");
646#endif
647 c.init(0.0);
648 c[0] = a.two_norm();
649}
650
651template <typename ordinal_type, typename value_type, typename node_type>
652void
653Stokhos::OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::
654abs(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
655 const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a)
656{
657#ifdef STOKHOS_TEUCHOS_TIME_MONITOR
658 TEUCHOS_FUNC_TIME_MONITOR("Stokhos::OrthogPolyExpansionBase::abs(OPA)");
659#endif
660 c.init(0.0);
661 c[0] = a.two_norm();
662}
663
664template <typename ordinal_type, typename value_type, typename node_type>
665void
666Stokhos::OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::
667max(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
668 const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a,
669 const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& b)
670{
671#ifdef STOKHOS_TEUCHOS_TIME_MONITOR
672 TEUCHOS_FUNC_TIME_MONITOR("Stokhos::OrthogPolyExpansionBase::max(OPA,OPA)");
673#endif
674 if (a.two_norm() >= b.two_norm())
675 c = a;
676 else
677 c = b;
678}
679
680template <typename ordinal_type, typename value_type, typename node_type>
681void
682Stokhos::OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::
683max(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
684 const value_type& a,
685 const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& b)
686{
687#ifdef STOKHOS_TEUCHOS_TIME_MONITOR
688 TEUCHOS_FUNC_TIME_MONITOR("Stokhos::OrthogPolyExpansionBase::max(const,OPA)");
689#endif
690 if (a >= b.two_norm()) {
691 c = OrthogPolyApprox<ordinal_type, value_type, node_type>(basis);
692 c[0] = a;
693 }
694 else
695 c = b;
696}
697
698template <typename ordinal_type, typename value_type, typename node_type>
699void
700Stokhos::OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::
701max(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
702 const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a,
703 const value_type& b)
704{
705#ifdef STOKHOS_TEUCHOS_TIME_MONITOR
706 TEUCHOS_FUNC_TIME_MONITOR("Stokhos::OrthogPolyExpansionBase::max(OPA,const)");
707#endif
708 if (a.two_norm() >= b)
709 c = a;
710 else {
711 c = OrthogPolyApprox<ordinal_type, value_type, node_type>(basis);
712 c[0] = b;
713 }
714}
715
716template <typename ordinal_type, typename value_type, typename node_type>
717void
718Stokhos::OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::
719min(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
720 const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a,
721 const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& b)
722{
723#ifdef STOKHOS_TEUCHOS_TIME_MONITOR
724 TEUCHOS_FUNC_TIME_MONITOR("Stokhos::OrthogPolyExpansionBase::min(OPA,OPA)");
725#endif
726 if (a.two_norm() <= b.two_norm())
727 c = a;
728 else
729 c = b;
730}
731
732template <typename ordinal_type, typename value_type, typename node_type>
733void
734Stokhos::OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::
735min(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
736 const value_type& a,
737 const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& b)
738{
739#ifdef STOKHOS_TEUCHOS_TIME_MONITOR
740 TEUCHOS_FUNC_TIME_MONITOR("Stokhos::OrthogPolyExpansionBase::min(const,OPA)");
741#endif
742 if (a <= b.two_norm()) {
743 c = OrthogPolyApprox<ordinal_type, value_type, node_type>(basis);
744 c[0] = a;
745 }
746 else
747 c = b;
748}
749
750template <typename ordinal_type, typename value_type, typename node_type>
751void
752Stokhos::OrthogPolyExpansionBase<ordinal_type, value_type, node_type>::
753min(Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& c,
754 const Stokhos::OrthogPolyApprox<ordinal_type, value_type, node_type>& a,
755 const value_type& b)
756{
757#ifdef STOKHOS_TEUCHOS_TIME_MONITOR
758 TEUCHOS_FUNC_TIME_MONITOR("Stokhos::OrthogPolyExpansionBase::min(OPA,const)");
759#endif
760 if (a.two_norm() <= b)
761 c = a;
762 else {
763 c = OrthogPolyApprox<ordinal_type, value_type, node_type>(basis);
764 c[0] = b;
765 }
766}
j
j
Definition Sacado_Fad_Exp_MP_Vector.hpp:527
val
expr val()
Stokhos_CGDivisionExpansionStrategy.hpp
Stokhos_DenseDirectDivisionExpansionStrategy.hpp
Stokhos_DynamicArrayTraits.hpp
Stokhos_GMRESDivisionExpansionStrategy.hpp
Stokhos_MeanBasedDivisionExpansionStrategy.hpp
Stokhos_SPDDenseDirectDivisionExpansionStrategy.hpp
Stokhos::CGDivisionExpansionStrategy
Strategy interface for computing PCE of a/b using only b[0].
Definition Stokhos_CGDivisionExpansionStrategy.hpp:73
Stokhos::DenseDirectDivisionExpansionStrategy
Strategy interface for computing PCE of a/b using only b[0].
Definition Stokhos_DenseDirectDivisionExpansionStrategy.hpp:64
Stokhos::GMRESDivisionExpansionStrategy
Strategy interface for computing PCE of a/b using only b[0].
Definition Stokhos_GMRESDivisionExpansionStrategy.hpp:70
Stokhos::MeanBasedDivisionExpansionStrategy
Strategy interface for computing PCE of a/b using only b[0].
Definition Stokhos_MeanBasedDivisionExpansionStrategy.hpp:58
Stokhos::OrthogPolyApprox
Class to store coefficients of a projection onto an orthogonal polynomial basis.
Definition Stokhos_OrthogPolyApprox.hpp:63
Stokhos::OrthogPolyApprox::init
void init(const value_type &v)
Initialize coefficients to value.
Definition Stokhos_OrthogPolyApproxImp.hpp:106
Stokhos::OrthogPolyApprox::two_norm
value_type two_norm() const
Compute the two-norm of expansion.
Definition Stokhos_OrthogPolyApproxImp.hpp:292
Stokhos::OrthogPolyApprox::resize
void resize(ordinal_type sz)
Resize coefficient array (coefficients are preserved)
Definition Stokhos_OrthogPolyApproxImp.hpp:152
Stokhos::OrthogPolyApprox::coeff
pointer coeff()
Return coefficient array.
Definition Stokhos_OrthogPolyApproxImp.hpp:168
Stokhos::OrthogPolyApprox::size
ordinal_type size() const
Return size.
Definition Stokhos_OrthogPolyApproxImp.hpp:160
Stokhos::OrthogPolyApprox::basis
Teuchos::RCP< const Stokhos::OrthogPolyBasis< ordinal_type, value_type > > basis() const
Return basis.
Definition Stokhos_OrthogPolyApproxImp.hpp:130
Stokhos::OrthogPolyBasis
Abstract base class for multivariate orthogonal polynomials.
Definition Stokhos_OrthogPolyBasis.hpp:74
Stokhos::Cijk_type::Cijk
Teuchos::RCP< const Stokhos::Sparse3Tensor< ordinal_type, value_type > > Cijk
Definition Stokhos_OrthogPolyExpansionBase.hpp:192
Stokhos::OrthogPolyExpansionBase::minusEqual
void minusEqual(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const value_type &x)
Definition Stokhos_OrthogPolyExpansionBaseImp.hpp:152
Stokhos::OrthogPolyExpansionBase::max
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)
Definition Stokhos_OrthogPolyExpansionBaseImp.hpp:667
Stokhos::OrthogPolyExpansionBase::OrthogPolyExpansionBase
OrthogPolyExpansionBase(const Teuchos::RCP< const OrthogPolyBasis< ordinal_type, value_type > > &basis, const Teuchos::RCP< const Stokhos::Sparse3Tensor< ordinal_type, value_type > > &Cijk, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Constructor.
Stokhos::OrthogPolyExpansionBase::fabs
void fabs(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
Definition Stokhos_OrthogPolyExpansionBaseImp.hpp:641
Stokhos::OrthogPolyExpansionBase::plus
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)
Definition Stokhos_OrthogPolyExpansionBaseImp.hpp:324
Stokhos::OrthogPolyExpansionBase::timesEqual
void timesEqual(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const value_type &x)
Definition Stokhos_OrthogPolyExpansionBaseImp.hpp:161
Stokhos::Cijk_type::params
Teuchos::RCP< Teuchos::ParameterList > params
Definition Stokhos_OrthogPolyExpansionBase.hpp:195
Stokhos::OrthogPolyExpansionBase::unaryMinus
void unaryMinus(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
Definition Stokhos_OrthogPolyExpansionBaseImp.hpp:121
Stokhos::OrthogPolyExpansionBase::divideEqual
void divideEqual(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const value_type &x)
Definition Stokhos_OrthogPolyExpansionBaseImp.hpp:176
Stokhos::OrthogPolyExpansionBase::divide
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)
Definition Stokhos_OrthogPolyExpansionBaseImp.hpp:599
Stokhos::OrthogPolyExpansionBase::minus
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)
Definition Stokhos_OrthogPolyExpansionBaseImp.hpp:402
Stokhos::Cijk_type::basis
Teuchos::RCP< const OrthogPolyBasis< ordinal_type, value_type > > basis
Definition Stokhos_OrthogPolyExpansionBase.hpp:186
Stokhos::Cijk_type::division_strategy
Teuchos::RCP< Stokhos::DivisionExpansionStrategy< ordinal_type, value_type, node_type > > division_strategy
Definition Stokhos_OrthogPolyExpansionBase.hpp:198
Stokhos::OrthogPolyExpansionBase::min
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)
Definition Stokhos_OrthogPolyExpansionBaseImp.hpp:719
Stokhos::OrthogPolyExpansionBase::times
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)
Definition Stokhos_OrthogPolyExpansionBaseImp.hpp:480
Stokhos::OrthogPolyExpansionBase::plusEqual
void plusEqual(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const value_type &x)
Definition Stokhos_OrthogPolyExpansionBaseImp.hpp:143
Stokhos::OrthogPolyExpansionBase::abs
void abs(OrthogPolyApprox< ordinal_type, value_type, node_type > &c, const OrthogPolyApprox< ordinal_type, value_type, node_type > &a)
Definition Stokhos_OrthogPolyExpansionBaseImp.hpp:654
Stokhos::SPDDenseDirectDivisionExpansionStrategy
Strategy interface for computing PCE of a/b using only b[0].
Definition Stokhos_SPDDenseDirectDivisionExpansionStrategy.hpp:64
Stokhos::Sparse3Tensor
Data structure storing a sparse 3-tensor C(i,j,k) in a a compressed format.
Definition Stokhos_Sparse3Tensor.hpp:56
Stokhos::Sparse3Tensor::ikj_iterator
j_sparse_array::const_iterator ikj_iterator
Iterator for looping over j entries given i and k.
Definition Stokhos_Sparse3Tensor.hpp:100
Stokhos::Sparse3Tensor::ik_iterator
kj_sparse_array::const_iterator ik_iterator
Iterator for looping over k entries given i.
Definition Stokhos_Sparse3Tensor.hpp:97
Stokhos::Sparse3Tensor::i_iterator
ikj_sparse_array::const_iterator i_iterator
Iterator for looping over i entries.
Definition Stokhos_Sparse3Tensor.hpp:91
Stokhos::Sparse3Tensor::value
SparseArrayIterator< index_iterator, value_iterator >::value_reference value(const SparseArrayIterator< index_iterator, value_iterator > &it)
Definition Stokhos_Sparse3Tensor.hpp:315
Stokhos::Sparse3Tensor::index
SparseArrayIterator< index_iterator, value_iterator >::value_type index(const SparseArrayIterator< index_iterator, value_iterator > &it)
Definition Stokhos_Sparse3Tensor.hpp:295
diag
@ diag
Definition experimental/linear2d_diffusion_pce.cpp:140
linear
@ linear
Definition experimental/linear2d_diffusion_pce.cpp:146
Stokhos::ds_array::destroy_and_release
static void destroy_and_release(T *m, int sz)
Destroy array elements and release memory.
Definition Stokhos_DynamicArrayTraits.hpp:115
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:82
Generated by doxygen 1.13.2