Stokhos Package Browser (Single Doxygen Collection) Version of the Day
Loading...
Searching...
No Matches
Stokhos_BasisInteractionGraph.cpp
Go to the documentation of this file.
1// $Id$
2// $Source$
3// @HEADER
4// ***********************************************************************
5//
6// Stokhos Package
7// Copyright (2009) Sandia Corporation
8//
9// Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
10// license for use of this work by or on behalf of the U.S. Government.
11//
12// Redistribution and use in source and binary forms, with or without
13// modification, are permitted provided that the following conditions are
14// met:
15//
16// 1. Redistributions of source code must retain the above copyright
17// notice, this list of conditions and the following disclaimer.
18//
19// 2. Redistributions in binary form must reproduce the above copyright
20// notice, this list of conditions and the following disclaimer in the
21// documentation and/or other materials provided with the distribution.
22//
23// 3. Neither the name of the Corporation nor the names of the
24// contributors may be used to endorse or promote products derived from
25// this software without specific prior written permission.
26//
27// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
28// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
29// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
30// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
31// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
32// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
33// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
34// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
35// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
36// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
37// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
38//
39// Questions? Contact Eric T. Phipps (etphipp@sandia.gov).
40//
41// ***********************************************************************
42// @HEADER
43
45
48
50 : vecLookup_(flm.vecLookup_), onlyUseLinear_(flm.onlyUseLinear_)
51{ }
52
54 : onlyUseLinear_(onlyUseLinear)
55{
56 using Teuchos::RCP;
57
58 RCP<const Stokhos::Sparse3Tensor<int,double> > Cijk;
59 if(porder<0)
60 Cijk = max_basis.computeTripleProductTensor();
61 else
62 Cijk = max_basis.computeLinearTripleProductTensor();
63 initialize(max_basis,*Cijk,porder);
64}
65
67 const Stokhos::ProductBasis<int,double> & rowBasis,
68 const Stokhos::ProductBasis<int,double> & colBasis,bool onlyUseLinear,int porder)
69 : onlyUseLinear_(onlyUseLinear)
70{
71 using Teuchos::RCP;
72
73 RCP<const Stokhos::Sparse3Tensor<int,double> > Cijk;
74 if(porder<0)
75 Cijk = masterBasis.computeTripleProductTensor();
76 else
77 Cijk = masterBasis.computeLinearTripleProductTensor();
78 initialize(masterBasis,*Cijk,rowBasis,colBasis,porder);
79}
80
83 bool onlyUseLinear,int porder)
84 : onlyUseLinear_(onlyUseLinear)
85{
86 initialize(max_basis,Cijk,porder);
87}
88
91 const Stokhos::ProductBasis<int,double> & rowBasis,
92 const Stokhos::ProductBasis<int,double> & colBasis,bool onlyUseLinear,int porder)
93 : onlyUseLinear_(onlyUseLinear)
94{
95 using Teuchos::RCP;
96
97 initialize(masterBasis,Cijk,rowBasis,colBasis,porder);
98}
99
102 int porder)
103{
104 using Teuchos::RCP;
105 typedef Stokhos::Sparse3Tensor<int,double> Cijk_type;
106
107 // // max it out if defaulted
108 // if(porder<0)
109 // porder = max_basis.size();
110
111 // RCP<Stokhos::Sparse3Tensor<int,double> > Cijk = max_basis.computeTripleProductTensor(porder);
112
113 Cijk_type::k_iterator k_end = Cijk.k_end();
114 if (onlyUseLinear_) {
115 int dim = max_basis.dimension();
116 k_end = Cijk.find_k(dim+1);
117 }
118
119 vecLookup_.resize(max_basis.size()); // defines number of rows
120 numCols_ = vecLookup_.size(); // set number of columns
121
122 // Loop over Cijk entries including a non-zero in the graph at
123 // indices (i,j) if there is any k for which Cijk is non-zero
124 for(Cijk_type::k_iterator k_it=Cijk.k_begin(); k_it!=k_end; ++k_it) {
125 for(Cijk_type::kj_iterator j_it = Cijk.j_begin(k_it); j_it != Cijk.j_end(k_it); ++j_it) {
126 int j = index(j_it);
127 for(Cijk_type::kji_iterator i_it = Cijk.i_begin(j_it); i_it != Cijk.i_end(j_it); ++i_it) {
128 int i = index(i_it);
129 vecLookup_[i].push_back(j);
130 }
131 }
132 }
133}
134
137 const Stokhos::ProductBasis<int,double> & rowBasis,
138 const Stokhos::ProductBasis<int,double> & colBasis,int porder)
139{
140 // for determining if their is an interaction or not
141 Stokhos::BasisInteractionGraph masterBig(masterBasis,Cijk,onlyUseLinear_,porder);
142
143 vecLookup_.resize(rowBasis.size()); // defines number of rows
144
145 // set number of columns
146 numCols_ = colBasis.size();
147
148 // build row basis terms
149 std::vector<int> rowIndexToMasterIndex(rowBasis.size());
150 for(int i=0;i<rowBasis.size();i++)
151 rowIndexToMasterIndex[i] = masterBasis.index(rowBasis.term(i));
152
153 // build column basis terms
154 std::vector<int> colIndexToMasterIndex(colBasis.size());
155 for(int i=0;i<colBasis.size();i++)
156 colIndexToMasterIndex[i] = masterBasis.index(colBasis.term(i));
157
158 // build graph by looking up sparsity in master basis
159 for(int r=0;r<rowBasis.size();r++) {
160 int masterRow = rowIndexToMasterIndex[r];
161 for(int c=0;c<colBasis.size();c++) {
162 int masterCol = colIndexToMasterIndex[c];
163
164 // is row and column active in master element
165 bool activeRC = masterBig(masterRow,masterCol);
166
167 // if active add to local graph
168 if(activeRC)
169 vecLookup_[r].push_back(c);
170 }
171 }
172}
173
174bool Stokhos::BasisInteractionGraph::operator()(std::size_t i,std::size_t j) const
175{
176 const std::vector<std::size_t> & indices = activeIndices(i);
177 return indices.end() != std::find(indices.begin(),indices.end(),j);
178}
179
181{
182 for(std::size_t r=0;r<rowCount();r++) {
183 for(std::size_t c=0;c<colCount();c++)
184 if(operator()(r,c)) os << " * ";
185 else os << " ";
186 os << std::endl;
187 }
188}
189
191{
192 std::size_t nnz = 0;
193 for(std::size_t r=0;r<rowCount();r++)
194 nnz += activeIndices(r).size();
195 return nnz;
196}
void initialize(const Stokhos::OrthogPolyBasis< int, double > &max_basis, const Stokhos::Sparse3Tensor< int, double > &Cijk, int porder=-1)
Setup the lookup graph.
bool operator()(std::size_t i, std::size_t j) const
Is there an entry for (i,j) in the graph.
std::size_t numNonZeros() const
How many non zeros are in this graph.
Abstract base class for multivariate orthogonal polynomials.
Abstract base class for multivariate orthogonal polynomials generated from tensor products of univari...
Data structure storing a sparse 3-tensor C(i,j,k) in a a compressed format.
kj_iterator j_end(const k_iterator &k) const
Iterator pointing to last j entry for given k.
k_iterator k_begin() const
Iterator pointing to first k entry.
kji_iterator i_begin(const kj_iterator &j) const
Iterator pointing to first i entry for given j and k.
kj_iterator j_begin(const k_iterator &k) const
Iterator pointing to first j entry for given k.
k_iterator find_k(ordinal_type k) const
Return k iterator for given index k.
kji_iterator i_end(const kj_iterator &j) const
Iterator pointing to last i entry for given j and k.
k_iterator k_end() const
Iterator pointing to last k entry.