Amesos Package Browser (Single Doxygen Collection) Development
amesos_cholmod_factor.c
Go to the documentation of this file.
00001 /* ========================================================================== */
00002 /* === Core/cholmod_factor ================================================== */
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 /* Core utility routines for the cholmod_factor object:
00015  *
00016  * The data structure for an LL' or LDL' factorization is too complex to
00017  * describe in one sentence.  This object can hold the symbolic analysis alone,
00018  * or in combination with a "simplicial" (similar to a sparse matrix) or
00019  * "supernodal" form of the numerical factorization.  Only the routine to free
00020  * a factor is primary, since a factor object is created by the factorization
00021  * routine (cholmod_factorize).  It must be freed with cholmod_free_factor.
00022  *
00023  * Primary routine:
00024  * ----------------
00025  * cholmod_free_factor    free a factor
00026  *
00027  * Secondary routines:
00028  * -------------------
00029  * cholmod_allocate_factor  allocate a symbolic factor (LL' or LDL')
00030  * cholmod_reallocate_factor  change the # entries in a factor 
00031  * cholmod_change_factor  change the type of factor (e.g., LDL' to LL')
00032  * cholmod_pack_factor    pack the columns of a factor
00033  * cholmod_reallocate_column  resize a single column of a factor
00034  * cholmod_factor_to_sparse create a sparse matrix copy of a factor
00035  * cholmod_copy_factor    create a copy of a factor
00036  *
00037  * Note that there is no cholmod_sparse_to_factor routine to create a factor
00038  * as a copy of a sparse matrix.  It could be done, after a fashion, but a
00039  * lower triangular sparse matrix would not necessarily have a chordal graph,
00040  * which would break the many CHOLMOD routines that rely on this property.
00041  *
00042  * The cholmod_factor_to_sparse routine is provided so that matrix operations
00043  * in the MatrixOps module may be applied to L.  Those operations operate on
00044  * cholmod_sparse objects, and they are not guaranteed to maintain the chordal
00045  * property of L.  Such a modified L cannot be safely converted back to a
00046  * cholmod_factor object.
00047  */
00048 
00049 #include "amesos_cholmod_internal.h"
00050 #include "amesos_cholmod_core.h"
00051 
00052 
00053 /* ========================================================================== */
00054 /* === cholmod_allocate_factor ============================================== */
00055 /* ========================================================================== */
00056 
00057 /* Allocate a simplicial symbolic factor, with L->Perm and L->ColCount allocated
00058  * and initialized to "empty" values (Perm [k] = k, and ColCount[k] = 1).
00059  * The integer and numerical parts of L are not allocated.  L->xtype is returned
00060  * as CHOLMOD_PATTERN and L->is_super are returned as FALSE.  L->is_ll is also
00061  * returned FALSE, but this may be modified when the matrix is factorized.
00062  *
00063  * This is sufficient (but far from ideal) for input to cholmod_factorize,
00064  * since the simplicial LL' or LDL' factorization (cholmod_rowfac) can
00065  * reallocate the columns of L as needed.  The primary purpose of this routine
00066  * is to allocate space for a symbolic factorization, for the "expert" user to
00067  * do his or her own symbolic analysis.  The typical user should use
00068  * cholmod_analyze instead of this routine.
00069  *
00070  * workspace: none
00071  */
00072 
00073 cholmod_factor *CHOLMOD(allocate_factor)
00074 (
00075     /* ---- input ---- */
00076     size_t n,   /* L is n-by-n */
00077     /* --------------- */
00078     cholmod_common *Common
00079 )
00080 {
00081     Int j ;
00082     Int *Perm, *ColCount ;
00083     cholmod_factor *L ;
00084     int ok = TRUE ;
00085 
00086     RETURN_IF_NULL_COMMON (FALSE) ;
00087     Common->status = CHOLMOD_OK ;
00088 
00089     /* ensure the dimension does not cause integer overflow */
00090     (void) CHOLMOD(add_size_t) (n, 2, &ok) ;
00091     if (!ok || n > Int_max)
00092     {
00093   ERROR (CHOLMOD_TOO_LARGE, "problem too large") ;
00094   return (NULL) ;
00095     }
00096 
00097     L = CHOLMOD(malloc) (sizeof (cholmod_factor), 1, Common) ;
00098     if (Common->status < CHOLMOD_OK)
00099     {
00100   return (NULL) ;     /* out of memory */
00101     }
00102     L->n = n ;
00103     L->is_ll = FALSE ;
00104     L->is_super = FALSE ;
00105     L->is_monotonic = TRUE ;
00106     L->itype = ITYPE ;
00107     L->xtype = CHOLMOD_PATTERN ;
00108     L->dtype = DTYPE ;
00109 
00110     /* allocate the purely symbolic part of L */
00111     L->ordering = CHOLMOD_NATURAL ;
00112     L->Perm = CHOLMOD(malloc) (n, sizeof (Int), Common) ;
00113     L->ColCount = CHOLMOD(malloc) (n, sizeof (Int), Common) ;
00114 
00115     /* simplicial part of L is empty */
00116     L->nzmax = 0 ;
00117     L->p = NULL ;
00118     L->i = NULL ;
00119     L->x = NULL ;
00120     L->z = NULL ;
00121     L->nz = NULL ;
00122     L->next = NULL ;
00123     L->prev = NULL ;
00124 
00125     /* supernodal part of L is also empty */
00126     L->nsuper = 0 ;
00127     L->ssize = 0 ;
00128     L->xsize = 0 ;
00129     L->maxesize = 0 ;
00130     L->maxcsize = 0 ;
00131     L->super = NULL ;
00132     L->pi = NULL ;
00133     L->px = NULL ;
00134     L->s = NULL ;
00135 
00136     /* L has not been factorized */
00137     L->minor = n ;
00138 
00139     if (Common->status < CHOLMOD_OK)
00140     {
00141   CHOLMOD(free_factor) (&L, Common) ;
00142   return (NULL) ;   /* out of memory */
00143     }
00144 
00145     /* initialize Perm and ColCount */
00146     Perm = L->Perm ;
00147     for (j = 0 ; j < ((Int) n) ; j++)
00148     {
00149   Perm [j] = j ;
00150     }
00151     ColCount = L->ColCount ;
00152     for (j = 0 ; j < ((Int) n) ; j++)
00153     {
00154   ColCount [j] = 1 ;
00155     }
00156 
00157     return (L) ;
00158 }
00159 
00160 
00161 /* ========================================================================== */
00162 /* === cholmod_free_factor ================================================== */
00163 /* ========================================================================== */
00164 
00165 /* Free a factor object.
00166  *
00167  * workspace: none
00168  */
00169 
00170 int CHOLMOD(free_factor)
00171 (
00172     /* ---- in/out --- */
00173     cholmod_factor **LHandle, /* factor to free, NULL on output */
00174     /* --------------- */
00175     cholmod_common *Common
00176 )
00177 {
00178     Int n, lnz, xs, ss, s ;
00179     cholmod_factor *L ;
00180 
00181     RETURN_IF_NULL_COMMON (FALSE) ;
00182 
00183     if (LHandle == NULL)
00184     {
00185   /* nothing to do */
00186   return (TRUE) ;
00187     }
00188     L = *LHandle ;
00189     if (L == NULL)
00190     {
00191   /* nothing to do */
00192   return (TRUE) ;
00193     }
00194 
00195     n = L->n ;
00196     lnz = L->nzmax ;
00197     s = L->nsuper + 1 ;
00198     xs = (L->is_super) ? ((Int) (L->xsize)) : (lnz) ;
00199     ss = L->ssize ;
00200 
00201     /* symbolic part of L */
00202     CHOLMOD(free) (n,   sizeof (Int), L->Perm,     Common) ;
00203     CHOLMOD(free) (n,   sizeof (Int), L->ColCount, Common) ;
00204 
00205     /* simplicial form of L */
00206     CHOLMOD(free) (n+1, sizeof (Int), L->p,        Common) ;
00207     CHOLMOD(free) (lnz, sizeof (Int), L->i,        Common) ;
00208     CHOLMOD(free) (n,   sizeof (Int), L->nz,       Common) ;
00209     CHOLMOD(free) (n+2, sizeof (Int), L->next,     Common) ;
00210     CHOLMOD(free) (n+2, sizeof (Int), L->prev,     Common) ;
00211 
00212     /* supernodal form of L */
00213     CHOLMOD(free) (s,   sizeof (Int), L->pi,       Common) ;
00214     CHOLMOD(free) (s,   sizeof (Int), L->px,       Common) ;
00215     CHOLMOD(free) (s,   sizeof (Int), L->super,    Common) ;
00216     CHOLMOD(free) (ss,  sizeof (Int), L->s,        Common) ;
00217 
00218     /* numerical values for both simplicial and supernodal L */
00219     if (L->xtype == CHOLMOD_REAL)
00220     {
00221   CHOLMOD(free) (xs, sizeof (double), L->x, Common) ;
00222     }
00223     else if (L->xtype == CHOLMOD_COMPLEX)
00224     {
00225   CHOLMOD(free) (xs, 2*sizeof (double), L->x, Common) ;
00226     }
00227     else if (L->xtype == CHOLMOD_ZOMPLEX)
00228     {
00229   CHOLMOD(free) (xs, sizeof (double), L->x, Common) ;
00230   CHOLMOD(free) (xs, sizeof (double), L->z, Common) ;
00231     }
00232 
00233     *LHandle = CHOLMOD(free) (1, sizeof (cholmod_factor), (*LHandle), Common) ;
00234     return (TRUE) ;
00235 }
00236 
00237 
00238 /* ========================================================================== */
00239 /* === cholmod_reallocate_factor ============================================ */
00240 /* ========================================================================== */
00241 
00242 /* Change the size of L->i and L->x, or allocate them if their current size
00243  * is zero.  L must be simplicial.
00244  *
00245  * workspace: none
00246  */
00247 
00248 int CHOLMOD(reallocate_factor)
00249 (
00250     /* ---- input ---- */
00251     size_t nznew, /* new # of entries in L */
00252     /* ---- in/out --- */
00253     cholmod_factor *L,  /* factor to modify */
00254     /* --------------- */
00255     cholmod_common *Common
00256 )
00257 {
00258     /* ---------------------------------------------------------------------- */
00259     /* get inputs */
00260     /* ---------------------------------------------------------------------- */
00261 
00262     RETURN_IF_NULL_COMMON (FALSE) ;
00263     RETURN_IF_NULL (L, FALSE) ;
00264     RETURN_IF_XTYPE_INVALID (L, CHOLMOD_REAL, CHOLMOD_ZOMPLEX, FALSE) ;
00265     PRINT1 (("realloc factor: xtype %d\n", L->xtype)) ;
00266     if (L->is_super)
00267     {
00268   /* L must be simplicial, and not symbolic */
00269   ERROR (CHOLMOD_INVALID, "L invalid") ;
00270   return (FALSE) ;
00271     }
00272     Common->status = CHOLMOD_OK ;
00273     PRINT1 (("realloc factor %g to %g\n", (double) L->nzmax, (double) nznew)) ;
00274 
00275     /* ---------------------------------------------------------------------- */
00276     /* resize (or allocate) the L->i and L->x components of the factor */
00277     /* ---------------------------------------------------------------------- */
00278 
00279     CHOLMOD(realloc_multiple) (nznew, 1, L->xtype, &(L->i), NULL,
00280       &(L->x), &(L->z), &(L->nzmax), Common) ;
00281     return (Common->status == CHOLMOD_OK) ;
00282 }
00283 
00284 
00285 /* ========================================================================== */
00286 /* === cholmod_reallocate_column =========================================== */
00287 /* ========================================================================== */
00288 
00289 /* Column j needs more space, reallocate it at the end of L->i and L->x.
00290  * If the reallocation fails, the factor is converted to a simplicial
00291  * symbolic factor (no pattern, just L->Perm and L->ColCount).
00292  *
00293  * workspace: none
00294  */
00295 
00296 int CHOLMOD(reallocate_column)
00297 (
00298     /* ---- input ---- */
00299     size_t j,   /* the column to reallocate */
00300     size_t need,  /* required size of column j */
00301     /* ---- in/out --- */
00302     cholmod_factor *L,  /* factor to modify */
00303     /* --------------- */
00304     cholmod_common *Common
00305 )
00306 {
00307     double xneed ;
00308     double *Lx, *Lz ;
00309     Int *Lp, *Lprev, *Lnext, *Li, *Lnz ;
00310     Int n, pold, pnew, len, k, tail ;
00311 
00312     /* ---------------------------------------------------------------------- */
00313     /* get inputs */
00314     /* ---------------------------------------------------------------------- */
00315 
00316     RETURN_IF_NULL_COMMON (FALSE) ;
00317     RETURN_IF_NULL (L, FALSE) ;
00318     RETURN_IF_XTYPE_INVALID (L, CHOLMOD_REAL, CHOLMOD_ZOMPLEX, FALSE) ;
00319     if (L->is_super)
00320     {
00321   ERROR (CHOLMOD_INVALID, "L must be simplicial") ;
00322   return (FALSE) ;
00323     }
00324     n = L->n ;
00325     if (j >= L->n || need == 0)
00326     {
00327   ERROR (CHOLMOD_INVALID, "j invalid") ;
00328   return (FALSE) ;      /* j out of range */
00329     }
00330     Common->status = CHOLMOD_OK ;
00331 
00332     DEBUG (CHOLMOD(dump_factor) (L, "start colrealloc", Common)) ;
00333 
00334     /* ---------------------------------------------------------------------- */
00335     /* increase the size of L if needed */
00336     /* ---------------------------------------------------------------------- */
00337 
00338     /* head = n+1 ; */
00339     tail = n ;
00340     Lp = L->p ;
00341     Lnz = L->nz ;
00342     Lprev = L->prev ;
00343     Lnext = L->next ;
00344 
00345     ASSERT (Lnz != NULL) ;
00346     ASSERT (Lnext != NULL && Lprev != NULL) ;
00347     PRINT1 (("col %g need %g\n", (double) j, (double) need)) ;
00348 
00349     /* column j cannot have more than n-j entries if all entries are present */
00350     need = MIN (need, n-j) ;
00351 
00352     /* compute need in double to avoid integer overflow */
00353     if (Common->grow1 >= 1.0)
00354     {
00355   xneed = (double) need ;
00356   xneed = Common->grow1 * xneed + Common->grow2 ;
00357   xneed = MIN (xneed, n-j) ;
00358   need = (Int) xneed ;
00359     }
00360     PRINT1 (("really new need %g current %g\n", (double) need,
00361       (double) (Lp [Lnext [j]] - Lp [j]))) ;
00362     ASSERT (need >= 1 && need <= n-j) ;
00363 
00364     if (Lp [Lnext [j]] - Lp [j] >= (Int) need)
00365     {
00366   /* no need to reallocate the column, it's already big enough */
00367   PRINT1 (("colrealloc: quick return %g %g\n",
00368       (double) (Lp [Lnext [j]] - Lp [j]), (double) need)) ;
00369   return (TRUE) ;
00370 
00371     }
00372 
00373     if (Lp [tail] + need > L->nzmax)
00374     {
00375   /* use double to avoid integer overflow */
00376   xneed = (double) need ;
00377   if (Common->grow0 < 1.2)      /* fl. pt. compare, false if NaN */
00378   {
00379       /* if grow0 is less than 1.2 or NaN, don't use it */
00380       xneed = 1.2 * (((double) L->nzmax) + xneed + 1) ;
00381   }
00382   else
00383   {
00384       xneed = Common->grow0 * (((double) L->nzmax) + xneed + 1) ;
00385   }
00386   if (xneed > Size_max ||
00387     !CHOLMOD(reallocate_factor) ((Int) xneed, L, Common))
00388   {
00389       /* out of memory, convert to simplicial symbolic */
00390       CHOLMOD(change_factor) (CHOLMOD_PATTERN, L->is_ll, FALSE, TRUE,
00391         TRUE, L, Common) ;
00392       ERROR (CHOLMOD_OUT_OF_MEMORY, "out of memory; L now symbolic") ;
00393       return (FALSE) ;      /* out of memory */
00394   }
00395   PRINT1 (("\n=== GROW L from %g to %g\n",
00396         (double) L->nzmax, (double) xneed)) ;
00397   /* pack all columns so that each column has at most grow2 free space */
00398   CHOLMOD(pack_factor) (L, Common) ;
00399   ASSERT (Common->status == CHOLMOD_OK) ;
00400   Common->nrealloc_factor++ ;
00401     }
00402 
00403     /* ---------------------------------------------------------------------- */
00404     /* reallocate the column */
00405     /* ---------------------------------------------------------------------- */
00406 
00407     Common->nrealloc_col++ ;
00408 
00409     Li = L->i ;
00410     Lx = L->x ;
00411     Lz = L->z ;
00412 
00413     /* remove j from its current position in the list */
00414     Lnext [Lprev [j]] = Lnext [j] ;
00415     Lprev [Lnext [j]] = Lprev [j] ;
00416 
00417     /* place j at the end of the list */
00418     Lnext [Lprev [tail]] = j ;
00419     Lprev [j] = Lprev [tail] ;
00420     Lnext [j] = n ;
00421     Lprev [tail] = j ;
00422 
00423     /* L is no longer monotonic; columns are out-of-order */
00424     L->is_monotonic = FALSE ;
00425 
00426     /* allocate space for column j */
00427     pold = Lp [j] ;
00428     pnew = Lp [tail] ;
00429     Lp [j] = pnew  ;
00430     Lp [tail] += need ;
00431 
00432     /* copy column j to the new space */
00433     len = Lnz [j] ;
00434     for (k = 0 ; k < len ; k++)
00435     {
00436   Li [pnew + k] = Li [pold + k] ;
00437     }
00438 
00439     if (L->xtype == CHOLMOD_REAL)
00440     {
00441   for (k = 0 ; k < len ; k++)
00442   {
00443       Lx [pnew + k] = Lx [pold + k] ;
00444   }
00445     }
00446     else if (L->xtype == CHOLMOD_COMPLEX)
00447     {
00448   for (k = 0 ; k < len ; k++)
00449   {
00450       Lx [2*(pnew + k)  ] = Lx [2*(pold + k)  ] ;
00451       Lx [2*(pnew + k)+1] = Lx [2*(pold + k)+1] ;
00452   }
00453     }
00454     else if (L->xtype == CHOLMOD_ZOMPLEX)
00455     {
00456   for (k = 0 ; k < len ; k++)
00457   {
00458       Lx [pnew + k] = Lx [pold + k] ;
00459       Lz [pnew + k] = Lz [pold + k] ;
00460   }
00461     }
00462 
00463     DEBUG (CHOLMOD(dump_factor) (L, "colrealloc done", Common)) ;
00464 
00465     /* successful reallocation of column j of L */
00466     return (TRUE) ;
00467 }
00468 
00469 
00470 /* ========================================================================== */
00471 /* === cholmod_pack_factor ================================================== */
00472 /* ========================================================================== */
00473 
00474 /* Pack the columns of a simplicial LDL' or LL' factor.  This can be followed
00475  * by a call to cholmod_reallocate_factor to reduce the size of L to the exact
00476  * size required by the factor, if desired.  Alternatively, you can leave the
00477  * size of L->i and L->x the same, to allow space for future updates/rowadds.
00478  *
00479  * Each column is reduced in size so that it has at most Common->grow2 free
00480  * space at the end of the column.
00481  *
00482  * Does nothing and returns silently if given any other type of factor.
00483  *
00484  * Does NOT force the columns of L to be monotonic.  It thus differs from
00485  * cholmod_change_factor (xtype, -, FALSE, TRUE, TRUE, L, Common), which
00486  * packs the columns and ensures that they appear in monotonic order.
00487  */
00488 
00489 int CHOLMOD(pack_factor)
00490 (
00491     /* ---- in/out --- */
00492     cholmod_factor *L,  /* factor to modify */
00493     /* --------------- */
00494     cholmod_common *Common
00495 )
00496 {
00497     double *Lx, *Lz ;
00498     Int *Lp, *Li, *Lnz, *Lnext ;
00499     Int pnew, j, k, pold, len, n, head, tail, grow2 ;
00500 
00501     /* ---------------------------------------------------------------------- */
00502     /* get inputs */
00503     /* ---------------------------------------------------------------------- */
00504 
00505     RETURN_IF_NULL_COMMON (FALSE) ;
00506     RETURN_IF_NULL (L, FALSE) ;
00507     RETURN_IF_XTYPE_INVALID (L, CHOLMOD_PATTERN, CHOLMOD_ZOMPLEX, FALSE) ;
00508     Common->status = CHOLMOD_OK ;
00509     DEBUG (CHOLMOD(dump_factor) (L, "start pack", Common)) ;
00510     PRINT1 (("PACK factor %d\n", L->is_super)) ;
00511 
00512     if (L->xtype == CHOLMOD_PATTERN || L->is_super)
00513     {
00514   /* nothing to do unless L is simplicial numeric */
00515   return (TRUE) ;
00516     }
00517 
00518     /* ---------------------------------------------------------------------- */
00519     /* pack */
00520     /* ---------------------------------------------------------------------- */
00521 
00522     grow2 = Common->grow2 ;
00523     PRINT1 (("\nPACK grow2 "ID"\n", grow2)) ;
00524 
00525     pnew = 0 ;
00526     n = L->n ;
00527     Lp = L->p ;
00528     Li = L->i ;
00529     Lx = L->x ;
00530     Lz = L->z ;
00531     Lnz = L->nz ;
00532     Lnext = L->next ;
00533 
00534     head = n+1 ;
00535     tail = n ;
00536 
00537     for (j = Lnext [head] ; j != tail ; j = Lnext [j])
00538     {
00539   /* pack column j */
00540   pold = Lp [j] ;
00541   len = Lnz [j] ;
00542   ASSERT (len > 0) ;
00543   PRINT2 (("col "ID" pnew "ID" pold "ID"\n", j, pnew, pold)) ;
00544   if (pnew < pold)
00545   {
00546       PRINT2 (("    pack this column\n")) ;
00547 
00548       for (k = 0 ; k < len ; k++)
00549       {
00550     Li [pnew + k] = Li [pold + k] ;
00551       }
00552 
00553       if (L->xtype == CHOLMOD_REAL)
00554       {
00555     for (k = 0 ; k < len ; k++)
00556     {
00557         Lx [pnew + k] = Lx [pold + k] ;
00558     }
00559       }
00560       else if (L->xtype == CHOLMOD_COMPLEX)
00561       {
00562     for (k = 0 ; k < len ; k++)
00563     {
00564         Lx [2*(pnew + k)  ] = Lx [2*(pold + k)  ] ;
00565         Lx [2*(pnew + k)+1] = Lx [2*(pold + k)+1] ;
00566     }
00567       }
00568       else if (L->xtype == CHOLMOD_ZOMPLEX)
00569       {
00570     for (k = 0 ; k < len ; k++)
00571     {
00572         Lx [pnew + k] = Lx [pold + k] ;
00573         Lz [pnew + k] = Lz [pold + k] ;
00574     }
00575       }
00576 
00577       Lp [j] = pnew ;
00578   }
00579   len = MIN (len + grow2, n - j) ;
00580   pnew = MIN (Lp [j] + len, Lp [Lnext [j]]) ;
00581     }
00582     PRINT2 (("final pnew = "ID"\n", pnew)) ;
00583     return (TRUE) ;
00584 }
00585 
00586 
00587 /* ========================================================================== */
00588 /* === cholmod_factor_to_sparse ============================================= */
00589 /* ========================================================================== */
00590 
00591 /* Constructs a column-oriented sparse matrix containing the pattern and values
00592  * of a simplicial or supernodal numerical factor, and then converts the factor
00593  * into a simplicial symbolic factor.  If L is already packed, monotonic,
00594  * and simplicial (which is the case when cholmod_factorize uses the simplicial
00595  * Cholesky factorization algorithm) then this routine requires only O(1)
00596  * memory and takes O(1) time.
00597  *
00598  * Only operates on numeric factors (real, complex, or zomplex).  Does not
00599  * change the numeric L->xtype (the resulting sparse matrix has the same xtype
00600  * as L).  If this routine fails, L is left unmodified.
00601  */
00602 
00603 cholmod_sparse *CHOLMOD(factor_to_sparse)
00604 (
00605     /* ---- in/out --- */
00606     cholmod_factor *L,  /* factor to copy, converted to symbolic on output */
00607     /* --------------- */
00608     cholmod_common *Common
00609 )
00610 {
00611     cholmod_sparse *Lsparse ;
00612 
00613     /* ---------------------------------------------------------------------- */
00614     /* get inputs */
00615     /* ---------------------------------------------------------------------- */
00616 
00617     RETURN_IF_NULL_COMMON (NULL) ;
00618     RETURN_IF_NULL (L, NULL) ;
00619     RETURN_IF_XTYPE_INVALID (L, CHOLMOD_REAL, CHOLMOD_ZOMPLEX, NULL) ;
00620     Common->status = CHOLMOD_OK ;
00621     DEBUG (CHOLMOD(dump_factor) (L, "start convert to matrix", Common)) ;
00622 
00623     /* ---------------------------------------------------------------------- */
00624     /* convert to packed, monotonic, simplicial, numeric */
00625     /* ---------------------------------------------------------------------- */
00626 
00627     /* leave as LL or LDL' */
00628     if (!CHOLMOD(change_factor) (L->xtype, L->is_ll, FALSE, TRUE, TRUE, L,
00629     Common))
00630     {
00631   ERROR (CHOLMOD_INVALID, "cannot convert L") ;
00632   return (NULL) ;
00633     }
00634 
00635     /* ---------------------------------------------------------------------- */
00636     /* create Lsparse */
00637     /* ---------------------------------------------------------------------- */
00638 
00639     /* allocate the header for Lsparse, the sparse matrix version of L */
00640     Lsparse = CHOLMOD(malloc) (sizeof (cholmod_sparse), 1, Common) ;
00641     if (Common->status < CHOLMOD_OK)
00642     {
00643   return (NULL) ;   /* out of memory */
00644     }
00645 
00646     /* transfer the contents from L to Lsparse */
00647     Lsparse->nrow = L->n ;
00648     Lsparse->ncol = L->n ;
00649     Lsparse->p = L->p ;
00650     Lsparse->i = L->i ;
00651     Lsparse->x = L->x ;
00652     Lsparse->z = L->z ;
00653     Lsparse->nz = NULL ;
00654     Lsparse->stype = 0 ;
00655     Lsparse->itype = L->itype ;
00656     Lsparse->xtype = L->xtype ;
00657     Lsparse->dtype = L->dtype ;
00658     Lsparse->sorted = TRUE ;
00659     Lsparse->packed = TRUE ;
00660     Lsparse->nzmax = L->nzmax ;
00661     ASSERT (CHOLMOD(dump_sparse) (Lsparse, "Lsparse", Common) >= 0) ;
00662 
00663     /* ---------------------------------------------------------------------- */
00664     /* convert L to symbolic, but do not free contents transfered to Lsparse */
00665     /* ---------------------------------------------------------------------- */
00666 
00667     L->p = NULL ;
00668     L->i = NULL ;
00669     L->x = NULL ;
00670     L->z = NULL ;
00671     L->xtype = CHOLMOD_PATTERN ;
00672     CHOLMOD(change_factor) (CHOLMOD_PATTERN, FALSE, FALSE, TRUE, TRUE, L,
00673       Common) ;
00674 
00675     return (Lsparse) ;
00676 }
00677 
00678 
00679 /* ========================================================================== */
00680 /* === cholmod_copy_factor ================================================== */
00681 /* ========================================================================== */
00682 
00683 /* Create an exact copy of a factor, with one exception:
00684  *
00685  * Entries in unused space are not copied (they might not be initialized,
00686  *  and copying them would cause program checkers such as purify and
00687  *  valgrind to complain).
00688  *
00689  * Note that a supernodal L cannot be zomplex.
00690  */
00691 
00692 cholmod_factor *CHOLMOD(copy_factor)
00693 (
00694     /* ---- input ---- */
00695     cholmod_factor *L,  /* factor to copy */
00696     /* --------------- */
00697     cholmod_common *Common
00698 )
00699 {
00700     cholmod_factor *L2 ;
00701     double *Lx, *L2x, *Lz, *L2z ;
00702     Int *Perm, *ColCount, *Lp, *Li, *Lnz, *Lnext, *Lprev, *Lsuper, *Lpi, *Lpx,
00703   *Ls, *Perm2, *ColCount2, *L2p, *L2i, *L2nz, *L2next, *L2prev, *L2super,
00704   *L2pi, *L2px, *L2s ;
00705     Int n, j, p, pend, s, xsize, ssize, nsuper ;
00706 
00707     /* ---------------------------------------------------------------------- */
00708     /* get inputs */
00709     /* ---------------------------------------------------------------------- */
00710 
00711     RETURN_IF_NULL_COMMON (NULL) ;
00712     RETURN_IF_NULL (L, NULL) ;
00713     RETURN_IF_XTYPE_INVALID (L, CHOLMOD_PATTERN, CHOLMOD_ZOMPLEX, NULL) ;
00714     Common->status = CHOLMOD_OK ;
00715     DEBUG (CHOLMOD(dump_factor) (L, "start copy", Common)) ;
00716 
00717     n = L->n ;
00718 
00719     /* ---------------------------------------------------------------------- */
00720     /* allocate a simplicial symbolic factor  */
00721     /* ---------------------------------------------------------------------- */
00722 
00723     /* allocates L2->Perm and L2->ColCount */
00724     L2 = CHOLMOD(allocate_factor) (n, Common) ;
00725     if (Common->status < CHOLMOD_OK)
00726     {
00727   return (NULL) ;     /* out of memory */
00728     }
00729     ASSERT (L2->xtype == CHOLMOD_PATTERN && !(L2->is_super)) ;
00730 
00731     Perm = L->Perm ;
00732     ColCount = L->ColCount ;
00733     Perm2 = L2->Perm ;
00734     ColCount2 = L2->ColCount ;
00735     L2->ordering = L->ordering ;
00736 
00737     for (j = 0 ; j < n ; j++)
00738     {
00739   Perm2 [j] = Perm [j] ;
00740     }
00741     for (j = 0 ; j < n ; j++)
00742     {
00743   ColCount2 [j] = ColCount [j] ;
00744     }
00745     L2->is_ll = L->is_ll ;
00746 
00747     /* ---------------------------------------------------------------------- */
00748     /* copy the rest of the factor */
00749     /* ---------------------------------------------------------------------- */
00750 
00751     if (L->xtype != CHOLMOD_PATTERN && !(L->super))
00752     {
00753 
00754   /* ------------------------------------------------------------------ */
00755   /* allocate a simplicial numeric factor */
00756   /* ------------------------------------------------------------------ */
00757 
00758   /* allocate L2->p, L2->nz, L2->prev, L2->next, L2->i, and L2->x.
00759    * packed = -1 so that cholmod_change_factor allocates space of
00760    * size L2->nzmax */
00761   L2->nzmax = L->nzmax ;
00762   if (!CHOLMOD(change_factor) (L->xtype, L->is_ll, FALSE, -1, TRUE,
00763         L2, Common))
00764   {
00765       CHOLMOD(free_factor) (&L2, Common) ;
00766       return (NULL) ; /* out of memory */
00767   }
00768   ASSERT (L->nzmax == L2->nzmax) ;
00769 
00770   /* ------------------------------------------------------------------ */
00771   /* copy the contents of a simplicial numeric factor */
00772   /* ------------------------------------------------------------------ */
00773 
00774   Lp = L->p ;
00775   Li = L->i ;
00776   Lx = L->x ;
00777   Lz = L->z ;
00778   Lnz = L->nz ;
00779   Lnext = L->next ;
00780   Lprev = L->prev ;
00781 
00782   L2p = L2->p ;
00783   L2i = L2->i ;
00784   L2x = L2->x ;
00785   L2z = L2->z ;
00786   L2nz = L2->nz ;
00787   L2next = L2->next ;
00788   L2prev = L2->prev ;
00789   L2->xtype = L->xtype ;
00790   L2->dtype = L->dtype ;
00791 
00792   for (j = 0 ; j <= n ; j++)
00793   {
00794       L2p [j] = Lp [j] ;
00795   }
00796 
00797   for (j = 0 ; j < n+2 ; j++)
00798   {
00799       L2prev [j] = Lprev [j] ;
00800   }
00801 
00802   for (j = 0 ; j < n+2 ; j++)
00803   {
00804       L2next [j] = Lnext [j] ;
00805   }
00806 
00807   for (j = 0 ; j < n ; j++)
00808   {
00809       L2nz [j] = Lnz [j] ;
00810   }
00811 
00812   for (j = 0 ; j < n ; j++)
00813   {
00814       p = Lp [j] ;
00815       pend = p + Lnz [j] ;
00816       for ( ; p < pend ; p++)
00817       {
00818     L2i [p] = Li [p] ;
00819       }
00820       p = Lp [j] ;
00821 
00822       if (L->xtype == CHOLMOD_REAL)
00823       {
00824     for ( ; p < pend ; p++)
00825     {
00826         L2x [p] = Lx [p] ;
00827     }
00828       }
00829       else if (L->xtype == CHOLMOD_COMPLEX)
00830       {
00831     for ( ; p < pend ; p++)
00832     {
00833         L2x [2*p  ] = Lx [2*p  ] ;
00834         L2x [2*p+1] = Lx [2*p+1] ;
00835     }
00836       }
00837       else if (L->xtype == CHOLMOD_ZOMPLEX)
00838       {
00839     for ( ; p < pend ; p++)
00840     {
00841         L2x [p] = Lx [p] ;
00842         L2z [p] = Lz [p] ;
00843     }
00844       }
00845 
00846   }
00847 
00848     }
00849     else if (L->is_super)
00850     {
00851 
00852   /* ------------------------------------------------------------------ */
00853   /* copy a supernodal factor */
00854   /* ------------------------------------------------------------------ */
00855 
00856   xsize = L->xsize ;
00857   ssize = L->ssize ;
00858   nsuper = L->nsuper ;
00859 
00860   L2->xsize = xsize ;
00861   L2->ssize = ssize ;
00862   L2->nsuper = nsuper ;
00863 
00864   /* allocate L2->super, L2->pi, L2->px, and L2->s.  Allocate L2->x if
00865    * L is numeric */
00866   if (!CHOLMOD(change_factor) (L->xtype, TRUE, TRUE, TRUE, TRUE, L2,
00867         Common))
00868   {
00869       CHOLMOD(free_factor) (&L2, Common) ;
00870       return (NULL) ; /* out of memory */
00871   }
00872 
00873   ASSERT (L2->s != NULL) ;
00874 
00875   /* ------------------------------------------------------------------ */
00876   /* copy the contents of a supernodal factor */
00877   /* ------------------------------------------------------------------ */
00878 
00879   Lsuper = L->super ;
00880   Lpi = L->pi ;
00881   Lpx = L->px ;
00882   Ls = L->s ;
00883   Lx = L->x ;
00884 
00885   L2super = L2->super ;
00886   L2pi = L2->pi ;
00887   L2px = L2->px ;
00888   L2s = L2->s ;
00889   L2x = L2->x ;
00890 
00891   L2->maxcsize = L->maxcsize ;
00892   L2->maxesize = L->maxesize ;
00893 
00894   for (s = 0 ; s <= nsuper ; s++)
00895   {
00896       L2super [s] = Lsuper [s] ;
00897   }
00898   for (s = 0 ; s <= nsuper ; s++)
00899   {
00900       L2pi [s] = Lpi [s] ;
00901   }
00902   for (s = 0 ; s <= nsuper ; s++)
00903   {
00904       L2px [s] = Lpx [s] ;
00905   }
00906 
00907   L2s [0] = 0 ;
00908   for (p = 0 ; p < ssize ; p++)
00909   {
00910       L2s [p] = Ls [p] ;
00911   }
00912 
00913   if (L->xtype == CHOLMOD_REAL)
00914   {
00915       for (p = 0 ; p < xsize ; p++)
00916       {
00917     L2x [p] = Lx [p] ;
00918       }
00919   }
00920   else if (L->xtype == CHOLMOD_COMPLEX)
00921   {
00922       for (p = 0 ; p < 2*xsize ; p++)
00923       {
00924     L2x [p] = Lx [p] ;
00925       }
00926   }
00927     }
00928 
00929     L2->minor = L->minor ;
00930     L2->is_monotonic = L->is_monotonic ;
00931 
00932     DEBUG (CHOLMOD(dump_factor) (L2, "L2 got copied", Common)) ;
00933     ASSERT (L2->xtype == L->xtype && L2->is_super == L->is_super) ;
00934     return (L2) ;
00935 }
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines