Amesos Package Browser (Single Doxygen Collection)  Development
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
CreateTridi.cpp
Go to the documentation of this file.
1 // @HEADER
2 // ***********************************************************************
3 //
4 // Amesos: Direct Sparse Solver 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., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301
23 // USA
24 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
25 //
26 // ***********************************************************************
27 // @HEADER
28 
29 //
30 // CreateTridi populates an empty EpetraCrsMatrix with a tridiagonal with
31 // -1 on the off-diagonals and 2 on the diagonal.
32 //
33 // CreateTridiPlus creates the same matrix as CreateTridi except that it adds
34 // -1 in the two off diagonal corners.
35 //
36 // This code was plaguerized from epetra/example/petra_power_method/cxx_main.cpp
37 // presumably written by Mike Heroux.
38 //
39 // Adapted by Ken Stanley - Aug 2003
40 //
41 #include "Epetra_Map.h"
42 #include "Epetra_CrsMatrix.h"
43 
45 {
46 
47  Epetra_Map Map = A.RowMap();
48  int NumMyElements = Map.NumMyElements();
49  int NumGlobalElements = Map.NumGlobalElements();
50 
51  int * MyGlobalElements = new int[NumMyElements];
52  Map.MyGlobalElements(MyGlobalElements);
53 
54  // Add rows one-at-a-time
55  // Need some vectors to help
56  // Off diagonal Values will always be -1
57 
58 
59  double *Values = new double[3];
60  int *Indices = new int[3];
61  int NumEntries;
62 
63  for (int i=0; i<NumMyElements; i++)
64  {
65  if (MyGlobalElements[i]==0)
66  {
67  Indices[0] = 0;
68  Indices[1] = 1;
69  Values[0] = 2.0;
70  Values[1] = -1.0;
71  NumEntries = 2;
72  }
73  else if (MyGlobalElements[i] == NumGlobalElements-1)
74  {
75  Indices[0] = NumGlobalElements-1;
76  Indices[1] = NumGlobalElements-2;
77  Values[0] = 2.0;
78  Values[1] = -1.0;
79  NumEntries = 2;
80  }
81  else
82  {
83  Indices[0] = MyGlobalElements[i]-1;
84  Indices[1] = MyGlobalElements[i];
85  Indices[2] = MyGlobalElements[i]+1;
86  Values[0] = -1.0;
87  Values[1] = 2.0;
88  Values[2] = -1.0;
89  NumEntries = 3;
90  }
91 
92  assert(A.InsertGlobalValues(MyGlobalElements[i], NumEntries, Values, Indices)==0);
93  // Put in the diagonal entry
94  // assert(A.InsertGlobalValues(MyGlobalElements[i], 1, &two, &MyGlobalElements[i])==0);
95  }
96 
97  // Finish up
98  assert(A.FillComplete()==0);
99 
100 
101  delete[] MyGlobalElements;
102  delete[] Values;
103  delete[] Indices;
104  return 0;
105 }
106 
108 {
109 
110  Epetra_Map Map = A.RowMap();
111  int NumMyElements = Map.NumMyElements();
112  int NumGlobalElements = Map.NumGlobalElements();
113 
114  int * MyGlobalElements = new int[NumMyElements];
115  Map.MyGlobalElements(MyGlobalElements);
116 
117  // Add rows one-at-a-time
118  // Need some vectors to help
119  // Off diagonal Values will always be -1
120 
121 
122  double *Values = new double[3];
123  int *Indices = new int[3];
124  int NumEntries;
125 
126  for (int i=0; i<NumMyElements; i++)
127  {
128  if (MyGlobalElements[i]==0)
129  {
130  Indices[0] = 0;
131  Indices[1] = 1;
132  Indices[2] = NumGlobalElements-1;
133  Values[0] = 2.0;
134  Values[1] = -1.0;
135  Values[2] = -0.5;
136  NumEntries = 3;
137  }
138  else if (MyGlobalElements[i] == NumGlobalElements-1)
139  {
140  Indices[0] = NumGlobalElements-1;
141  Indices[1] = NumGlobalElements-2;
142  Indices[2] = 0;
143  Values[0] = 2.0;
144  Values[1] = -1.0;
145  Values[2] = -0.5;
146  NumEntries = 3;
147  }
148  else
149  {
150  Indices[0] = MyGlobalElements[i]-1;
151  Indices[1] = MyGlobalElements[i];
152  Indices[2] = MyGlobalElements[i]+1;
153  Values[0] = -1.0;
154  Values[1] = 2.0;
155  Values[2] = -1.0;
156  NumEntries = 3;
157  }
158 
159  assert(A.InsertGlobalValues(MyGlobalElements[i], NumEntries, Values, Indices)==0);
160  // Put in the diagonal entry
161  // assert(A.InsertGlobalValues(MyGlobalElements[i], 1, &two, &MyGlobalElements[i])==0);
162  }
163 
164  // Finish up
165  assert(A.FillComplete()==0);
166 
167 
168  delete[] MyGlobalElements;
169  delete[] Values;
170  delete[] Indices;
171  return 0;
172 }
173 
int NumGlobalElements() const
int MyGlobalElements(int *MyGlobalElementList) const
virtual int InsertGlobalValues(int GlobalRow, int NumEntries, const double *Values, const int *Indices)
int CreateTridi(Epetra_CrsMatrix &A)
Definition: CreateTridi.cpp:44
int FillComplete(bool OptimizeDataStorage=true)
const Epetra_Map & RowMap() const
int NumMyElements() const
int CreateTridiPlus(Epetra_CrsMatrix &A)