Epetra_C_wrappers.cpp

Go to the documentation of this file.
00001 
00002 //@HEADER
00003 // ************************************************************************
00004 // 
00005 //               Epetra: Linear Algebra Services Package 
00006 //                 Copyright (2001) Sandia Corporation
00007 // 
00008 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
00009 // license for use of this work by or on behalf of the U.S. Government.
00010 // 
00011 // This library is free software; you can redistribute it and/or modify
00012 // it under the terms of the GNU Lesser General Public License as
00013 // published by the Free Software Foundation; either version 2.1 of the
00014 // License, or (at your option) any later version.
00015 //  
00016 // This library is distributed in the hope that it will be useful, but
00017 // WITHOUT ANY WARRANTY; without even the implied warranty of
00018 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00019 // Lesser General Public License for more details.
00020 //  
00021 // You should have received a copy of the GNU Lesser General Public
00022 // License along with this library; if not, write to the Free Software
00023 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00024 // USA
00025 // Questions? Contact Michael A. Heroux (maherou@sandia.gov) 
00026 // 
00027 // ************************************************************************
00028 //@HEADER
00029 
00030 #ifdef EPETRA_MPI
00031 #include <mpi.h>
00032 #endif
00033 
00034 #include "Epetra_Object.h"
00035 #include "Epetra_Comm.h"
00036 #include "Epetra_SerialComm.h"
00037 #include "Epetra_Map.h"
00038 #include "Epetra_LocalMap.h"
00039 #include "Epetra_BlockMap.h"
00040 #include "Epetra_MultiVector.h"
00041 #include "Epetra_Vector.h"
00042 #include "Epetra_VbrMatrix.h"
00043 #include "Epetra_CrsMatrix.h"
00044 #include "Epetra_C_wrappers.h"
00045 #ifdef EPETRA_MPI
00046 #include "Epetra_MpiComm.h"
00047 #endif
00048 
00049 #ifdef __cplusplus
00050 extern "C" {
00051 #endif
00052 
00054 //                  Epetra_Comm                    //
00056 
00057 #ifdef EPETRA_MPI
00058   EPETRA_OBJECT_PTR MANGLE(epetra_mpicomm_create1)() {
00059     Epetra_Comm *comm_ = new Epetra_MpiComm(MPI_COMM_WORLD);
00060     return((EPETRA_OBJECT_PTR ) comm_);
00061   }
00062   EPETRA_OBJECT_PTR MANGLE(epetra_mpicomm_create2)(MPI_Comm * comm) {
00063     Epetra_Comm *comm_ = new Epetra_MpiComm(*comm);
00064     return((EPETRA_OBJECT_PTR ) comm_);
00065   }
00066 #endif
00067 
00068   EPETRA_OBJECT_PTR MANGLE(epetra_serialcomm_create)() {
00069     Epetra_Comm *comm = new Epetra_SerialComm();
00070     return((EPETRA_OBJECT_PTR ) comm);
00071   }
00072 
00073   int MANGLE(epetra_comm_mypid)(EPETRA_OBJECT_REF comm) {
00074     Epetra_Comm *comm_ = (Epetra_Comm *) comm;
00075     return(comm_->MyPID());
00076   
00077   }
00078   int MANGLE(epetra_comm_numproc)(EPETRA_OBJECT_REF comm) {
00079     Epetra_Comm* comm_ = (Epetra_Comm *) comm;
00080     return(comm_->NumProc());
00081   
00082   }
00083 
00084   void MANGLE(epetra_comm_barrier)(EPETRA_OBJECT_REF comm) {
00085     Epetra_Comm* comm_ = (Epetra_Comm *) comm;
00086     comm_->Barrier();
00087   
00088   }
00089 
00090   void MANGLE(epetra_comm_destroy)(EPETRA_OBJECT_REF comm) {
00091     delete (Epetra_Comm *) comm;
00092   }
00093 
00095   //                  Epetra_Map                     //
00097 
00098   EPETRA_OBJECT_PTR MANGLE(epetra_map_create1)(EPETRA_INT numGlobalElements,
00099                  EPETRA_INT indexBase,
00100                  EPETRA_OBJECT_REF comm) {
00101     Epetra_Comm& comm_ = *(Epetra_Comm *) comm;
00102     Epetra_Map *map = new Epetra_Map(EPETRA_DEREF(numGlobalElements), EPETRA_DEREF(indexBase), comm_);
00103     return((EPETRA_OBJECT_PTR ) map);
00104   }
00105 
00106   EPETRA_OBJECT_PTR MANGLE(epetra_map_create2)(EPETRA_INT numGlobalElements,
00107                  EPETRA_INT numMyElements,
00108                  EPETRA_INT indexBase,
00109                  EPETRA_OBJECT_REF comm) {
00110     Epetra_Comm& comm_ = *(Epetra_Comm *) comm;
00111     Epetra_Map *map = new Epetra_Map(EPETRA_DEREF(numGlobalElements), EPETRA_DEREF(numMyElements), 
00112              EPETRA_DEREF(indexBase), comm_);
00113     return((EPETRA_OBJECT_PTR ) map);
00114   }
00115 
00116   EPETRA_OBJECT_PTR MANGLE(epetra_map_create3)(EPETRA_INT numGlobalElements,
00117                  EPETRA_INT numLocalElements,
00118                  int *updateList, 
00119                  EPETRA_INT indexBase,
00120                  EPETRA_OBJECT_REF comm) {
00121     Epetra_Comm& comm_ = *(Epetra_Comm *) comm;
00122     Epetra_Map *map = new Epetra_Map(EPETRA_DEREF(numGlobalElements), EPETRA_DEREF(numLocalElements),
00123              updateList, EPETRA_DEREF(indexBase), comm_);
00124     return((EPETRA_OBJECT_PTR ) map);
00125   }
00126 
00127   int MANGLE(epetra_map_nummyelements)(EPETRA_OBJECT_REF map) {
00128     Epetra_Map * map_ = (Epetra_Map *) map;
00129     return(map_->NumMyElements());
00130   }
00131 
00132   int MANGLE(epetra_map_numglobalelements)(EPETRA_OBJECT_REF map) {
00133     Epetra_Map * map_ = (Epetra_Map *) map;
00134     return(map_->NumGlobalElements());
00135   }
00136 #ifndef EPETRA_FORTRAN  /* Fortran cannot receive a pointer to int */
00137   int * MANGLE(epetra_map_myglobalelements)(EPETRA_OBJECT_REF map) {
00138     Epetra_Map * map_ = (Epetra_Map *) map;
00139     return(map_->MyGlobalElements());
00140   }
00141 #endif
00142   EPETRA_OBJECT_PTR MANGLE(epetra_map_comm)(EPETRA_OBJECT_REF map) {
00143     Epetra_Map * map_ = (Epetra_Map *) map;
00144     return((EPETRA_OBJECT_PTR) &(map_->Comm()));
00145   }
00146 
00147   void MANGLE(epetra_map_destroy)(EPETRA_OBJECT_REF map)
00148   {
00149     delete (Epetra_Map *) map;
00150   }
00151 
00153   //                  Epetra_Vector                  //
00155 
00156   EPETRA_OBJECT_PTR MANGLE(epetra_vector_create1)(EPETRA_OBJECT_REF map) {
00157     Epetra_Map& map_ = *(Epetra_Map *) map;
00158     Epetra_Vector *vector = new Epetra_Vector(map_);
00159     return((EPETRA_OBJECT_PTR ) vector);
00160   }
00161 
00162   EPETRA_OBJECT_PTR MANGLE(epetra_vector_create2)(EPETRA_INT CopyValues, EPETRA_OBJECT_REF map,
00163               double * V) {
00164     Epetra_DataAccess CV = View;
00165     if (EPETRA_DEREF(CopyValues)==1) CV = Copy;
00166     Epetra_Map& map_ = *(Epetra_Map *) map;
00167     Epetra_Vector *vector = new Epetra_Vector(CV, map_, V);
00168     return((EPETRA_OBJECT_PTR ) vector);
00169   }
00170 
00171   int MANGLE(epetra_vector_putscalar)(EPETRA_OBJECT_REF x, EPETRA_DOUBLE scalar) {
00172     Epetra_Vector *x_ = (Epetra_Vector *) x;
00173     return(x_->PutScalar(EPETRA_DEREF(scalar)));
00174   }
00175 
00176   int MANGLE(epetra_vector_norm1)(EPETRA_OBJECT_REF x, double *scalar) {
00177     Epetra_Vector *x_ = (Epetra_Vector *) x;
00178     return(x_->Norm1(scalar));
00179   }
00180 
00181   int MANGLE(epetra_vector_norm2)(EPETRA_OBJECT_REF x, double *scalar) {
00182     Epetra_Vector *x_ = (Epetra_Vector *) x;
00183     return(x_->Norm2(scalar));
00184   }
00185 
00186   int MANGLE(epetra_vector_random)(EPETRA_OBJECT_REF x) {
00187     Epetra_Vector *x_ = (Epetra_Vector *) x;
00188     return(x_->Random());
00189   }
00190 
00191   int MANGLE(epetra_vector_update)(EPETRA_OBJECT_REF x, EPETRA_DOUBLE scalara, EPETRA_OBJECT_REF a, 
00192            EPETRA_DOUBLE scalarb, EPETRA_OBJECT_REF b, EPETRA_DOUBLE scalarx) {
00193     Epetra_Vector *x_ = (Epetra_Vector *) x;
00194     Epetra_Vector& a_ = *(Epetra_Vector *) a;
00195     Epetra_Vector& b_ = *(Epetra_Vector *) b;
00196     return(x_->Update(EPETRA_DEREF(scalara), a_, EPETRA_DEREF(scalarb), b_, EPETRA_DEREF(scalarx)));
00197   }
00198 
00199   void MANGLE(epetra_vector_print)(EPETRA_OBJECT_REF x) {
00200     cout << *(Epetra_Vector *) x;
00201   }
00202 
00203   void MANGLE(epetra_vector_destroy)(EPETRA_OBJECT_REF x) {
00204     delete (Epetra_Vector *) x;
00205   }
00206 
00208 #ifdef SKIP4NOW /* Comment this out for now */
00210   //                  Epetra_DVBR_Matrix         //
00212 
00213 
00214   EPETRA_OBJECT_PTR MANGLE(epetra_rdp_dvbr_matrix_create)
00215     (EPETRA_MAP rowmap)
00216   {
00217     Epetra_BlockMap& rowmap_ = *(Epetra_BlockMap *) rowmap;
00218     Epetra_DVBR_Matrix *B = new Epetra_DVBR_Matrix(rowmap_);
00219     return((EPETRA_OBJECT_PTR) B);
00220   }
00221 
00222   int MANGLE(epetra_rdp_dvbr_matrix_allocate)
00223     (EPETRA_MATRIX A, int* numNzBlks, int* blkColInds)
00224   {
00225     Epetra_DVBR_Matrix *B = (Epetra_DVBR_Matrix *) A;
00226     return(B->allocate(numNzBlks, blkColInds));
00227   }
00228   int MANGLE(epetra_rdp_dvbr_matrix_putblockrow)
00229     (EPETRA_MATRIX A, EPETRA_INT blk_row, EPETRA_INT num_nz_blocks, 
00230      double* vals, int* blk_col_inds)
00231   {
00232     Epetra_DVBR_Matrix *B = (Epetra_DVBR_Matrix *) A;
00233     return(B->putBlockRow( EPETRA_DEREF(blk_row), EPETRA_DEREF(num_nz_blocks), vals, 
00234          blk_col_inds));
00235   }
00236 
00237   int MANGLE(epetra_rdp_dvbr_matrix_fillcomplete)(EPETRA_MATRIX A)
00238   {
00239     Epetra_DVBR_Matrix *B = (Epetra_DVBR_Matrix *) A;
00240     return(B->fillComplete());
00241   }
00242 
00243   int MANGLE(epetra_rdp_dvbr_matrix_matvec)(EPETRA_MATRIX A, EPETRA_VECTOR x,
00244               EPETRA_VECTOR y)
00245   {
00246     Epetra_DVBR_Matrix *B = (Epetra_DVBR_Matrix *) A;
00247     const Epetra_Vector& x1 = *(Epetra_Vector *) x;
00248     Epetra_Vector& y1 = *(Epetra_Vector *) y;
00249     return(B->matvec(x1, y1));
00250   }
00251 
00252   int MANGLE(epetra_rdp_dvbr_matrix_matmultivec)(EPETRA_MATRIX A,
00253              EPETRA_MULTIVECTOR x,
00254              EPETRA_MULTIVECTOR y)
00255   {
00256     Epetra_DVBR_Matrix *B = (Epetra_DVBR_Matrix *) A;
00257     const Epetra_MultiVector& x1 = *(Epetra_MultiVector *) x;
00258     Epetra_MultiVector& y1 = *(Epetra_MultiVector *) y;
00259     return(B->matvec(x1, y1));
00260   }
00261 
00262   void MANGLE(epetra_rdp_dvbr_matrix_destroy)(EPETRA_MATRIX A)
00263   {
00264     delete (Epetra_DVBR_Matrix *) A;
00265   }
00266 
00268   //                  Epetra_DCRS_Matrix         //
00270 
00271 
00272   EPETRA_OBJECT_PTR MANGLE(epetra_rdp_dcrs_matrix_create) (EPETRA_MAP rowmap)
00273   {
00274     Epetra_Map& rowmap_ = *(Epetra_Map *) rowmap;
00275     Epetra_DCRS_Matrix *B = new Epetra_DCRS_Matrix(rowmap_);
00276     return((EPETRA_OBJECT_PTR) B);
00277   }
00278 
00279   int MANGLE(epetra_rdp_dcrs_matrix_allocate)
00280     (EPETRA_MATRIX A, int* rowLengths)
00281   {
00282     Epetra_DCRS_Matrix *B = (Epetra_DCRS_Matrix *) A;
00283     return(B->allocate(rowLengths));
00284   }
00285   int MANGLE(epetra_rdp_dcrs_matrix_putrow)(EPETRA_MATRIX A, EPETRA_INT row,
00286               EPETRA_INT num_nz, 
00287               double* vals, int* col_inds)
00288   {
00289     Epetra_DCRS_Matrix *B = (Epetra_DCRS_Matrix *) A;
00290     return(B->putRow( EPETRA_DEREF(row), EPETRA_DEREF(num_nz), vals, col_inds));
00291   }
00292 
00293   int MANGLE(epetra_rdp_dcrs_matrix_sumintodiagonal)
00294     (EPETRA_MATRIX A, double* diagonal)
00295   {
00296     Epetra_DCRS_Matrix *B = (Epetra_DCRS_Matrix *) A;
00297     return(B->sumIntoDiagonal( diagonal));
00298   }
00299 
00300   int MANGLE(epetra_rdp_dcrs_matrix_fillcomplete)(EPETRA_MATRIX A)
00301   {
00302     Epetra_DCRS_Matrix *B = (Epetra_DCRS_Matrix *) A;
00303     return(B->fillComplete());
00304   }
00305 
00306   int MANGLE(epetra_rdp_dcrs_matrix_matvec)(EPETRA_MATRIX A, EPETRA_VECTOR x,
00307               EPETRA_VECTOR y)
00308   {
00309     Epetra_DCRS_Matrix *B = (Epetra_DCRS_Matrix *) A;
00310     const Epetra_Vector& x1 = *(Epetra_Vector *) x;
00311     Epetra_Vector& y1 = *(Epetra_Vector *) y;
00312     return(B->matvec(x1, y1));
00313   }
00314 
00315   int MANGLE(epetra_rdp_dcrs_matrix_matmultivec)(EPETRA_MATRIX A,
00316              EPETRA_MULTIVECTOR x,
00317              EPETRA_MULTIVECTOR y)
00318   {
00319     Epetra_DCRS_Matrix *B = (Epetra_DCRS_Matrix *) A;
00320     const Epetra_MultiVector& x1 = *(Epetra_MultiVector *) x;
00321     Epetra_MultiVector& y1 = *(Epetra_MultiVector *) y;
00322     return(B->matvec(x1, y1));
00323   }
00324 
00325 
00326   void MANGLE(epetra_rdp_dcrs_matrix_destroy)(EPETRA_MATRIX A)
00327   {
00328     delete (Epetra_DCRS_Matrix *) A;
00329   }
00330 
00331   //                  Epetra_MultiVector         //
00333 
00334   EPETRA_OBJECT_PTR MANGLE(epetra_rdp_multivector_create)()
00335   {
00336     Epetra_MultiVector *vector = new Epetra_MultiVector();
00337     return((EPETRA_OBJECT_PTR) vector);
00338   }
00339 
00340   EPETRA_OBJECT_PTR MANGLE(epetra_rdp_multivector_create1)
00341     (EPETRA_MAP map, EPETRA_INT numVectors)
00342   {
00343     Epetra_Map& map_ = *(Epetra_Map *) map;
00344     Epetra_MultiVector *vector = new Epetra_MultiVector(map_, EPETRA_DEREF(numVectors));
00345     return((EPETRA_OBJECT_PTR) vector);
00346   }
00347 
00348   EPETRA_OBJECT_PTR MANGLE(epetra_rdp_multivector_create2)(EPETRA_MAP map, 
00349                     double *A, EPETRA_INT lda, EPETRA_INT numVectors)
00350   {
00351     Epetra_Map& map_ = *(Epetra_Map *) map;
00352     Epetra_MultiVector *vector = new Epetra_MultiVector(map_, A, EPETRA_DEREF(lda),
00353               EPETRA_DEREF(numVectors));
00354     return((EPETRA_OBJECT_PTR) vector);
00355   }
00356 
00357   EPETRA_OBJECT_PTR MANGLE(epetra_rdp_multivector_create3)(EPETRA_MAP map, 
00358                     double **in_multiVector, EPETRA_INT numVectors)
00359   {
00360     Epetra_Map& map_ = *(Epetra_Map *) map;
00361     Epetra_MultiVector *vector = new Epetra_MultiVector(map_, in_multiVector,
00362               EPETRA_DEREF(numVectors));
00363     return((EPETRA_OBJECT_PTR) vector);
00364   }
00365 
00366   EPETRA_OBJECT_PTR MANGLE(epetra_rdp_multivector_create4)
00367     (EPETRA_MULTIVECTOR in_multiVector)
00368   {
00369     Epetra_MultiVector & in_multiVector_ = *(Epetra_MultiVector *) in_multiVector;
00370     Epetra_MultiVector *vector = new Epetra_MultiVector(in_multiVector_);
00371     return((EPETRA_OBJECT_PTR) vector);
00372   }
00373 
00374   EPETRA_OBJECT_PTR MANGLE(epetra_rdp_multivector_create5)(EPETRA_MULTIVECTOR 
00375                     in_multiVector, EPETRA_INT numVectors, int *vecIndices)
00376   {
00377     Epetra_MultiVector & in_multiVector_ = *(Epetra_MultiVector *) in_multiVector;
00378     Epetra_MultiVector *vector = new Epetra_MultiVector(in_multiVector_,
00379               EPETRA_DEREF(numVectors), vecIndices));
00380   return((EPETRA_OBJECT_PTR) vector);
00381 }
00382 
00383 EPETRA_OBJECT_PTR MANGLE(epetra_rdp_multivector_create6)(EPETRA_MULTIVECTOR
00384                   in_multiVector, EPETRA_INT startindex, EPETRA_INT numvectors)
00385 {
00386   Epetra_MultiVector & in_multiVector_ = *(Epetra_MultiVector *) in_multiVector;
00387   Epetra_MultiVector *vector = new Epetra_MultiVector(in_multiVector_, EPETRA_DEREF(startindex),
00388                   EPETRA_DEREF(numvectors));
00389   return((EPETRA_OBJECT_PTR) vector);
00390 }
00391 
00392 int MANGLE(epetra_rdp_multivector_putmultivector)
00393   (EPETRA_MULTIVECTOR multiVector, 
00394    double **in_multiVector)
00395 {
00396   Epetra_MultiVector *multiVector_ = (Epetra_MultiVector *) multiVector;
00397   const double ** t = (const double **) in_multiVector;
00398   return(multiVector_->putMultiVector(t));
00399 }
00400 
00401 int MANGLE(epetra_rdp_multivector_allocate)(EPETRA_MULTIVECTOR multiVector, 
00402               EPETRA_MAP map, EPETRA_INT numVectors)
00403 {
00404   Epetra_Map& map_ = *(Epetra_Map *) map;
00405   Epetra_MultiVector *multiVector_ = (Epetra_MultiVector *) multiVector;
00406   return(multiVector_->allocate(map_, EPETRA_DEREF(numVectors)));
00407 }
00408 
00409 int MANGLE(epetra_rdp_multivector_putscalar)
00410   (EPETRA_MULTIVECTOR multiVector, EPETRA_DOUBLE scalar)
00411 {
00412   Epetra_MultiVector *multiVector_ = (Epetra_MultiVector *) multiVector;
00413   return(multiVector_->putScalar(EPETRA_DEREF(scalar)));
00414 }
00415 
00416 int MANGLE(epetra_rdp_multivector_scale)
00417   (EPETRA_MULTIVECTOR multiVector, EPETRA_DOUBLE scalar)
00418 {
00419   Epetra_MultiVector *multiVector_ = (Epetra_MultiVector *) multiVector;
00420   return(multiVector_->scale(EPETRA_DEREF(scalar)));
00421 }
00422 
00423 int MANGLE(epetra_rdp_multivector_scalecopy)
00424   (EPETRA_MULTIVECTOR multiVector, EPETRA_MULTIVECTOR multiVector_in,
00425    EPETRA_DOUBLE scalar)
00426 {
00427   Epetra_MultiVector *multiVector_ = (Epetra_MultiVector *) multiVector;
00428   Epetra_MultiVector& multiVector_in_ = *(Epetra_MultiVector *) multiVector_in;
00429   return(multiVector_->scaleCopy(multiVector_in_, EPETRA_DEREF(scalar)));
00430 }
00431 
00432 int MANGLE(epetra_rdp_multivector_dotprod)
00433   (EPETRA_MULTIVECTOR multiVector, EPETRA_MULTIVECTOR multiVector_in,
00434    double *scalar)
00435 {
00436   Epetra_MultiVector *multiVector_ = (Epetra_MultiVector *) multiVector;
00437   Epetra_MultiVector& multiVector_in_ = *(Epetra_MultiVector *) multiVector_in;
00438   return(multiVector_->dotProd(multiVector_in_, scalar));
00439 }
00440 
00441 int MANGLE(epetra_rdp_multivector_addvec)
00442   (EPETRA_MULTIVECTOR multiVector, EPETRA_DOUBLE scalar, 
00443    EPETRA_MULTIVECTOR multiVector_in)
00444 {
00445   Epetra_MultiVector *multiVector_ = (Epetra_MultiVector *) multiVector;
00446   Epetra_MultiVector& multiVector_in_ = *(Epetra_MultiVector *) multiVector_in;
00447   return(multiVector_->addVec(EPETRA_DEREF(scalar), multiVector_in_)));
00448 }
00449 
00450 int MANGLE(epetra_rdp_multivector_norm1)
00451   (EPETRA_MULTIVECTOR multiVector, double *scalar)
00452 {
00453   Epetra_MultiVector *multiVector_ = (Epetra_MultiVector *) multiVector;
00454   return(multiVector_->norm1(scalar));
00455 }
00456 
00457 int MANGLE(epetra_rdp_multivector_norm2)
00458   (EPETRA_MULTIVECTOR multiVector, double *scalar)
00459 {
00460   Epetra_MultiVector *multiVector_ = (Epetra_MultiVector *) multiVector;
00461   return(multiVector_->norm2(scalar));
00462 }
00463 
00464 int MANGLE(epetra_rdp_multivector_lincomb)(EPETRA_MULTIVECTOR multiVector,
00465              EPETRA_MULTIVECTOR b, 
00466              EPETRA_DOUBLE scalar, EPETRA_MULTIVECTOR c)
00467 {
00468   Epetra_MultiVector *multiVector_ = (Epetra_MultiVector *) multiVector;
00469   Epetra_MultiVector& b_ = *(Epetra_MultiVector *) b;
00470   Epetra_MultiVector& c_ = *(Epetra_MultiVector *) c;
00471   return(multiVector_->linComb(b_,EPETRA_DEREF(scalar,c_)));
00472 }
00473 
00474 int MANGLE(epetra_rdp_multivector_random)
00475   (EPETRA_MULTIVECTOR multiVector)
00476 {
00477   Epetra_MultiVector *multiVector_ = (Epetra_MultiVector *) multiVector;
00478   return(multiVector_->random());
00479 }
00480 
00481 int MANGLE(epetra_rdp_multivector_reduce)(EPETRA_MULTIVECTOR multiVector)
00482 {
00483   Epetra_MultiVector *multiVector_ = (Epetra_MultiVector *) multiVector;
00484   return(multiVector_->reduce());
00485 }
00486 
00487 int MANGLE(epetra_rdp_multivector_numvectors)(EPETRA_MULTIVECTOR multiVector)
00488 {
00489   Epetra_MultiVector *multiVector_ = (Epetra_MultiVector *) multiVector;
00490   return(multiVector_->numVectors());
00491 }
00492 
00493 int MANGLE(epetra_rdp_multivector_gemm)(EPETRA_MULTIVECTOR multiVector,
00494           EPETRA_INT transa, EPETRA_INT transb, EPETRA_DOUBLE alpha,
00495           EPETRA_MULTIVECTOR A, EPETRA_MULTIVECTOR B,
00496           EPETRA_DOUBLE beta )
00497 {
00498   Epetra_MultiVector *multiVector_ = (Epetra_MultiVector *) multiVector;
00499   Epetra_MultiVector& A_ = *(Epetra_MultiVector *) A;
00500   Epetra_MultiVector& B_ = *(Epetra_MultiVector *) B;
00501   bool transa_ = !(EPETRA_DEREF(transa==0));
00502   bool transb_ = !(EPETRA_DEREF(transb==0));
00503   return(multiVector_->GEMM(transa_, transb_, EPETRA_DEREF(alpha), A_, B_, EPETRA_DEREF(beta)));
00504 }
00505 
00506 void MANGLE(epetra_rdp_multivector_destroy)(EPETRA_MULTIVECTOR multiVector)
00507 {
00508   delete (Epetra_MultiVector *) multiVector;
00509 }
00510 
00512 //                  Epetra_BlockMap                //
00514 
00515 EPETRA_OBJECT_PTR MANGLE(epetra_blockmap_create1)(
00516               EPETRA_INT numGlobalElements, EPETRA_INT numLocalElements, int *updateList,
00517               EPETRA_INT numGlobalBlocks, EPETRA_INT numLocalBlocks, 
00518               int *blockUpdateList,
00519               int* blockSizes, EPETRA_INT indexBase, EPETRA_COMM comm)
00520 {
00521   Epetra_Comm& comm_ = *(Epetra_Comm *) comm;
00522   Epetra_BlockMap *blockmap = new Epetra_BlockMap(EPETRA_DEREF(numGlobalElements),
00523               EPETRA_DEREF(numLocalElements), updateList,
00524               EPETRA_DEREF(numGlobalBlocks), EPETRA_DEREF(numLocalBlocks),
00525               blockUpdateList,
00526               blockSizes, EPETRA_DEREF(indexBase), comm_);
00527   return((EPETRA_OBJECT_PTR ) blockmap);
00528 }
00529 
00530 EPETRA_OBJECT_PTR MANGLE(epetra_blockmap_create2)(
00531               EPETRA_INT numGlobalBlocks, EPETRA_INT numLocalBlocks, 
00532               int *blockUpdateList,
00533               int* blockSizes, EPETRA_INT indexBase, EPETRA_COMM comm)
00534 {
00535   Epetra_Comm& comm_ = *(Epetra_Comm *) comm;
00536   Epetra_BlockMap *blockmap = new Epetra_BlockMap(
00537               EPETRA_DEREF(numGlobalBlocks), EPETRA_DEREF(numLocalBlocks), 
00538               blockUpdateList,
00539               blockSizes, EPETRA_DEREF(indexBase), comm_);
00540   return((EPETRA_OBJECT_PTR ) blockmap);
00541 }
00542 
00543 void MANGLE(epetra_blockmap_destroy)(EPETRA_BLOCKMAP blockmap)
00544 {
00545   delete (Epetra_BlockMap *) blockmap;
00546 }
00547 
00549 //                  Epetra_LocalMap                //
00551 
00552 EPETRA_OBJECT_PTR MANGLE(epetra_localmap_create)(EPETRA_INT numLocalElements,
00553                    EPETRA_INT indexBase, EPETRA_COMM comm)
00554 {
00555   Epetra_Comm& comm_ = *(Epetra_Comm *) comm;
00556   Epetra_LocalMap *localmap = new Epetra_LocalMap(EPETRA_DEREF(numLocalElements),
00557               EPETRA_DEREF(indexBase), comm_);
00558   return((EPETRA_OBJECT_PTR ) localmap);
00559 }
00560 
00561 void MANGLE(epetra_localmap_destroy)(EPETRA_LOCALMAP localmap)
00562 {
00563   delete (Epetra_LocalMap *) localmap;
00564 }
00565 
00567 //                  Epetra_LocalBlockMap           //
00569 
00570 EPETRA_OBJECT_PTR MANGLE(epetra_localblockmap_create1)(
00571                   EPETRA_INT numLocalElements,
00572                   EPETRA_INT numLocalBlocks,
00573                   int* blockSizes,
00574                   EPETRA_INT indexBase, EPETRA_COMM comm)
00575 {
00576   Epetra_Comm& comm_ = *(Epetra_Comm *) comm;
00577   Epetra_LocalBlockMap *localblockmap = new
00578     Epetra_LocalBlockMap(EPETRA_DEREF(numLocalElements),
00579        EPETRA_DEREF(numLocalBlocks),
00580        blockSizes,
00581        EPETRA_DEREF(indexBase), comm_);
00582   return((EPETRA_OBJECT_PTR ) localblockmap);
00583 }
00584 
00585 EPETRA_OBJECT_PTR MANGLE(epetra_localblockmap_create2)(
00586                   EPETRA_INT numLocalBlocks,
00587                   int* blockSizes,
00588                   EPETRA_INT indexBase, EPETRA_COMM comm)
00589 {
00590   Epetra_Comm& comm_ = *(Epetra_Comm *) comm;
00591   Epetra_LocalBlockMap *localblockmap = new
00592     Epetra_LocalBlockMap(EPETRA_DEREF(numLocalBlocks),
00593        blockSizes,
00594        EPETRA_DEREF(indexBase), comm_);
00595   return((EPETRA_OBJECT_PTR ) localblockmap);
00596 }
00597 
00598 void MANGLE(epetra_localblockmap_destroy)(EPETRA_LOCALBLOCKMAP localblockmap)
00599 {
00600   delete (Epetra_LocalBlockMap *) localblockmap;
00601 }
00602 
00603 
00604 #endif /* 0 */
00605 
00606 #ifdef __cplusplus
00607 }
00608 #endif

Generated on Wed May 12 21:41:05 2010 for Epetra Package Browser (Single Doxygen Collection) by  doxygen 1.4.7