Amesos Package Browser (Single Doxygen Collection) Development
CreateTridi.cpp
Go to the documentation of this file.
00001 // @HEADER
00002 // ***********************************************************************
00003 // 
00004 //                Amesos: Direct Sparse Solver Package
00005 //                 Copyright (2004) Sandia Corporation
00006 // 
00007 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
00008 // license for use of this work by or on behalf of the U.S. Government.
00009 // 
00010 // This library is free software; you can redistribute it and/or modify
00011 // it under the terms of the GNU Lesser General Public License as
00012 // published by the Free Software Foundation; either version 2.1 of the
00013 // License, or (at your option) any later version.
00014 //  
00015 // This library is distributed in the hope that it will be useful, but
00016 // WITHOUT ANY WARRANTY; without even the implied warranty of
00017 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00018 // Lesser General Public License for more details.
00019 //  
00020 // You should have received a copy of the GNU Lesser General Public
00021 // License along with this library; if not, write to the Free Software
00022 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00023 // USA
00024 // Questions? Contact Michael A. Heroux (maherou@sandia.gov) 
00025 // 
00026 // ***********************************************************************
00027 // @HEADER
00028 
00029 //
00030 //  CreateTridi populates an empty EpetraCrsMatrix with a tridiagonal with 
00031 //  -1 on the off-diagonals and 2 on the diagonal.  
00032 //
00033 //  CreateTridiPlus creates the same matrix as CreateTridi except that it adds
00034 //  -1 in the two off diagonal corners.
00035 //
00036 //  This code was plaguerized from epetra/example/petra_power_method/cxx_main.cpp
00037 //  presumably written by Mike Heroux.
00038 //
00039 //  Adapted by Ken Stanley - Aug 2003 
00040 //
00041 #include "Epetra_Map.h"
00042 #include "Epetra_CrsMatrix.h"
00043 
00044 int CreateTridi(Epetra_CrsMatrix& A)
00045 {
00046 
00047   Epetra_Map Map = A.RowMap();
00048   int NumMyElements = Map.NumMyElements();
00049   int NumGlobalElements = Map.NumGlobalElements();
00050 
00051   int * MyGlobalElements = new int[NumMyElements];
00052     Map.MyGlobalElements(MyGlobalElements);
00053 
00054   // Add  rows one-at-a-time
00055   // Need some vectors to help
00056   // Off diagonal Values will always be -1
00057 
00058 
00059   double *Values = new double[3];
00060   int *Indices = new int[3];
00061   int NumEntries;
00062   
00063   for (int i=0; i<NumMyElements; i++)
00064     {
00065     if (MyGlobalElements[i]==0)
00066       {
00067   Indices[0] = 0;
00068   Indices[1] = 1;
00069   Values[0] = 2.0;
00070   Values[1] = -1.0;
00071   NumEntries = 2;
00072       }
00073     else if (MyGlobalElements[i] == NumGlobalElements-1)
00074       {
00075   Indices[0] = NumGlobalElements-1;
00076   Indices[1] = NumGlobalElements-2;
00077   Values[0] = 2.0;
00078   Values[1] = -1.0;
00079   NumEntries = 2;
00080       }
00081     else
00082       {
00083   Indices[0] = MyGlobalElements[i]-1;
00084   Indices[1] = MyGlobalElements[i];
00085   Indices[2] = MyGlobalElements[i]+1;
00086   Values[0] = -1.0; 
00087   Values[1] = 2.0;
00088   Values[2] = -1.0;
00089   NumEntries = 3;
00090       }
00091     
00092     assert(A.InsertGlobalValues(MyGlobalElements[i], NumEntries, Values, Indices)==0);
00093      // Put in the diagonal entry
00094      //     assert(A.InsertGlobalValues(MyGlobalElements[i], 1, &two, &MyGlobalElements[i])==0);
00095     }
00096   
00097   // Finish up
00098   assert(A.FillComplete()==0);
00099 
00100 
00101   delete[] MyGlobalElements;
00102   delete[] Values;
00103   delete[] Indices;
00104   return 0;
00105 }
00106 
00107 int CreateTridiPlus(Epetra_CrsMatrix& A)
00108 {
00109 
00110   Epetra_Map Map = A.RowMap();
00111   int NumMyElements = Map.NumMyElements();
00112   int NumGlobalElements = Map.NumGlobalElements();
00113 
00114   int * MyGlobalElements = new int[NumMyElements];
00115     Map.MyGlobalElements(MyGlobalElements);
00116 
00117   // Add  rows one-at-a-time
00118   // Need some vectors to help
00119   // Off diagonal Values will always be -1
00120 
00121 
00122   double *Values = new double[3];
00123   int *Indices = new int[3];
00124   int NumEntries;
00125   
00126   for (int i=0; i<NumMyElements; i++)
00127     {
00128     if (MyGlobalElements[i]==0)
00129       {
00130   Indices[0] = 0;
00131   Indices[1] = 1;
00132   Indices[2] = NumGlobalElements-1;
00133   Values[0] = 2.0;
00134   Values[1] = -1.0;
00135   Values[2] = -0.5;
00136   NumEntries = 3;
00137       }
00138     else if (MyGlobalElements[i] == NumGlobalElements-1)
00139       {
00140   Indices[0] = NumGlobalElements-1;
00141   Indices[1] = NumGlobalElements-2;
00142   Indices[2] = 0;
00143   Values[0] = 2.0;
00144   Values[1] = -1.0;
00145   Values[2] = -0.5;
00146   NumEntries = 3;
00147       }
00148     else
00149       {
00150   Indices[0] = MyGlobalElements[i]-1;
00151   Indices[1] = MyGlobalElements[i];
00152   Indices[2] = MyGlobalElements[i]+1;
00153   Values[0] = -1.0; 
00154   Values[1] = 2.0;
00155   Values[2] = -1.0;
00156   NumEntries = 3;
00157       }
00158     
00159     assert(A.InsertGlobalValues(MyGlobalElements[i], NumEntries, Values, Indices)==0);
00160      // Put in the diagonal entry
00161      //     assert(A.InsertGlobalValues(MyGlobalElements[i], 1, &two, &MyGlobalElements[i])==0);
00162     }
00163   
00164   // Finish up
00165   assert(A.FillComplete()==0);
00166 
00167 
00168   delete[] MyGlobalElements;
00169   delete[] Values;
00170   delete[] Indices;
00171   return 0;
00172 }
00173 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines