Epetra Package Browser (Single Doxygen Collection) Development
Loading...
Searching...
No Matches
test/CrsRectMatrix/cxx_main.cpp
Go to the documentation of this file.
1//@HEADER
2// ************************************************************************
3//
4// Epetra: Linear Algebra Services Package
5// Copyright 2011 Sandia Corporation
6//
7// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
8// the U.S. Government retains certain rights in this software.
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? Contact Michael A. Heroux (maherou@sandia.gov)
38//
39// ************************************************************************
40//@HEADER
41
42
43#include "Epetra_LocalMap.h"
44#include "Epetra_Map.h"
45#include "Epetra_Time.h"
46#include "Epetra_Vector.h"
47#include "Epetra_CrsMatrix.h"
48#include "Epetra_Flops.h"
49#include "Epetra_Export.h"
50#include "Epetra_Import.h"
53#ifdef EPETRA_MPI
54#include "Epetra_MpiComm.h"
55#include "mpi.h"
56#else
57#include "Epetra_SerialComm.h"
58#endif
59#include "../epetra_test_err.h"
60#include "Epetra_Version.h"
61
62int main(int argc, char *argv[])
63{
64 int ierr = 0, i, j, forierr = 0;
65 int ntrials = 1;
66
67#ifdef EPETRA_MPI
68
69 // Initialize MPI
70
71 MPI_Init(&argc,&argv);
72 Epetra_MpiComm Comm( MPI_COMM_WORLD );
73
74#else
75
77
78#endif
79
80 bool verbose = false;
81
82 // Check if we should print results to standard out
83 if (argc>1) if (argv[1][0]=='-' && argv[1][1]=='v') verbose = true;
84
85
86
87 // char tmp; if (Comm.MyPID()==0) { cout << "Press any key to continue..."<< endl; cin >> tmp;} Comm.Barrier();
88
89 Comm.SetTracebackMode(0); // This should shut down any error traceback reporting
90 int MyPID = Comm.MyPID();
91 int NumProc = Comm.NumProc();
92
93 if (verbose && Comm.MyPID()==0)
94 cout << Epetra_Version() << endl << endl;
95
96 if (verbose) cout << "Processor "<<MyPID<<" of "<< NumProc
97 << " is alive."<<endl;
98
99 bool verbose1 = verbose;
100
101 // Redefine verbose to only print on PE 0
102 if (verbose && Comm.MyPID()!=0) verbose = false;
103
104 int NumMyEquations = 10000;
105
106 int NumGlobalEquations = NumMyEquations*NumProc;
107 int NumGlobalVariables = 2 * NumGlobalEquations+1;
108
109 // Construct a Map that puts approximately the same Number of equations on each processor
110
111 Epetra_Map RowMap(NumGlobalEquations, 0, Comm);
112 Epetra_Map XMap(NumGlobalVariables, 0, Comm);
113 Epetra_Map& YMap = RowMap;
114
115 // Get update list and number of local equations from newly created Map
116 int * MyGlobalElements = new int[RowMap.NumMyElements()];
117 RowMap.MyGlobalElements(MyGlobalElements);
118
119 // Get update list and number of local equations from newly created XMap
120 int * XGlobalElements = new int[XMap.NumMyElements()];
121 XMap.MyGlobalElements(XGlobalElements);
122
123 // Get update list and number of local variables from newly created YMap
124 int * YGlobalElements = new int[YMap.NumMyElements()];
125 YMap.MyGlobalElements(YGlobalElements);
126
127 // We need vectors to compute:
128 // X = A^T*Y
129 // AATY = A*A^T*Y = A*X
130 // and
131 // BY = B*Y
132
133 Epetra_Vector Y(YMap);
134 Epetra_Vector X(XMap);
135 Epetra_Vector AATY(YMap);
136 Epetra_Vector BY(YMap);
137
138
139 // Fill Y Vector
140 Y.Random();
141 //Y.PutScalar(1.0);
142
143 // To create A^T explicitly we need an assembly map that is two elements longer than
144 // the XMap, because each processor will be making contributions to two rows beyond what
145 // it will own.
146 int ATAssemblyNumMyElements = 2*MyGlobalElements[NumMyEquations-1] + 2 - 2*MyGlobalElements[0] + 1;
147 int * ATAssemblyGlobalElements = new int[ATAssemblyNumMyElements];
148
149 for (i=0; i<ATAssemblyNumMyElements; i++) ATAssemblyGlobalElements[i] = 2*MyGlobalElements[0] + i;
150 Epetra_Map ATAssemblyMap(-1, ATAssemblyNumMyElements, ATAssemblyGlobalElements, 0, Comm);
151
152 // Create a Epetra_Matrix with the values of A
153 // A is a simple 1D weighted average operator that mimics a restriction operator
154 // that might be found in a multigrid code.
155 // Also create A^T explicitly
156
157 Epetra_CrsMatrix A(Copy, RowMap, 3);
158 Epetra_CrsMatrix ATAssembly(Copy, ATAssemblyMap, 2);
159 Epetra_CrsMatrix AT(Copy, XMap, 2);
160
161 //cout << "ATAssemblyMap = "<< endl<< ATAssemblyMap << endl
162 // << "XMap = " << endl << XMap << endl
163 // << "RowMap = " << endl << RowMap << endl;
164 // Add rows one-at-a-time
165 // Need some vectors to help
166 // Off diagonal Values will always be -1
167
168
169 double *Values = new double[3];
170 int *Indices = new int[3];
171 int NumEntries;
172 /*
173 Values[0] = 0.25;
174 Values[1] = 0.5;
175 Values[2] = 0.25;
176 */
177 Values[0] = 0.5;
178 Values[1] = 0.25;
179 Values[2] = 0.25;
180 forierr = 0;
181 for (i=0; i<NumMyEquations; i++)
182 {
183 /*
184 Indices[0] = 2*MyGlobalElements[i];
185 Indices[1] = 2*MyGlobalElements[i]+1;
186 Indices[2] = 2*MyGlobalElements[i]+2;
187 */
188 Indices[0] = 2*MyGlobalElements[i]+1;
189 Indices[1] = 2*MyGlobalElements[i]+2;
190 Indices[2] = 2*MyGlobalElements[i];
191 NumEntries = 3;
192 forierr += !(A.InsertGlobalValues(MyGlobalElements[i], NumEntries, Values, Indices)==0);
193 for (j=0; j<3; j++)
194 forierr += !(ATAssembly.InsertGlobalValues(Indices[j],1, &(Values[j]), &(MyGlobalElements[i]))>=0);
195 }
196 EPETRA_TEST_ERR(forierr,ierr);
197
198
199 EPETRA_TEST_ERR(!(ATAssembly.FillComplete()==0),ierr);
200 // Gather AT values from ATAssembly matrix
201 Epetra_Export Exporter(ATAssemblyMap, XMap);
202 EPETRA_TEST_ERR(!(AT.Export(ATAssembly, Exporter, Add)==0),ierr);
203
204 // Finish up
205 EPETRA_TEST_ERR(!(A.FillComplete(XMap, YMap)==0),ierr);
206 EPETRA_TEST_ERR(!(AT.FillComplete(YMap, XMap)==0),ierr);
207
208
209 if (verbose1 && NumGlobalEquations<20) {
210 if (verbose) cout << "\n\n Matrix A\n" << endl;
211 cout << A << endl;
212 if (verbose) cout << " \n\n Matrix A Transpose\n" << endl;
213 cout << AT << endl;
214 }
215
216
217 // Create a Epetra_Matrix containing B = A*A^T.
218 // This matrix will be a square tridiagonal matrix. We will use it to compare the results
219 // of A*(A^T*X) using two methods: (1) with two calls to Multiply using A^T and then A and
220 // (2) using B directly.
221
222 Epetra_CrsMatrix B(Copy, RowMap, 3);
223
224 Values[0] = 1.0/16.0;
225 Values[1] = 3.0/8.0;
226 Values[2] = 1.0/16.0;
227 int Valstart;
228 forierr = 0;
229 for (i=0; i<NumMyEquations; i++)
230 {
231 if (MyGlobalElements[i] == 0) {
232 Indices[0] = MyGlobalElements[i];
233 Indices[1] = MyGlobalElements[i]+1;
234 NumEntries = 2;
235 Valstart = 1;
236 }
237 else {
238 Indices[0] = MyGlobalElements[i]-1;
239 Indices[1] = MyGlobalElements[i];
240 Indices[2] = MyGlobalElements[i]+1;
241 NumEntries = 3;
242 Valstart = 0;
243 }
244 if (MyGlobalElements[i] == NumGlobalEquations-1) NumEntries--;
245 forierr += !(B.InsertGlobalValues(MyGlobalElements[i], NumEntries, Values+Valstart, Indices)==0);
246 }
247 EPETRA_TEST_ERR(forierr,ierr);
248
249 // Finish up
250 EPETRA_TEST_ERR(!(B.FillComplete()==0),ierr);
251 if (verbose && NumGlobalEquations<20) cout << "\n\nMatrix B \n" << endl;
252 if (verbose1 && NumGlobalEquations<20) cout << B << endl;
253
254
255 Epetra_Flops counter;
256 A.SetFlopCounter(counter);
257 B.SetFlopCounter(A);
258 Epetra_Time timer(Comm);
259 for (i=0; i<ntrials; i++) {
260 EPETRA_TEST_ERR(!(B.Multiply(false, Y, BY)==0),ierr); // Compute BY = B*Y
261 }
262 double elapsed_time = timer.ElapsedTime();
263 double total_flops = B.Flops();
264 counter.ResetFlops();
265
266 double MFLOPs = total_flops/elapsed_time/1000000.0;
267
268 if (verbose) cout << "\n\nTotal MFLOPs for B*Y = " << MFLOPs << endl<< endl;
269 if (verbose && NumGlobalEquations<20) cout << "\n\nVector Z = B*Y \n";
270 if (verbose1 && NumGlobalEquations<20) cout << BY << endl;
271
272
274
275 timer.ResetStartTime();
276 for (i=0; i<ntrials; i++) {
277 EPETRA_TEST_ERR(!(A.Multiply(true, Y, X)==0),ierr); // Compute X = A^T*Y
278 }
279 elapsed_time = timer.ElapsedTime();
280 total_flops = A.Flops();
281 counter.ResetFlops();
282 MFLOPs = total_flops/elapsed_time/1000000.0;
283
284 if (verbose) cout << "\n\nTotal MFLOPs for A^T*Y using A and trans=true = " << MFLOPs << endl<< endl;
285 if (verbose && NumGlobalEquations<20) cout << "\n\nVector Z = AT*Y \n";
286 if (verbose1 && NumGlobalEquations<20) cout << X << endl;
287
289
290 timer.ResetStartTime();
291 EPETRA_TEST_ERR(!(A.Multiply(false, X, AATY)==0),ierr); // Compute AATY = A*X
292 elapsed_time = timer.ElapsedTime();
293 total_flops = A.Flops();
294 MFLOPs = total_flops/elapsed_time/1000000.0;
295 counter.ResetFlops();
296 Epetra_Vector resid(YMap);
297 resid.Update(1.0, BY, -1.0, AATY, 0.0);
298 double residual;
299 resid.Norm2(&residual);
300
301 if (verbose) cout << "\n\nTotal MFLOPs for A*X using A and trans=false = " << MFLOPs << endl<< endl;
302 if (verbose) cout << "Residual = " << residual << endl<< endl;
303 if (verbose && NumGlobalEquations<20) cout << "\n\nVector Z = A*ATY \n";
304 if (verbose1 && NumGlobalEquations<20) cout << AATY << endl;
305
307
308 AT.SetFlopCounter(counter);
309 timer.ResetStartTime();
310 for (i=0; i<ntrials; i++) {
311 EPETRA_TEST_ERR(!(AT.Multiply(false, Y, X)==0),ierr); // Compute X = A^T*Y
312 }
313 elapsed_time = timer.ElapsedTime();
314 total_flops = AT.Flops();
315 counter.ResetFlops();
316 MFLOPs = total_flops/elapsed_time/1000000.0;
317
318 if (verbose) cout << "\n\nTotal MFLOPs for A^T*Y using AT and trans=false = " << MFLOPs << endl<< endl;
319 if (verbose && NumGlobalEquations<20) cout << "\n\nVector Z = AT*Y \n";
320 if (verbose1 && NumGlobalEquations<20) cout << X << endl;
321
323
324 timer.ResetStartTime();
325 for (i=0; i<ntrials; i++) {
326 EPETRA_TEST_ERR(!(AT.Multiply(true, X, AATY)==0),ierr); // Compute AATY = A*X
327 }
328 elapsed_time = timer.ElapsedTime();
329 total_flops = AT.Flops();
330 MFLOPs = total_flops/elapsed_time/1000000.0;
331 counter.ResetFlops();
332 resid.Update(1.0, BY, -1.0, AATY, 0.0);
333 resid.Norm2(&residual);
334
335 if (verbose) cout << "\n\nTotal MFLOPs for A*X using AT and trans=true = " << MFLOPs << endl<< endl;
336 if (verbose) cout << "Residual = " << residual << endl<< endl;
337 if (verbose && NumGlobalEquations<20) cout << "\n\nVector Z = A*ATY \n";
338 if (verbose1 && NumGlobalEquations<20) cout <<AATY << endl;
339
340
342 // Now test use of Epetra_LocalMap vectors: First test case of local replicated range vector
343
344 {
345 Epetra_CrsMatrix AL(Copy, RowMap, 3);
346 for (i=0; i<NumMyEquations; i++)
347 {
348 forierr += !(A.ExtractGlobalRowCopy(MyGlobalElements[i], 3, NumEntries, Values, Indices)==0);
349 forierr += !(AL.InsertGlobalValues(MyGlobalElements[i], NumEntries, Values, Indices)==0);
350 }
351 EPETRA_TEST_ERR(forierr,ierr);
352
353 Epetra_LocalMap YLMap(NumGlobalEquations, 0, Comm);
354 EPETRA_TEST_ERR(!(AL.FillComplete(XMap, YLMap)==0),ierr);
355 AL.SetFlopCounter(A);
356 Epetra_Vector YL(YLMap);
357 Epetra_Vector ALX(YLMap);
358
359 timer.ResetStartTime();
360 for (i=0; i<ntrials; i++) {
361 EPETRA_TEST_ERR(!(A.Multiply(false, X, Y)==0),ierr); // Compute Y= A*X
362 }
363 elapsed_time = timer.ElapsedTime();
364 total_flops = A.Flops();
365 counter.ResetFlops();
366 MFLOPs = total_flops/elapsed_time/1000000.0;
367
368
369
370 if (verbose) cout << "\n\nTotal MFLOPs for Y=A*X using global distributed Y = " << MFLOPs << endl<< endl;
371 if (verbose && NumGlobalEquations<20) cout << "\n\nVector Y = A*X using distributed Y \n";
372 if (verbose1 && NumGlobalEquations<20) cout << Y << endl;
373 if (verbose && NumGlobalEquations<20) cout << "\n\nA using dist Y range map\n";
374 if (verbose1 && NumGlobalEquations<20) cout << A << endl;
375
376 timer.ResetStartTime();
377 for (i=0; i<ntrials; i++) {
378 EPETRA_TEST_ERR(!(AL.Multiply(false, X, ALX)==0),ierr); // Compute YL= A*X
379 }
380 elapsed_time = timer.ElapsedTime();
381 total_flops = AL.Flops();
382 counter.ResetFlops();
383 MFLOPs = total_flops/elapsed_time/1000000.0;
384
385 if (verbose) cout << "\n\nTotal MFLOPs for Y=A*X using Local replicated Y = " << MFLOPs << endl<< endl;
386 if (verbose && NumGlobalEquations<20) cout << "\n\nVector YL = AL*X using local replicated Y \n";
387 if (verbose1 && NumGlobalEquations<20) cout << ALX << endl;
388 if (verbose && NumGlobalEquations<20) cout << "\n\nA using local Y range map\n";
389 if (verbose1 && NumGlobalEquations<20) cout << AL << endl;
390
391 // Now gather Y values from the distributed Y and compare them to the local replicated Y values
392 Epetra_Import g2limporter(YLMap, YMap); // Import from distributed Y map to local Y map
393 EPETRA_TEST_ERR(!(YL.Import(Y, g2limporter, Insert)==0),ierr);
394 if (verbose && NumGlobalEquations<20) cout << "\n\nVector YL = imported from distributed Y \n";
395 if (verbose1 && NumGlobalEquations<20) cout << YL << endl;
396 EPETRA_TEST_ERR(!(YL.Update(-1.0, ALX, 1.0)==0),ierr);
397 EPETRA_TEST_ERR(!(YL.Norm2(&residual)==0),ierr);
398 if (verbose) cout << "Residual = " << residual << endl<< endl;
399
400
401 //
402 // Multiply by transpose
403 //
404
405 timer.ResetStartTime();
406 for (i=0; i<ntrials; i++) {
407 EPETRA_TEST_ERR(!(A.Multiply(true, Y, X)==0),ierr); // Compute X = A^TY
408 }
409 elapsed_time = timer.ElapsedTime();
410 total_flops = A.Flops();
411 counter.ResetFlops();
412 MFLOPs = total_flops/elapsed_time/1000000.0;
413
414
415
416 if (verbose) cout << "\n\nTotal MFLOPs for X=A^TY using global distributed Y = " << MFLOPs << endl<< endl;
417 if (verbose && NumGlobalEquations<20) cout << "\n\nVector X using distributed Y \n";
418 if (verbose1 && NumGlobalEquations<20) cout << Y << endl;
419 if (verbose && NumGlobalEquations<20) cout << "\n\nA using dist Y range map\n";
420 if (verbose1 && NumGlobalEquations<20) cout << A << endl;
421
422 Epetra_Vector X1(XMap);
423
424 timer.ResetStartTime();
425 for (i=0; i<ntrials; i++) {
426 EPETRA_TEST_ERR(!(AL.Multiply(true, ALX, X1)==0),ierr); // Compute X1 = AL^T*Y
427 }
428 elapsed_time = timer.ElapsedTime();
429 total_flops = AL.Flops();
430 counter.ResetFlops();
431 MFLOPs = total_flops/elapsed_time/1000000.0;
432
433 if (verbose) cout << "\n\nTotal MFLOPs for X1=A^T*Y using Local replicated Y = " << MFLOPs << endl<< endl;
434 if (verbose && NumGlobalEquations<20) cout << "\n\nVector X1 using local replicated Y \n";
435 if (verbose1 && NumGlobalEquations<20) cout << X1 << endl;
436
437 EPETRA_TEST_ERR(!(X1.Update(-1.0, X, 1.0)==0),ierr);
438 EPETRA_TEST_ERR(!(X1.Norm2(&residual)==0),ierr);
439 if (verbose) cout << "Residual = " << residual << endl<< endl;
440 }
441
443 // Finally test use of Epetra_LocalMap vectors using local replicated domain vector
444
445 {
446 Epetra_CrsMatrix AL(Copy, RowMap, 3);
447 for (i=0; i<NumMyEquations; i++)
448 {
449 forierr += !(A.ExtractGlobalRowCopy(MyGlobalElements[i], 3, NumEntries, Values, Indices)==0);
450 forierr += !(AL.InsertGlobalValues(MyGlobalElements[i], NumEntries, Values, Indices)==0);
451 }
452 EPETRA_TEST_ERR(forierr,ierr);
453
454 Epetra_LocalMap XLMap(NumGlobalVariables, 0, Comm);
455 EPETRA_TEST_ERR(!(AL.FillComplete(XLMap, YMap)==0),ierr);
456 AL.SetFlopCounter(A);
457 Epetra_Vector XL(XLMap);
458 Epetra_Vector ALX(XLMap);
459
460 timer.ResetStartTime();
461 for (i=0; i<ntrials; i++) {
462 EPETRA_TEST_ERR(!(A.Multiply(false, X, Y)==0),ierr); // Compute Y= A*X
463 }
464 elapsed_time = timer.ElapsedTime();
465 total_flops = A.Flops();
466 counter.ResetFlops();
467 MFLOPs = total_flops/elapsed_time/1000000.0;
468
469
470
471 if (verbose) cout << "\n\nTotal MFLOPs for Y=A*X using global distributed X = " << MFLOPs << endl<< endl;
472 if (verbose && NumGlobalEquations<20) cout << "\n\nVector Y = A*X using distributed X \n";
473 if (verbose1 && NumGlobalEquations<20) cout << Y << endl;
474 //if (verbose && NumGlobalEquations<20) cout << "\n\nA using dist X range map\n";
475 //if (verbose1 && NumGlobalEquations<20) cout << A << endl;
476
477 // Now gather X values from the distributed X
478 Epetra_Import g2limporter(XLMap, XMap); // Import from distributed X map to local X map
479 EPETRA_TEST_ERR(!(XL.Import(X, g2limporter, Insert)==0),ierr);
480 if (verbose && NumGlobalEquations<20) cout << "\n\nVector XL = imported from distributed X \n";
481 if (verbose1 && NumGlobalEquations<20) cout << XL << endl;
482 Epetra_Vector Y1(Y);
483 timer.ResetStartTime();
484 for (i=0; i<ntrials; i++) {
485 EPETRA_TEST_ERR(!(AL.Multiply(false, XL, Y1)==0),ierr); // Compute Y1= AL*XL
486 }
487 elapsed_time = timer.ElapsedTime();
488 total_flops = AL.Flops();
489 counter.ResetFlops();
490 MFLOPs = total_flops/elapsed_time/1000000.0;
491
492 if (verbose) cout << "\n\nTotal MFLOPs for Y1=A*XL using Local replicated X = " << MFLOPs << endl<< endl;
493 if (verbose && NumGlobalEquations<20) cout << "\n\nVector Y1 = AL*XL using local replicated X \n";
494 if (verbose1 && NumGlobalEquations<20) cout << Y1 << endl;
495 //if (verbose && NumGlobalEquations<20) cout << "\n\nA using local X domain map\n";
496 //if (verbose1 && NumGlobalEquations<20) cout << AL << endl;
497
498 EPETRA_TEST_ERR(!(Y1.Update(-1.0, Y, 1.0)==0),ierr);
499 EPETRA_TEST_ERR(!(Y1.Norm2(&residual)==0),ierr);
500 if (verbose) cout << "Residual = " << residual << endl<< endl;
501
502
503 //
504 // Multiply by transpose
505 //
506
507 timer.ResetStartTime();
508 for (i=0; i<ntrials; i++) {
509 EPETRA_TEST_ERR(!(A.Multiply(true, Y, X)==0),ierr); // Compute X = A^TY
510 }
511 elapsed_time = timer.ElapsedTime();
512 total_flops = A.Flops();
513 counter.ResetFlops();
514 MFLOPs = total_flops/elapsed_time/1000000.0;
515
516
517
518 if (verbose) cout << "\n\nTotal MFLOPs for X=A^TY using global distributed X = " << MFLOPs << endl<< endl;
519 if (verbose && NumGlobalEquations<20) cout << "\n\nVector X using distributed X \n";
520 if (verbose1 && NumGlobalEquations<20) cout << X << endl;
521 //if (verbose && NumGlobalEquations<20) cout << "\n\nA using dist X domain map\n";
522 //if (verbose1 && NumGlobalEquations<20) cout << A << endl;
523
524 timer.ResetStartTime();
525 for (i=0; i<ntrials; i++) {
526 EPETRA_TEST_ERR(!(AL.Multiply(true, Y, XL)==0),ierr); // Compute XL = AL^T*Y
527 }
528 elapsed_time = timer.ElapsedTime();
529 total_flops = AL.Flops();
530 counter.ResetFlops();
531 MFLOPs = total_flops/elapsed_time/1000000.0;
532
533 if (verbose) cout << "\n\nTotal MFLOPs for XL=A^T*Y1 using Local replicated X = " << MFLOPs << endl<< endl;
534 if (verbose && NumGlobalEquations<20) cout << "\n\nVector XL using local replicated X \n";
535 if (verbose1 && NumGlobalEquations<20) cout << XL << endl;
536
537 Epetra_Vector XL1(XLMap);
538 EPETRA_TEST_ERR(!(XL1.Import(X, g2limporter, Insert)==0),ierr);
539 EPETRA_TEST_ERR(!(XL1.Update(-1.0, XL, 1.0)==0),ierr);
540 EPETRA_TEST_ERR(!(XL1.Norm2(&residual)==0),ierr);
541 if (verbose) cout << "Residual = " << residual << endl<< endl;
542 }
543 // Release all objects
544 delete [] Values;
545 delete [] Indices;
546 delete [] MyGlobalElements;
547 delete [] XGlobalElements;
548 delete [] YGlobalElements;
549 delete [] ATAssemblyGlobalElements;
550
551
552#ifdef EPETRA_MPI
553 MPI_Finalize() ;
554#endif
555
556/* end main
557*/
558return ierr ;
559}
560
std::string Epetra_Version()
int MyGlobalElements(int *MyGlobalElementList) const
Puts list of global elements on this processor into the user-provided array.
int NumMyElements() const
Number of elements on the calling processor.
void SetFlopCounter(const Epetra_Flops &FlopCounter_in)
Set the internal Epetra_Flops() pointer.
double Flops() const
Returns the number of floating point operations with this multi-vector.
Epetra_CrsMatrix: A class for constructing and using real-valued double-precision sparse compressed r...
int FillComplete(bool OptimizeDataStorage=true)
Signal that data entry is complete. Perform transformations to local index space.
int ExtractGlobalRowCopy(int GlobalRow, int Length, int &NumEntries, double *Values, int *Indices) const
Returns a copy of the specified global row in user-provided arrays.
int Multiply(bool TransA, const Epetra_Vector &x, Epetra_Vector &y) const
Returns the result of a Epetra_CrsMatrix multiplied by a Epetra_Vector x in y.
virtual int InsertGlobalValues(int GlobalRow, int NumEntries, const double *Values, const int *Indices)
Insert a list of elements in a given global row of the matrix.
int Import(const Epetra_SrcDistObject &A, const Epetra_Import &Importer, Epetra_CombineMode CombineMode, const Epetra_OffsetIndex *Indexor=0)
Imports an Epetra_DistObject using the Epetra_Import object.
int Export(const Epetra_SrcDistObject &A, const Epetra_Import &Importer, Epetra_CombineMode CombineMode, const Epetra_OffsetIndex *Indexor=0)
Exports an Epetra_DistObject using the Epetra_Import object.
Epetra_Export: This class builds an export object for efficient exporting of off-processor elements.
Epetra_Flops: The Epetra Floating Point Operations Class.
void ResetFlops()
Resets the number of floating point operations to zero for this multi-vector.
Epetra_Import: This class builds an import object for efficient importing of off-processor elements.
Epetra_LocalMap: A class for replicating vectors and matrices across multiple processors.
Epetra_Map: A class for partitioning vectors and matrices.
Definition Epetra_Map.h:119
Epetra_MpiComm: The Epetra MPI Communication Class.
int NumProc() const
Returns total number of processes.
int MyPID() const
Return my process ID.
int Update(double ScalarA, const Epetra_MultiVector &A, double ScalarThis)
Update multi-vector values with scaled values of A, this = ScalarThis*this + ScalarA*A.
int Random()
Set multi-vector values to random numbers.
int Norm2(double *Result) const
Compute 2-norm of each vector in multi-vector.
static void SetTracebackMode(int TracebackModeValue)
Set the value of the Epetra_Object error traceback report mode.
Epetra_SerialComm: The Epetra Serial Communication Class.
Epetra_Time: The Epetra Timing Class.
Definition Epetra_Time.h:75
void ResetStartTime(void)
Epetra_Time function to reset the start time for a timer object.
double ElapsedTime(void) const
Epetra_Time elapsed time function.
Epetra_Vector: A class for constructing and using dense vectors on a parallel computer.
#define EPETRA_TEST_ERR(a, b)
int main(int argc, char *argv[])