~chaffra/+junk/trilinos

« back to all changes in this revision

Viewing changes to packages/kokkos/test/OskiMatrix/GenerateOskiProblem.hpp

  • Committer: Bazaar Package Importer
  • Author(s): Christophe Prud'homme, Christophe Prud'homme, Johannes Ring
  • Date: 2009-12-13 12:53:22 UTC
  • mfrom: (5.1.2 sid)
  • Revision ID: james.westby@ubuntu.com-20091213125322-in0nrdjc55deqsw9
Tags: 10.0.3.dfsg-1
[Christophe Prud'homme]
* New upstream release

[Johannes Ring]
* debian/patches/libname.patch: Add prefix 'libtrilinos_' to all
  libraries. 
* debian/patches/soname.patch: Add soversion to libraries.
* debian/watch: Update download URL.
* debian/control:
  - Remove python-numeric from Build-Depends (virtual package).
  - Remove automake and autotools from Build-Depends and add cmake to
    reflect switch to CMake.
  - Add python-support to Build-Depends.
* debian/rules: 
  - Cleanup and updates for switch to CMake.

Show diffs side-by-side

added added

removed removed

Lines of Context:
1
 
//@HEADER
2
 
// ************************************************************************
3
 
// 
4
 
//                 Kokkos: A Fast Kernel Package
5
 
//              Copyright (2004) 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
 
// This library is free software; you can redistribute it and/or modify
11
 
// it under the terms of the GNU Lesser General Public License as
12
 
// published by the Free Software Foundation; either version 2.1 of the
13
 
// License, or (at your option) any later version.
14
 
//  
15
 
// This library is distributed in the hope that it will be useful, but
16
 
// WITHOUT ANY WARRANTY; without even the implied warranty of
17
 
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
18
 
// Lesser General Public License for more details.
19
 
//  
20
 
// You should have received a copy of the GNU Lesser General Public
21
 
// License along with this library; if not, write to the Free Software
22
 
// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
23
 
// USA
24
 
// Questions? Contact Michael A. Heroux (maherou@sandia.gov) 
25
 
// 
26
 
// ************************************************************************
27
 
//@HEADER
28
 
 
29
 
#ifndef KOKKOS_TEST_GENERATEOSKIPROBLEM_H
30
 
#define KOKKOS_TEST_GENERATEOSKIPROBLEM_H
31
 
#include "Kokkos_OskiVector.hpp"
32
 
#include "Kokkos_OskiMultiVector.hpp"
33
 
#include "Kokkos_OskiMatrix.hpp"
34
 
 
35
 
namespace KokkosTest {
36
 
 
37
 
 
38
 
//! KokkosTest::GenerateOskiProblem: Constructs a 2D PDE finite difference matrix using the list of x and y offsets.
39
 
  /*!
40
 
 
41
 
  nx      (In) - number of grid points in x direction
42
 
 
43
 
  ny      (In) - number of grid points in y direction
44
 
 
45
 
  The total number of equations will be nx*ny ordered such that the x direction changes
46
 
  most rapidly: 
47
 
  First equation is at point (0,0)
48
 
  Second at                  (1,0)
49
 
  ...
50
 
  nx equation at             (nx-1,0)
51
 
  nx+1st equation at         (0,1)
52
 
 
53
 
  npoints (In) - number of points in finite difference stencil
54
 
  xoff    (In) - stencil offsets in x direction (of length npoints)
55
 
  yoff    (In) - stencil offsets in y direction (of length npoints)
56
 
  A standard 5-point finite difference stencil would be described as:
57
 
  npoints = 5
58
 
  xoff = [-1, 1, 0,  0, 0]
59
 
  yoff = [ 0, 0, 0, -1, 1]
60
 
 
61
 
  nrhs - Number of rhs to generate. (First interface produces vectors, so nrhs is not needed
62
 
 
63
 
  A      (Out) - Kokkos::OskiMatrix constructed for nx by ny grid using prescribed stencil
64
 
  Off-diagonal values are random between 0 and 1.  If diagonal is part of stencil,
65
 
  diagonal will be slightly diag dominant.
66
 
  x      (Out) - Initial guess vector set to zero.
67
 
  b      (Out) - Generated RHS.  Values satisfy b = A*xexact
68
 
  xexact (Out) - Generated exact solution to Ax = b.
69
 
 
70
 
  Note: Caller of this function is responsible for deleting all output objects.
71
 
  */
72
 
 
73
 
  template<typename OrdinalType, typename ScalarType>
74
 
  class GenerateOskiProblem {
75
 
  public:
76
 
 
77
 
    //@{ \name Constructors/Destructor.
78
 
 
79
 
    //! Single RHS constuctor.
80
 
/* don't need this because Kokkos::OskiVector inherits from Kokkos::OskiMultivector
81
 
    GenerateOskiProblem(bool isRowOriented,
82
 
                      OrdinalType nx, OrdinalType ny, OrdinalType npoints, 
83
 
                      OrdinalType * xoff, OrdinalType * yoff,
84
 
                      Kokkos::CisMatrix<OrdinalType, ScalarType> *& A, 
85
 
                      Kokkos::Vector<OrdinalType, ScalarType> *& x, 
86
 
                      Kokkos::Vector<OrdinalType, ScalarType> *& b,
87
 
                      Kokkos::Vector<OrdinalType, ScalarType> *&xexact,
88
 
                      OrdinalType & numEntries);
89
 
*/
90
 
    //! Multi RHS constuctor.
91
 
    GenerateOskiProblem(bool isRowOriented,
92
 
                      OrdinalType nx, OrdinalType ny, OrdinalType npoints, 
93
 
                      OrdinalType * xoff, OrdinalType * yoff, OrdinalType nrhs,
94
 
                      Kokkos::CisMatrix<OrdinalType, ScalarType> *& A, 
95
 
                      Kokkos::MultiVector<OrdinalType, ScalarType> *& x, 
96
 
                      Kokkos::MultiVector<OrdinalType, ScalarType> *& b,
97
 
                      Kokkos::MultiVector<OrdinalType, ScalarType> *&xexact,
98
 
                      OrdinalType & numEntries);
99
 
 
100
 
  public:
101
 
 
102
 
    //! OskiMatrix Destructor
103
 
    virtual ~GenerateOskiProblem();
104
 
    //@}
105
 
  
106
 
   private:
107
 
 
108
 
    //! Copy constructor (Not implemented).
109
 
    GenerateOskiProblem(const GenerateOskiProblem<OrdinalType, ScalarType>& source){};
110
 
 
111
 
/*   void GenerateProblem(bool isRowOriented,
112
 
                           OrdinalType nx, OrdinalType ny, OrdinalType npoints, 
113
 
                           OrdinalType * xoff, OrdinalType * yoff,
114
 
                           Kokkos::CisMatrix<OrdinalType, ScalarType> *& A, 
115
 
                           Kokkos::Vector<OrdinalType, ScalarType> *& x, 
116
 
                           Kokkos::Vector<OrdinalType, ScalarType> *& b,
117
 
                           Kokkos::Vector<OrdinalType, ScalarType> *&xexact,
118
 
                           OrdinalType & numEntries);
119
 
*/
120
 
      void GenerateProblem(bool isRowOriented,
121
 
                       OrdinalType nx, OrdinalType ny, OrdinalType npoints, 
122
 
                       OrdinalType * xoff, OrdinalType * yoff, OrdinalType nrhs,
123
 
                       Kokkos::CisMatrix<OrdinalType, ScalarType> *& A, 
124
 
                       Kokkos::MultiVector<OrdinalType, ScalarType> *& x, 
125
 
                       Kokkos::MultiVector<OrdinalType, ScalarType> *& b,
126
 
                       Kokkos::MultiVector<OrdinalType, ScalarType> *&xexact,
127
 
                             OrdinalType & numEntries);
128
 
 
129
 
  private:
130
 
 
131
 
    Kokkos::OskiVector<OrdinalType, ScalarType> * xd1;
132
 
    Kokkos::OskiVector<OrdinalType, ScalarType> *bd1;
133
 
    Kokkos::OskiVector<OrdinalType, ScalarType> *xexactd1;
134
 
    Kokkos::OskiMultiVector<OrdinalType, ScalarType> * xd;
135
 
    Kokkos::OskiMultiVector<OrdinalType, ScalarType> * bd;
136
 
    Kokkos::OskiMultiVector<OrdinalType, ScalarType> * xexactd;
137
 
    Kokkos::OskiMatrix<OrdinalType, ScalarType> * Ad; // Construct matrix
138
 
    ScalarType ** xv;
139
 
    ScalarType ** bv;
140
 
    ScalarType ** xexactv;
141
 
    OrdinalType numEquations;
142
 
    OrdinalType nrhsv;
143
 
 
144
 
  OrdinalType ** indices;
145
 
  //ScalarType ** values;
146
 
  OrdinalType * pointers;
147
 
  OrdinalType * allIndices;
148
 
  ScalarType * allValues;
149
 
  //OrdinalType * profiles;
150
 
  };
151
 
} // namespace KokkosTest
152
 
 
153
 
using namespace KokkosTest;
154
 
//==============================================================================
155
 
/*template<typename OrdinalType, typename ScalarType>
156
 
GenerateOskiProblem<OrdinalType, ScalarType>::
157
 
GenerateOskiProblem(bool isRowOriented,
158
 
                      OrdinalType nx, OrdinalType ny, OrdinalType npoints, 
159
 
                      OrdinalType * xoff, OrdinalType * yoff,
160
 
                      Kokkos::CisMatrix<OrdinalType, ScalarType> *& A, 
161
 
                      Kokkos::Vector<OrdinalType, ScalarType> *& x, 
162
 
                      Kokkos::Vector<OrdinalType, ScalarType> *& b,
163
 
                      Kokkos::Vector<OrdinalType, ScalarType> *&xexact,
164
 
                      OrdinalType & numEntries):
165
 
  xd1(0),
166
 
  bd1(0),
167
 
  xexactd1(0),
168
 
  xd(0),
169
 
  bd(0),
170
 
  xexactd(0),
171
 
  Ad(0),
172
 
  xv(0),
173
 
  bv(0),
174
 
  xexactv(0),
175
 
  numEquations(0),
176
 
  nrhsv(0),
177
 
  indices(0),
178
 
  values(0),
179
 
  pointers(0),
180
 
  allIndices(0),
181
 
  allValues(0),
182
 
  profiles(0) {
183
 
  GenerateProblem(isRowOriented, nx, ny, npoints, xoff, yoff, A, x, b, xexact, numEntries);
184
 
}
185
 
*/  
186
 
//==============================================================================
187
 
template<typename OrdinalType, typename ScalarType>
188
 
GenerateOskiProblem<OrdinalType, ScalarType>::
189
 
GenerateOskiProblem(bool isRowOriented,
190
 
                      OrdinalType nx, OrdinalType ny, OrdinalType npoints, 
191
 
                      OrdinalType * xoff, OrdinalType * yoff, OrdinalType nrhs,
192
 
                      Kokkos::CisMatrix<OrdinalType, ScalarType> *& A, 
193
 
                      Kokkos::MultiVector<OrdinalType, ScalarType> *& x, 
194
 
                      Kokkos::MultiVector<OrdinalType, ScalarType> *& b,
195
 
                      Kokkos::MultiVector<OrdinalType, ScalarType> *&xexact,
196
 
                      OrdinalType & numEntries):
197
 
  xd1(0),
198
 
  bd1(0),
199
 
  xexactd1(0),
200
 
  xd(0),
201
 
  bd(0),
202
 
  xexactd(0),
203
 
  Ad(0),
204
 
  xv(0),
205
 
  bv(0),
206
 
  xexactv(0),
207
 
  numEquations(0),
208
 
  nrhsv(0),
209
 
  indices(0),
210
 
 // values(0),
211
 
  pointers(0),
212
 
  allIndices(0),
213
 
  allValues(0) {
214
 
  //profiles(0) {
215
 
  GenerateProblem(isRowOriented, nx, ny, npoints, xoff, yoff, nrhs, A, x, b, xexact, numEntries);
216
 
}
217
 
 /* 
218
 
//==============================================================================
219
 
template<typename OrdinalType, typename ScalarType>
220
 
void GenerateOskiProblem<OrdinalType, ScalarType>::
221
 
GenerateProblem(bool isRowOriented,
222
 
                OrdinalType nx, OrdinalType ny, OrdinalType npoints, 
223
 
                OrdinalType * xoff, OrdinalType * yoff,
224
 
                Kokkos::CisMatrix<OrdinalType, ScalarType> *& A, 
225
 
                Kokkos::Vector<OrdinalType, ScalarType> *& x, 
226
 
                Kokkos::Vector<OrdinalType, ScalarType> *& b,
227
 
                Kokkos::Vector<OrdinalType, ScalarType> *&xexact,
228
 
                OrdinalType & numEntries) {
229
 
 
230
 
  Kokkos::MultiVector<OrdinalType, ScalarType> * x1, * b1, * xexact1;
231
 
        
232
 
  GenerateProblem(isRowOriented, nx, ny, npoints, xoff, yoff, 1, A, x1, b1, xexact1, numEntries);
233
 
 
234
 
  xd1 = new Kokkos::OskiVector<OrdinalType, ScalarType>();
235
 
  bd1 = new Kokkos::OskiVector<OrdinalType, ScalarType>();
236
 
  xexactd1 = new Kokkos::OskiVector<OrdinalType, ScalarType>();
237
 
 
238
 
  xd1->initializeValues(x1->getNumRows(), x1->getValues(0), x1->getColInc());
239
 
  bd1->initializeValues(b1->getNumRows(), b1->getValues(0), b1->getColInc());
240
 
  xexactd1->initializeValues(xexact1->getNumRows(), xexact1->getValues(0), xexact1->getColInc());
241
 
 
242
 
  x = dynamic_cast<Kokkos::Vector<OrdinalType, ScalarType> *>(xd1);
243
 
  b = dynamic_cast<Kokkos::Vector<OrdinalType, ScalarType> *>(bd1);
244
 
  xexact = dynamic_cast<Kokkos::Vector<OrdinalType, ScalarType> *>(xexactd1);
245
 
 
246
 
  return;
247
 
}
248
 
*/
249
 
template<typename OrdinalType, typename ScalarType>
250
 
void GenerateOskiProblem<OrdinalType, ScalarType>::
251
 
GenerateProblem(bool isRowOriented,
252
 
                       OrdinalType nx, OrdinalType ny, OrdinalType npoints, 
253
 
                       OrdinalType * xoff, OrdinalType * yoff, OrdinalType nrhs,
254
 
                       Kokkos::CisMatrix<OrdinalType, ScalarType> *& A, 
255
 
                       Kokkos::MultiVector<OrdinalType, ScalarType> *& x, 
256
 
                       Kokkos::MultiVector<OrdinalType, ScalarType> *& b,
257
 
                       Kokkos::MultiVector<OrdinalType, ScalarType> *&xexact,
258
 
                       OrdinalType & numEntries) {
259
 
 
260
 
  // Number of equations is nx*ny.
261
 
  numEquations = nx*ny;
262
 
  nrhsv = nrhs;
263
 
  OrdinalType i, j, k;
264
 
 
265
 
  xv = new ScalarType*[nrhs];
266
 
  bv = new ScalarType*[nrhs];
267
 
  xexactv = new ScalarType*[nrhs];
268
 
 
269
 
  for (i=0; i<nrhs; i++) {
270
 
    xv[i] = new ScalarType[numEquations];
271
 
    bv[i] = new ScalarType[numEquations];
272
 
    xexactv[i] = new ScalarType[numEquations];
273
 
  }
274
 
 
275
 
  for (i=0; i<nrhs; i++)
276
 
    for (j=0; j<numEquations; j++) {
277
 
      xexactv[i][j] = ((ScalarType) rand())/ ((ScalarType) RAND_MAX);
278
 
      xv[i][j] = 0.0;
279
 
      bv[i][j] = 0.0;
280
 
    }
281
 
  indices = 0;
282
 
  //values = 0;
283
 
  pointers = 0;
284
 
  allIndices = 0;
285
 
  allValues = 0;
286
 
  //profiles = 0;
287
 
 
288
 
  allIndices = new OrdinalType[numEquations*npoints];
289
 
  allValues = new ScalarType[numEquations*npoints];
290
 
  pointers = new OrdinalType[numEquations+1];
291
 
 
292
 
  OrdinalType * curIndices;
293
 
  ScalarType * curValues;
294
 
 
295
 
  ScalarType dnpoints = (ScalarType) npoints;
296
 
  numEntries = 0;
297
 
  pointers[0] = 0;
298
 
  for (i=0; i<numEquations; i++) {
299
 
 
300
 
    OrdinalType rowID = i;
301
 
    OrdinalType numIndices = 0;
302
 
    
303
 
   curIndices = allIndices+numEntries;
304
 
   curValues = allValues+numEntries;
305
 
 
306
 
    for (j=0; j<npoints; j++) {
307
 
      OrdinalType colID = rowID + xoff[j] + nx*yoff[j]; // Compute column ID based on stencil offsets
308
 
      if (colID>-1 && colID<numEquations) {
309
 
        curIndices[numIndices] = colID;
310
 
        ScalarType value = - ((ScalarType) rand())/ ((ScalarType) RAND_MAX);
311
 
        if (colID==rowID)
312
 
          curValues[numIndices] = dnpoints - value; // Make diagonal dominant
313
 
        else
314
 
          curValues[numIndices] = -value;
315
 
        if (isRowOriented)
316
 
          for (k=0; k<nrhs; k++)
317
 
            bv[k][i] += curValues[numIndices]*xexactv[k][curIndices[numIndices]];
318
 
        else
319
 
          for (k=0; k<nrhs; k++)
320
 
            bv[k][curIndices[numIndices]] += curValues[numIndices]*xexactv[k][i];
321
 
        numEntries++;
322
 
        numIndices++;
323
 
      }
324
 
    }
325
 
  pointers[i+1] = pointers[i]+numIndices;
326
 
  }
327
 
  // Now construct Kokkos objects and register pointers
328
 
  xd = new Kokkos::OskiMultiVector<OrdinalType, ScalarType>();
329
 
  bd = new Kokkos::OskiMultiVector<OrdinalType, ScalarType>();
330
 
  xexactd = new Kokkos::OskiMultiVector<OrdinalType, ScalarType>();
331
 
  Ad = new Kokkos::OskiMatrix<OrdinalType, ScalarType>(); // Construct matrix
332
 
 
333
 
  xd->initializeValues(numEquations, nrhs, xv[0], numEquations, 1);
334
 
  bd->initializeValues(numEquations, nrhs, bv[0], numEquations, 1);
335
 
  xexactd->initializeValues(numEquations, nrhs, xexactv[0], numEquations, 1);
336
 
 
337
 
  Ad->initializeStructure(numEquations, numEquations, isRowOriented, pointers, allIndices);
338
 
  Ad->initializeValues(allValues);
339
 
 
340
 
  x = dynamic_cast<Kokkos::MultiVector<OrdinalType, ScalarType> *>(xd);
341
 
  b = dynamic_cast<Kokkos::MultiVector<OrdinalType, ScalarType> *>(bd);
342
 
  xexact = dynamic_cast<Kokkos::MultiVector<OrdinalType, ScalarType> *>(xexactd);
343
 
  A = dynamic_cast<Kokkos::CisMatrix<OrdinalType, ScalarType> *>(Ad);
344
 
  return;
345
 
}
346
 
 
347
 
 
348
 
template<typename OrdinalType, typename ScalarType>
349
 
GenerateOskiProblem<OrdinalType, ScalarType>::
350
 
~GenerateOskiProblem() {
351
 
 
352
 
  OrdinalType i;
353
 
  if (xd1!=0) delete xd1;
354
 
  if (bd1!=0) delete bd1;
355
 
  if (xexactd1!=0) delete xexactd1;
356
 
  if (xd!=0) delete xd;
357
 
  if (bd!=0) delete bd;
358
 
  if (xexactd!=0) delete xexactd;
359
 
  if (Ad!=0) delete Ad;
360
 
 
361
 
  if (xv!=0) {
362
 
    for (i=0; i<nrhsv; i++) if (xv[i]!=0) delete [] xv[i];
363
 
    delete [] xv;
364
 
  }
365
 
  if (bv!=0) {
366
 
    for (i=0; i<nrhsv; i++) if (bv[i]!=0) delete [] bv[i];
367
 
    delete [] bv;
368
 
  }
369
 
  if (xexactv!=0) {
370
 
    for (i=0; i<nrhsv; i++) if (xexactv[i]!=0) delete [] xexactv[i];
371
 
    delete [] xexactv;
372
 
  }
373
 
  if (indices!=0) {
374
 
    for (i=0; i<numEquations; i++) if (indices[i]!=0) delete [] indices[i];
375
 
    delete [] indices;
376
 
  }
377
 
/*  if (values!=0) {
378
 
    for (i=0; i<numEquations; i++) if (values[i]!=0) delete [] values[i];
379
 
    delete [] values;
380
 
  }*/
381
 
  if (pointers!=0) delete [] pointers;
382
 
  if (allIndices!=0) delete [] allIndices;
383
 
  if (allValues!=0) delete [] allValues;
384
 
}
385
 
#endif /* KOKKOS_TEST_GENERATEOSKIPROBLEM_H */