kusano 7d535a
kusano 7d535a
<meta content="text/html;charset=UTF-8" http-equiv="Content-Type">
kusano 7d535a
<title>SuperLU: SRC/dgstrf.c File Reference</title>
kusano 7d535a
<link href="doxygen.css" rel="stylesheet" type="text/css">
kusano 7d535a
<link href="tabs.css" rel="stylesheet" type="text/css">
kusano 7d535a
kusano 7d535a
kusano 7d535a
kusano 7d535a
  
kusano 7d535a
    
    kusano 7d535a
          
  • Main Page
  • kusano 7d535a
          
  • Data Structures
  • kusano 7d535a
          
  • Files
  • kusano 7d535a
        
    kusano 7d535a
      
    kusano 7d535a
    kusano 7d535a
    kusano 7d535a

    SRC/dgstrf.c File Reference

    Computes an LU factorization of a general sparse matrix. More...
    kusano 7d535a

    kusano 7d535a
    #include "slu_ddefs.h"
    kusano 7d535a
    kusano 7d535a
    kusano 7d535a

    Functions

    kusano 7d535a
    void dgstrf (superlu_options_t *options, SuperMatrix *A, int relax, int panel_size, int *etree, void *work, int lwork, int *perm_c, int *perm_r, SuperMatrix *L, SuperMatrix *U, SuperLUStat_t *stat, int *info)
    kusano 7d535a
    kusano 7d535a
    kusano 7d535a

    Detailed Description

    kusano 7d535a
    kusano 7d535a
     -- SuperLU routine (version 3.0) --
    kusano 7d535a
     Univ. of California Berkeley, Xerox Palo Alto Research Center,
    kusano 7d535a
     and Lawrence Berkeley National Lab.
    kusano 7d535a
     October 15, 2003

    kusano 7d535a
     Copyright (c) 1994 by Xerox Corporation.  All rights reserved.

    kusano 7d535a
     THIS MATERIAL IS PROVIDED AS IS, WITH ABSOLUTELY NO WARRANTY
    kusano 7d535a
     EXPRESSED OR IMPLIED.  ANY USE IS AT YOUR OWN RISK.

    kusano 7d535a
     Permission is hereby granted to use or copy this program for any
    kusano 7d535a
     purpose, provided the above notices are retained on all copies.
    kusano 7d535a
     Permission to modify the code and to distribute modified code is
    kusano 7d535a
     granted, provided the above notices are retained, and a notice that
    kusano 7d535a
     the code was modified is included with the above copyright notice.
    kusano 7d535a
      

    Function Documentation

    kusano 7d535a
    kusano 7d535a
    kusano 7d535a
    kusano 7d535a
          
    kusano 7d535a
            
    kusano 7d535a
              void dgstrf           
    kusano 7d535a
              (
    kusano 7d535a
              superlu_options_t
    kusano 7d535a
               options, 
    kusano 7d535a
            
    kusano 7d535a
            
    kusano 7d535a
              
    kusano 7d535a
              
    kusano 7d535a
              SuperMatrix
    kusano 7d535a
               A, 
    kusano 7d535a
            
    kusano 7d535a
            
    kusano 7d535a
              
    kusano 7d535a
              
    kusano 7d535a
              int 
    kusano 7d535a
               relax, 
    kusano 7d535a
            
    kusano 7d535a
            
    kusano 7d535a
              
    kusano 7d535a
              
    kusano 7d535a
              int 
    kusano 7d535a
               panel_size, 
    kusano 7d535a
            
    kusano 7d535a
            
    kusano 7d535a
              
    kusano 7d535a
              
    kusano 7d535a
              int * 
    kusano 7d535a
               etree, 
    kusano 7d535a
            
    kusano 7d535a
            
    kusano 7d535a
              
    kusano 7d535a
              
    kusano 7d535a
              void * 
    kusano 7d535a
               work, 
    kusano 7d535a
            
    kusano 7d535a
            
    kusano 7d535a
              
    kusano 7d535a
              
    kusano 7d535a
              int 
    kusano 7d535a
               lwork, 
    kusano 7d535a
            
    kusano 7d535a
            
    kusano 7d535a
              
    kusano 7d535a
              
    kusano 7d535a
              int * 
    kusano 7d535a
               perm_c, 
    kusano 7d535a
            
    kusano 7d535a
            
    kusano 7d535a
              
    kusano 7d535a
              
    kusano 7d535a
              int * 
    kusano 7d535a
               perm_r, 
    kusano 7d535a
            
    kusano 7d535a
            
    kusano 7d535a
              
    kusano 7d535a
              
    kusano 7d535a
              SuperMatrix
    kusano 7d535a
               L, 
    kusano 7d535a
            
    kusano 7d535a
            
    kusano 7d535a
              
    kusano 7d535a
              
    kusano 7d535a
              SuperMatrix
    kusano 7d535a
               U, 
    kusano 7d535a
            
    kusano 7d535a
            
    kusano 7d535a
              
    kusano 7d535a
              
    kusano 7d535a
              SuperLUStat_t
    kusano 7d535a
               stat, 
    kusano 7d535a
            
    kusano 7d535a
            
    kusano 7d535a
              
    kusano 7d535a
              
    kusano 7d535a
              int * 
    kusano 7d535a
               info 
    kusano 7d535a
            
    kusano 7d535a
            
    kusano 7d535a
              
    kusano 7d535a
              )
    kusano 7d535a
              
    kusano 7d535a
            
    kusano 7d535a
          
    kusano 7d535a
    kusano 7d535a
    kusano 7d535a
    kusano 7d535a

    kusano 7d535a
    kusano 7d535a
     Purpose
    kusano 7d535a
     =======

    kusano 7d535a
     DGSTRF computes an LU factorization of a general sparse m-by-n
    kusano 7d535a
     matrix A using partial pivoting with row interchanges.
    kusano 7d535a
     The factorization has the form
    kusano 7d535a
         Pr * A = L * U
    kusano 7d535a
     where Pr is a row permutation matrix, L is lower triangular with unit
    kusano 7d535a
     diagonal elements (lower trapezoidal if A->nrow > A->ncol), and U is upper 
    kusano 7d535a
     triangular (upper trapezoidal if A->nrow < A->ncol).

    kusano 7d535a
     See supermatrix.h for the definition of 'SuperMatrix' structure.

    kusano 7d535a
     Arguments
    kusano 7d535a
     =========

    kusano 7d535a
     options (input) superlu_options_t*
    kusano 7d535a
             The structure defines the input parameters to control
    kusano 7d535a
             how the LU decomposition will be performed.

    kusano 7d535a
     A        (input) SuperMatrix*
    kusano 7d535a
    	    Original matrix A, permuted by columns, of dimension
    kusano 7d535a
              (A->nrow, A->ncol). The type of A can be:
    kusano 7d535a
              Stype = SLU_NCP; Dtype = SLU_D; Mtype = SLU_GE.

    kusano 7d535a
     relax    (input) int
    kusano 7d535a
              To control degree of relaxing supernodes. If the number
    kusano 7d535a
              of nodes (columns) in a subtree of the elimination tree is less
    kusano 7d535a
              than relax, this subtree is considered as one supernode,
    kusano 7d535a
              regardless of the row structures of those columns.

    kusano 7d535a
     panel_size (input) int
    kusano 7d535a
              A panel consists of at most panel_size consecutive columns.

    kusano 7d535a
     etree    (input) int*, dimension (A->ncol)
    kusano 7d535a
              Elimination tree of A'*A.
    kusano 7d535a
              Note: etree is a vector of parent pointers for a forest whose
    kusano 7d535a
              vertices are the integers 0 to A->ncol-1; etree[root]==A->ncol.
    kusano 7d535a
              On input, the columns of A should be permuted so that the
    kusano 7d535a
              etree is in a certain postorder.

    kusano 7d535a
     work     (input/output) void*, size (lwork) (in bytes)
    kusano 7d535a
              User-supplied work space and space for the output data structures.
    kusano 7d535a
              Not referenced if lwork = 0;

    kusano 7d535a
     lwork   (input) int
    kusano 7d535a
             Specifies the size of work array in bytes.
    kusano 7d535a
             = 0:  allocate space internally by system malloc;
    kusano 7d535a
             > 0:  use user-supplied work array of length lwork in bytes,
    kusano 7d535a
                   returns error if space runs out.
    kusano 7d535a
             = -1: the routine guesses the amount of space needed without
    kusano 7d535a
                   performing the factorization, and returns it in
    kusano 7d535a
                   *info; no other side effects.

    kusano 7d535a
     perm_c   (input) int*, dimension (A->ncol)
    kusano 7d535a
    	    Column permutation vector, which defines the 
    kusano 7d535a
              permutation matrix Pc; perm_c[i] = j means column i of A is 
    kusano 7d535a
              in position j in A*Pc.
    kusano 7d535a
              When searching for diagonal, perm_c[*] is applied to the
    kusano 7d535a
              row subscripts of A, so that diagonal threshold pivoting
    kusano 7d535a
              can find the diagonal of A, rather than that of A*Pc.

    kusano 7d535a
     perm_r   (input/output) int*, dimension (A->nrow)
    kusano 7d535a
              Row permutation vector which defines the permutation matrix Pr,
    kusano 7d535a
              perm_r[i] = j means row i of A is in position j in Pr*A.
    kusano 7d535a
              If options->Fact = SamePattern_SameRowPerm, the pivoting routine
    kusano 7d535a
                 will try to use the input perm_r, unless a certain threshold
    kusano 7d535a
                 criterion is violated. In that case, perm_r is overwritten by
    kusano 7d535a
                 a new permutation determined by partial pivoting or diagonal
    kusano 7d535a
                 threshold pivoting.
    kusano 7d535a
              Otherwise, perm_r is output argument;

    kusano 7d535a
     L        (output) SuperMatrix*
    kusano 7d535a
              The factor L from the factorization Pr*A=L*U; use compressed row 
    kusano 7d535a
              subscripts storage for supernodes, i.e., L has type: 
    kusano 7d535a
              Stype = SLU_SC, Dtype = SLU_D, Mtype = SLU_TRLU.

    kusano 7d535a
     U        (output) SuperMatrix*
    kusano 7d535a
    	    The factor U from the factorization Pr*A*Pc=L*U. Use column-wise
    kusano 7d535a
              storage scheme, i.e., U has types: Stype = SLU_NC, 
    kusano 7d535a
              Dtype = SLU_D, Mtype = SLU_TRU.

    kusano 7d535a
     stat     (output) SuperLUStat_t*
    kusano 7d535a
              Record the statistics on runtime and floating-point operation count.
    kusano 7d535a
              See slu_util.h for the definition of 'SuperLUStat_t'.

    kusano 7d535a
     info     (output) int*
    kusano 7d535a
              = 0: successful exit
    kusano 7d535a
              < 0: if info = -i, the i-th argument had an illegal value
    kusano 7d535a
              > 0: if info = i, and i is
    kusano 7d535a
                 <= A->ncol: U(i,i) is exactly zero. The factorization has
    kusano 7d535a
                    been completed, but the factor U is exactly singular,
    kusano 7d535a
                    and division by zero will occur if it is used to solve a
    kusano 7d535a
                    system of equations.
    kusano 7d535a
                 > A->ncol: number of bytes allocated when memory allocation
    kusano 7d535a
                    failure occurred, plus A->ncol. If lwork = -1, it is
    kusano 7d535a
                    the estimated amount of space needed, plus A->ncol.

    kusano 7d535a
     ======================================================================

    kusano 7d535a
     Local Working Arrays: 
    kusano 7d535a
     ======================
    kusano 7d535a
       m = number of rows in the matrix
    kusano 7d535a
       n = number of columns in the matrix

    kusano 7d535a
       xprune[0:n-1]: xprune[*] points to locations in subscript 
    kusano 7d535a
    	vector lsub[*]. For column i, xprune[i] denotes the point where 
    kusano 7d535a
    	structural pruning begins. I.e. only xlsub[i],..,xprune[i]-1 need 
    kusano 7d535a
    	to be traversed for symbolic factorization.

    kusano 7d535a
       marker[0:3*m-1]: marker[i] = j means that node i has been 
    kusano 7d535a
    	reached when working on column j.
    kusano 7d535a
    	Storage: relative to original row subscripts
    kusano 7d535a
    	NOTE: There are 3 of them: marker/marker1 are used for panel dfs, 
    kusano 7d535a
    	      see dpanel_dfs.c; marker2 is used for inner-factorization,
    kusano 7d535a
                see dcolumn_dfs.c.

    kusano 7d535a
       parent[0:m-1]: parent vector used during dfs
    kusano 7d535a
          Storage: relative to new row subscripts

    kusano 7d535a
       xplore[0:m-1]: xplore[i] gives the location of the next (dfs) 
    kusano 7d535a
    	unexplored neighbor of i in lsub[*]

    kusano 7d535a
       segrep[0:nseg-1]: contains the list of supernodal representatives
    kusano 7d535a
    	in topological order of the dfs. A supernode representative is the 
    kusano 7d535a
    	last column of a supernode.
    kusano 7d535a
          The maximum size of segrep[] is n.

    kusano 7d535a
       repfnz[0:W*m-1]: for a nonzero segment U[*,j] that ends at a 
    kusano 7d535a
    	supernodal representative r, repfnz[r] is the location of the first 
    kusano 7d535a
    	nonzero in this segment.  It is also used during the dfs: repfnz[r]>0
    kusano 7d535a
    	indicates the supernode r has been explored.
    kusano 7d535a
    	NOTE: There are W of them, each used for one column of a panel.

    kusano 7d535a
       panel_lsub[0:W*m-1]: temporary for the nonzeros row indices below 
    kusano 7d535a
          the panel diagonal. These are filled in during dpanel_dfs(), and are
    kusano 7d535a
          used later in the inner LU factorization within the panel.
    kusano 7d535a
    	panel_lsub[]/dense[] pair forms the SPA data structure.
    kusano 7d535a
    	NOTE: There are W of them.

    kusano 7d535a
       dense[0:W*m-1]: sparse accumulating (SPA) vector for intermediate values;
    kusano 7d535a
    	    	   NOTE: there are W of them.

    kusano 7d535a
       tempv[0:*]: real temporary used for dense numeric kernels;
    kusano 7d535a
    	The size of this array is defined by NUM_TEMPV() in slu_ddefs.h.
    kusano 7d535a
      
    kusano 7d535a
    kusano 7d535a

    kusano 7d535a
    kusano 7d535a

    <address style="text-align: right;"><small>Generated on Mon Nov 22 10:23:47 2010 for SuperLU by </small></address>
    kusano 7d535a
    kusano 7d535a
    doxygen 1.5.5 
    kusano 7d535a
    kusano 7d535a