Xpetra Version of the Day
Loading...
Searching...
No Matches
Xpetra_TpetraBlockCrsMatrix_def.hpp
Go to the documentation of this file.
1// @HEADER
2//
3// ***********************************************************************
4//
5// Xpetra: A linear algebra interface package
6// Copyright 2012 Sandia Corporation
7//
8// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
9// the U.S. Government retains certain rights in this software.
10//
11// Redistribution and use in source and binary forms, with or without
12// modification, are permitted provided that the following conditions are
13// met:
14//
15// 1. Redistributions of source code must retain the above copyright
16// notice, this list of conditions and the following disclaimer.
17//
18// 2. Redistributions in binary form must reproduce the above copyright
19// notice, this list of conditions and the following disclaimer in the
20// documentation and/or other materials provided with the distribution.
21//
22// 3. Neither the name of the Corporation nor the names of the
23// contributors may be used to endorse or promote products derived from
24// this software without specific prior written permission.
25//
26// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
27// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
28// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
29// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
30// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
31// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
32// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
33// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
34// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
35// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
36// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
37//
38// Questions? Contact
39// Jonathan Hu (jhu@sandia.gov)
40// Andrey Prokopenko (aprokop@sandia.gov)
41// Ray Tuminaro (rstumin@sandia.gov)
42//
43// ***********************************************************************
44//
45// @HEADER
46#ifndef XPETRA_TPETRABLOCKCRSMATRIX_DEF_HPP
47#define XPETRA_TPETRABLOCKCRSMATRIX_DEF_HPP
48
49
51#include "Xpetra_TpetraCrsGraph.hpp"
52
53namespace Xpetra {
54
55
57 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
60 size_t maxNumEntriesPerRow,
62 {
63 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix function not implemented in"+std::string(__FILE__)+":"+std::to_string(__LINE__));
64 }
65
66
68 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
71 const ArrayRCP< const size_t > &NumEntriesPerRowToAlloc,
73 {
74 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix function not implemented in"+std::string(__FILE__)+":"+std::to_string(__LINE__));
75 }
76
77
79 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
83 size_t maxNumEntriesPerRow,
85 {
86 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix function not implemented in"+std::string(__FILE__)+":"+std::to_string(__LINE__));
87 }
88
89
91 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
95 const ArrayRCP< const size_t > &NumEntriesPerRowToAlloc,
97 {
98 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix function not implemented in"+std::string(__FILE__)+":"+std::to_string(__LINE__));
99 }
100
101
103 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
107 // : mtx_(Teuchos::rcp(new Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >(toTpetra(graph), params)))
108 // * there is no Tpetra::BlockCrsMatrix(graph, params) c'tor. We throw anyways here so no need to set mtx_.
109 {
110 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix function not implemented in"+std::string(__FILE__)+":"+std::to_string(__LINE__));
111 }
112
113
115 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
118 const LocalOrdinal blockSize)
119 : mtx_(Teuchos::rcp(new Tpetra::BlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>(*toTpetra(graph), blockSize)))
120 { }
121
122
124 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
127 const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& pointDomainMap,
128 const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& pointRangeMap,
129 const LocalOrdinal blockSize)
130 : mtx_(Teuchos::rcp(new Tpetra::BlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>(*toTpetra(graph), *toTpetra(pointDomainMap), *toTpetra(pointRangeMap),blockSize)))
131 { }
132
133
135 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
137 TpetraBlockCrsMatrix(const Teuchos::RCP<const Tpetra::BlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >& sourceMatrix,
139 const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& domainMap,
142 {
143 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix function not implemented in"+std::string(__FILE__)+":"+std::to_string(__LINE__));
144 }
145
146
148 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
150 TpetraBlockCrsMatrix(const Teuchos::RCP<const Tpetra::BlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >& sourceMatrix,
152 const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& domainMap,
155 {
156 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix function not implemented in"+std::string(__FILE__)+":"+std::to_string(__LINE__));
157 }
158
159
161 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
163 TpetraBlockCrsMatrix(const Teuchos::RCP<const Tpetra::BlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >& sourceMatrix,
164 const Import<LocalOrdinal,GlobalOrdinal,Node> & RowImporter,
165 const Teuchos::RCP<const Import<LocalOrdinal,GlobalOrdinal,Node> > DomainImporter,
166 const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& domainMap,
169 {
170 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix function not implemented in"+std::string(__FILE__)+":"+std::to_string(__LINE__));
171 }
172
173
175 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
177 TpetraBlockCrsMatrix(const Teuchos::RCP<const Tpetra::BlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >& sourceMatrix,
178 const Export<LocalOrdinal,GlobalOrdinal,Node> & RowExporter,
179 const Teuchos::RCP<const Export<LocalOrdinal,GlobalOrdinal,Node> > DomainExporter,
180 const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& domainMap,
183 {
184 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix function not implemented in"+std::string(__FILE__)+":"+std::to_string(__LINE__));
185 }
186
187
189 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
192
193
195
196
198 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
199 void
201 insertGlobalValues(GlobalOrdinal globalRow,
203 const ArrayView< const Scalar > &vals)
204 {
205 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix function not implemented in"+std::string(__FILE__)+":"+std::to_string(__LINE__));
206 }
207
208
210 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
211 void
213 insertLocalValues(LocalOrdinal localRow,
215 const ArrayView< const Scalar > &vals)
216 {
217 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix function not implemented in"+std::string(__FILE__)+":"+std::to_string(__LINE__));
218 }
219
220
222 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
223 void
225 replaceGlobalValues(GlobalOrdinal globalRow,
227 const ArrayView< const Scalar > &vals)
228 {
229 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix function not implemented in"+std::string(__FILE__)+":"+std::to_string(__LINE__));
230 }
231
232
234 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
235 void
237 replaceLocalValues (LocalOrdinal localRow,const ArrayView<const LocalOrdinal> &cols,const ArrayView<const Scalar> &vals)
238 {
239 XPETRA_MONITOR("TpetraBlockCrsMatrix::replaceLocalValues");
240 mtx_->replaceLocalValues(localRow,cols.getRawPtr(),vals.getRawPtr(),cols.size());
241 }
242
243
245 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
246 void
248 setAllToScalar(const Scalar &alpha)
249 {
250 XPETRA_MONITOR("TpetraBlockCrsMatrix::setAllToScalar"); mtx_->setAllToScalar(alpha);
251 }
252
253
255 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
256 void
258 scale(const Scalar &alpha)
259 {
260 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix function not implemented in "+std::string(__FILE__)+":"+std::to_string(__LINE__));
261 }
262
263
265 //** \warning This is an expert-only routine and should not be called from user code. (not implemented)
266 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
267 void
269 allocateAllValues(size_t numNonZeros,ArrayRCP<size_t> & rowptr, ArrayRCP<LocalOrdinal> & colind, ArrayRCP<Scalar> & values)
270 {
271 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix function not implemented in "+std::string(__FILE__)+":"+std::to_string(__LINE__));
272 }
273
274
276 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
277 void
279 setAllValues(const ArrayRCP<size_t> & rowptr, const ArrayRCP<LocalOrdinal> & colind, const ArrayRCP<Scalar> & values)
280 {
281 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix function not implemented in "+std::string(__FILE__)+":"+std::to_string(__LINE__));
282 }
283
284
286 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
287 void
291 ArrayRCP<const Scalar>& values) const
292 {
293 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix function not implemented in "+std::string(__FILE__)+":"+std::to_string(__LINE__));
294 }
295
297 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
298 void
301 {
302 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix function not implemented in "+std::string(__FILE__)+":"+std::to_string(__LINE__));
303 }
304
305
306
308
309 // Transformational Methods
311
312
313 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
314 void
320
321
322 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
323 void
331
332
333 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
334 void
340
341
342 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
343 void
347 {
348 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix function not implemented in "+std::string(__FILE__)+":"+std::to_string(__LINE__));
349 }
350
351
352 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
353 void
356 const RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> > & rangeMap,
357 const RCP<const Import<LocalOrdinal,GlobalOrdinal,Node> > &importer,
358 const RCP<const Export<LocalOrdinal,GlobalOrdinal,Node> > &exporter,
359 const RCP<ParameterList> &params)
360 {
361 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix function not implemented in "+std::string(__FILE__)+":"+std::to_string(__LINE__));
362 }
363
365
366
368
369
370 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
373 getRowMap() const
374 {
375 XPETRA_MONITOR("TpetraBlockCrsMatrix::getRowMap"); return toXpetra(mtx_->getRowMap());
376 }
377
378
379 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
382 getColMap() const
383 {
384 XPETRA_MONITOR("TpetraBlockCrsMatrix::getColMap"); return toXpetra(mtx_->getColMap());
385 }
386
387
388 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
391 getCrsGraph() const
392 {
393 XPETRA_MONITOR("TpetraBlockCrsMatrix::getCrsGraph");
394 using G_t = Tpetra::CrsGraph<LocalOrdinal,GlobalOrdinal,Node>;
396 RCP<G_t> t_graph = Teuchos::rcp_const_cast<G_t>(Teuchos::rcpFromRef(mtx_->getCrsGraph()));
397 RCP<const G_x> x_graph = rcp(new G_x(t_graph));
398 return x_graph;
399 }
400
401
402 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
405 getGlobalNumRows() const
406 { XPETRA_MONITOR("TpetraBlockCrsMatrix::getGlobalNumRows"); return mtx_->getGlobalNumRows(); }
407
408
409 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
412 getGlobalNumCols() const
413 { XPETRA_MONITOR("TpetraBlockCrsMatrix::getGlobalNumCols"); return mtx_->getGlobalNumCols(); }
414
415
416 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
417 size_t
419 getLocalNumRows() const
420 { XPETRA_MONITOR("TpetraBlockCrsMatrix::getLocalNumRows"); return mtx_->getLocalNumRows(); }
421
422
423 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
424 size_t
426 getLocalNumCols() const
427 { XPETRA_MONITOR("TpetraBlockCrsMatrix::getLocalNumCols"); return mtx_->getLocalNumCols(); }
428
429
430 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
434 { XPETRA_MONITOR("TpetraBlockCrsMatrix::getGlobalNumEntries"); return mtx_->getGlobalNumEntries(); }
435
436
437 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
438 size_t
441 { XPETRA_MONITOR("TpetraBlockCrsMatrix::getLocalNumEntries"); return mtx_->getLocalNumEntries(); }
442
443
444 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
445 size_t
447 getNumEntriesInLocalRow(LocalOrdinal localRow) const
448 { XPETRA_MONITOR("TpetraBlockCrsMatrix::getNumEntriesInLocalRow"); return mtx_->getNumEntriesInLocalRow(localRow); }
449
450
451 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
452 size_t
454 getNumEntriesInGlobalRow(GlobalOrdinal globalRow) const
455 { XPETRA_MONITOR("TpetraBlockCrsMatrix::getNumEntriesInGlobalRow"); return mtx_->getNumEntriesInGlobalRow(globalRow); }
456
457
458 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
460 { XPETRA_MONITOR("TpetraBlockCrsMatrix::getGlobalMaxNumRowEntries"); return mtx_->getGlobalMaxNumRowEntries(); }
461
462
463 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
465 { XPETRA_MONITOR("TpetraBlockCrsMatrix::getLocalMaxNumRowEntries"); return mtx_->getLocalMaxNumRowEntries(); }
466
467
468 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
470 { XPETRA_MONITOR("TpetraBlockCrsMatrix::isLocallyIndexed"); return mtx_->isLocallyIndexed(); }
471
472
473 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
475 { XPETRA_MONITOR("TpetraBlockCrsMatrix::isGloballyIndexed"); return mtx_->isGloballyIndexed(); }
476
477
478 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
480 { XPETRA_MONITOR("TpetraBlockCrsMatrix::isFillComplete"); return mtx_->isFillComplete(); }
481
482
483 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
485 { XPETRA_MONITOR("TpetraBlockCrsMatrix::isFillActive"); return false; }
486
487
488 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
490 { XPETRA_MONITOR("TpetraBlockCrsMatrix::getFrobeniusNorm"); return mtx_->getFrobeniusNorm(); }
491
492
493 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
495 { XPETRA_MONITOR("TpetraBlockCrsMatrix::supportsRowViews"); return mtx_->supportsRowViews(); }
496
497
498 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
499 void
501 getLocalRowCopy(LocalOrdinal LocalRow,
502 const ArrayView< LocalOrdinal > &Indices,
503 const ArrayView< Scalar > &Values,
504 size_t &NumEntries) const
505 {
506 XPETRA_MONITOR("TpetraBlockCrsMatrix::getLocalRowCopy");
507 typename Tpetra::CrsGraph<LocalOrdinal,GlobalOrdinal,Node>::nonconst_local_inds_host_view_type indices("indices",Indices.size());
508 typename Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::nonconst_values_host_view_type values("values", Values.size());
509
510 mtx_->getLocalRowCopy(LocalRow, indices, values, NumEntries);
511 for (size_t i=0;i<NumEntries; ++i) {
512 Indices[i]=indices(i);
513 Values[i]=values(i);
514 }
515 }
516 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
517 void
519 getLocalRowView(LocalOrdinal LocalRow, ArrayView< const LocalOrdinal > &Indices,
520 ArrayView< const Scalar > &Values) const
521 {
522 XPETRA_MONITOR("TpetraBlockCrsMatrix::getLocalRowView");
523 typename Tpetra::CrsGraph<LocalOrdinal,GlobalOrdinal,Node>::local_inds_host_view_type indices;
524 typename Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::values_host_view_type values;
525
526 mtx_->getLocalRowView(LocalRow, indices, values);
527 Indices = ArrayView<const LocalOrdinal> (indices.data(), indices.extent(0));
528 Values = ArrayView<const Scalar> (reinterpret_cast<const Scalar*>(values.data()), values.extent(0));
529 }
530
531
532 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
533 void
535 getGlobalRowView(GlobalOrdinal GlobalRow,
537 ArrayView< const Scalar > &Values) const
538 {
539 XPETRA_MONITOR("TpetraBlockCrsMatrix::getGlobalRowView");
540 typename Tpetra::CrsGraph<LocalOrdinal,GlobalOrdinal,Node>::global_inds_host_view_type indices;
541 typename Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::values_host_view_type values;
542
543 mtx_->getGlobalRowView(GlobalRow, indices, values);
544 Indices = ArrayView<const GlobalOrdinal> (indices.data(), indices.extent(0));
545 Values = ArrayView<const Scalar> (reinterpret_cast<const Scalar*>(values.data()), values.extent(0));
546 }
547
548
549 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
550 void
552 getGlobalRowCopy(GlobalOrdinal GlobalRow,
553 const ArrayView< GlobalOrdinal > &Indices,
554 const ArrayView< Scalar > &Values,
555 size_t &NumEntries) const
556 {
557 XPETRA_MONITOR("TpetraBlockCrsMatrix::getGlobalRowCopy");
558 typename Tpetra::CrsGraph<LocalOrdinal,GlobalOrdinal,Node>::nonconst_global_inds_host_view_type indices("indices",Indices.size());
559 typename Tpetra::CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>::nonconst_values_host_view_type values("values", Values.size());
560
561 mtx_->getGlobalRowCopy(GlobalRow, indices, values, NumEntries);
562 for (size_t i=0;i<NumEntries; ++i) {
563 Indices[i]=indices(i);
564 Values[i]=values(i);
565 }
566 }
567
568
569 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
570 bool
574
575
576 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
580 Teuchos::ETransp mode,
581 Scalar alpha,
582 Scalar beta) const
583 { XPETRA_MONITOR("TpetraBlockCrsMatrix::apply"); mtx_->apply(toTpetra(X), toTpetra(Y), mode, alpha, beta); }
584
585 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
589 Teuchos::ETransp mode,
590 Scalar alpha,
591 Scalar beta,
592 bool sumInterfaceValues,
593 const RCP<Import<LocalOrdinal, GlobalOrdinal, Node> >& regionInterfaceImporter,
594 const Teuchos::ArrayRCP<LocalOrdinal>& regionInterfaceLIDs) const
595 { }
596
597
598 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
601 getDomainMap() const
602 { XPETRA_MONITOR("TpetraBlockCrsMatrix::getDomainMap"); return toXpetra(mtx_->getDomainMap()); }
603
604
605 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
608 getRangeMap() const
609 { XPETRA_MONITOR("TpetraBlockCrsMatrix::getRangeMap"); return toXpetra(mtx_->getRangeMap()); }
610
611
612 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
613 std::string
615 description() const
616 { XPETRA_MONITOR("TpetraBlockCrsMatrix::description"); return mtx_->description(); }
617
618
619 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
620 void
623 const Teuchos::EVerbosityLevel verbLevel) const
624 {
625 XPETRA_MONITOR("TpetraBlockCrsMatrix::describe");
626 mtx_->describe(out, verbLevel);
627 }
628
629
630 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
631 void
633 setObjectLabel( const std::string &objectLabel )
634 {
635 XPETRA_MONITOR("TpetraCrsMatrix::setObjectLabel");
637 mtx_->setObjectLabel(objectLabel);
638 }
639
640
641
642 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
643 void
646 {
647 XPETRA_MONITOR("TpetraBlockCrsMatrix::getLocalDiagCopy");
649 diag,
650 tDiag,
651 "Xpetra::TpetraBlockCrsMatrix.getLocalDiagCopy() only accept Xpetra::TpetraVector as input arguments.");
652 mtx_->getLocalDiagCopy(*tDiag.getTpetra_Vector());
653 }
654
655
657 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
658 void
661 const Teuchos::ArrayView<const size_t> &offsets) const
662 {
663 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix function not implemented in "+std::string(__FILE__)+":"+std::to_string(__LINE__));
664 }
665
666
667
668 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
669 void
672 {
673 XPETRA_MONITOR("TpetraBlockCrsMatrix::getLocalDiagOffsets");
674
675 const size_t lclNumRows = mtx_->getGraph()->getLocalNumRows();
676 if (static_cast<size_t>(offsets.size()) < lclNumRows)
677 {
678 offsets.resize(lclNumRows);
679 }
680
681 // The input ArrayRCP must always be a host pointer. Thus, if
682 // device_type::memory_space is Kokkos::HostSpace, it's OK for us
683 // to write to that allocation directly as a Kokkos::View.
684 typedef typename Node::device_type device_type;
685 typedef typename device_type::memory_space memory_space;
686 if (std::is_same<memory_space, Kokkos::HostSpace>::value)
687 {
688 // It is always syntactically correct to assign a raw host
689 // pointer to a device View, so this code will compile correctly
690 // even if this branch never runs.
691 typedef Kokkos::View<size_t*, device_type, Kokkos::MemoryUnmanaged> output_type;
692 output_type offsetsOut (offsets.getRawPtr(), offsets.size());
693 mtx_->getLocalDiagOffsets(offsetsOut);
694 }
695 else
696 {
697 Kokkos::View<size_t*, device_type> offsetsTmp ("diagOffsets", offsets.size());
698 mtx_->getLocalDiagOffsets(offsetsTmp);
699 typedef Kokkos::View<size_t*, Kokkos::HostSpace, Kokkos::MemoryUnmanaged> output_type;
700 output_type offsetsOut(offsets.getRawPtr(), offsets.size());
701 Kokkos::deep_copy(offsetsOut, offsetsTmp);
702 }
703 }
704
705
706 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
707 void
710 {
711 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix::replaceDiag: function not implemented in "+std::string(__FILE__)+":"+std::to_string(__LINE__));
712 }
713
714
715 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
716 void
719 {
720 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix function not implemented in "+std::string(__FILE__)+":"+std::to_string(__LINE__));
721 }
722
723
724 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
725 void
728 {
729 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix function not implemented in "+std::string(__FILE__)+":"+std::to_string(__LINE__));
730 }
731
732
733 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
736 getMap() const
737 {
738 XPETRA_MONITOR("TpetraBlockCrsMatrix::getMap"); return rcp( new TpetraMap< LocalOrdinal, GlobalOrdinal, Node >(mtx_->getMap()) );
739 }
740
741
743 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
744 void
748 {
749 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix function not implemented in "+std::string(__FILE__)+":"+std::to_string(__LINE__));
750 }
751
752
754 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
755 void
759 {
760 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix function not implemented in "+std::string(__FILE__)+":"+std::to_string(__LINE__));
761 }
762
763
765 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
766 void
770 {
771 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix function not implemented in "+std::string(__FILE__)+":"+std::to_string(__LINE__));
772 }
773
774
776 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
777 void
781 {
782 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix function not implemented in "+std::string(__FILE__)+":"+std::to_string(__LINE__));
783 }
784
785
786 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
787 void
790 {
791 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix function not implemented in "+std::string(__FILE__)+":"+std::to_string(__LINE__));
792 }
793
794
795template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
796bool
798hasMatrix() const
799{
800 return !mtx_.is_null();
801}
802
803
804template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
806TpetraBlockCrsMatrix(const Teuchos::RCP<Tpetra::BlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > &mtx)
807: mtx_(mtx)
808{ }
809
810
811template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
818
819
820// TODO: remove
821template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
828
829#ifdef HAVE_XPETRA_TPETRA
830
831// was: typedef typename Xpetra::CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::local_matrix_type local_matrix_type;
832//using local_matrix_type = typename CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node>::local_matrix_type;
833
834template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
838{
839 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix does not support getLocalMatrix due to missing Kokkos::CrsMatrix in Tpetra's experimental implementation in "+std::string(__FILE__)+":"+std::to_string(__LINE__));
840
841#ifndef __NVCC__
843#endif // __NVCC__
844
846}
847
848template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
851getLocalMatrixHost () const
852{
853 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix does not support getLocalMatrix due to missing Kokkos::CrsMatrix in Tpetra's experimental implementation in "+std::string(__FILE__)+":"+std::to_string(__LINE__));
854
855#ifndef __NVCC__
856 typename local_matrix_type::HostMirror ret;
857#endif // __NVCC__
858
860}
861
862
863template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
864void
866setAllValues (const typename local_matrix_type::row_map_type& ptr,
867 const typename local_matrix_type::StaticCrsGraphType::entries_type::non_const_type& ind,
868 const typename local_matrix_type::values_type& val)
869{
870 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix does not support setAllValues due to missing Kokkos::CrsMatrix in Tpetra's experimental implementation in "+std::string(__FILE__)+":"+std::to_string(__LINE__));
871}
872
873#endif // HAVE_XPETRA_TPETRA
874
875
876#ifdef HAVE_XPETRA_EPETRA
877
878#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT))) || \
879 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_INT))))
880
881
882 // specialization of TpetraBlockCrsMatrix for GO=LO=int and Node=EpetraNode
883 template <class Scalar>
884 class TpetraBlockCrsMatrix<Scalar,int,int,EpetraNode>
885 : public CrsMatrix<Scalar,int,int,EpetraNode>//, public TpetraRowMatrix<Scalar,int,int,Node>
886 {
887
888 // The following typedef are used by the XPETRA_DYNAMIC_CAST() macro.
889 typedef int LocalOrdinal;
890 typedef int GlobalOrdinal;
896
897 public:
898
900
905
910
915
920
925
930
932 TpetraBlockCrsMatrix(const Teuchos::RCP<const Tpetra::BlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >& sourceMatrix,
934 const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& domainMap = Teuchos::null,
935 const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& rangeMap = Teuchos::null,
936 const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null)
938
940 TpetraBlockCrsMatrix(const Teuchos::RCP<const Tpetra::BlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >& sourceMatrix,
942 const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& domainMap = Teuchos::null,
943 const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& rangeMap = Teuchos::null,
944 const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null)
946
948 TpetraBlockCrsMatrix(const Teuchos::RCP<const Tpetra::BlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >& sourceMatrix,
949 const Import<LocalOrdinal,GlobalOrdinal,Node> & RowImporter,
950 const Teuchos::RCP<const Import<LocalOrdinal,GlobalOrdinal,Node> > DomainImporter,
951 const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& domainMap,
955
957 TpetraBlockCrsMatrix(const Teuchos::RCP<const Tpetra::BlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >& sourceMatrix,
958 const Export<LocalOrdinal,GlobalOrdinal,Node> & RowExporter,
959 const Teuchos::RCP<const Export<LocalOrdinal,GlobalOrdinal,Node> > DomainExporter,
960 const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& domainMap,
963 {
965 }
966
969
970
972
976
980
984
988
990 void setAllToScalar(const Scalar &alpha) {}
991
993 void scale(const Scalar &alpha)
994 {}
995
997 //** \warning This is an expert-only routine and should not be called from user code. (not implemented)
998 void allocateAllValues(size_t numNonZeros,ArrayRCP<size_t> & rowptr, ArrayRCP<LocalOrdinal> & colind, ArrayRCP<Scalar> & values)
999 {}
1000
1002 void setAllValues(const ArrayRCP<size_t> & rowptr, const ArrayRCP<LocalOrdinal> & colind, const ArrayRCP<Scalar> & values)
1003 {}
1004
1008
1009
1012 {}
1013
1015
1017 void resumeFill(const RCP< ParameterList > &params=null) { /*noop*/ }
1018
1020 void fillComplete(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const RCP< ParameterList > &params=null) { /*noop*/ }
1021
1023 void fillComplete(const RCP< ParameterList > &params=null) { /*noop*/ }
1024
1025
1029
1032 const RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> > & rangeMap,
1033 const RCP<const Import<LocalOrdinal,GlobalOrdinal,Node> > &importer=Teuchos::null,
1034 const RCP<const Export<LocalOrdinal,GlobalOrdinal,Node> > &exporter=Teuchos::null,
1035 const RCP<ParameterList> &params=Teuchos::null)
1036 {}
1037
1038
1040
1042 const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getRowMap() const { return Teuchos::null; }
1043
1045 const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getColMap() const { return Teuchos::null; }
1046
1050
1052 global_size_t getGlobalNumRows() const { return 0; }
1053
1055 global_size_t getGlobalNumCols() const { return 0; }
1056
1058 size_t getLocalNumRows() const { return 0; }
1059
1061 size_t getLocalNumCols() const { return 0; }
1062
1065
1067 size_t getLocalNumEntries() const { return 0; }
1068
1070 size_t getNumEntriesInLocalRow(LocalOrdinal localRow) const { return 0; }
1071
1073 size_t getNumEntriesInGlobalRow(GlobalOrdinal globalRow) const { return 0; }
1074
1076 size_t getGlobalMaxNumRowEntries() const { return 0; }
1077
1079 size_t getLocalMaxNumRowEntries() const { return 0; }
1080
1082 bool isLocallyIndexed() const { return false; }
1083
1085 bool isGloballyIndexed() const { return false; }
1086
1088 bool isFillComplete() const { return false; }
1089
1091 bool isFillActive() const { return false; }
1092
1095
1097 bool supportsRowViews() const { return false; }
1098
1100 void getLocalRowCopy(LocalOrdinal LocalRow, const ArrayView< LocalOrdinal > &Indices, const ArrayView< Scalar > &Values, size_t &NumEntries) const { }
1101
1104
1106 void getGlobalRowCopy(GlobalOrdinal GlobalRow, const ArrayView< GlobalOrdinal > &indices, const ArrayView< Scalar > &values, size_t &numEntries) const { }
1107
1110
1112 bool haveGlobalConstants() const {return false;}
1113
1114
1116
1119
1122
1125
1126
1128
1130 std::string description() const { return std::string(""); }
1131
1134
1135
1138
1141
1144
1148
1150
1153
1154
1156
1159
1164
1169
1174
1179
1182
1183
1184
1186
1188 bool hasMatrix() const { return false; }
1189
1191 TpetraBlockCrsMatrix(const Teuchos::RCP<Tpetra::BlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > &mtx) {
1193 }
1194
1197
1200
1201#ifdef HAVE_XPETRA_TPETRA
1203
1205 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix does not support getLocalMatrix due to missing Kokkos::CrsMatrix in Tpetra's experimental implementation in "+std::string(__FILE__)+":"+std::to_string(__LINE__));
1207 return ret; // make compiler happy
1208 }
1209
1210 void setAllValues (const typename local_matrix_type::row_map_type& ptr,
1211 const typename local_matrix_type::StaticCrsGraphType::entries_type::non_const_type& ind,
1212 const typename local_matrix_type::values_type& val)
1213 {
1214 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix does not support setAllValues due to missing Kokkos::CrsMatrix in Tpetra's experimental implementation in "+std::string(__FILE__)+":"+std::to_string(__LINE__));
1215 }
1216#endif // HAVE_XPETRA_TPETRA
1217
1218 }; // TpetraBlockCrsMatrix class
1219
1220
1221#endif // #if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_INT)))
1222
1223
1224
1225
1226#if ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_LONG_LONG))) || \
1227 (!defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_SERIAL) || !defined(HAVE_TPETRA_INST_INT_LONG_LONG))))
1228
1229 // specialization of TpetraBlockCrsMatrix for GO=long long and Node=EpetraNode
1230 template <class Scalar>
1231 class TpetraBlockCrsMatrix<Scalar,int,long long,EpetraNode>
1232 : public CrsMatrix<Scalar,int,long long,EpetraNode>//, public TpetraRowMatrix<Scalar,int,int,Node>
1233 {
1234
1235 // The following typedef are used by the XPETRA_DYNAMIC_CAST() macro.
1236 typedef int LocalOrdinal;
1237 typedef long long GlobalOrdinal;
1243
1244 public:
1245
1247
1251
1255
1259
1263
1267
1271
1272
1273
1274
1276 TpetraBlockCrsMatrix(const Teuchos::RCP<const Tpetra::BlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >& sourceMatrix,
1278 const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& domainMap = Teuchos::null,
1279 const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& rangeMap = Teuchos::null,
1280 const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null)
1282
1284 TpetraBlockCrsMatrix(const Teuchos::RCP<const Tpetra::BlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >& sourceMatrix,
1286 const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& domainMap = Teuchos::null,
1287 const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& rangeMap = Teuchos::null,
1288 const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null)
1290
1292 TpetraBlockCrsMatrix(const Teuchos::RCP<const Tpetra::BlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >& sourceMatrix,
1293 const Import<LocalOrdinal,GlobalOrdinal,Node> & RowImporter,
1294 const Teuchos::RCP<const Import<LocalOrdinal,GlobalOrdinal,Node> > DomainImporter,
1295 const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& domainMap,
1296 const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& rangeMap,
1299
1301 TpetraBlockCrsMatrix(const Teuchos::RCP<const Tpetra::BlockCrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >& sourceMatrix,
1302 const Export<LocalOrdinal,GlobalOrdinal,Node> & RowExporter,
1303 const Teuchos::RCP<const Export<LocalOrdinal,GlobalOrdinal,Node> > DomainExporter,
1304 const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& domainMap,
1305 const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& rangeMap,
1308
1311
1312
1314
1318
1322
1326
1330
1332 void setAllToScalar(const Scalar &alpha) {}
1333
1335 void scale(const Scalar &alpha)
1336 {}
1337
1339 //** \warning This is an expert-only routine and should not be called from user code. (not implemented)
1340 void allocateAllValues(size_t numNonZeros,ArrayRCP<size_t> & rowptr, ArrayRCP<LocalOrdinal> & colind, ArrayRCP<Scalar> & values)
1341 {}
1342
1344 void setAllValues(const ArrayRCP<size_t> & rowptr, const ArrayRCP<LocalOrdinal> & colind, const ArrayRCP<Scalar> & values)
1345 {}
1346
1350
1351
1354 {}
1355
1356
1358
1360 void resumeFill(const RCP< ParameterList > &params=null) { /*noop*/ }
1361
1363 void fillComplete(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const RCP< ParameterList > &params=null) { /*noop*/ }
1364
1366 void fillComplete(const RCP< ParameterList > &params=null) { /*noop*/ }
1367
1368
1372
1375 const RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> > & rangeMap,
1376 const RCP<const Import<LocalOrdinal,GlobalOrdinal,Node> > &importer=Teuchos::null,
1377 const RCP<const Export<LocalOrdinal,GlobalOrdinal,Node> > &exporter=Teuchos::null,
1378 const RCP<ParameterList> &params=Teuchos::null)
1379 {}
1380
1381
1383
1385 const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getRowMap() const { return Teuchos::null; }
1386
1388 const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getColMap() const { return Teuchos::null; }
1389
1393
1395 global_size_t getGlobalNumRows() const { return 0; }
1396
1398 global_size_t getGlobalNumCols() const { return 0; }
1399
1401 size_t getLocalNumRows() const { return 0; }
1402
1404 size_t getLocalNumCols() const { return 0; }
1405
1408
1410 size_t getLocalNumEntries() const { return 0; }
1411
1413 size_t getNumEntriesInLocalRow(LocalOrdinal localRow) const { return 0; }
1414
1416 size_t getNumEntriesInGlobalRow(GlobalOrdinal globalRow) const { return 0; }
1417
1419 size_t getGlobalMaxNumRowEntries() const { return 0; }
1420
1422 size_t getLocalMaxNumRowEntries() const { return 0; }
1423
1425 bool isLocallyIndexed() const { return false; }
1426
1428 bool isGloballyIndexed() const { return false; }
1429
1431 bool isFillComplete() const { return false; }
1432
1434 bool isFillActive() const { return false; }
1435
1438
1440 bool supportsRowViews() const { return false; }
1441
1443 void getLocalRowCopy(LocalOrdinal LocalRow, const ArrayView< LocalOrdinal > &Indices, const ArrayView< Scalar > &Values, size_t &NumEntries) const { }
1444
1447
1449 void getGlobalRowCopy(GlobalOrdinal GlobalRow, const ArrayView< GlobalOrdinal > &indices, const ArrayView< Scalar > &values, size_t &numEntries) const { }
1450
1453
1455 bool haveGlobalConstants() const {return true;}
1456
1457
1459
1462
1465
1468
1469
1471
1473 std::string description() const { return std::string(""); }
1474
1477
1480
1483
1486
1490
1492
1495
1497
1500
1505
1510
1515
1520
1523
1524
1525
1527
1529 bool hasMatrix() const { return false; }
1530
1532 TpetraBlockCrsMatrix(const Teuchos::RCP<Tpetra::BlockCrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> > &mtx) {
1534 }
1535
1538
1541
1542#ifdef HAVE_XPETRA_TPETRA
1544
1546 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix does not support getLocalMatrix due to missing Kokkos::CrsMatrix in Tpetra's experimental implementation in "+std::string(__FILE__)+":"+std::to_string(__LINE__));
1549 }
1550
1551 void setAllValues (const typename local_matrix_type::row_map_type& ptr,
1552 const typename local_matrix_type::StaticCrsGraphType::entries_type::non_const_type& ind,
1553 const typename local_matrix_type::values_type& val)
1554 {
1555 throw std::runtime_error("Xpetra::TpetraBlockCrsMatrix does not support setAllValues due to missing Kokkos::CrsMatrix in Tpetra's experimental implementation in "+std::string(__FILE__)+":"+std::to_string(__LINE__));
1556 }
1557#endif // HAVE_XPETRA_TPETRA
1558
1559 }; // TpetraBlockCrsMatrix class
1560
1561
1562#endif // IF ((defined(EPETRA_HAVE_OMP) && (!defined(HAVE_TPETRA_INST_OPENMP) || !defined(HAVE_TPETRA_INST_INT_LONG_LONG)))
1563
1564
1565
1566
1567#endif // HAVE_XPETRA_EPETRA
1568
1569
1570} // Xpetra namespace
1571
1572#endif // XPETRA_TPETRABLOCKCRSMATRIX_DEF_HPP
1573
1574
#define XPETRA_MONITOR(funcName)
#define XPETRA_DYNAMIC_CAST(type, obj, newObj, exceptionMsg)
#define XPETRA_TPETRA_ETI_EXCEPTION(cl, obj, go, node)
size_type size() const
T * getRawPtr() const
void resize(const size_type n, const T &val=T())
size_type size() const
T * getRawPtr() const
static const EVerbosityLevel verbLevel_default
virtual void setObjectLabel(const std::string &objectLabel)
KokkosSparse::CrsMatrix< impl_scalar_type, LocalOrdinal, execution_space, void, typename local_graph_type::size_type > local_matrix_type
The specialization of Kokkos::CrsMatrix that represents the part of the sparse matrix on each MPI pro...
void setAllValues(const typename local_matrix_type::row_map_type &ptr, const typename local_matrix_type::StaticCrsGraphType::entries_type::non_const_type &ind, const typename local_matrix_type::values_type &val)
Xpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::local_matrix_type local_matrix_type
void setAllToScalar(const Scalar &alpha)
Set all matrix entries equal to scalarThis.
RCP< Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_BlockCrsMatrixNonConst() const
Get the underlying Tpetra matrix.
void getLocalRowCopy(LocalOrdinal LocalRow, const ArrayView< LocalOrdinal > &Indices, const ArrayView< Scalar > &Values, size_t &NumEntries) const
Extract a list of entries in a specified local row of the matrix. Put into storage allocated by calli...
size_t getLocalNumCols() const
Returns the number of columns connected to the locally owned rows of this matrix.
bool isFillActive() const
Returns true if the matrix is in edit mode.
TpetraBlockCrsMatrix(const Teuchos::RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > &graph, const LocalOrdinal blockSize)
Constructor specifying a previously constructed graph & blocksize.
size_t getNumEntriesInLocalRow(LocalOrdinal localRow) const
Returns the current number of entries on this node in the specified local row.
TpetraBlockCrsMatrix(const Teuchos::RCP< const Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &sourceMatrix, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap=Teuchos::null, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap=Teuchos::null, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Constructor for a fused export (not implemented(.
void getGlobalRowCopy(GlobalOrdinal GlobalRow, const ArrayView< GlobalOrdinal > &indices, const ArrayView< Scalar > &values, size_t &numEntries) const
Extract a list of entries in a specified global row of this matrix. Put into pre-allocated storage.
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getRangeMap() const
Returns the Map associated with the range of this operator, which must be compatible with Y....
void replaceGlobalValues(GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Replace matrix entries, using global IDs (not implemented)
size_t getLocalMaxNumRowEntries() const
Returns the maximum number of entries across all rows/columns on this node.
void fillComplete(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const RCP< ParameterList > &params=null)
Signal that data entry is complete, specifying domain and range maps.
void setAllValues(const ArrayRCP< size_t > &rowptr, const ArrayRCP< LocalOrdinal > &colind, const ArrayRCP< Scalar > &values)
Sets the 1D pointer arrays of the graph (not impelmented)
TpetraBlockCrsMatrix(const Teuchos::RCP< const Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &sourceMatrix, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap=Teuchos::null, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap=Teuchos::null, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Constructor for a fused import ( not implemented )
void apply(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, Scalar alpha=ScalarTraits< Scalar >::one(), Scalar beta=ScalarTraits< Scalar >::zero()) const
Computes the sparse matrix-multivector multiplication.
std::string description() const
A simple one-line description of this object.
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getRowMap() const
Returns the Map that describes the row distribution in this matrix.
TpetraBlockCrsMatrix(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, size_t maxNumEntriesPerRow, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Constructor specifying fixed number of entries for each row (not implemented)
void fillComplete(const RCP< ParameterList > &params=null)
Signal that data entry is complete.
size_t getLocalNumRows() const
Returns the number of matrix rows owned on the calling node.
void expertStaticFillComplete(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const RCP< const Import< LocalOrdinal, GlobalOrdinal, Node > > &importer=Teuchos::null, const RCP< const Export< LocalOrdinal, GlobalOrdinal, Node > > &exporter=Teuchos::null, const RCP< ParameterList > &params=Teuchos::null)
Expert static fill complete.
TpetraBlockCrsMatrix(const Teuchos::RCP< const Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &sourceMatrix, const Import< LocalOrdinal, GlobalOrdinal, Node > &RowImporter, const Teuchos::RCP< const Import< LocalOrdinal, GlobalOrdinal, Node > > DomainImporter, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const Teuchos::RCP< Teuchos::ParameterList > &params)
Constructor for a fused import ( not implemented )
global_size_t getGlobalNumEntries() const
Returns the global number of entries in this matrix.
void removeEmptyProcessesInPlace(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &newMap)
void rightScale(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &x)
Right scale matrix using the given vector entries.
void replaceLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Replace matrix entries, using local IDs.
bool isFillComplete() const
Returns true if the matrix is in compute mode, i.e. if fillComplete() has been called.
bool supportsRowViews() const
Returns true if getLocalRowView() and getGlobalRowView() are valid for this class.
void getLocalDiagCopy(Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &diag, const Teuchos::ArrayView< const size_t > &offsets) const
Get a copy of the diagonal entries owned by this node, with local row indices.
TpetraBlockCrsMatrix(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &colMap, size_t maxNumEntriesPerRow, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Constructor specifying column Map and fixed number of entries for each row (not implemented)
bool isLocallyIndexed() const
If matrix indices are in the local range, this function returns true. Otherwise, this function return...
TpetraBlockCrsMatrix(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, const ArrayRCP< const size_t > &NumEntriesPerRowToAlloc, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Constructor specifying (possibly different) number of entries in each row (not implemented)
void replaceDomainMapAndImporter(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &newDomainMap, Teuchos::RCP< const Import< LocalOrdinal, GlobalOrdinal, Node > > &newImporter)
Replaces the current domainMap and importer with the user-specified objects.
ScalarTraits< Scalar >::magnitudeType getFrobeniusNorm() const
Returns the Frobenius norm of the matrix.
void allocateAllValues(size_t numNonZeros, ArrayRCP< size_t > &rowptr, ArrayRCP< LocalOrdinal > &colind, ArrayRCP< Scalar > &values)
Allocates and returns ArrayRCPs of the Crs arrays — This is an Xpetra-only routine.
void getLocalRowView(LocalOrdinal LocalRow, ArrayView< const LocalOrdinal > &indices, ArrayView< const Scalar > &values) const
Extract a const, non-persisting view of local indices in a specified row of the matrix.
size_t getNumEntriesInGlobalRow(GlobalOrdinal globalRow) const
Returns the current number of entries in the (locally owned) global row.
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getColMap() const
Returns the Map that describes the column distribution in this matrix.
void getLocalDiagOffsets(Teuchos::ArrayRCP< size_t > &offsets) const
Get offsets of the diagonal entries in the matrix.
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const
Print the object with some verbosity level to an FancyOStream object.
RCP< const Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_BlockCrsMatrix() const
Get the underlying Tpetra matrix.
TpetraBlockCrsMatrix(const TpetraBlockCrsMatrix &matrix)
Deep copy constructor.
void getGlobalRowView(GlobalOrdinal GlobalRow, ArrayView< const GlobalOrdinal > &indices, ArrayView< const Scalar > &values) const
Extract a const, non-persisting view of global indices in a specified row of the matrix.
void doImport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &source, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM)
Import (using an Exporter).
void getAllValues(ArrayRCP< Scalar > &values)
Gets the 1D pointer arrays of the graph (not implemented)
void doImport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &source, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Import.
RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > getCrsGraph() const
Returns the CrsGraph associated with this matrix.
void replaceDiag(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &diag)
Replace the diagonal entries of the matrix.
bool isGloballyIndexed() const
If matrix indices are in the global range, this function returns true. Otherwise, this function retur...
TpetraBlockCrsMatrix(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &colMap, const ArrayRCP< const size_t > &NumEntriesPerRowToAlloc, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Constructor specifying column Map and number of entries in each row (not implemented)
bool haveGlobalConstants() const
Returns true if globalConstants have been computed; false otherwise.
void doExport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &dest, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM)
Export (using an Importer).
TpetraBlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > TpetraBlockCrsMatrixClass
void getLocalDiagCopy(Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &diag) const
Get a copy of the diagonal entries owned by this node, with local row idices.
TpetraBlockCrsMatrix(const Teuchos::RCP< Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &mtx)
TpetraBlockCrsMatrix constructor to wrap a Tpetra::BlockCrsMatrix object.
TpetraBlockCrsMatrix(const Teuchos::RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > &graph, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Constructor specifying a previously constructed graph ( not implemented )
size_t getGlobalMaxNumRowEntries() const
Returns the maximum number of entries across all rows/columns on all nodes.
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getDomainMap() const
Returns the Map associated with the domain of this operator. This will be null until fillComplete() i...
TpetraBlockCrsMatrix(const Teuchos::RCP< const Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &sourceMatrix, const Export< LocalOrdinal, GlobalOrdinal, Node > &RowExporter, const Teuchos::RCP< const Export< LocalOrdinal, GlobalOrdinal, Node > > DomainExporter, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const Teuchos::RCP< Teuchos::ParameterList > &params)
Constructor for a fused export (not implemented(.
void leftScale(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &x)
Left scale matrix using the given vector entries.
global_size_t getGlobalNumRows() const
Number of global elements in the row map of this matrix.
void getAllValues(ArrayRCP< const size_t > &rowptr, ArrayRCP< const LocalOrdinal > &colind, ArrayRCP< const Scalar > &values) const
Gets the 1D pointer arrays of the graph (not implemented)
global_size_t getGlobalNumCols() const
Number of global columns in the matrix.
void doExport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &dest, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Export.
void insertLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Insert matrix entries, using local IDs (not implemented)
TpetraVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > TpetraVectorClass
size_t getLocalNumEntries() const
Returns the local number of entries in this matrix.
void scale(const Scalar &alpha)
Scale the current values of a matrix, this = alpha*this (not implemented)
void insertGlobalValues(GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Insert matrix entries, using global IDs (not implemented)
Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getMap() const
Implements DistObject interface.
TpetraBlockCrsMatrix(const Teuchos::RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > &graph, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Constructor specifying a previously constructed graph ( not implemented )
void replaceLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Replace matrix entries, using local IDs.
void rightScale(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &x)
Right scale matrix using the given vector entries.
void setAllValues(const typename local_matrix_type::row_map_type &ptr, const typename local_matrix_type::StaticCrsGraphType::entries_type::non_const_type &ind, const typename local_matrix_type::values_type &val)
TpetraBlockCrsMatrix(const Teuchos::RCP< const Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &sourceMatrix, const Import< LocalOrdinal, GlobalOrdinal, Node > &RowImporter, const Teuchos::RCP< const Import< LocalOrdinal, GlobalOrdinal, Node > > DomainImporter, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const Teuchos::RCP< Teuchos::ParameterList > &params)
Constructor for a fused import (not implemented)
void getAllValues(ArrayRCP< Scalar > &values)
Gets the 1D pointer arrays of the graph (not implemented)
void replaceGlobalValues(GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Replace matrix entries, using global IDs (not implemented)
size_t getNumEntriesInLocalRow(LocalOrdinal localRow) const
Returns the current number of entries on this node in the specified local row.
TpetraBlockCrsMatrix(const Teuchos::RCP< Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &mtx)
TpetraBlockCrsMatrix constructor to wrap a Tpetra::BlockCrsMatrix object.
void doExport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &dest, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM)
Export (using an Importer).
void removeEmptyProcessesInPlace(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &newMap)
bool isGloballyIndexed() const
If matrix indices are in the global range, this function returns true. Otherwise, this function retur...
void getLocalRowCopy(LocalOrdinal LocalRow, const ArrayView< LocalOrdinal > &Indices, const ArrayView< Scalar > &Values, size_t &NumEntries) const
Extract a list of entries in a specified local row of the matrix. Put into storage allocated by calli...
TpetraBlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > TpetraBlockCrsMatrixClass
TpetraBlockCrsMatrix(const Teuchos::RCP< const Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &sourceMatrix, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap=Teuchos::null, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap=Teuchos::null, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Constructor for a fused import (not implemented)
TpetraBlockCrsMatrix(const Teuchos::RCP< const Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &sourceMatrix, const Export< LocalOrdinal, GlobalOrdinal, Node > &RowExporter, const Teuchos::RCP< const Export< LocalOrdinal, GlobalOrdinal, Node > > DomainExporter, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const Teuchos::RCP< Teuchos::ParameterList > &params)
Constructor for a fused export (not implemented)
RCP< Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_BlockCrsMatrixNonConst() const
Get the underlying Tpetra matrix.
TpetraBlockCrsMatrix(const Teuchos::RCP< const Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > &sourceMatrix, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap=Teuchos::null, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap=Teuchos::null, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Constructor for a fused export (not implemented)
bool isFillComplete() const
Returns true if the matrix is in compute mode, i.e. if fillComplete() has been called.
void doImport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &source, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Import.
size_t getNumEntriesInGlobalRow(GlobalOrdinal globalRow) const
Returns the current number of entries in the (locally owned) global row.
size_t getLocalMaxNumRowEntries() const
Returns the maximum number of entries across all rows/columns on this node.
global_size_t getGlobalNumRows() const
Number of global elements in the row map of this matrix.
size_t getLocalNumEntries() const
Returns the local number of entries in this matrix.
void getLocalDiagCopy(Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &diag, const Teuchos::ArrayView< const size_t > &offsets) const
Get a copy of the diagonal entries owned by this node, with local row indices.
TpetraBlockCrsMatrix(const TpetraBlockCrsMatrix &matrix)
Deep copy constructor.
size_t getLocalNumRows() const
Returns the number of matrix rows owned on the calling node.
void fillComplete(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const RCP< ParameterList > &params=null)
Signal that data entry is complete, specifying domain and range maps.
void insertGlobalValues(GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Insert matrix entries, using global IDs (not implemented)
size_t getLocalNumCols() const
Returns the number of columns connected to the locally owned rows of this matrix.
void fillComplete(const RCP< ParameterList > &params=null)
Signal that data entry is complete.
Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getMap() const
Implements DistObject interface.
void insertLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Insert matrix entries, using local IDs (not implemented)
TpetraBlockCrsMatrix(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &colMap, size_t maxNumEntriesPerRow, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Constructor specifying column Map and fixed number of entries for each row (not implemented)
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const
Print the object with some verbosity level to an FancyOStream object.
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getDomainMap() const
Returns the Map associated with the domain of this operator. This will be null until fillComplete() i...
void getGlobalRowView(GlobalOrdinal GlobalRow, ArrayView< const GlobalOrdinal > &indices, ArrayView< const Scalar > &values) const
Extract a const, non-persisting view of global indices in a specified row of the matrix.
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getRangeMap() const
Returns the Map associated with the range of this operator, which must be compatible with Y....
void expertStaticFillComplete(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const RCP< const Import< LocalOrdinal, GlobalOrdinal, Node > > &importer=Teuchos::null, const RCP< const Export< LocalOrdinal, GlobalOrdinal, Node > > &exporter=Teuchos::null, const RCP< ParameterList > &params=Teuchos::null)
Expert static fill complete.
void replaceDomainMapAndImporter(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &newDomainMap, Teuchos::RCP< const Import< LocalOrdinal, GlobalOrdinal, Node > > &newImporter)
Replaces the current domainMap and importer with the user-specified objects.
void allocateAllValues(size_t numNonZeros, ArrayRCP< size_t > &rowptr, ArrayRCP< LocalOrdinal > &colind, ArrayRCP< Scalar > &values)
Allocates and returns ArrayRCPs of the Crs arrays — This is an Xpetra-only routine.
TpetraBlockCrsMatrix(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &colMap, const ArrayRCP< const size_t > &NumEntriesPerRowToAlloc, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Constructor specifying column Map and number of entries in each row (not implemented)
Xpetra::CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::local_matrix_type local_matrix_type
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getColMap() const
Returns the Map that describes the column distribution in this matrix.
void replaceDiag(Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &diag) const
ScalarTraits< Scalar >::magnitudeType getFrobeniusNorm() const
Returns the Frobenius norm of the matrix.
global_size_t getGlobalNumEntries() const
Returns the global number of entries in this matrix.
RCP< const Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_BlockCrsMatrix() const
Get the underlying Tpetra matrix.
void getAllValues(ArrayRCP< const size_t > &rowptr, ArrayRCP< const LocalOrdinal > &colind, ArrayRCP< const Scalar > &values) const
Gets the 1D pointer arrays of the graph (not implemented)
void scale(const Scalar &alpha)
Scale the current values of a matrix, this = alpha*this (not implemented)
size_t getGlobalMaxNumRowEntries() const
Returns the maximum number of entries across all rows/columns on all nodes.
void setAllValues(const ArrayRCP< size_t > &rowptr, const ArrayRCP< LocalOrdinal > &colind, const ArrayRCP< Scalar > &values)
Sets the 1D pointer arrays of the graph (not impelmented)
void setAllToScalar(const Scalar &alpha)
Set all matrix entries equal to scalarThis.
bool haveGlobalConstants() const
Returns true if globalConstants have been computed; false otherwise.
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getRowMap() const
Returns the Map that describes the row distribution in this matrix.
TpetraBlockCrsMatrix(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, size_t maxNumEntriesPerRow, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Constructor specifying fixed number of entries for each row (not implemented)
TpetraBlockCrsMatrix(const Teuchos::RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > &graph, const LocalOrdinal blockSize)
Constructor specifying a previously constructed graph & blocksize.
RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > getCrsGraph() const
Returns the CrsGraph associated with this matrix.
void getLocalRowView(LocalOrdinal LocalRow, ArrayView< const LocalOrdinal > &indices, ArrayView< const Scalar > &values) const
Extract a const, non-persisting view of local indices in a specified row of the matrix.
void apply(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, Scalar alpha=ScalarTraits< Scalar >::one(), Scalar beta=ScalarTraits< Scalar >::zero()) const
Computes the sparse matrix-multivector multiplication.
std::string description() const
A simple one-line description of this object.
void getGlobalRowCopy(GlobalOrdinal GlobalRow, const ArrayView< GlobalOrdinal > &indices, const ArrayView< Scalar > &values, size_t &numEntries) const
Extract a list of entries in a specified global row of this matrix. Put into pre-allocated storage.
void doExport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &dest, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Export.
void doImport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &source, const Export< LocalOrdinal, GlobalOrdinal, Node > &exporter, CombineMode CM)
Import (using an Exporter).
global_size_t getGlobalNumCols() const
Number of global columns in the matrix.
bool isLocallyIndexed() const
If matrix indices are in the local range, this function returns true. Otherwise, this function return...
bool supportsRowViews() const
Returns true if getLocalRowView() and getGlobalRowView() are valid for this class.
void getLocalDiagOffsets(Teuchos::ArrayRCP< size_t > &offsets) const
Get offsets of the diagonal entries in the matrix.
TpetraBlockCrsMatrix(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, const ArrayRCP< const size_t > &NumEntriesPerRowToAlloc, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Constructor specifying (possibly different) number of entries in each row (not implemented)
void leftScale(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &x)
Left scale matrix using the given vector entries.
TpetraVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > TpetraVectorClass
void getLocalDiagCopy(Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &diag) const
Get a copy of the diagonal entries owned by this node, with local row idices.
void insertGlobalValues(GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Insert matrix entries, using global IDs (not implemented)
TpetraBlockCrsMatrix(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rowMap, size_t maxNumEntriesPerRow, const Teuchos::RCP< Teuchos::ParameterList > &params=Teuchos::null)
Constructor specifying fixed number of entries for each row (not implemented)
void setObjectLabel(const std::string &objectLabel)
local_matrix_type::HostMirror getLocalMatrixHost() const
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getRowMap() const
Returns the Map that describes the row distribution in this matrix.
void setAllToScalar(const Scalar &alpha)
Set all matrix entries equal to scalarThis.
std::string description() const
A simple one-line description of this object.
virtual bool haveGlobalConstants() const
Returns true if globalConstants have been computed; false otherwise.
size_t getLocalNumCols() const
Returns the number of columns connected to the locally owned rows of this matrix.
void allocateAllValues(size_t numNonZeros, ArrayRCP< size_t > &rowptr, ArrayRCP< LocalOrdinal > &colind, ArrayRCP< Scalar > &values)
Allocates and returns ArrayRCPs of the Crs arrays — This is an Xpetra-only routine.
void getLocalRowCopy(LocalOrdinal LocalRow, const ArrayView< LocalOrdinal > &Indices, const ArrayView< Scalar > &Values, size_t &NumEntries) const
Extract a list of entries in a specified local row of the matrix. Put into storage allocated by calli...
size_t getGlobalMaxNumRowEntries() const
Returns the maximum number of entries across all rows/columns on all nodes.
bool supportsRowViews() const
Returns true if getLocalRowView() and getGlobalRowView() are valid for this class.
bool isLocallyIndexed() const
If matrix indices are in the local range, this function returns true. Otherwise, this function return...
global_size_t getGlobalNumRows() const
Number of global elements in the row map of this matrix.
size_t getNumEntriesInLocalRow(LocalOrdinal localRow) const
Returns the current number of entries on this node in the specified local row.
void leftScale(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &x)
Left scale operator with given vector values.
void getAllValues(ArrayRCP< const size_t > &rowptr, ArrayRCP< const LocalOrdinal > &colind, ArrayRCP< const Scalar > &values) const
Gets the 1D pointer arrays of the graph (not implemented)
size_t getLocalNumRows() const
Returns the number of matrix rows owned on the calling node.
void getLocalDiagCopy(Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &diag) const
Get a copy of the diagonal entries owned by this node, with local row idices.
bool hasMatrix() const
Does this have an underlying matrix.
void getLocalDiagOffsets(Teuchos::ArrayRCP< size_t > &offsets) const
Get offsets of the diagonal entries in the matrix.
void setAllValues(const ArrayRCP< size_t > &rowptr, const ArrayRCP< LocalOrdinal > &colind, const ArrayRCP< Scalar > &values)
Sets the 1D pointer arrays of the graph (not impelmented)
void replaceDiag(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &diag)
Replace the diagonal entries of the matrix.
global_size_t getGlobalNumEntries() const
Returns the global number of entries in this matrix.
bool isGloballyIndexed() const
If matrix indices are in the global range, this function returns true. Otherwise, this function retur...
void rightScale(const Vector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &x)
Right scale operator with given vector values.
RCP< const Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_BlockCrsMatrix() const
Get the underlying Tpetra matrix.
void scale(const Scalar &alpha)
Scale the current values of a matrix, this = alpha*this (not implemented)
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getColMap() const
Returns the Map that describes the column distribution in this matrix.
void insertLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Insert matrix entries, using local IDs (not implemented)
typename CrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node >::local_matrix_type local_matrix_type
void apply(const MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &X, MultiVector< Scalar, LocalOrdinal, GlobalOrdinal, Node > &Y, Teuchos::ETransp mode=Teuchos::NO_TRANS, Scalar alpha=ScalarTraits< Scalar >::one(), Scalar beta=ScalarTraits< Scalar >::zero()) const
Computes the sparse matrix-multivector multiplication.
RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > getCrsGraph() const
Returns the CrsGraph associated with this matrix.
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getDomainMap() const
Returns the Map associated with the domain of this operator. This will be null until fillComplete() i...
const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getRangeMap() const
Returns the Map associated with the range of this operator, which must be compatible with Y....
size_t getLocalNumEntries() const
Returns the local number of entries in this matrix.
void fillComplete(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const RCP< ParameterList > &params=null)
Signal that data entry is complete, specifying domain and range maps.
ScalarTraits< Scalar >::magnitudeType getFrobeniusNorm() const
Returns the Frobenius norm of the matrix.
void resumeFill(const RCP< ParameterList > &params=null)
void doImport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &source, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Import.
global_size_t getGlobalNumCols() const
Number of global columns in the matrix.
size_t getNumEntriesInGlobalRow(GlobalOrdinal globalRow) const
Returns the current number of entries in the (locally owned) global row.
bool isFillComplete() const
Returns true if the matrix is in compute mode, i.e. if fillComplete() has been called.
void replaceLocalValues(LocalOrdinal localRow, const ArrayView< const LocalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Replace matrix entries, using local IDs.
void getGlobalRowCopy(GlobalOrdinal GlobalRow, const ArrayView< GlobalOrdinal > &indices, const ArrayView< Scalar > &values, size_t &numEntries) const
Extract a list of entries in a specified global row of this matrix. Put into pre-allocated storage.
void expertStaticFillComplete(const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &domainMap, const RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &rangeMap, const RCP< const Import< LocalOrdinal, GlobalOrdinal, Node > > &importer=Teuchos::null, const RCP< const Export< LocalOrdinal, GlobalOrdinal, Node > > &exporter=Teuchos::null, const RCP< ParameterList > &params=Teuchos::null)
Expert static fill complete.
void doExport(const DistObject< char, LocalOrdinal, GlobalOrdinal, Node > &dest, const Import< LocalOrdinal, GlobalOrdinal, Node > &importer, CombineMode CM)
Export.
void replaceGlobalValues(GlobalOrdinal globalRow, const ArrayView< const GlobalOrdinal > &cols, const ArrayView< const Scalar > &vals)
Replace matrix entries, using global IDs (not implemented)
bool isFillActive() const
Returns true if the matrix is in edit mode.
void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const
Print the object with some verbosity level to an FancyOStream object.
void replaceDomainMapAndImporter(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &newDomainMap, Teuchos::RCP< const Import< LocalOrdinal, GlobalOrdinal, Node > > &newImporter)
Replaces the current domainMap and importer with the user-specified objects.
size_t getLocalMaxNumRowEntries() const
Returns the maximum number of entries across all rows/columns on this node.
RCP< Tpetra::BlockCrsMatrix< Scalar, LocalOrdinal, GlobalOrdinal, Node > > getTpetra_BlockCrsMatrixNonConst() const
Get the underlying Tpetra matrix.
void getGlobalRowView(GlobalOrdinal GlobalRow, ArrayView< const GlobalOrdinal > &indices, ArrayView< const Scalar > &values) const
Extract a const, non-persisting view of global indices in a specified row of the matrix.
void getLocalRowView(LocalOrdinal LocalRow, ArrayView< const LocalOrdinal > &indices, ArrayView< const Scalar > &values) const
Extract a const, non-persisting view of local indices in a specified row of the matrix.
Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > getMap() const
Implements DistObject interface.
void removeEmptyProcessesInPlace(const Teuchos::RCP< const Map< LocalOrdinal, GlobalOrdinal, Node > > &newMap)
#define TEUCHOS_UNREACHABLE_RETURN(dummyReturnVal)
Xpetra namespace
size_t global_size_t
Global size_t object.
RCP< const CrsGraph< int, GlobalOrdinal, Node > > toXpetra(const Epetra_CrsGraph &g)
RCP< const Tpetra::CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > toTpetra(const RCP< const CrsGraph< LocalOrdinal, GlobalOrdinal, Node > > &graph)
CombineMode
Xpetra::Combine Mode enumerable type.
static magnitudeType magnitude(T a)