Stokhos Package Browser (Single Doxygen Collection) Version of the Day
Loading...
Searching...
No Matches
Stokhos_Sparse3TensorImp.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_ConfigDefs.h"
43#include "Teuchos_Assert.hpp"
44
45template <typename ordinal_type, typename value_type>
51
52template <typename ordinal_type, typename value_type>
53void
55add_term(ordinal_type i, ordinal_type j, ordinal_type k, const value_type& c)
56{
57#ifdef STOKHOS_DEBUG
58 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == true, std::logic_error,
59 "You can't call add_term() after calling fillComplete()!");
60#endif
61
62 kji_data[k][j][i] = c;
63 ikj_data[i][k][j] = c;
64}
65
66template <typename ordinal_type, typename value_type>
67void
69sum_term(ordinal_type i, ordinal_type j, ordinal_type k, const value_type& c)
70{
71#ifdef STOKHOS_DEBUG
72 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == true, std::logic_error,
73 "You can't call sum_term() after calling fillComplete()!");
74#endif
75
76 kji_data[k][j][i] += c;
77 ikj_data[i][k][j] += c;
78}
79
80template <typename ordinal_type, typename value_type>
81void
84{
85 if (fill_completed)
86 return;
87
88 ikj_array.resize(ikj_data.size());
89 ordinal_type n = 0;
90 for (typename ikj_map::const_iterator i_it = ikj_data.begin();
91 i_it != ikj_data.end(); ++i_it) {
92 ikj_array.indices[n] = i_it->first;
93 ikj_array.values[n].resize(i_it->second.size());
94 ordinal_type l = 0;
95 for (typename kj_map::const_iterator k_it = i_it->second.begin();
96 k_it != i_it->second.end(); ++k_it) {
97 ikj_array.values[n].indices[l] = k_it->first;
98 ikj_array.values[n].values[l].resize(k_it->second.size());
99 ordinal_type m = 0;
100 for (typename j_map::const_iterator j_it = k_it->second.begin();
101 j_it != k_it->second.end(); ++j_it) {
102 ikj_array.values[n].values[l].indices[m] = j_it->first;
103 ikj_array.values[n].values[l].values[m] = j_it->second;
104 m++;
105 }
106 l++;
107 }
108 n++;
110
111 kji_array.resize(kji_data.size());
112 n = 0;
113 for (typename kji_map::const_iterator k_it = kji_data.begin();
114 k_it != kji_data.end(); ++k_it) {
115 kji_array.indices[n] = k_it->first;
116 kji_array.values[n].resize(k_it->second.size());
117 ordinal_type l = 0;
118 for (typename ji_map::const_iterator j_it = k_it->second.begin();
119 j_it != k_it->second.end(); ++j_it) {
120 kji_array.values[n].indices[l] = j_it->first;
121 kji_array.values[n].values[l].resize(j_it->second.size());
122 ordinal_type m = 0;
123 for (typename i_map::const_iterator i_it = j_it->second.begin();
124 i_it != j_it->second.end(); ++i_it) {
125 kji_array.values[n].values[l].indices[m] = i_it->first;
126 kji_array.values[n].values[l].values[m] = i_it->second;
127 m++;
129 l++;
130 }
131 n++;
132 }
133
134 ikj_data.clear();
135 kji_data.clear();
136
137 fill_completed = true;
138}
139
140template <typename ordinal_type, typename value_type>
141bool
144{
145 return fill_completed;
147
148template <typename ordinal_type, typename value_type>
149void
151print(std::ostream& os) const
153#ifdef STOKHOS_DEBUG
154 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
155 "You must call fillComplete() before calling print()!");
156#endif
157
158 for (k_iterator k=k_begin(); k!=k_end(); ++k)
159 for (kj_iterator j=j_begin(k); j!=j_end(k); ++j)
160 for (kji_iterator i=i_begin(j); i!=i_end(j); ++i)
161 os << "k = " << index(k)
162 << ", j = " << index(j)
163 << ", i = " << index(i)
164 << ", Cijk = " << value(i) << std::endl;
165}
166
167template <typename ordinal_type, typename value_type>
168value_type
170getValue(ordinal_type i, ordinal_type j, ordinal_type k) const
171{
172#ifdef STOKHOS_DEBUG
173 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
174 "You must call fillComplete() before calling getValue()!");
175#endif
177 k_iterator k_it = find_k(k);
178 if (k_it == k_end())
179 return value_type(0);
180
181 kj_iterator j_it = find_j(k_it, j);
182 if (j_it == j_end(k_it))
183 return value_type(0);
184
185 kji_iterator i_it = find_i(j_it, i);
186 if (i_it == i_end(j_it))
187 return value_type(0);
188
189 return i_it.value();
190}
191
192template <typename ordinal_type, typename value_type>
193ordinal_type
195num_entries() const
197#ifdef STOKHOS_DEBUG
198 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
199 "You must call fillComplete() before calling num_entries()!");
200#endif
201
202 ordinal_type num = 0;
203 for (k_iterator k = k_begin(); k != k_end(); ++k)
204 for (kj_iterator j = j_begin(k); j != j_end(k); ++j)
205 for (kji_iterator i = i_begin(j); i != i_end(j); ++i)
206 ++num;
207
208 return num;
209}
210
211template <typename ordinal_type, typename value_type>
212ordinal_type
214num_k() const
215{
216#ifdef STOKHOS_DEBUG
217 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
218 "You must call fillComplete() before calling num_k()!");
219#endif
220 return kji_array.size();
221}
222
223template <typename ordinal_type, typename value_type>
224ordinal_type
226num_j(const k_iterator& k) const
227{
228#ifdef STOKHOS_DEBUG
229 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
230 "You must call fillComplete() before calling num_j()!");
231#endif
233 return k.value().size();
234}
236template <typename ordinal_type, typename value_type>
237ordinal_type
239num_i(const kj_iterator& j) const
240{
241#ifdef STOKHOS_DEBUG
242 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
243 "You must call fillComplete() before calling num_i()!");
244#endif
245 return j.value().size();
246}
247
248template <typename ordinal_type, typename value_type>
251find_k(ordinal_type k) const
252{
253#ifdef STOKHOS_DEBUG
254 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
255 "You must call fillComplete() before calling find_k()!");
256#endif
257 return kji_array.find(k);
258}
259
260template <typename ordinal_type, typename value_type>
263find_j(const k_iterator& k, ordinal_type j) const
264{
265#ifdef STOKHOS_DEBUG
266 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
267 "You must call fillComplete() before calling find_j()!");
268#endif
269 return k.value().find(j);
270}
271
272template <typename ordinal_type, typename value_type>
275find_i(const kj_iterator& j, ordinal_type i) const
276{
277#ifdef STOKHOS_DEBUG
278 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
279 "You must call fillComplete() before calling find_i()!");
280#endif
281 return j.value().find(i);
282}
283
284template <typename ordinal_type, typename value_type>
287k_begin() const
288{
289#ifdef STOKHOS_DEBUG
290 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
291 "You must call fillComplete() before calling k_begin()!");
292#endif
293 return kji_array.begin();
294}
295
296template <typename ordinal_type, typename value_type>
299k_end() const
300{
301#ifdef STOKHOS_DEBUG
302 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
303 "You must call fillComplete() before calling k_end()!");
304#endif
305 return kji_array.end();
306}
307
308template <typename ordinal_type, typename value_type>
311k_rbegin() const
312{
313#ifdef STOKHOS_DEBUG
314 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
315 "You must call fillComplete() before calling k_rbegin()!");
316#endif
317 return kji_array.rbegin();
318}
319
320template <typename ordinal_type, typename value_type>
323k_rend() const
324{
325#ifdef STOKHOS_DEBUG
326 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
327 "You must call fillComplete() before calling k_rend()!");
328#endif
329 return kji_array.rend();
330}
331
332template <typename ordinal_type, typename value_type>
335j_begin(const k_iterator& k) const
336{
337#ifdef STOKHOS_DEBUG
338 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
339 "You must call fillComplete() before calling j_begin()!");
340#endif
341 return k.value().begin();
342}
343
344template <typename ordinal_type, typename value_type>
347j_end(const k_iterator& k) const
348{
349#ifdef STOKHOS_DEBUG
350 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
351 "You must call fillComplete() before calling j_end()!");
352#endif
353 return k.value().end();
354}
355
356template <typename ordinal_type, typename value_type>
359j_begin(const k_reverse_iterator& k) const
360{
361#ifdef STOKHOS_DEBUG
362 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
363 "You must call fillComplete() before calling j_begin()!");
364#endif
365 return k.value().begin();
366}
367
368template <typename ordinal_type, typename value_type>
371j_end(const k_reverse_iterator& k) const
372{
373#ifdef STOKHOS_DEBUG
374 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
375 "You must call fillComplete() before calling j_end()!");
376#endif
377 return k.value().end();
378}
379
380template <typename ordinal_type, typename value_type>
383i_begin(const kj_iterator& j) const
384{
385#ifdef STOKHOS_DEBUG
386 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
387 "You must call fillComplete() before calling i_begin()!");
388#endif
389 return j.value().begin();
390}
391
392template <typename ordinal_type, typename value_type>
395i_end(const kj_iterator& j) const
396{
397#ifdef STOKHOS_DEBUG
398 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
399 "You must call fillComplete() before calling i_end()!");
400#endif
401 return j.value().end();
402}
403
404template <typename ordinal_type, typename value_type>
405ordinal_type
407num_i() const
408{
409#ifdef STOKHOS_DEBUG
410 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
411 "You must call fillComplete() before calling num_i()!");
412#endif
413 return ikj_array.size();
414}
415
416template <typename ordinal_type, typename value_type>
417ordinal_type
419num_k(const i_iterator& i) const
420{
421#ifdef STOKHOS_DEBUG
422 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
423 "You must call fillComplete() before calling num_k()!");
424#endif
425 return i.value().size();
426}
427
428template <typename ordinal_type, typename value_type>
429ordinal_type
431num_j(const ik_iterator& k) const
432{
433#ifdef STOKHOS_DEBUG
434 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
435 "You must call fillComplete() before calling num_j()!");
436#endif
437 return k.value().size();
438}
439
440template <typename ordinal_type, typename value_type>
443find_i(ordinal_type i) const
444{
445#ifdef STOKHOS_DEBUG
446 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
447 "You must call fillComplete() before calling find_i()!");
448#endif
449 return ikj_array.find(i);
450}
451
452template <typename ordinal_type, typename value_type>
455find_k(const i_iterator& i, ordinal_type k) const
456{
457#ifdef STOKHOS_DEBUG
458 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
459 "You must call fillComplete() before calling find_k()!");
460#endif
461 return i.value().find(k);
462}
463
464template <typename ordinal_type, typename value_type>
467find_j(const ik_iterator& k, ordinal_type j) const
468{
469#ifdef STOKHOS_DEBUG
470 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
471 "You must call fillComplete() before calling find_j()!");
472#endif
473 return k.value().find(j);
474}
475
476template <typename ordinal_type, typename value_type>
479i_begin() const
480{
481#ifdef STOKHOS_DEBUG
482 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
483 "You must call fillComplete() before calling i_begin()!");
484#endif
485 return ikj_array.begin();
486}
487
488template <typename ordinal_type, typename value_type>
491i_end() const
492{
493#ifdef STOKHOS_DEBUG
494 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
495 "You must call fillComplete() before calling i_end()!");
496#endif
497 return ikj_array.end();
498}
499
500template <typename ordinal_type, typename value_type>
503i_rbegin() const
504{
505#ifdef STOKHOS_DEBUG
506 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
507 "You must call fillComplete() before calling i_rbegin()!");
508#endif
509 return ikj_array.rbegin();
510}
511
512template <typename ordinal_type, typename value_type>
515i_rend() const
516{
517#ifdef STOKHOS_DEBUG
518 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
519 "You must call fillComplete() before calling i_rend()!");
520#endif
521 return ikj_array.rend();
522}
523
524template <typename ordinal_type, typename value_type>
527k_begin(const i_iterator& i) const
528{
529#ifdef STOKHOS_DEBUG
530 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
531 "You must call fillComplete() before calling k_begin()!");
532#endif
533 return i.value().begin();
534}
535
536template <typename ordinal_type, typename value_type>
539k_end(const i_iterator& i) const
540{
541#ifdef STOKHOS_DEBUG
542 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
543 "You must call fillComplete() before calling k_end()!");
544#endif
545 return i.value().end();
546}
547
548template <typename ordinal_type, typename value_type>
551k_begin(const i_reverse_iterator& i) const
552{
553#ifdef STOKHOS_DEBUG
554 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
555 "You must call fillComplete() before calling k_begin()!");
556#endif
557 return i.value().begin();
558}
559
560template <typename ordinal_type, typename value_type>
563k_end(const i_reverse_iterator& i) const
564{
565#ifdef STOKHOS_DEBUG
566 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
567 "You must call fillComplete() before calling k_end()!");
568#endif
569 return i.value().end();
570}
571
572template <typename ordinal_type, typename value_type>
575j_begin(const ik_iterator& k) const
576{
577#ifdef STOKHOS_DEBUG
578 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
579 "You must call fillComplete() before calling j_begin()!");
580#endif
581 return k.value().begin();
582}
583
584template <typename ordinal_type, typename value_type>
587j_end(const ik_iterator& k) const
588{
589#ifdef STOKHOS_DEBUG
590 TEUCHOS_TEST_FOR_EXCEPTION(fill_completed == false, std::logic_error,
591 "You must call fillComplete() before calling j_end()!");
592#endif
593 return k.value().end();
594}
ordinal_type num_j(const k_iterator &k) const
Number of j entries in C(i,j,k) for given k.
i_iterator i_end() const
Iterator pointing to last k entry.
void sum_term(ordinal_type i, ordinal_type j, ordinal_type k, const value_type &c)
Add new term for given (i,j,k) and sum in if already there.
ikj_sparse_array::const_reverse_iterator i_reverse_iterator
Iterator for looping over i entries in reverse.
k_reverse_iterator k_rbegin() const
Reverse iterator pointing to last k entry.
i_reverse_iterator i_rbegin() const
Reverse iterator pointing to last k entry.
k_reverse_iterator k_rend() const
Reverse iterator pointing to first k entry.
j_sparse_array::const_iterator ikj_iterator
Iterator for looping over j entries given i and k.
kj_iterator j_end(const k_iterator &k) const
Iterator pointing to last j entry for given k.
void print(std::ostream &os) const
Print tensor.
ordinal_type num_i() const
Number of i entries in C(i,j,k)
k_iterator k_begin() const
Iterator pointing to first k entry.
kji_sparse_array::const_iterator k_iterator
Iterator for looping over k entries.
value_type getValue(ordinal_type i, ordinal_type j, ordinal_type k) const
Get Cijk value for a given i, j, k indices.
i_reverse_iterator i_rend() const
Reverse iterator pointing to first k entry.
ikj_sparse_array::const_iterator i_iterator
Iterator for looping over i entries.
void fillComplete()
Signal all terms have been added.
kj_iterator j_begin(const k_iterator &k) const
Iterator pointing to first j entry for given k.
k_iterator find_k(ordinal_type k) const
Return k iterator for given index k.
ordinal_type num_k() const
Number of k entries in C(i,j,k)
i_iterator i_begin() const
Iterator pointing to first k entry.
kj_iterator find_j(const k_iterator &k, ordinal_type j) const
Return j iterator given k iterator and index j.
kji_sparse_array::const_reverse_iterator k_reverse_iterator
Iterator for looping over k entries in reverse.
ordinal_type num_entries() const
Return number of non-zero entries.
j_sparse_array::const_iterator kji_iterator
Iterator for looping over i entries given k and j.
ji_sparse_array::const_iterator kj_iterator
Iterator for looping over j entries given k.
kji_iterator find_i(const kj_iterator &j, ordinal_type i) const
Return i iterator given j iterator and index i.
k_iterator k_end() const
Iterator pointing to last k entry.
void add_term(ordinal_type i, ordinal_type j, ordinal_type k, const value_type &c)
Add new term for given (i,j,k)
bool fillCompleted() const
Return whether fillComplete() has been called.
Bi-directional iterator for traversing a sparse array.
value_reference value() const
Return value associated with iterator.