Epetra Package Browser (Single Doxygen Collection) Development
Loading...
Searching...
No Matches
FEVbrMatrix/ExecuteTestProblems.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_BLAS.h"
44#include "ExecuteTestProblems.h"
45#include "Epetra_Comm.h"
46#include "Epetra_Vector.h"
47#include "Epetra_Map.h"
48#include "Epetra_FEVbrMatrix.h"
49#include <Epetra_FEVector.h>
50
51#include "../epetra_test_err.h"
52
53
54int quad1(const Epetra_Map& map, bool verbose)
55{
56 const Epetra_Comm & Comm = map.Comm();
57 int numProcs = Comm.NumProc();
58 int localProc = Comm.MyPID();
59
60 Comm.Barrier();
61 if (verbose && localProc == 0) {
62 cout << "====================== quad1 =============================="<<endl;
63 }
64
65 //Set up a simple finite-element mesh containing 2-D quad elements, 1 per proc.
66 //
67 // *-----*-----*-----*
68 // 0| 2| 4| 6|
69 // | 0 | 1 | np-1|
70 // | | | |
71 // *-----*-----*-----*
72 // 1 3 5 7
73 //
74 //In the above drawing, 'np' means num-procs. node-numbers are to the
75 //lower-left of each node (*).
76 //
77 //Each processor owns element 'localProc', and each processor owns
78 //nodes 'localProc*2' and 'localProc*2+1' except for the last processor,
79 //which also owns the last two nodes.
80 //
81 //There will be 3 degrees-of-freedom per node, so each element-matrix is
82 //of size 12x12. (4 nodes per element, X 3 dof per node)
83 //
84
85 int myFirstNode = localProc*2;
86 int myLastNode = localProc*2+1;
87 if (localProc == numProcs-1) {
88 myLastNode += 2;
89 }
90
91 int numMyNodes = myLastNode - myFirstNode + 1;
92 int* myNodes = new int[numMyNodes];
93 int i, j;
94 int ierr = 0;
95 for(i=0; i<numMyNodes; ++i) {
96 myNodes[i] = myFirstNode + i;
97 }
98
99 int dofPerNode = 3; //degrees-of-freedom per node
100 int indexBase = 0;
101
102 Epetra_BlockMap blkMap(-1, numMyNodes, myNodes, dofPerNode,
103 indexBase, Comm);
104
105 int rowLengths = 3; //each element-matrix will have 4 block-columns.
106 //the rows of the assembled matrix will be longer than
107 //this, but we don't need to worry about that because the
108 //VbrMatrix will add memory as needed. For a real
109 //application where efficiency is a concern, better
110 //performance would be obtained by giving more accurate
111 //row-lengths here.
112
113 Epetra_FEVbrMatrix A(Copy, blkMap, rowLengths);
114
115 int nodesPerElem = 4;
116 int* elemNodes = new int[nodesPerElem];
117 for(i=0; i<nodesPerElem; ++i) elemNodes[i] = myFirstNode+i;
118
119 int elemMatrixDim = nodesPerElem*dofPerNode;
120 int len = elemMatrixDim*elemMatrixDim;
121 double* elemMatrix = new double[len];
122
123 //In an actual finite-element problem, we would calculate and fill
124 //meaningful element stiffness matrices. But for this simple matrix assembly
125 //test, we're just going to fill our element matrix with 1.0's. This will
126 //make it easy to see whether the matrix is correct after it's assembled.
127
128 for(i=0; i<len; ++i) elemMatrix[i] = 1.0;
129
130 //For filling in the matrix block-entries, we would ordinarily have to
131 //carefully copy, or set pointers to, appropriate sections of the
132 //element-matrix. But for this simple case we know that the element-matrix
133 //is all 1's, so we'll just set our block-entry pointer to point to the
134 //beginning of the element-matrix and leave it at that.
135 //Note that the matrix class will refer to dofPerNode X dofPerNode (==9)
136 //positions in the memory pointed to by 'blockEntry'.
137
138 double* blockEntry = elemMatrix;
139
140 //The element-matrix is a 4x4 (nodesPerElem X nodesPerElem) matrix of
141 //3x3 block-entries. We'll now load our element-matrix into the global
142 //matrix by looping over it and loading block-entries individually.
143
144 for(i=0; i<nodesPerElem; ++i) {
145 int blkrow = myFirstNode+i;
146 EPETRA_TEST_ERR( A.BeginInsertGlobalValues(blkrow, nodesPerElem, elemNodes),
147 ierr);
148
149 for(j=0; j<nodesPerElem; ++j) {
150 for(int ii=0; ii<dofPerNode*dofPerNode; ii++) {
151 blockEntry[ii] = blkrow+elemNodes[j];
152 }
153 EPETRA_TEST_ERR( A.SubmitBlockEntry( blockEntry, dofPerNode,
154 dofPerNode, dofPerNode), ierr);
155 }
156
157 int err = A.EndSubmitEntries();
158 if (err < 0) {
159 cout << "quad1: error in A.EndSubmitEntries: "<<err<<endl;
160 return(-1);
161 }
162 }
163
165
166 if (verbose && localProc==0) {
167 cout << "after globalAssemble"<<endl;
168 }
169 if (verbose) {
170 A.Print(cout);
171 }
172
173 int numMyRows = A.NumMyRows();
174 int correct_numMyRows = dofPerNode*numMyNodes;
175
176 if (numMyRows != correct_numMyRows) {
177 cout << "proc " << localProc << ", numMyRows("<<numMyRows<<") doesn't match"
178 << " correct_numMyRows("<<correct_numMyRows<<")."<<endl;
179 return(-1);
180 }
181
182 int numMyNonzeros = A.NumMyNonzeros();
183 int correct_numMyNonzeros = nodesPerElem*nodesPerElem*dofPerNode*dofPerNode;
184
185 if (numProcs > 1) {
186 if (localProc == numProcs-1) {
187 correct_numMyNonzeros += dofPerNode*dofPerNode*4;
188 }
189 else if (localProc > 0) {
190 correct_numMyNonzeros -= dofPerNode*dofPerNode*4;
191 }
192 else {
193 //localProc==0 && numProcs > 1
194 correct_numMyNonzeros -= dofPerNode*dofPerNode*8;
195 }
196 }
197
198 if (numMyNonzeros != correct_numMyNonzeros) {
199 cout << "proc " << localProc << ", numMyNonzeros(" << numMyNonzeros
200 <<") != correct_numMyNonzeros("<<correct_numMyNonzeros<<")"<<endl;
201 return(-1);
202 }
203
204
205 delete [] elemMatrix;
206 delete [] myNodes;
207 delete [] elemNodes;
208
209 Comm.Barrier();
210
211 return(0);
212}
213
214int quad2(const Epetra_Map& map, bool verbose)
215{
216 const Epetra_Comm & Comm = map.Comm();
217 int numProcs = Comm.NumProc();
218 int localProc = Comm.MyPID();
219 Comm.Barrier();
220 if (verbose && localProc == 0) {
221 cout << "====================== quad2 =============================="<<endl;
222 }
223
224 //Set up a simple finite-element mesh containing 2-D quad elements, 2 per proc.
225 //(This test is similar to quad1() above, except for having 2 elements-per-proc
226 // rather than 1.)
227 //
228 // *-----*-----*-----*-------*
229 // 0| 2| 4| 6| 8|
230 // | 0 | 1 | 2 | 2*np-1|
231 // | | | | |
232 // *-----*-----*-----*-------*
233 // 1 3 5 7 9
234 //
235 //In the above drawing, 'np' means num-procs. node-numbers are to the
236 //lower-left of each node (*).
237 //
238 //Each processor owns element 'localProc' and 'localProc+1', and each processor
239 //owns nodes 'localProc*4' through 'localProc*4+3' except for the last
240 //processor, which also owns the last two nodes in the mesh.
241 //
242 //There will be 3 degrees-of-freedom per node, so each element-matrix is
243 //of size 12x12.
244 //
245
246 int myFirstNode = localProc*4;
247 int myLastNode = localProc*4+3;
248 if (localProc == numProcs-1) {
249 myLastNode += 2;
250 }
251
252 int numMyElems = 2;
253
254 int numMyNodes = myLastNode - myFirstNode + 1;
255 int* myNodes = new int[numMyNodes];
256 int i, j;
257 int ierr = 0;
258 for(i=0; i<numMyNodes; ++i) {
259 myNodes[i] = myFirstNode + i;
260 }
261
262 int dofPerNode = 3; //degrees-of-freedom per node
263 int indexBase = 0;
264
265 Epetra_BlockMap blkMap(-1, numMyNodes, myNodes, dofPerNode,
266 indexBase, Comm);
267
268 int rowLengths = 4; //each element-matrix will have 4 block-columns.
269 //the rows of the assembled matrix will be longer than
270 //this, but we don't need to worry about that because the
271 //VbrMatrix will add memory as needed. For a real
272 //application where efficiency is a concern, better
273 //performance would be obtained by giving a more accurate
274 //row-length here.
275
276 Epetra_FEVbrMatrix A(Copy, blkMap, rowLengths);
277
278 int nodesPerElem = 4;
279 int* elemNodes = new int[nodesPerElem];
280
281 int elemMatrixDim = nodesPerElem*dofPerNode;
282 int len = elemMatrixDim*elemMatrixDim;
283 double* elemMatrix = new double[len];
284
285 //In an actual finite-element problem, we would calculate and fill
286 //meaningful element stiffness matrices. But for this simple matrix assembly
287 //test, we're just going to fill our element matrix with 1.0's. This will
288 //make it easy to see whether the matrix is correct after it's assembled.
289
290 for(i=0; i<len; ++i) elemMatrix[i] = 1.0;
291
292 //For filling in the matrix block-entries, we would ordinarily have to
293 //carefully copy, or set pointers to, appropriate sections of the
294 //element-matrix. But for this simple case we know that the element-matrix
295 //is all 1's, so we'll just set our block-entry pointer to point to the
296 //beginning of the element-matrix and leave it at that.
297 //Note that the matrix class will refer to dofPerNode X dofPerNode (==9)
298 //positions in the memory pointed to by 'blockEntry'.
299
300 double* blockEntry = elemMatrix;
301
302 //Each element-matrix is a 4x4 (nodesPerElem X nodesPerElem) matrix of
303 //3x3 block-entries. We'll now load our element-matrices into the global
304 //matrix by looping over them and loading block-entries individually.
305
306 int firstNode = myFirstNode;
307 for(int el=0; el<numMyElems; ++el) {
308 for(i=0; i<nodesPerElem; ++i) {
309
310 for(int n=0; n<nodesPerElem; ++n) elemNodes[n] = firstNode+n;
311
312 int blkrow = firstNode+i;
313 EPETRA_TEST_ERR( A.BeginInsertGlobalValues(blkrow, nodesPerElem, elemNodes),
314 ierr);
315
316 for(j=0; j<nodesPerElem; ++j) {
317 EPETRA_TEST_ERR( A.SubmitBlockEntry( blockEntry, dofPerNode,
318 dofPerNode, dofPerNode), ierr);
319 }
320
321 int this_err = A.EndSubmitEntries();
322 if (this_err < 0) {
323 cerr << "error in quad2, A.EndSubmitEntries(): " << this_err << endl;
324 return(this_err);
325 }
326 }
327
328 firstNode += 2;
329 }
330
332
333 if (verbose && localProc==0) {
334 cout << "after globalAssemble"<<endl;
335 }
336
337 if (verbose) {
338 A.Print(cout);
339 }
340
341 //now let's make sure that we can perform a matvec...
342 Epetra_FEVector x(blkMap, 1), y(blkMap, 1);
343 x.PutScalar(1.0);
344
345 EPETRA_TEST_ERR( A.Multiply(false, x, y), ierr);
346
347 if (verbose && localProc==0) {
348 cout << "quad2, y:"<<endl;
349 }
350 if (verbose) {
351 y.Print(cout);
352 }
353
354 delete [] elemMatrix;
355 delete [] myNodes;
356 delete [] elemNodes;
357
358 return(0);
359}
360
361int MultiVectorTests(const Epetra_Map & Map, int NumVectors, bool verbose)
362{
363 const Epetra_Comm & Comm = Map.Comm();
364 int ierr = 0, i, j;
365
366 /* get number of processors and the name of this processor */
367
368 int MyPID = Comm.MyPID();
369
370 // Construct FEVbrMatrix
371
372 if (verbose && MyPID==0) cout << "constructing Epetra_FEVbrMatrix" << endl;
373
374 //
375 //we'll set up a tri-diagonal matrix.
376 //
377
378 int numGlobalRows = Map.NumGlobalElements();
379 int minLocalRow = Map.MinMyGID();
380 int rowLengths = 3;
381
382 Epetra_FEVbrMatrix A(Copy, Map, rowLengths);
383
384 if (verbose && MyPID==0) {
385 cout << "calling A.InsertGlobalValues with 1-D data array"<<endl;
386 }
387
388 int numCols = 3;
389 int* ptIndices = new int[numCols];
390 for(int k=0; k<numCols; ++k) {
391 ptIndices[k] = minLocalRow+k;
392 }
393
394 double* values_1d = new double[numCols*numCols];
395 for(j=0; j<numCols*numCols; ++j) {
396 values_1d[j] = 3.0;
397 }
398
399 //For an extreme test, we'll have all processors sum into all rows.
400
401 int minGID = Map.MinAllGID();
402
403 //For now we're going to assume that there's just one point associated with
404 //each GID (element).
405
406 double* ptCoefs = new double[3];
407
408 {for(i=0; i<numGlobalRows; ++i) {
409 if (i>0 && i<numGlobalRows-1) {
410 ptIndices[0] = minGID+i-1;
411 ptIndices[1] = minGID+i;
412 ptIndices[2] = minGID+i+1;
413 ptCoefs[0] = -1.0;
414 ptCoefs[1] = 2.0;
415 ptCoefs[2] = -1.0;
416 numCols = 3;
417 }
418 else if (i == 0) {
419 ptIndices[0] = minGID+i;
420 ptIndices[1] = minGID+i+1;
421 ptIndices[2] = minGID+i+2;
422 ptCoefs[0] = 2.0;
423 ptCoefs[1] = -1.0;
424 ptCoefs[2] = -1.0;
425 numCols = 3;
426 }
427 else {
428 ptIndices[0] = minGID+i-2;
429 ptIndices[1] = minGID+i-1;
430 ptIndices[2] = minGID+i;
431 ptCoefs[0] = -1.0;
432 ptCoefs[1] = -1.0;
433 ptCoefs[2] = 2.0;
434 numCols = 3;
435 }
436
437 int row = minGID+i;
438
439 EPETRA_TEST_ERR( A.BeginInsertGlobalValues(row, rowLengths, ptIndices), ierr);
440
441 for(j=0; j<rowLengths; ++j) {
442 EPETRA_TEST_ERR( A.SubmitBlockEntry(&(ptCoefs[j]), 1, 1, 1), ierr);
443 }
444
446
447 }}
448
449 if (verbose&&MyPID==0) {
450 cout << "calling A.GlobalAssemble()" << endl;
451 }
452
453 EPETRA_TEST_ERR( A.GlobalAssemble(), ierr );
454
455 if (verbose&&MyPID==0) {
456 cout << "after globalAssemble"<<endl;
457 }
458 if (verbose) {
459 A.Print(cout);
460 }
461
462 delete [] values_1d;
463 delete [] ptIndices;
464 delete [] ptCoefs;
465
466 return(ierr);
467}
468
469int four_quads(const Epetra_Comm& Comm, bool preconstruct_graph, bool verbose)
470{
471 if (verbose) {
472 cout << "******************* four_quads ***********************"<<endl;
473 }
474
475 //This function assembles a matrix representing a finite-element mesh
476 //of four 2-D quad elements. There are 9 nodes in the problem. The
477 //same problem is assembled no matter how many processors are being used
478 //(within reason). It may not work if more than 9 processors are used.
479 //
480 // *------*------*
481 // 6| 7| 8|
482 // | E2 | E3 |
483 // *------*------*
484 // 3| 4| 5|
485 // | E0 | E1 |
486 // *------*------*
487 // 0 1 2
488 //
489 //Nodes are denoted by * with node-numbers below and left of each node.
490 //E0, E1 and so on are element-numbers.
491 //
492 //Each processor will contribute a sub-matrix of size 4x4, filled with 1's,
493 //for each element. Thus, the coefficient value at position 0,0 should end up
494 //being 1.0*numProcs, the value at position 4,4 should be 1.0*4*numProcs, etc.
495 //
496 //Depending on the number of processors being used, the locations of the
497 //specific matrix positions (in terms of which processor owns them) will vary.
498 //
499
500 int numProcs = Comm.NumProc();
501
502 int numNodes = 9;
503 int numElems = 4;
504 int numNodesPerElem = 4;
505
506 int blockSize = 1;
507 int indexBase = 0;
508
509 //Create a map using epetra-defined linear distribution.
510 Epetra_BlockMap map(numNodes, blockSize, indexBase, Comm);
511
512 Epetra_CrsGraph* graph = NULL;
513
514 int* nodes = new int[numNodesPerElem];
515 int i, j, k, err = 0;
516
517 if (preconstruct_graph) {
518 graph = new Epetra_CrsGraph(Copy, map, 1);
519
520 //we're going to fill the graph with indices, but remember it will only
521 //accept indices in rows for which map.MyGID(row) is true.
522
523 for(i=0; i<numElems; ++i) {
524 switch(i) {
525 case 0:
526 nodes[0] = 0; nodes[1] = 1; nodes[2] = 4; nodes[3] = 3;
527 break;
528 case 1:
529 nodes[0] = 1; nodes[1] = 2; nodes[2] = 5; nodes[3] = 4;
530 break;
531 case 2:
532 nodes[0] = 3; nodes[1] = 4; nodes[2] = 7; nodes[3] = 6;
533 break;
534 case 3:
535 nodes[0] = 4; nodes[1] = 5; nodes[2] = 8; nodes[3] = 7;
536 break;
537 }
538
539 for(j=0; j<numNodesPerElem; ++j) {
540 if (map.MyGID(nodes[j])) {
541 err = graph->InsertGlobalIndices(nodes[j], numNodesPerElem,
542 nodes);
543 if (err<0) return(err);
544 }
545 }
546 }
547
548 EPETRA_CHK_ERR( graph->FillComplete() );
549 }
550
551 Epetra_FEVbrMatrix* A = NULL;
552
553 if (preconstruct_graph) {
554 A = new Epetra_FEVbrMatrix(Copy, *graph);
555 }
556 else {
557 A = new Epetra_FEVbrMatrix(Copy, map, 1);
558 }
559
560 //EPETRA_CHK_ERR( A->PutScalar(0.0) );
561
562 double* values_1d = new double[numNodesPerElem*numNodesPerElem];
563 double** values_2d = new double*[numNodesPerElem];
564
565 for(i=0; i<numNodesPerElem*numNodesPerElem; ++i) values_1d[i] = 1.0;
566
567 int offset = 0;
568 for(i=0; i<numNodesPerElem; ++i) {
569 values_2d[i] = &(values_1d[offset]);
570 offset += numNodesPerElem;
571 }
572
573 for(i=0; i<numElems; ++i) {
574 switch(i) {
575 case 0:
576 nodes[0] = 0; nodes[1] = 1; nodes[2] = 4; nodes[3] = 3;
577 break;
578
579 case 1:
580 nodes[0] = 1; nodes[1] = 2; nodes[2] = 5; nodes[3] = 4;
581 break;
582
583 case 2:
584 nodes[0] = 3; nodes[1] = 4; nodes[2] = 7; nodes[3] = 6;
585 break;
586
587 case 3:
588 nodes[0] = 4; nodes[1] = 5; nodes[2] = 8; nodes[3] = 7;
589 break;
590 }
591
592 for(j=0; j<numNodesPerElem; ++j) {
593 if (preconstruct_graph) {
594 err = A->BeginSumIntoGlobalValues(nodes[j], numNodesPerElem, nodes);
595 if (err<0) return(err);
596 }
597 else {
598 err = A->BeginInsertGlobalValues(nodes[j], numNodesPerElem, nodes);
599 if (err<0) return(err);
600 }
601
602 for(k=0; k<numNodesPerElem; ++k) {
603 err = A->SubmitBlockEntry(values_1d, blockSize, blockSize, blockSize);
604 if (err<0) return(err);
605 }
606
607 err = A->EndSubmitEntries();
608 if (err<0) return(err);
609 }
610 }
611
613
614 Epetra_FEVbrMatrix* Acopy = new Epetra_FEVbrMatrix(*A);
615
616 if (verbose) {
617 cout << "A:"<<*A << endl;
618 cout << "Acopy:"<<*Acopy<<endl;
619 }
620
621 Epetra_Vector x(A->RowMap()), y(A->RowMap());
622
623 x.PutScalar(1.0); y.PutScalar(0.0);
624
625 Epetra_Vector x2(Acopy->RowMap()), y2(Acopy->RowMap());
626
627 x2.PutScalar(1.0); y2.PutScalar(0.0);
628
629 A->Multiply(false, x, y);
630
631 Acopy->Multiply(false, x2, y2);
632
633 double ynorm2, y2norm2;
634
635 y.Norm2(&ynorm2);
636 y2.Norm2(&y2norm2);
637 if (ynorm2 != y2norm2) {
638 cerr << "norm2(A*ones) != norm2(*Acopy*ones)"<<endl;
639 return(-99);
640 }
641
642 Epetra_FEVbrMatrix* Acopy2 =
643 new Epetra_FEVbrMatrix(Copy, A->RowMap(), A->ColMap(), 1);
644
645 *Acopy2 = *Acopy;
646
647 Epetra_Vector x3(Acopy->RowMap()), y3(Acopy->RowMap());
648
649 x3.PutScalar(1.0); y3.PutScalar(0.0);
650
651 Acopy2->Multiply(false, x3, y3);
652
653 double y3norm2;
654 y3.Norm2(&y3norm2);
655
656 if (y3norm2 != y2norm2) {
657 cerr << "norm2(Acopy*ones) != norm2(Acopy2*ones)"<<endl;
658 return(-999);
659 }
660
661 int len = 20;
662 int* indices = new int[len];
663 double* values = new double[len];
664 int numIndices;
665
666 if (map.MyGID(0)) {
667 int lid = map.LID(0);
668 EPETRA_CHK_ERR( A->ExtractMyRowCopy(lid, len, numIndices,
669 values, indices) );
670 if (numIndices != 4) {
671 return(-1);
672 }
673 if (indices[0] != lid) {
674 return(-2);
675 }
676
677 if (values[0] != 1.0*numProcs) {
678 cout << "ERROR: values[0] ("<<values[0]<<") should be "<<numProcs<<endl;
679 return(-3);
680 }
681 }
682
683 if (map.MyGID(4)) {
684 int lid = map.LID(4);
685 EPETRA_CHK_ERR( A->ExtractMyRowCopy(lid, len, numIndices,
686 values, indices) );
687
688 if (numIndices != 9) {
689 return(-4);
690 }
691 int lcid = A->LCID(4);
692// if (indices[lcid] != 4) {
693// cout << "ERROR: indices[4] ("<<indices[4]<<") should be "
694// <<A->LCID(4)<<endl;
695// return(-5);
696// }
697 if (values[lcid] != 4.0*numProcs) {
698 cout << "ERROR: values["<<lcid<<"] ("<<values[lcid]<<") should be "
699 <<4*numProcs<<endl;
700 return(-6);
701 }
702 }
703
704 delete [] values_2d;
705 delete [] values_1d;
706 delete [] nodes;
707 delete [] indices;
708 delete [] values;
709
710 delete A;
711 delete Acopy2;
712 delete Acopy;
713 delete graph;
714
715 return(0);
716}
#define EPETRA_CHK_ERR(a)
@ Copy
int quad1(const Epetra_Map &map, bool verbose)
int quad2(const Epetra_Map &map, bool verbose)
int MultiVectorTests(const Epetra_Map &Map, int NumVectors, bool verbose)
int four_quads(const Epetra_Comm &Comm, bool preconstruct_graph, bool verbose)
Epetra_BlockMap: A class for partitioning block element vectors and matrices.
int MinAllGID() const
Returns the minimum global ID across the entire map.
int MinMyGID() const
Returns the minimum global ID owned by this processor.
int LID(int GID) const
Returns local ID of global ID, return -1 if not found on this processor.
int NumGlobalElements() const
Number of elements across all processors.
const Epetra_Comm & Comm() const
Access function for Epetra_Comm communicator.
bool MyGID(int GID_in) const
Returns true if the GID passed in belongs to the calling processor in this map, otherwise returns fal...
Epetra_Comm: The Epetra Communication Abstract Base Class.
Definition: Epetra_Comm.h:73
virtual int NumProc() const =0
Returns total number of processes.
virtual int MyPID() const =0
Return my process ID.
virtual void Barrier() const =0
Epetra_Comm Barrier function.
Epetra_CrsGraph: A class for constructing and using sparse compressed row graphs.
int InsertGlobalIndices(int GlobalRow, int NumIndices, int *Indices)
Enter a list of elements in a specified global row of the graph.
int FillComplete()
Tranform to local index space. Perform other operations to allow optimal matrix operations.
Epetra Finite-Element VbrMatrix.
int BeginInsertGlobalValues(int BlockRow, int NumBlockEntries, int *BlockIndices)
Initiate insertion of a list of elements in a given global row of the matrix, values are inserted via...
int GlobalAssemble(bool callFillComplete=true)
int BeginSumIntoGlobalValues(int BlockRow, int NumBlockEntries, int *BlockIndices)
Initiate summing into current values with this list of entries for a given global row of the matrix,...
int EndSubmitEntries()
Completes processing of all data passed in for the current block row.
int SubmitBlockEntry(double *Values, int LDA, int NumRows, int NumCols)
Submit a block entry to the indicated block row and column specified in the Begin routine.
Epetra Finite-Element Vector.
Epetra_Map: A class for partitioning vectors and matrices.
Definition: Epetra_Map.h:119
virtual void Print(std::ostream &os) const
Print method.
int Norm2(double *Result) const
Compute 2-norm of each vector in multi-vector.
int PutScalar(double ScalarConstant)
Initialize all values in a multi-vector with constant value.
virtual void Print(std::ostream &os) const
Print method.
int NumMyNonzeros() const
Returns the number of nonzero entriesowned by the calling processor .
int LCID(int GCID_in) const
Returns the local column index for given global column index, returns -1 if no local column for this ...
int NumMyRows() const
Returns the number of matrix rows owned by the calling processor.
const Epetra_BlockMap & ColMap() const
Returns the ColMap as an Epetra_BlockMap (the Epetra_Map base class) needed for implementing Epetra_R...
int Multiply(bool TransA, const Epetra_MultiVector &X, Epetra_MultiVector &Y) const
Returns the result of a Epetra_VbrMatrix multiplied by a Epetra_MultiVector X in Y.
int ExtractMyRowCopy(int MyRow, int Length, int &NumEntries, double *Values, int *Indices) const
Returns a copy of the specified local row in user-provided arrays.
const Epetra_BlockMap & RowMap() const
Returns the RowMap object as an Epetra_BlockMap (the Epetra_Map base class) needed for implementing E...
Epetra_Vector: A class for constructing and using dense vectors on a parallel computer.
#define EPETRA_TEST_ERR(a, b)