EpetraExt Development
EpetraExt_MultiVectorIn.cpp
Go to the documentation of this file.
00001 //@HEADER
00002 // ***********************************************************************
00003 // 
00004 //     EpetraExt: Epetra Extended - Linear Algebra Services Package
00005 //                 Copyright (2001) 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 #include "EpetraExt_MultiVectorIn.h"
00029 #include "Epetra_Comm.h"
00030 #include "Epetra_MultiVector.h"
00031 #include "Epetra_Vector.h"
00032 #include "Epetra_BlockMap.h"
00033 
00034 using namespace EpetraExt;
00035 namespace EpetraExt {
00036 
00037 int MatrixMarketFileToMultiVector( const char *filename, const Epetra_BlockMap & map, Epetra_MultiVector * & A) {
00038 
00039   const int lineLength = 1025;
00040   const int tokenLength = 35;
00041   char line[lineLength];
00042   char token1[tokenLength];
00043   char token2[tokenLength];
00044   char token3[tokenLength];
00045   char token4[tokenLength];
00046   char token5[tokenLength];
00047   int M, N;
00048 
00049   FILE * handle = 0;
00050 
00051   handle = fopen(filename,"r");  // Open file
00052   if (handle == 0)
00053     EPETRA_CHK_ERR(-1); // file not found
00054 
00055   // Check first line, which should be "%%MatrixMarket matrix coordinate real general" (without quotes)
00056   if(fgets(line, lineLength, handle)==0) return(-1);
00057   if(sscanf(line, "%s %s %s %s %s", token1, token2, token3, token4, token5 )==0) return(-1);
00058   if (strcmp(token1, "%%MatrixMarket") ||
00059       strcmp(token2, "matrix") ||
00060       strcmp(token3, "array") ||
00061       strcmp(token4, "real") ||
00062       strcmp(token5, "general")) return(-1);
00063 
00064   // Next, strip off header lines (which start with "%")
00065   do {
00066     if(fgets(line, lineLength, handle)==0) return(-1);
00067   } while (line[0] == '%');
00068 
00069   // Next get problem dimensions: M, N
00070   if(sscanf(line, "%d %d", &M, &N)==0) return(-1);
00071 
00072   // Compute the offset for each processor for when it should start storing values
00073   int numMyPoints = map.NumMyPoints();
00074   int offset;
00075   map.Comm().ScanSum(&numMyPoints, &offset, 1); // ScanSum will compute offsets for us
00076   offset -= numMyPoints; // readjust for my PE
00077 
00078   // Now construct vector/multivector
00079   if (N==1)
00080     A = new Epetra_Vector(map);
00081   else
00082     A = new Epetra_MultiVector(map, N);
00083 
00084   double ** Ap = A->Pointers();
00085 
00086   for (int j=0; j<N; j++) {
00087     double * v = Ap[j];
00088 
00089     // Now read in lines that we will discard
00090     for (int i=0; i<offset; i++)
00091       if(fgets(line, lineLength, handle)==0) return(-1);
00092     
00093     // Now read in each value and store to the local portion of the the  if the row is owned.
00094     double V;
00095     for (int i=0; i<numMyPoints; i++) {
00096       if(fgets(line, lineLength, handle)==0) return(-1);
00097       if(sscanf(line, "%lg\n", &V)==0) return(-1);
00098       v[i] = V;
00099     }
00100     // Now read in the rest of the lines to discard
00101     for (int i=0; i < M-numMyPoints-offset; i++) {
00102       if(fgets(line, lineLength, handle)==0) return(-1);
00103     }
00104   }
00105 
00106   if (fclose(handle)) return(-1);
00107   
00108   return(0);
00109 }
00110 
00111 } // namespace EpetraExt
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines