ifp_c_wrappers.cpp

00001 //@HEADER
00002 // ***********************************************************************
00003 // 
00004 //       Ifpack: Object-Oriented Algebraic Preconditioner Package
00005 //                 Copyright (2002) 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 #include "stdio.h"
00030 #include "ifp_c_wrappers.h"
00031 #include "ifp_BlockMat.h"
00032 #include "ifp_brelax.h"
00033 #include "ifp_biluk.h"
00034 #include "az_aztec.h"
00035 #include "ifp_ifpack.h"
00036 
00037 #ifdef __cplusplus
00038 using namespace std;
00039 extern "C" {
00040 #endif
00041 
00042 // prototypes for this file
00043 static void set_localprecon(ifp_GlobalPrecon *M, int local, 
00044   double lparam1, double lparam2);
00045 
00046 // catch memory errors
00047 //static void freeStoreException()
00048 //{
00049 //    cout << endl; // flush standard output
00050 //    ifp_error("ifp_Fortran: free store exhausted", 0);
00051 //}
00052 
00053 // void ifp_initialize()
00054 // {
00055 //     void (*_new_handler)();
00056 //     _new_handler = freeStoreException;
00057 // }
00058 
00059 void az2ifp_blockmatrix (void **bmat, AZ_MATRIX *Amat)
00060 {
00061 
00062   double *val    = Amat->val;
00063   int *bindx  = Amat->bindx;
00064   int *indx = Amat->indx;
00065   int *bpntr = Amat->bpntr;
00066   int *rpntr = Amat->rpntr;
00067   int *cpntr = Amat->cpntr;
00068   int *data_org = Amat->data_org;
00069   int nrow = data_org[AZ_N_int_blk] + data_org[AZ_N_bord_blk];
00070   int ncol = nrow; /* Assume single processor for now; */
00071   int nnz = bpntr[nrow]; /* number of block nonzeros */
00072   int nnzs = indx[nnz]; /* number of scalar nonzeros */
00073 
00074   *bmat = (void *) new ifp_BlockMat(val, indx, bindx, rpntr, cpntr,
00075            bpntr, nrow, ncol, nnz, nnzs);
00076 }
00077 
00078 
00079 void ifp_freeblockmatrix(void *bmat)
00080 {
00081     delete (ifp_BlockMat *) bmat;
00082 }
00083 
00084 void ifp_preconditioner(
00085         void    **precon,
00086   const void    *bmat,
00087   const int     global,
00088   const double  gparam1,
00089   const double  gparam2,
00090   const int     local,
00091   const double  lparam1,
00092   const double  lparam2)
00093 {
00094   ifp_GlobalPrecon *M;
00095   ifp_BlockMat *B = (ifp_BlockMat *) bmat;
00096 
00097   switch (global)
00098     {
00099     case 0:
00100       M = new ifp_None;
00101       break;
00102     case 1:
00103       M = new ifp_BJacobi;
00104       set_localprecon(M, (int)local, lparam1, lparam2);
00105       ((ifp_BJacobi *)M)->setup(*B);
00106       
00107       break;
00108     case 2:
00109       M = new ifp_BSOR;
00110       set_localprecon(M, (int)local, lparam1, lparam2);
00111       ((ifp_BSOR *)M)->setup(*B, (double)gparam1, (int)gparam2);
00112       break;
00113     case 3:
00114       M = new ifp_BSSOR;
00115       set_localprecon(M, (int)local, lparam1, lparam2);
00116       ((ifp_BSSOR *)M)->setup(*B, (double)gparam1, (int)gparam2);
00117       break;
00118     case 4:
00119       M = new ifp_biluk;
00120       set_localprecon(M, (int)local, lparam1, lparam2);
00121       ((ifp_biluk *)M)->setup(*B, (int)gparam1);
00122       
00123       break;
00124     default:
00125       ifp_error("ifp_c_wrapper: no such global preconditioner", 0);
00126     }
00127   *precon = (void *) M;
00128 }
00129 
00130 void ifp_freebiluk(void *precon)
00131 {
00132     ifp_biluk * bilukPrec = (ifp_biluk *) precon;
00133     ifp_BlockMat * bilukMat = bilukPrec->A();
00134     delete bilukPrec;
00135     delete bilukMat;
00136 }
00137 void ifp_freepreconditioner(void *precon)
00138 {
00139     delete (ifp_GlobalPrecon *) precon;
00140 }
00141 
00142 // static functions
00143 
00144 #define INT(d) ((int)(d+0.5)) // round to integer
00145 
00146 static void set_localprecon(ifp_GlobalPrecon *M, int local_, 
00147   double lparam1, double lparam2)
00148 {
00149     LocalPreconName local = (LocalPreconName) local_;
00150 
00151     switch (local)
00152     {
00153     case LP_LU:
00154     M->localprecon(local);
00155     break;
00156     case LP_INVERSE:
00157     M->localprecon(local);
00158     break;
00159     case LP_SVD:
00160     M->localprecon(local, lparam1, lparam2);
00161     break;
00162     case LP_DIAG:
00163     M->localprecon(local);
00164     break;
00165     case LP_SOR:
00166     M->localprecon(local, lparam1, INT(lparam2));
00167     break;
00168     case LP_SSOR:
00169     M->localprecon(local, lparam1, INT(lparam2));
00170     break;
00171     default:
00172         ifp_error("ifp_Fortran: no such local preconditioner", 0);
00173     }
00174 }
00175 
00176 void ifp_matvec(
00177   void    *bmat,
00178   int     nr,
00179   int     nc,
00180   const double *u,
00181   int     ldu,
00182   double *v,
00183   int     ldv)
00184 {
00185     ifp_BlockMat *B = (ifp_BlockMat *) bmat;
00186     B->mult(nr, nc, u, ldu, v, ldv);
00187 }
00188 
00189 void ifp_apply(
00190   void    *prec,
00191   int     nr,
00192   int     nc,
00193   const double *u,
00194   int     ldu,
00195   double *v,
00196   int     ldv)
00197 {
00198     ifp_GlobalPrecon *M = (ifp_GlobalPrecon *) prec;
00199     M->apply(nr, nc, u, ldu, v, ldv);
00200 }
00201 
00202 void ifp_BJacobi_condest(void *M)
00203 {
00204       double cond_number = ((ifp_BJacobi *)M)->condest();
00205       printf ("Condition number estimate = %.4e\n",cond_number);
00206 }
00207 
00208 void ifp_biluk_condest(void *M)
00209 {
00210       double cond_number = ((ifp_biluk *)M)->condest();
00211       printf ("Condition number estimate = %.4e\n",cond_number);
00212 }
00213 
00214 void ifp_biluk_stats(void *M)
00215 {
00216   ifp_biluk * Mp = (ifp_biluk * )M;
00217   ifp_BlockMat* Ap = Mp->A();
00218   int NumEntries_M = Mp->NumEntries();
00219   int NumNonzeros_M = Mp->NumNonzeros();
00220   int NumEntries_A = Ap->NumEntries();
00221   int NumNonzeros_A = Ap->NumNonzeros();
00222   double Entries_Ratio = ((double) NumEntries_M) / ((double) NumEntries_A);
00223   double Nonzeros_Ratio = ((double) NumNonzeros_M) / ((double) NumNonzeros_A);
00224   printf ("\n");
00225   printf ("********************************************************************\n");
00226   printf ("***** Number of Block Entries in User Matrix          = %d\n",NumEntries_A);
00227   printf ("***** Number of Nonzero Values in User Matrix         = %d\n",NumNonzeros_A);
00228   printf ("***** Number of Block Entries in IFPACK BILU Factors  = %d\n",NumEntries_M);
00229   printf ("***** Number of Nonzero Values in IFPACK BILU Factors = %d\n",NumNonzeros_M);
00230   printf ("***** Ratio of Block Entries in M to A                = %.4e\n",Entries_Ratio);
00231   printf ("***** Ratio of Nonzero Values in M to A               = %.4e\n",Nonzeros_Ratio);
00232   printf ("********************************************************************\n");
00233   printf ("\n");
00234 }
00235 
00236 #ifdef __cplusplus
00237 }
00238 #endif

Generated on Tue Oct 20 12:48:53 2009 for IFPACK by doxygen 1.4.7