Stokhos Package Browser (Single Doxygen Collection)  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
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 
45 template <typename ordinal_type, typename value_type>
48  fill_completed(false)
49 {
50 }
51 
52 template <typename ordinal_type, typename value_type>
53 void
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 
66 template <typename ordinal_type, typename value_type>
67 void
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 
80 template <typename ordinal_type, typename value_type>
81 void
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++;
109  }
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++;
128  }
129  l++;
130  }
131  n++;
132  }
133 
134  ikj_data.clear();
135  kji_data.clear();
136 
137  fill_completed = true;
138 }
139 
140 template <typename ordinal_type, typename value_type>
141 bool
144 {
145  return fill_completed;
146 }
147 
148 template <typename ordinal_type, typename value_type>
149 void
151 print(std::ostream& os) const
152 {
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 
167 template <typename ordinal_type, typename value_type>
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
176 
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 
192 template <typename ordinal_type, typename value_type>
195 num_entries() const
196 {
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 
211 template <typename ordinal_type, typename value_type>
214 num_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 
223 template <typename ordinal_type, typename value_type>
226 num_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
232 
233  return k.value().size();
234 }
235 
236 template <typename ordinal_type, typename value_type>
239 num_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 
248 template <typename ordinal_type, typename value_type>
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 
260 template <typename ordinal_type, typename value_type>
263 find_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 
272 template <typename ordinal_type, typename value_type>
275 find_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 
284 template <typename ordinal_type, typename value_type>
287 k_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 
296 template <typename ordinal_type, typename value_type>
299 k_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 
308 template <typename ordinal_type, typename value_type>
311 k_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 
320 template <typename ordinal_type, typename value_type>
323 k_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 
332 template <typename ordinal_type, typename value_type>
335 j_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 
344 template <typename ordinal_type, typename value_type>
347 j_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 
356 template <typename ordinal_type, typename value_type>
359 j_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 
368 template <typename ordinal_type, typename value_type>
371 j_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 
380 template <typename ordinal_type, typename value_type>
383 i_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 
392 template <typename ordinal_type, typename value_type>
395 i_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 
404 template <typename ordinal_type, typename value_type>
407 num_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 
416 template <typename ordinal_type, typename value_type>
419 num_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 
428 template <typename ordinal_type, typename value_type>
431 num_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 
440 template <typename ordinal_type, typename value_type>
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 
452 template <typename ordinal_type, typename value_type>
455 find_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 
464 template <typename ordinal_type, typename value_type>
467 find_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 
476 template <typename ordinal_type, typename value_type>
479 i_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 
488 template <typename ordinal_type, typename value_type>
491 i_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 
500 template <typename ordinal_type, typename value_type>
503 i_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 
512 template <typename ordinal_type, typename value_type>
515 i_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 
524 template <typename ordinal_type, typename value_type>
527 k_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 
536 template <typename ordinal_type, typename value_type>
539 k_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 
548 template <typename ordinal_type, typename value_type>
551 k_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 
560 template <typename ordinal_type, typename value_type>
563 k_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 
572 template <typename ordinal_type, typename value_type>
575 j_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 
584 template <typename ordinal_type, typename value_type>
587 j_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 }
k_iterator k_begin() const
Iterator pointing to first k entry.
ordinal_type num_j(const k_iterator &k) const
Number of j entries in C(i,j,k) for given k.
ordinal_type num_k() const
Number of k entries in C(i,j,k)
kj_iterator j_begin(const k_iterator &k) const
Iterator pointing to first j entry for given k.
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
void fillComplete()
Signal all terms have been added.
i_reverse_iterator i_rbegin() const
Reverse iterator pointing to last k entry.
kj_iterator j_end(const k_iterator &k) const
Iterator pointing to last j entry for given k.
Bi-directional iterator for traversing a sparse array.
Bi-directional reverse iterator for traversing a sparse array.
void add_term(ordinal_type i, ordinal_type j, ordinal_type k, const value_type &c)
Add new term for given (i,j,k)
value_reference value() const
Return value associated with iterator.
i_reverse_iterator i_rend() const
Reverse iterator pointing to first k entry.
k_iterator k_end() const
Iterator pointing to last k entry.
k_iterator find_k(ordinal_type k) const
Return k iterator for given index k.
i_iterator i_begin() const
Iterator pointing to first k entry.
bool fillCompleted() const
Return whether fillComplete() has been called.
kji_iterator find_i(const kj_iterator &j, ordinal_type i) const
Return i iterator given j iterator and index i.
ordinal_type num_i() const
Number of i entries in C(i,j,k)
void print(std::ostream &os) const
Print tensor.
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.
k_reverse_iterator k_rbegin() const
Reverse iterator pointing to last k entry.
value_type getValue(ordinal_type i, ordinal_type j, ordinal_type k) const
Get Cijk value for a given i, j, k indices.
value_reference value() const
Return value associated with iterator.
k_reverse_iterator k_rend() const
Reverse iterator pointing to first k entry.
int n
i_iterator i_end() const
Iterator pointing to last k entry.
kj_iterator find_j(const k_iterator &k, ordinal_type j) const
Return j iterator given k iterator and index j.
ordinal_type num_entries() const
Return number of non-zero entries.