Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
GlobalMatrix.cpp
Go to the documentation of this file.
1 ///////////////////////////////////////////////////////////////////////////////
2 //
3 // File: GlobalMatrix.cpp
4 //
5 // For more information, please see: http://www.nektar.info
6 //
7 // The MIT License
8 //
9 // Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
10 // Department of Aeronautics, Imperial College London (UK), and Scientific
11 // Computing and Imaging Institute, University of Utah (USA).
12 //
13 // License for the specific language governing rights and limitations under
14 // Permission is hereby granted, free of charge, to any person obtaining a
15 // copy of this software and associated documentation files (the "Software"),
16 // to deal in the Software without restriction, including without limitation
17 // the rights to use, copy, modify, merge, publish, distribute, sublicense,
18 // and/or sell copies of the Software, and to permit persons to whom the
19 // Software is furnished to do so, subject to the following conditions:
20 //
21 // The above copyright notice and this permission notice shall be included
22 // in all copies or substantial portions of the Software.
23 //
24 // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
25 // OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
26 // FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
27 // THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
28 // LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
29 // FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
30 // DEALINGS IN THE SOFTWARE.
31 //
32 // Description: GlobalMatrix definition
33 //
34 ///////////////////////////////////////////////////////////////////////////////
35 
40 
41 #include <iomanip>
42 #include <fstream>
43 
44 using namespace std;
45 
46 namespace Nektar
47 {
48  namespace MultiRegions
49  {
50  std::string GlobalMatrix::def = LibUtilities::SessionReader::
51  RegisterDefaultSolverInfo("GlobalMatrixStorageType","SmvBSR");
52  std::string GlobalMatrix::lookupIds[3] = {
53  LibUtilities::SessionReader::RegisterEnumValue(
54  "GlobalMatrixStorageType", "SmvBSR", MultiRegions::eSmvBSR)
55  };
56 
57 
58  /**
59  * @class GlobalMatrix
60  * This matrix is essentially a wrapper around a DNekSparseMat.
61  */
62 
63  /**
64  * Allocates a new DNekSparseMat object from the given specification.
65  * @param rows Number of rows in matrix.
66  * @param columns Number of columns in matrix.
67  * @param cooMat ?
68  */
69  GlobalMatrix::GlobalMatrix(
71  unsigned int rows,
72  unsigned int columns,
73  const COOMatType &cooMat,
74  const MatrixStorage& matStorage):
75  m_smvbsrmatrix(),
76  m_rows(rows),
77  m_mulCallsCounter(0)
78  {
79  MatrixStorageType storageType = pSession->
80  GetSolverInfoAsEnum<MatrixStorageType>("GlobalMatrixStorageType");
81 
82  unsigned int brows, bcols;
83 
84  // Size of dense matrix sub-blocks
85  int block_size = 1;
86 
87  BCOMatType bcoMat;
88 
89  // assuming current sparse format allows
90  // block-sparse data representation
91 
92  if(pSession->DefinesParameter("SparseBlockSize"))
93  {
94  pSession->LoadParameter("SparseBlockSize", block_size);
95  ASSERTL1(block_size > 0,"SparseBlockSize parameter must to be positive");
96  }
97 
98  brows = rows / block_size + (rows % block_size > 0);
99  bcols = columns / block_size + (columns % block_size > 0);
100 
101  if (rows % block_size > 0) m_copyOp = true;
102 
103  if (m_copyOp)
104  {
105  m_tmpin = Array<OneD, NekDouble> (brows*block_size, 0.0);
106  m_tmpout = Array<OneD, NekDouble> (brows*block_size, 0.0);
107  }
108 
109  convertCooToBco(brows, bcols, block_size, cooMat, bcoMat);
110 
111  size_t matBytes;
112  switch(storageType)
113  {
114  case eSmvBSR:
115  {
116 
117  // Create zero-based Smv-multiply BSR sparse storage holder
121  brows, bcols, block_size, bcoMat, matStorage );
122 
123  // Create sparse matrix
125  AllocateSharedPtr( sparseStorage );
126 
127  matBytes = m_smvbsrmatrix->GetMemoryFootprint();
128 
129  }
130  break;
131 
132  default:
133  NEKERROR(ErrorUtil::efatal,"Unsupported sparse storage type chosen");
134  }
135 
136  cout << "Global matrix storage type: "
137  << MatrixStorageTypeMap[storageType] << endl;
138  std::cout << "Global matrix memory, bytes = " << matBytes;
139  if (matBytes/(1024*1024) > 0)
140  {
141  std::cout << " ("<< matBytes/(1024*1024) <<" MB)" << std::endl;
142  }
143  else
144  {
145  std::cout << " ("<< matBytes/1024 <<" KB)" << std::endl;
146  }
147  std::cout << "Sparse storage block size = " << block_size << std::endl;
148  }
149 
150  /**
151  * Performs a matrix-vector multiply using the sparse format-specific
152  * multiply routine.
153  * @param in Input vector.
154  * @param out Output vector.
155  */
158  {
159  if (!m_copyOp)
160  {
161  if (m_smvbsrmatrix) m_smvbsrmatrix->Multiply(in,out);
162  }
163  else
164  {
165  // if block size makes the last row/column bigger, one needs
166  // using temporary storage for rhs and result vectors.
167  Vmath::Vcopy(m_rows, &in[0], 1, &m_tmpin[0], 1);
168 
170 
171  Vmath::Vcopy(m_rows, &m_tmpout[0], 1, &out[0], 1);
172  }
173 
175  }
176 
177  const unsigned long GlobalMatrix::GetMulCallsCounter() const
178  {
179  if (m_smvbsrmatrix) return m_smvbsrmatrix->GetMulCallsCounter();
180  return -1;
181  }
182 
183  const unsigned int GlobalMatrix::GetNumNonZeroEntries() const
184  {
185  if (m_smvbsrmatrix) return m_smvbsrmatrix->GetNumNonZeroEntries();
186  return -1;
187  }
188 
189 
190  } //end of namespace
191 } //end of namespace
192 
boost::shared_ptr< SparseStorageType > SparseStorageSharedPtr
#define NEKERROR(type, msg)
Assert Level 0 – Fundamental assert which is used whether in FULLDEBUG, DEBUG or OPT compilation mod...
Definition: ErrorUtil.hpp:191
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
std::map< CoordType, NekDouble > COOMatType
Array< OneD, NekDouble > m_tmpout
Definition: GlobalMatrix.h:78
const unsigned int GetNumNonZeroEntries() const
const char *const MatrixStorageTypeMap[]
STL namespace.
const unsigned long GetMulCallsCounter() const
boost::shared_ptr< SessionReader > SessionReaderSharedPtr
Definition: MeshPartition.h:51
void Multiply(const Array< OneD, const NekDouble > &in, Array< OneD, NekDouble > &out)
Perform a matrix-vector multiply.
Array< OneD, NekDouble > m_tmpin
Definition: GlobalMatrix.h:77
std::map< CoordType, BCOEntryType > BCOMatType
DNekSmvBsrMatSharedPtr m_smvbsrmatrix
Pointer to a double-precision Nektar++ sparse matrix.
Definition: GlobalMatrix.h:74
void convertCooToBco(const unsigned int blkRows, const unsigned int blkColumns, const unsigned int blkDim, const COOMatType &cooMat, BCOMatType &bcoMat)
Definition: SparseUtils.cpp:48
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:228
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1061