Ifpack Package Browser (Single Doxygen Collection) Development
Vec_dh.c
Go to the documentation of this file.
00001 /*@HEADER
00002 // ***********************************************************************
00003 //
00004 //       Ifpack: Object-Oriented Algebraic Preconditioner Package
00005 //                 Copyright (2009) 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 #include <stdlib.h>
00031 #include "Vec_dh.h"
00032 #include "Mem_dh.h"
00033 #include "SubdomainGraph_dh.h"
00034 #include "io_dh.h"
00035 
00036 #undef __FUNC__
00037 #define __FUNC__ "Vec_dhCreate"
00038 void
00039 Vec_dhCreate (Vec_dh * v)
00040 {
00041   START_FUNC_DH
00042     struct _vec_dh *tmp =
00043     (struct _vec_dh *) MALLOC_DH (sizeof (struct _vec_dh));
00044   CHECK_V_ERROR;
00045   *v = tmp;
00046   tmp->n = 0;
00047   tmp->vals = NULL;
00048 END_FUNC_DH}
00049 
00050 #undef __FUNC__
00051 #define __FUNC__ "Vec_dhDestroy"
00052 void
00053 Vec_dhDestroy (Vec_dh v)
00054 {
00055   START_FUNC_DH if (v->vals != NULL)
00056     FREE_DH (v->vals);
00057   CHECK_V_ERROR;
00058   FREE_DH (v);
00059   CHECK_V_ERROR;
00060 END_FUNC_DH}
00061 
00062 
00063 #undef __FUNC__
00064 #define __FUNC__ "Vec_dhInit"
00065 void
00066 Vec_dhInit (Vec_dh v, int size)
00067 {
00068   START_FUNC_DH v->n = size;
00069   v->vals = (double *) MALLOC_DH (size * sizeof (double));
00070   CHECK_V_ERROR;
00071 END_FUNC_DH}
00072 
00073 #undef __FUNC__
00074 #define __FUNC__ "Vec_dhCopy"
00075 void
00076 Vec_dhCopy (Vec_dh x, Vec_dh y)
00077 {
00078   START_FUNC_DH if (x->vals == NULL)
00079     SET_V_ERROR ("x->vals is NULL");
00080   if (y->vals == NULL)
00081     SET_V_ERROR ("y->vals is NULL");
00082   if (x->n != y->n)
00083     SET_V_ERROR ("x and y are different lengths");
00084   memcpy (y->vals, x->vals, x->n * sizeof (double));
00085 END_FUNC_DH}
00086 
00087 
00088 #undef __FUNC__
00089 #define __FUNC__ "Vec_dhDuplicate"
00090 void
00091 Vec_dhDuplicate (Vec_dh v, Vec_dh * out)
00092 {
00093   START_FUNC_DH Vec_dh tmp;
00094   int size = v->n;
00095   if (v->vals == NULL)
00096     SET_V_ERROR ("v->vals is NULL");
00097   Vec_dhCreate (out);
00098   CHECK_V_ERROR;
00099   tmp = *out;
00100   tmp->n = size;
00101   tmp->vals = (double *) MALLOC_DH (size * sizeof (double));
00102   CHECK_V_ERROR;
00103 END_FUNC_DH}
00104 
00105 #undef __FUNC__
00106 #define __FUNC__ "Vec_dhSet"
00107 void
00108 Vec_dhSet (Vec_dh v, double value)
00109 {
00110   START_FUNC_DH int i, m = v->n;
00111   double *vals = v->vals;
00112   if (v->vals == NULL)
00113     SET_V_ERROR ("v->vals is NULL");
00114   for (i = 0; i < m; ++i)
00115     vals[i] = value;
00116 END_FUNC_DH}
00117 
00118 #undef __FUNC__
00119 #define __FUNC__ "Vec_dhSetRand"
00120 void
00121 Vec_dhSetRand (Vec_dh v)
00122 {
00123   START_FUNC_DH int i, m = v->n;
00124   double max = 0.0;
00125   double *vals = v->vals;
00126 
00127   if (v->vals == NULL)
00128     SET_V_ERROR ("v->vals is NULL");
00129 
00130 #ifdef WIN32
00131   for (i = 0; i < m; ++i)
00132     vals[i] = rand ();
00133 #else
00134   for (i = 0; i < m; ++i)
00135     vals[i] = rand ();
00136 #endif
00137 
00138   /* find largest value in vector, and scale vector,
00139    * so all values are in [0.0,1.0]
00140    */
00141   for (i = 0; i < m; ++i)
00142     max = MAX (max, vals[i]);
00143   for (i = 0; i < m; ++i)
00144     vals[i] = vals[i] / max;
00145 END_FUNC_DH}
00146 
00147 
00148 #undef __FUNC__
00149 #define __FUNC__ "Vec_dhPrint"
00150 void
00151 Vec_dhPrint (Vec_dh v, SubdomainGraph_dh sg, char *filename)
00152 {
00153   START_FUNC_DH double *vals = v->vals;
00154   int pe, i, m = v->n;
00155   FILE *fp;
00156 
00157   if (v->vals == NULL)
00158     SET_V_ERROR ("v->vals is NULL");
00159 
00160   /*--------------------------------------------------------
00161    * case 1: no permutation information
00162    *--------------------------------------------------------*/
00163   if (sg == NULL)
00164     {
00165       for (pe = 0; pe < np_dh; ++pe)
00166   {
00167     MPI_Barrier (comm_dh);
00168     if (pe == myid_dh)
00169       {
00170         if (pe == 0)
00171     {
00172       fp = openFile_dh (filename, "w");
00173       CHECK_V_ERROR;
00174     }
00175         else
00176     {
00177       fp = openFile_dh (filename, "a");
00178       CHECK_V_ERROR;
00179     }
00180 
00181         for (i = 0; i < m; ++i)
00182     fprintf (fp, "%g\n", vals[i]);
00183 
00184         closeFile_dh (fp);
00185         CHECK_V_ERROR;
00186       }
00187   }
00188     }
00189 
00190   /*--------------------------------------------------------
00191    * case 2: single mpi task, multiple subdomains
00192    *--------------------------------------------------------*/
00193   else if (np_dh == 1)
00194     {
00195       int i, j;
00196 
00197       fp = openFile_dh (filename, "w");
00198       CHECK_V_ERROR;
00199 
00200       for (i = 0; i < sg->blocks; ++i)
00201   {
00202     int oldBlock = sg->n2o_sub[i];
00203     int beg_row = sg->beg_rowP[oldBlock];
00204     int end_row = beg_row + sg->row_count[oldBlock];
00205 
00206     printf ("seq: block= %i  beg= %i  end= %i\n", oldBlock, beg_row,
00207       end_row);
00208 
00209 
00210     for (j = beg_row; j < end_row; ++j)
00211       {
00212         fprintf (fp, "%g\n", vals[j]);
00213       }
00214   }
00215     }
00216 
00217   /*--------------------------------------------------------
00218    * case 3: multiple mpi tasks, one subdomain per task
00219    *--------------------------------------------------------*/
00220   else
00221     {
00222       int id = sg->o2n_sub[myid_dh];
00223       for (pe = 0; pe < np_dh; ++pe)
00224   {
00225     MPI_Barrier (comm_dh);
00226     if (id == pe)
00227       {
00228         if (pe == 0)
00229     {
00230       fp = openFile_dh (filename, "w");
00231       CHECK_V_ERROR;
00232     }
00233         else
00234     {
00235       fp = openFile_dh (filename, "a");
00236       CHECK_V_ERROR;
00237     }
00238 
00239         fprintf (stderr, "par: block= %i\n", id);
00240 
00241         for (i = 0; i < m; ++i)
00242     {
00243       fprintf (fp, "%g\n", vals[i]);
00244     }
00245 
00246         closeFile_dh (fp);
00247         CHECK_V_ERROR;
00248       }
00249   }
00250     }
00251 END_FUNC_DH}
00252 
00253 
00254 #undef __FUNC__
00255 #define __FUNC__ "Vec_dhPrintBIN"
00256 void
00257 Vec_dhPrintBIN (Vec_dh v, SubdomainGraph_dh sg, char *filename)
00258 {
00259   START_FUNC_DH if (np_dh > 1)
00260     {
00261       SET_V_ERROR ("only implemented for a single MPI task");
00262     }
00263   if (sg != NULL)
00264     {
00265       SET_V_ERROR ("not implemented for reordered vector; ensure sg=NULL");
00266     }
00267 
00268   io_dh_print_ebin_vec_private (v->n, 0, v->vals, NULL, NULL, NULL, filename);
00269   CHECK_V_ERROR;
00270 END_FUNC_DH}
00271 
00272 #define MAX_JUNK 200
00273 
00274 #undef __FUNC__
00275 #define __FUNC__ "Vec_dhRead"
00276 void
00277 Vec_dhRead (Vec_dh * vout, int ignore, char *filename)
00278 {
00279   START_FUNC_DH Vec_dh tmp;
00280   FILE *fp;
00281   int items, n, i;
00282   double *v, w;
00283   char junk[MAX_JUNK];
00284 
00285   Vec_dhCreate (&tmp);
00286   CHECK_V_ERROR;
00287   *vout = tmp;
00288 
00289   if (np_dh > 1)
00290     {
00291       SET_V_ERROR ("only implemented for a single MPI task");
00292     }
00293 
00294   fp = openFile_dh (filename, "w");
00295   CHECK_V_ERROR;
00296 
00297   /* skip over file lines */
00298   if (ignore)
00299     {
00300       printf ("Vec_dhRead:: ignoring following header lines:\n");
00301       printf
00302   ("--------------------------------------------------------------\n");
00303       for (i = 0; i < ignore; ++i)
00304   {
00305     fgets (junk, MAX_JUNK, fp);
00306     printf ("%s", junk);
00307   }
00308       printf
00309   ("--------------------------------------------------------------\n");
00310     }
00311 
00312   /* count floating point entries in file */
00313   n = 0;
00314   while (!feof (fp))
00315     {
00316       items = fscanf (fp, "%lg", &w);
00317       if (items != 1)
00318   {
00319     break;
00320   }
00321       ++n;
00322     }
00323 
00324   printf ("Vec_dhRead:: n= %i\n", n);
00325 
00326   /* allocate storage */
00327   tmp->n = n;
00328   v = tmp->vals = (double *) MALLOC_DH (n * sizeof (double));
00329   CHECK_V_ERROR;
00330 
00331   /* reset file, and skip over header again */
00332   rewind (fp);
00333   rewind (fp);
00334   for (i = 0; i < ignore; ++i)
00335     {
00336       fgets (junk, MAX_JUNK, fp);
00337     }
00338 
00339   /* read values */
00340   for (i = 0; i < n; ++i)
00341     {
00342       items = fscanf (fp, "%lg", v + i);
00343       if (items != 1)
00344   {
00345     sprintf (msgBuf_dh, "failed to read value %i of %i", i + 1, n);
00346   }
00347     }
00348 
00349   closeFile_dh (fp);
00350   CHECK_V_ERROR;
00351 END_FUNC_DH}
00352 
00353 #undef __FUNC__
00354 #define __FUNC__ "Vec_dhReadBIN"
00355 extern void
00356 Vec_dhReadBIN (Vec_dh * vout, char *filename)
00357 {
00358   START_FUNC_DH Vec_dh tmp;
00359 
00360   Vec_dhCreate (&tmp);
00361   CHECK_V_ERROR;
00362   *vout = tmp;
00363   io_dh_read_ebin_vec_private (&tmp->n, &tmp->vals, filename);
00364   CHECK_V_ERROR;
00365 END_FUNC_DH}
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines