Amesos Package Browser (Single Doxygen Collection) Development
amesos_t_cholmod_triplet.c
Go to the documentation of this file.
00001 /* ========================================================================== */
00002 /* === Core/t_cholmod_triplet =============================================== */
00003 /* ========================================================================== */
00004 
00005 /* -----------------------------------------------------------------------------
00006  * CHOLMOD/Core Module.  Copyright (C) 2005-2006,
00007  * Univ. of Florida.  Author: Timothy A. Davis
00008  * The CHOLMOD/Core Module is licensed under Version 2.1 of the GNU
00009  * Lesser General Public License.  See lesser.txt for a text of the license.
00010  * CHOLMOD is also available under other licenses; contact authors for details.
00011  * http://www.cise.ufl.edu/research/sparse
00012  * -------------------------------------------------------------------------- */
00013 
00014 /* Template routine for cholmod_triplet.  All xtypes supported */
00015 
00016 #include "amesos_cholmod_template.h"
00017 
00018 /* ========================================================================== */
00019 /* === t_cholmod_triplet_to_sparse ========================================== */
00020 /* ========================================================================== */
00021 
00022 static size_t TEMPLATE (cholmod_triplet_to_sparse)
00023 (
00024     /* ---- input ---- */
00025     cholmod_triplet *T, /* matrix to copy */
00026     /* ---- in/out --- */
00027     cholmod_sparse *R,  /* output matrix */
00028     /* --------------- */
00029     cholmod_common *Common
00030 )
00031 {
00032     double *Rx, *Rz, *Tx, *Tz ;
00033     Int *Wj, *Rp, *Ri, *Rnz, *Ti, *Tj  ;
00034     Int i, j, p, p1, p2, pdest, pj, k, stype, nrow, ncol, nz ;
00035     size_t anz ;
00036 
00037     /* ---------------------------------------------------------------------- */
00038     /* get inputs */
00039     /* ---------------------------------------------------------------------- */
00040 
00041     /* Wj contains a copy of Rp on input [ */
00042     Wj = Common->Iwork ;  /* size MAX (nrow,ncol). (i/l/l) */
00043 
00044     Rp = R->p ;
00045     Ri = R->i ;
00046     Rnz = R->nz ;
00047     Rx = R->x ;
00048     Rz = R->z ;
00049 
00050     Ti = T->i ;
00051     Tj = T->j ;
00052     Tx = T->x ;
00053     Tz = T->z ;
00054     nz = T->nnz ;
00055     nrow = T->nrow ;
00056     ncol = T->ncol ;
00057     stype = SIGN (T->stype) ;
00058 
00059     /* ---------------------------------------------------------------------- */
00060     /* construct the row form */
00061     /* ---------------------------------------------------------------------- */
00062 
00063     /* if Ti is jumbled, this part dominates the run time */
00064 
00065     if (stype > 0)
00066     {
00067   for (k = 0 ; k < nz ; k++)
00068   {
00069       i = Ti [k] ;
00070       j = Tj [k] ;
00071       if (i < j)
00072       {
00073     /* place triplet (j,i,x) in column i of R */
00074     p = Wj [i]++ ;
00075     Ri [p] = j ;
00076       }
00077       else
00078       {
00079     /* place triplet (i,j,x) in column j of R */
00080     p = Wj [j]++ ;
00081     Ri [p] = i ;
00082       }
00083       ASSIGN (Rx, Rz, p, Tx, Tz, k) ;
00084   }
00085     }
00086     else if (stype < 0)
00087     {
00088   for (k = 0 ; k < nz ; k++)
00089   {
00090       i = Ti [k] ;
00091       j = Tj [k] ;
00092       if (i > j)
00093       {
00094     /* place triplet (j,i,x) in column i of R */
00095     p = Wj [i]++ ;
00096     Ri [p] = j ;
00097       }
00098       else
00099       {
00100     /* place triplet (i,j,x) in column j of R */
00101     p = Wj [j]++ ;
00102     Ri [p] = i ;
00103       }
00104       ASSIGN (Rx, Rz, p, Tx, Tz, k) ;
00105   }
00106     }
00107     else
00108     {
00109   for (k = 0 ; k < nz ; k++)
00110   {
00111       /* place triplet (i,j,x) in column i of R */
00112       p = Wj [Ti [k]]++ ;
00113       Ri [p] = Tj [k] ;
00114       ASSIGN (Rx, Rz, p, Tx, Tz, k) ;
00115   }
00116     }
00117 
00118     /* done using Wj (i/l/l) as temporary row pointers ] */
00119 
00120     /* ---------------------------------------------------------------------- */
00121     /* sum up duplicates */
00122     /* ---------------------------------------------------------------------- */
00123 
00124     /* use Wj (i/l/l) of size ncol to keep track of duplicates in each row [ */
00125     for (j = 0 ; j < ncol ; j++)
00126     {
00127   Wj [j] = EMPTY ;
00128     }
00129 
00130     anz = 0 ;
00131     for (i = 0 ; i < nrow ; i++)
00132     {
00133   p1 = Rp [i] ;
00134   p2 = Rp [i+1] ;
00135   pdest = p1 ;
00136   /* at this point Wj [j] < p1 holds true for all columns j, because
00137    * Ri/Rx is stored in row oriented manner */
00138   for (p = p1 ; p < p2 ; p++)
00139   {
00140       j = Ri [p] ;
00141       pj = Wj [j] ;
00142       if (pj >= p1)
00143       {
00144     /* this column index j is already in row i at position pj;
00145      * sum up the duplicate entry */
00146     /* Rx [pj] += Rx [p] ; */
00147     ASSEMBLE (Rx, Rz, pj, Rx, Rz, p) ;
00148       }
00149       else
00150       {
00151     /* keep the entry and keep track in Wj [j] for case above */
00152     Wj [j] = pdest ;
00153     if (pdest != p)
00154     {
00155         Ri [pdest] = j ;
00156         ASSIGN (Rx, Rz, pdest, Rx, Rz, p) ;
00157     }
00158     pdest++ ;
00159       }
00160   }
00161   Rnz [i] = pdest - p1 ;
00162   anz += (pdest - p1) ;
00163     }
00164     /* done using Wj to keep track of duplicate entries in each row ] */
00165 
00166     /* ---------------------------------------------------------------------- */
00167     /* return number of entries after summing up duplicates */
00168     /* ---------------------------------------------------------------------- */
00169 
00170     return (anz) ;
00171 }
00172 
00173 #undef PATTERN
00174 #undef REAL
00175 #undef COMPLEX
00176 #undef ZOMPLEX
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines