Galeri Package Browser (Single Doxygen Collection)
Version of the Day
Loading...
Searching...
No Matches
example
VbrMatrix.cpp
Go to the documentation of this file.
1
// @HEADER
2
// ************************************************************************
3
//
4
// Galeri: Finite Element and Matrix Generation Package
5
// Copyright (2006) ETHZ/Sandia Corporation
6
//
7
// Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8
// license for use of this work by or on behalf of the U.S. Government.
9
//
10
// Redistribution and use in source and binary forms, with or without
11
// modification, are permitted provided that the following conditions are
12
// met:
13
//
14
// 1. Redistributions of source code must retain the above copyright
15
// notice, this list of conditions and the following disclaimer.
16
//
17
// 2. Redistributions in binary form must reproduce the above copyright
18
// notice, this list of conditions and the following disclaimer in the
19
// documentation and/or other materials provided with the distribution.
20
//
21
// 3. Neither the name of the Corporation nor the names of the
22
// contributors may be used to endorse or promote products derived from
23
// this software without specific prior written permission.
24
//
25
// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26
// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27
// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28
// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29
// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30
// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31
// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32
// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33
// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34
// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35
// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36
//
37
// Questions about Galeri? Contact Marzio Sala (marzio.sala _AT_ gmail.com)
38
//
39
// ************************************************************************
40
// @HEADER
41
42
#include "Galeri_Maps.h"
43
#include "Galeri_CrsMatrices.h"
44
#include "Galeri_VbrMatrices.h"
45
#ifdef HAVE_MPI
46
#include "
Epetra_MpiComm.h
"
47
#include "mpi.h"
48
#else
49
#include "
Epetra_SerialComm.h
"
50
#endif
51
#include "
Epetra_Map.h
"
52
#include "
Epetra_CrsMatrix.h
"
53
#include "
Epetra_VbrMatrix.h
"
54
#include "Teuchos_ParameterList.hpp"
55
56
using namespace
Galeri;
57
58
// =========== //
59
// main driver //
60
// =========== //
61
62
int
main
(
int
argc,
char
* argv[])
63
{
64
#ifdef HAVE_MPI
65
MPI_Init(&argc, &argv);
66
Epetra_MpiComm
Comm(MPI_COMM_WORLD);
67
#else
68
Epetra_SerialComm
Comm;
69
#endif
70
71
// pointer to the map to be created
72
Epetra_Map
* Map;
73
// pointer to the CrsMatrix to be created
74
Epetra_CrsMatrix
* CrsMatrix;
75
// pointer to the VbrMatrix to be created
76
Epetra_VbrMatrix
* VbrMatrix;
77
// container for parameters
78
Teuchos::ParameterList GaleriList;
79
// here we specify the global dimension of the problem
80
GaleriList.set(
"n"
, 2 * Comm.NumProc());
81
82
try
83
{
84
// Creates a simple linear map; for more details on the map creation
85
// refer to the documentation
86
Map = CreateMap(
"Linear"
, Comm, GaleriList);
87
88
// Creates a diagonal matrix with 1's on the diagonal
89
CrsMatrix = CreateCrsMatrix(
"Diag"
, Map, GaleriList);
90
91
// at this point we can create the VbrMatrix. Each block of the Vbr will
92
// be 2 x 2. (Any positive value is accepted, like 1.)
93
VbrMatrix = CreateVbrMatrix(CrsMatrix, 2);
94
95
// print out the matrix
96
cout << *VbrMatrix;
97
98
// To created objects must be free'd using delete
99
delete
Map;
100
delete
CrsMatrix;
101
delete
VbrMatrix;
102
}
103
catch
(Exception& rhs)
104
{
105
if
(Comm.MyPID() == 0)
106
{
107
cerr <<
"Caught exception: "
;
108
rhs.
Print
();
109
}
110
}
111
112
#ifdef HAVE_MPI
113
MPI_Finalize();
114
#endif
115
116
return
(EXIT_SUCCESS);
117
}
Epetra_CrsMatrix.h
Epetra_Map.h
Epetra_MpiComm.h
Epetra_SerialComm.h
Epetra_VbrMatrix.h
main
int main(int argc, char *argv[])
Definition
VbrMatrix.cpp:62
Epetra_CrsMatrix
Epetra_CrsMatrix::Print
virtual void Print(std::ostream &os) const
Epetra_Map
Epetra_MpiComm
Epetra_SerialComm
Epetra_VbrMatrix
Generated by
1.10.0