Epetra Package Browser (Single Doxygen Collection) Development
FEVector/ExecuteTestProblems.cpp
Go to the documentation of this file.
00001 //@HEADER
00002 // ************************************************************************
00003 // 
00004 //               Epetra: Linear Algebra Services Package 
00005 //                 Copyright 2011 Sandia Corporation
00006 // 
00007 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
00008 // the U.S. Government retains certain rights in this software.
00009 //
00010 // Redistribution and use in source and binary forms, with or without
00011 // modification, are permitted provided that the following conditions are
00012 // met:
00013 //
00014 // 1. Redistributions of source code must retain the above copyright
00015 // notice, this list of conditions and the following disclaimer.
00016 //
00017 // 2. Redistributions in binary form must reproduce the above copyright
00018 // notice, this list of conditions and the following disclaimer in the
00019 // documentation and/or other materials provided with the distribution.
00020 //
00021 // 3. Neither the name of the Corporation nor the names of the
00022 // contributors may be used to endorse or promote products derived from
00023 // this software without specific prior written permission.
00024 //
00025 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
00026 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
00027 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
00028 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
00029 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
00030 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
00031 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
00032 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
00033 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
00034 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
00035 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00036 //
00037 // Questions? Contact Michael A. Heroux (maherou@sandia.gov) 
00038 // 
00039 // ************************************************************************
00040 //@HEADER
00041 
00042 
00043 #include "Epetra_BLAS.h"
00044 #include "ExecuteTestProblems.h"
00045 #include "Epetra_Comm.h"
00046 #include "Epetra_Vector.h"
00047 #include "Epetra_IntSerialDenseVector.h"
00048 #include "Epetra_SerialDenseVector.h"
00049 #include <vector>
00050 
00051 int MultiVectorTests(const Epetra_BlockMap & Map, int NumVectors, bool verbose)
00052 {
00053   (void)NumVectors;
00054   const Epetra_Comm & Comm = Map.Comm();
00055   int ierr = 0;
00056   
00057   /* get number of processors and the name of this processor */
00058   
00059   // int NumProc = Comm.getNumProc();
00060   int MyPID   = Comm.MyPID();
00061   
00062   // Construct FEVector
00063   
00064   if (verbose&&MyPID==0) cout << "constructing Epetra_FEVector" << endl;
00065 
00066   Epetra_FEVector A(Map, 1);
00067  
00068   //For an extreme test, we'll have each processor sum-in a 1.0 for All
00069   //global ids.
00070 
00071   int minGID = Map.MinAllGID();
00072   int numGlobalIDs = Map.NumGlobalElements();
00073 
00074   //For now we're going to have just one point associated with
00075   //each GID (element).
00076 
00077   int* ptIndices = new int[numGlobalIDs];
00078   double* ptCoefs = new double[numGlobalIDs];
00079 
00080   Epetra_IntSerialDenseVector epetra_indices(View, ptIndices, numGlobalIDs);
00081   Epetra_SerialDenseVector epetra_coefs(View, ptCoefs, numGlobalIDs);
00082 
00083   {for(int i=0; i<numGlobalIDs; ++i) {
00084     ptIndices[i] = minGID+i;
00085     ptCoefs[i] = 1.0;
00086   }}
00087 
00088   if (verbose&&MyPID==0) {
00089     cout << "calling A.SumIntoGlobalValues with " << numGlobalIDs << " values"<<endl;
00090   }
00091   EPETRA_TEST_ERR( A.SumIntoGlobalValues(numGlobalIDs, ptIndices, ptCoefs), ierr);
00092 
00093   if (verbose&&MyPID==0) {
00094     cout << "calling A.SumIntoGlobalValues with " << numGlobalIDs << " values"<<endl;
00095   }
00096   EPETRA_TEST_ERR( A.SumIntoGlobalValues(epetra_indices, epetra_coefs), ierr);
00097 
00098   if (verbose&&MyPID==0) {
00099     cout << "calling A.GlobalAssemble()" << endl;
00100   }
00101 
00102   EPETRA_TEST_ERR( A.GlobalAssemble(), ierr );
00103 
00104   if (verbose&&MyPID==0) {
00105   cout << "after globalAssemble"<<endl;
00106   }
00107   if (verbose) {
00108   A.Print(cout);
00109   }
00110 
00111   //now do a quick test of the copy constructor
00112   Epetra_FEVector B(A);
00113 
00114   double nrm2a, nrm2b;
00115   A.Norm2(&nrm2a);
00116   B.Norm2(&nrm2b);
00117 
00118   if (nrm2a != nrm2b) {
00119     cerr << "copy-constructor test failed, norm of copy doesn't equal"
00120          << " norm of original."<<endl;
00121     return(-1);
00122   }
00123 
00124   delete [] ptIndices;
00125   delete [] ptCoefs;
00126 
00127   return(ierr);
00128 }
00129 
00130 int fevec0(Epetra_Comm& Comm, bool verbose)
00131 {
00132   int ierr = 0;
00133   int NumGlobalRows = 4;
00134   int indexBase = 0;
00135   Epetra_Map Map(NumGlobalRows, indexBase, Comm);
00136 
00137   int Numprocs = Comm.NumProc();
00138   int MyPID   = Comm.MyPID();
00139 
00140   if (Numprocs != 2) return(0);
00141 
00142 
00143   int NumCols = 3;
00144   int* Indices = new int[NumCols];
00145 
00146   double* Values = new double[NumCols];
00147 
00148 // Create vectors
00149 
00150   Epetra_FEVector b(Map, 1);
00151   Epetra_FEVector x0(Map, 1);
00152 
00153 // source terms
00154   NumCols = 2;
00155 
00156   if(MyPID==0)  // indices corresponding to element 0 on processor 0
00157   {
00158     Indices[0] = 0;
00159     Indices[1] = 3;
00160 
00161     Values[0] = 1./2.;
00162     Values[1] = 1./2.;
00163 
00164    }
00165    else
00166    {
00167     Indices[0] = 1;
00168     Indices[1] = 2;
00169 
00170     Values[0] = 0;
00171     Values[1] = 0;
00172    }
00173 
00174   EPETRA_TEST_ERR( b.SumIntoGlobalValues(NumCols, Indices, Values),
00175                    ierr);
00176 
00177   EPETRA_TEST_ERR( b.GlobalAssemble(), ierr);
00178 
00179   if (verbose&&MyPID==0) {
00180     cout << "b:"<<endl;
00181   }
00182 
00183   if (verbose) {
00184   b.Print(cout);
00185   }
00186 
00187   x0 = b;
00188 
00189   if (verbose&&MyPID==0) {
00190     cout << "x:"<<endl;
00191   }
00192 
00193   if (verbose) {
00194   x0.Print(cout);
00195   }
00196 
00197   delete [] Values;
00198   delete [] Indices;
00199 
00200   return(0);
00201 }
00202 
00203 int fevec1(Epetra_Comm& Comm, bool verbose)
00204 {
00205   int Numprocs = Comm.NumProc();
00206 
00207   if (Numprocs != 2) return(0);
00208   int MyPID   = Comm.MyPID();
00209 
00210   int ierr = 0;
00211   int NumGlobalRows = 6;
00212   const int NumVectors = 4;
00213   int indexBase = 0;
00214   Epetra_Map Map(NumGlobalRows, indexBase, Comm);
00215 
00216   const int Num = 4;
00217   int Indices[Num];
00218  
00219   double Values[Num];
00220  
00221 // Create vectors
00222 
00223   Epetra_FEVector b(Map, NumVectors);
00224   Epetra_FEVector x0(Map, NumVectors);
00225  
00226 // source terms
00227  
00228   if(MyPID==0)  // indices corresponding to element 0 on processor 0
00229   {
00230     Indices[0] = 0;
00231     Indices[1] = 1;
00232     Indices[2] = 4;
00233     Indices[3] = 5;
00234  
00235     Values[0] = 1./2.;
00236     Values[1] = 1./2.;
00237     Values[2] = 1./2.;
00238     Values[3] = 1./2.;
00239 
00240    }
00241    else
00242    {
00243     Indices[0] = 1;
00244     Indices[1] = 2;
00245     Indices[2] = 3;
00246     Indices[3] = 4;
00247  
00248     Values[0] = 0;
00249     Values[1] = 0;
00250     Values[2] = 0;
00251     Values[3] = 0;
00252    }
00253 
00254   for(int i=0; i<NumVectors; ++i) {
00255     EPETRA_TEST_ERR( b.SumIntoGlobalValues(Num, Indices, Values, i),
00256        ierr);
00257   }
00258 
00259   EPETRA_TEST_ERR( b.GlobalAssemble(), ierr);
00260 
00261   double nrm2[NumVectors];
00262 
00263   b.Norm2(nrm2);
00264 
00265   for(int i=1; i<NumVectors; ++i) {
00266     if (fabs(nrm2[i]-nrm2[0]) > 1.e-12) {
00267       EPETRA_TEST_ERR(-1, ierr);
00268       return(-1);
00269     }
00270   }
00271 
00272 
00273   //now sum-in again, to make sure the previous call to GlobalAssemble
00274   //didn't do something nasty to internal non-local data structures.
00275   //(This is a specific case that has bitten me. Hence this test...)
00276   for(int i=0; i<NumVectors; ++i) {
00277     EPETRA_TEST_ERR( b.SumIntoGlobalValues(Num, Indices, Values, i),
00278                    ierr);
00279   }
00280 
00281   //and now GlobalAssemble again...
00282   EPETRA_TEST_ERR( b.GlobalAssemble(), ierr);
00283 
00284 
00285   if (verbose&&MyPID==0) {
00286     cout << "b:"<<endl;
00287   }
00288 
00289   if (verbose) {
00290     b.Print(cout);
00291   }
00292 
00293   x0 = b;
00294 
00295   if (verbose&&MyPID==0) {
00296     cout << "x:"<<endl;
00297   }
00298 
00299   if (verbose) {
00300     x0.Print(cout);
00301   }
00302 
00303   return(0);
00304 }
00305 
00306 int fevec2(Epetra_Comm& Comm, bool verbose)
00307 {
00308   int ierr = 0;
00309   int NumGlobalElems = 4;
00310   int elemSize = 3;
00311   int indexBase = 0;
00312   Epetra_BlockMap Map(NumGlobalElems, elemSize, indexBase, Comm);
00313 
00314   int Numprocs = Comm.NumProc();
00315   int MyPID   = Comm.MyPID();
00316 
00317   if (Numprocs != 2) return(0);
00318 
00319   int NumCols = 3;
00320   int* Indices = new int[NumCols];
00321   int* numValuesPerID = new int[NumCols];
00322   for(int i=0; i<NumCols; ++i) {
00323     numValuesPerID[i] = elemSize;
00324   }
00325  
00326   double* Values = new double[NumCols*elemSize];
00327  
00328 // Create vectors
00329 
00330   Epetra_FEVector b(Map, 1);
00331   Epetra_FEVector x0(Map, 1);
00332  
00333 // source terms
00334   NumCols = 2;
00335  
00336   if(MyPID==0)  // indices corresponding to element 0 on processor 0
00337   {
00338     Indices[0] = 0;
00339     Indices[1] = 3;
00340  
00341     Values[0] = 1./2.;
00342     Values[1] = 1./2.;
00343     Values[2] = 1./2.;
00344     Values[3] = 1./2.;
00345     Values[4] = 1./2.;
00346     Values[5] = 1./2.;
00347 
00348    }
00349    else
00350    {
00351     Indices[0] = 1;
00352     Indices[1] = 2;
00353  
00354     Values[0] = 0;
00355     Values[1] = 0;
00356     Values[2] = 0;
00357     Values[3] = 0;
00358     Values[4] = 0;
00359     Values[5] = 0;
00360    }
00361 
00362   EPETRA_TEST_ERR( b.SumIntoGlobalValues(NumCols, Indices,
00363            numValuesPerID, Values),
00364        ierr);
00365 
00366   EPETRA_TEST_ERR( b.GlobalAssemble(), ierr);
00367 
00368   if (verbose&&MyPID==0) {
00369     cout << "b:"<<endl;
00370   }
00371 
00372   if (verbose) {
00373   b.Print(cout);
00374   }
00375 
00376   x0 = b;
00377 
00378   if (verbose&&MyPID==0) {
00379     cout << "x:"<<endl;
00380   }
00381 
00382   if (verbose) {
00383   x0.Print(cout);
00384   }
00385 
00386   delete [] Values;
00387   delete [] Indices;
00388   delete [] numValuesPerID;
00389 
00390   return(0);
00391 }
00392 
00393 int fevec3(Epetra_Comm& Comm, bool verbose)
00394 {
00395   int ierr = 0;
00396   int NumGlobalElems = 4;
00397   int elemSize = 40;
00398   int indexBase = 0;
00399   Epetra_BlockMap Map(NumGlobalElems, elemSize, indexBase, Comm);
00400 
00401   int Numprocs = Comm.NumProc();
00402   int MyPID   = Comm.MyPID();
00403 
00404   if (Numprocs != 2) return(0);
00405 
00406   int NumCols = 3;
00407   int* Indices = new int[NumCols];
00408   int* numValuesPerID = new int[NumCols];
00409   for(int i=0; i<NumCols; ++i) {
00410     numValuesPerID[i] = elemSize;
00411   }
00412  
00413   double* Values = new double[NumCols*elemSize];
00414  
00415 // Create vectors
00416 
00417   Epetra_FEVector b(Map, 1);
00418   Epetra_FEVector x0(Map, 1);
00419  
00420 // source terms
00421   NumCols = 2;
00422  
00423   if(MyPID==0)  // indices corresponding to element 0 on processor 0
00424   {
00425     Indices[0] = 0;
00426     Indices[1] = 3;
00427  
00428     for(int ii=0; ii<NumCols*elemSize; ++ii) {
00429       Values[ii] = 1./2.;
00430     }
00431 
00432   }
00433   else
00434   {
00435     Indices[0] = 1;
00436     Indices[1] = 2;
00437  
00438     for(int ii=0; ii<NumCols*elemSize; ++ii) {
00439       Values[ii] = 0.;
00440     }
00441 
00442   }
00443 
00444   EPETRA_TEST_ERR( b.SumIntoGlobalValues(NumCols, Indices,
00445            numValuesPerID, Values),
00446        ierr);
00447 
00448   EPETRA_TEST_ERR( b.GlobalAssemble(), ierr);
00449 
00450   if (verbose&&MyPID==0) {
00451     cout << "b:"<<endl;
00452   }
00453 
00454   if (verbose) {
00455   b.Print(cout);
00456   }
00457 
00458   x0 = b;
00459 
00460   if (verbose&&MyPID==0) {
00461     cout << "x:"<<endl;
00462   }
00463 
00464   if (verbose) {
00465   x0.Print(cout);
00466   }
00467 
00468   delete [] Values;
00469   delete [] Indices;
00470   delete [] numValuesPerID;
00471 
00472   return(0);
00473 }
00474 
00475 int fevec4(Epetra_Comm& Comm, bool verbose)
00476 {
00477   int NumElements = 4;
00478   Epetra_Map     Map(NumElements, 0, Comm);
00479   Epetra_FEVector x1(Map);
00480   const double value = 1.;
00481   x1.PutScalar (value);
00482         // replace one element by itself. processor 0
00483         // does not own this element
00484   const int GID = 3;
00485   x1.ReplaceGlobalValues(1, &GID, &value);
00486   x1.GlobalAssemble (Insert);
00487 
00488   if (Map.MyGID(3)) {
00489     //insist that the value for GID==3 is 1:
00490     if (std::abs(x1.Values()[Map.LID(3)] - 1) > 1.e-9) return -1;
00491   }
00492 
00493   std::cout << x1;
00494 
00495   Comm.Barrier();
00496 
00497         // re-apply GlobalAssemble. Nothing should
00498         // happen
00499   x1.GlobalAssemble (Insert);
00500   std::cout << x1;
00501   if (Map.MyGID(3)) {
00502     //insist that the value for GID==3 is 1:
00503     if (std::abs(x1.Values()[Map.LID(3)] - 1) > 1.e-9) return -1;
00504   }
00505 
00506   return 0;
00507 }
00508 
00509 int fevec5(Epetra_Comm& Comm, bool verbose)
00510 {
00511   int NumElements = 4;
00512   Epetra_Map     Map(NumElements, 0, Comm);
00513   Epetra_FEVector x1(Map);
00514   x1.PutScalar (0);
00515 
00516         // let all processors set global entry 0 to 1
00517   const int GID = 0;
00518   const double value = 1;
00519   x1.ReplaceGlobalValues(1, &GID, &value);
00520   x1.GlobalAssemble (Insert);
00521   if (Comm.MyPID()==0)
00522     std::cout << "Entry " << GID << " after construct & set: " 
00523         << x1[0][0] << std::endl;
00524 
00525         // copy vector
00526   Epetra_FEVector x2 (x1);
00527 
00528   x2.PutScalar(0);
00529 
00530         // re-apply 1 to the vector, but only on the
00531         // owning processor. should be enough to set
00532         // the value (as non-local data in x1 should
00533         // have been eliminated after calling
00534         // GlobalAssemble).
00535   if (Comm.MyPID()==0)
00536     x2.ReplaceGlobalValues(1, &GID, &value);
00537   x2.GlobalAssemble (Insert);
00538 
00539   if (Comm.MyPID()==0)
00540     std::cout << "Entry " << GID << " after copy & set:      " 
00541         << x2[0][0] << std::endl;
00542 
00543   return 0;
00544 }
00545 
00546 int fevec6(Epetra_Comm& Comm, bool verbose)
00547 {
00548   int NumElements = 4;
00549   Epetra_Map     Map(NumElements, 0, Comm);
00550   Epetra_FEVector x1(Map);
00551   x1.PutScalar (0);
00552 
00553         // let all processors set global entry 0 to 1
00554   const int GID = 0;
00555   const double value = 1;
00556   x1.ReplaceGlobalValues(1, &GID, &value);
00557   x1.GlobalAssemble (Insert);
00558   if (Comm.MyPID()==0)
00559     std::cout << "Entry " << GID << " after construct & set: " 
00560         << x1[0][0] << std::endl;
00561 
00562   x1.PutScalar(0);
00563 
00564         // re-apply 1 to the vector, but only on the
00565         // owning processor. should be enough to set
00566         // the value (as non-local data in x1 should
00567         // have been eliminated after calling
00568         // GlobalAssemble).
00569   if (Comm.MyPID()==0)
00570     x1.ReplaceGlobalValues(1, &GID, &value);
00571   x1.GlobalAssemble (Insert);
00572 
00573   if (Comm.MyPID()==0) {
00574     std::cout << "Entry " << GID << " after PutScalar & set:      " 
00575         << x1[0][0] << std::endl;
00576     if (x1[0][0] != value) return -1;
00577   }
00578 
00579   return 0;
00580 }
00581 
00582 int fevec7(Epetra_Comm& Comm, bool verbose)
00583 {
00584   const int NumVectors = 4;
00585   const int NumElements = 4;
00586   Epetra_Map     Map(NumElements, 0, Comm);
00587   std::vector<double> mydata(NumElements*NumVectors, 1.0);
00588   Epetra_FEVector x1(View, Map, &mydata[0], NumElements, NumVectors);
00589 
00590   x1.PutScalar (0);
00591 
00592         // let all processors set global entry 0 to 1
00593   const int GID = 0;
00594   const double value = 1;
00595   x1.ReplaceGlobalValues(1, &GID, &value);
00596   x1.GlobalAssemble (Insert);
00597 
00598   if (Comm.MyPID()==0 && x1[0][0] != value) return -1;
00599   return 0;
00600 }
00601 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines