home *** CD-ROM | disk | FTP | other *** search
/ SGI Developer Toolbox 6.1 / SGI Developer Toolbox 6.1 - Disc 4.iso / lib / mathlib / skyline / README
Encoding:
Text File  |  1994-08-02  |  14.2 KB  |  384 lines
  1.                                  SKYLIB
  2.  
  3. SKYLIB is a library of Fortran 77 routines for solving linear systems
  4. of equations where the matrix is banded, with the bands varying in
  5. width, stored in a packed form known as skyline. These routines handle
  6. two types of matrices: symmetric (LDL**T factorization) and symmetric
  7. in shape, but not in values (LU factorization). Further details can
  8. be found below.
  9.  
  10. The current library comes in four formats:
  11.  
  12. o mips1: compilation flags -O2,
  13.         will run on all systems,
  14.         This is NOT a parallel version.
  15.  
  16. o mips2: compiled with the "-O2 -mips2" flags,
  17.         will run only on systems with R4000 CPUs,
  18.         This is NOT a parallel version.
  19.  
  20. o mp_mips1: compiled "-O2 -mp"
  21.         will run all systems.
  22.         Requires Fortran compilers since it uses the Fortran MP library.
  23.         This version will run in parallel on Power Series systems.
  24.  
  25. o mp_mips2: compiled "-O2 -mips2 -mp"
  26.         will run all Challenge and Onyx systems only.
  27.         Requires Fortran compilers since it uses the Fortran MP library.
  28.         This version will run in parallel.
  29.  
  30.  
  31. The example directory contains example routines that use both solvers. 
  32. The makefile in this directory will create single and double precision
  33. test programs for both solvers. You have to modify the Makefile to point 
  34. to the appropriate skylib for your machine
  35.  
  36.  
  37.  
  38.                          LDL SKYLINE SOLVER
  39.                  ( for symmetric matrices )
  40.  
  41. The routines for the skyline solver operate on matrices that are 
  42. banded and symmetric in form and in values. The operations performed 
  43. are the triangular factorization LDL(t) WITHOUT pivoting of 
  44. matrices and the solution of simultaneous linear equations.
  45.  
  46. The system parameter routine is in: dparam.f
  47.  
  48. The factorization routines are in the files: dskydc.f
  49.                                              dspslv.f
  50.                                              dmpslv.f
  51.                                              dskyul.f
  52. The solve routines are in the file: djesol.f
  53.                                     drjsol.f
  54.  
  55.  
  56. ( Note: The routines should be linked to the blas library by including
  57.   ----  -lblas or -lblas_mp and the appropriate link flags -mips2 -mp
  58.   on the link command.  Also, the system parameters need to be initialized 
  59.   before the factorization routines are called. )
  60.  
  61.  
  62.  
  63. DSKYDC computes a LDL(t) factorization without pivoting of a
  64.        real*8 matrix V stored in either Jenning's or reverse
  65.        Jenning's storage scheme.
  66.  
  67.  
  68.    Arguments
  69.    ---------
  70.  
  71.       ON ENTRY
  72.  
  73.          V       DOUBLE PRECISION ( PD(N) )
  74.                  The matrix to be factored stored as a 1-dim
  75.                  vector. The storage scheme is described in detail
  76.                  below.
  77.  
  78.          N       INTEGER
  79.                  The order of the matrix (number of columns).
  80.  
  81.          PD      INTEGER (N), or INTEGER (N+1)
  82.                  An integer vector pointing to the position of the
  83.                  diagonal elements in V.
  84.                  V( PD(1) ) will have the first diagonal element,
  85.                  V( PD(2) ) will have the second diagonal element,
  86.                  .
  87.                  .
  88.                  V( PD(N) ) will have the last diagonal element,
  89.  
  90.                  For the Jenning's method,
  91.                  PD(N) will also have the number of elements in V.
  92.                  For the reverse Jenning's method,
  93.                  PD(N+1) will have the number of elements in V plus 1.
  94.  
  95.          ITYP    INTEGER
  96.                  The type of storage schemes
  97.             .EQ.1 == Jennings 
  98.                     .EQ.2 == reverse Jenning's (NOT AVAILABLE YET)
  99.  
  100.  
  101.       ON RETURN
  102.  
  103.          V       DOUBLE PRECISION ( PD(N) )
  104.                  The resultant L matrix as a 1-dim vector, except that
  105.                  the V(PD(i)) term corresponds to the i_th term in D(-1).
  106.  
  107.          INFO    INTEGER
  108.                  = 0  Successful exit.
  109.                  = K  If  V(pd(k)) .EQ. 0.0 .  The factorization has
  110.                       not been completed due to this singularity.
  111.  
  112.  
  113. DRJSOL  solves the double precision system
  114.         V * X = B using the factorized matrix
  115.         stored in reverse Jenning's format
  116.         (This routine has NOT been implemented yet)
  117.  
  118.  
  119. DJESOL  solves the double precision system
  120.         V * X = B  using the factorized matrix 
  121.         stored in Jenning's format
  122.  
  123.  
  124.  
  125.    Arguments
  126.    ---------
  127.  
  128.       ON ENTRY
  129.  
  130.          V       DOUBLE PRECISION ( PD(N) )
  131.                  The output from DSKYDC stored in Jenning's skyline 
  132.                  form. The storage scheme is described in DSKYDC.
  133.  
  134.          N       INTEGER
  135.                  The order of the matrix (number of rows/columns).
  136.  
  137.  
  138.          PD      INTEGER (N)
  139.                  An integer vector pointing to the position of the
  140.                  diagonal elements in V.
  141.                  V( PD(1) ) will have the first diagonal element,
  142.                  V( PD(2) ) will have the second diagonal element,
  143.                  .
  144.                  .
  145.                  V( PD(N) ) will have the last diagonal element,
  146.                  PD(N) will also have the number of elements in V.
  147.  
  148.  
  149.          B       DOUBLE PRECISION ( N )
  150.                  The right hand side vector.
  151.  
  152.  
  153.       ON RETURN
  154.  
  155.          B       THE SOLUTION VECTOR  X .
  156.  
  157.  
  158.  
  159.  
  160.    Jenning's Storage Method
  161.    ------------------------
  162.  
  163.       The variable band or skyline storage scheme is illustrated
  164.       by the following example, where:
  165.  
  166.          N = 7
  167.          PD = ( 1, 3, 6, 9, 13, 16, 21 )
  168.          V = ( 11, 12, 22, 13, 23, 33, 24, 0, 44, 25,  
  169.                35, 0, 55, 46, 56, 66, 37, 0, 57, 0, 77 )
  170.  
  171.  
  172.   The elements of V                      correponding to the elements of A
  173.  
  174.   ------------------------------------   ------------------------------------
  175.   |  1 |  2 |  4 |  * |  * |  * |  * |   | 11 | 12 | 13 |  0 |  0 |  0 |  0 |
  176.   |----|    |    |    |    |    |    |   |----|    |    |    |    |    |    |
  177.        |  3 |  5 |  7 | 10 |  * |  * |   | 12 | 22 | 23 | 24 | 25 |  0 |  0 |
  178.        -----|    |    |    |    |    |   |---------|    |    |    |    |    |
  179.             |  6 |  8 | 11 |  * | 17 |   | 13 | 23 | 33 |  0 | 35 |  0 | 37 |
  180.             -----|    |    |    |    |   |--------------|    |    |    |    |
  181.                  |  9 | 12 | 14 | 18 |   |  0   24    0 | 44 |  0 | 46 |  0 |
  182.                  -----|    |    |    |   |-------------------|    |    |    |
  183.                       | 13 | 15 | 19 |   |  0   25   35    0 | 55 | 56 | 57 |
  184.                       -----|    |    |   |------------------------|    |    |
  185.                            | 16 | 20 |   |  0    0    0   46   56 | 66 |  0 |
  186.                            -----|    |   |-----------------------------|    |
  187.                                 | 21 |   |  0    0   37    0   57    0 | 77 |
  188.                                 ------   ------------------------------------
  189.  
  190.       Due to symmetry, only the upper triangular part of the array needs
  191.       to be stored.  The columns are stored in a top-down fashion, i.e. 
  192.       from the first nonzero entry to the diagonal, into V.  The zero
  193.       entries outside the skyline profile are not stored.  The book-keeping
  194.       array PD keeps track of the locations of the diagonal elements in V.
  195.       
  196.  
  197.    Reverse Jenning's Method
  198.    ------------------------
  199.  
  200.       The variable band or skyline storage scheme is illustrated
  201.       by the following example, where:
  202.  
  203.          N = 7
  204.          PD = ( 1, 2, 4, 7, 10, 14, 17, 22 )
  205.          V = ( 11, 22, 12, 33, 23, 13, 44, 0, 24, 55,  
  206.                0, 35, 25, 66, 56, 46, 77, 0, 57, 0, 37 )
  207.  
  208.  
  209.   The elements of V                      correponding to the elements of A
  210.  
  211.   ------------------------------------   ------------------------------------
  212.   |  1 |  3 |  6 |  * |  * |  * |  * |   | 11 | 12 | 13 |  0 |  0 |  0 |  0 |
  213.   |----|    |    |    |    |    |    |   |----|    |    |    |    |    |    |
  214.        |  2 |  5 |  9 | 13 |  * |  * |   | 12 | 22 | 23 | 24 | 25 |  0 |  0 |
  215.        -----|    |    |    |    |    |   |---------|    |    |    |    |    |
  216.             |  4 |  8 | 12 |  * | 21 |   | 13 | 23 | 33 |  0 | 35 |  0 | 37 |
  217.             -----|    |    |    |    |   |--------------|    |    |    |    |
  218.                  |  7 | 11 | 16 | 20 |   |  0   24    0 | 44 |  0 | 46 |  0 |
  219.                  -----|    |    |    |   |-------------------|    |    |    |
  220.                       | 10 | 15 | 19 |   |  0   25   35    0 | 55 | 56 | 57 |
  221.                       -----|    |    |   |------------------------|    |    |
  222.                            | 14 | 18 |   |  0    0    0   46   56 | 66 |  0 |
  223.                            -----|    |   |-----------------------------|    |
  224.                                 | 17 |   |  0    0   37    0   57    0 | 77 |
  225.                                 ------   ------------------------------------
  226.  
  227.       Due to symmetry, only the upper triangular part of the array needs
  228.       to be stored.  The columns are stored in a bottom-up fashion, i.e. 
  229.       from the diagonal to the first nonzero entry, into V.  The zero
  230.       entries outside the skyline profile are not stored.  The book-keeping
  231.       array PD keeps track of the locations of the diagonal elements in V.
  232.       Note that there are N+1 terms in PD.  The last term provides the
  233.       information of the total number of terms in V plus 1.
  234.       
  235.  
  236.  
  237.  
  238.  
  239.                            LU SKYLINE SOLVER
  240.            (for symmetrical shaped matrices)
  241.  
  242. The routines for the skyline solver operate on matrices that are 
  243. variable band and symmetric in form but not in values. The operations
  244. performed are the triangular factorization (LU) WITHOUT pivoting of 
  245. matrices and the solution of simultaneous linear equations.
  246.  
  247. The factorization routines are in the files: dskyfa.f
  248.                                              dskyf1.f
  249.                                              dskyf2.f
  250. The solve routines are in the file: dskysl.f
  251.                                     dskyf1.f
  252.  
  253.  
  254. ( Note: The routines should be linked to the blas library by including
  255.   ----  -lblas on the link command. )
  256.  
  257.  
  258.  
  259. DSKYFA computes an LU factorization without pivoting of a
  260.        real n-by-n variable band (skyline) matrix V stored in a
  261.        fishbone packed form.
  262.  
  263.     Arguments
  264.     ---------
  265.  
  266.      ON ENTRY
  267.  
  268.         V       DOUBLE PRECISION ( PD(N) )
  269.                 The matrix to be factored stored as a 1-dim
  270.                 vector. The storage scheme is described in detail
  271.                 below.
  272.  
  273.         N       INTEGER
  274.                 The order of the matrix (number of columns).
  275.  
  276.         PD      INTEGER (N)
  277.                 An integer vector pointing to the position of the
  278.                 diagonal elements in V.
  279.                 V( PD(1) ) will have the first diagonal element,
  280.                 V( PD(2) ) will have the second diagonal element,
  281.                 .
  282.                 .
  283.                 V( PD(N) ) will have the last diagonal element,
  284.                 PD(N) will also have the number of elements in V.
  285.  
  286.  
  287.      ON RETURN
  288.  
  289.         V       An upper triangular matrix stored as a 1-dim vector
  290.                 and the multipliers which were used to obtain it.
  291.  
  292.         INFO    INTEGER
  293.                 = 0  Successful exit.
  294.                 = K  If  V(pd(k)) .EQ. 0.0 .  The factorization has
  295.                      been completed, but the factor U is exactly
  296.                      singular. Division by zero will occur if it is
  297.                      used to solve a system of equations in DSKYSL.
  298.  
  299.  
  300.  
  301. DGESL solves the double precision system
  302.       V * X = B  OR  TRANS(V) * X = B
  303.       using the factors computed by DSKYFA.
  304.  
  305.    Arguments
  306.    ---------
  307.  
  308.      ON ENTRY
  309.  
  310.         V       DOUBLE PRECISION ( PD(N) )
  311.                 The output from DSKYSL stored in skyline
  312.                 packed form. The storage scheme is described
  313.                 in more detail below.
  314.  
  315.         N       INTEGER
  316.                 The order of the matrix (number of rows/columns).
  317.  
  318.  
  319.         PD      INTEGER (N)
  320.                 An integer vector pointing to the position of the
  321.                 diagonal elements in V.
  322.                 V( PD(1) ) will have the first diagonal element,
  323.                 V( PD(2) ) will have the second diagonal element,
  324.                 .
  325.                 .
  326.                 V( PD(N) ) will have the last diagonal element,
  327.                 PD(N) will also have the number of elements in V.
  328.  
  329.  
  330.         B       DOUBLE PRECISION ( N )
  331.                 The right hand side vector.
  332.  
  333.         TRANS   CHARACTER*1
  334.                 The operation applied to V.
  335.                 = 'N':  Solve  V * X = B  (No transpose)
  336.                 = 'T':  Solve  V'* X = B  (Transpose)
  337.  
  338.  
  339.      ON RETURN
  340.  
  341.         B       THE SOLUTION VECTOR  X .
  342.  
  343.         INFO    INTEGER
  344.                 = 0:  successful exit
  345.                 < 0: if INFO = -k, the k-th argument had an illegal value
  346.  
  347.  
  348. Storage
  349. -------
  350. The variable band or skyline storage scheme is illustrated
  351. by the following example, where:
  352.  
  353.         N = 7
  354.         PD = ( 1, 4, 9, 14, 21, 26, 35 )
  355.  
  356.  
  357.  The elements of V                      correponding to the elements of A
  358.  
  359.  ------------------------------------   ------------------------------------
  360.  |  1 |  3 |  7 |  * |  * |  * |  * |   | 11 | 12 | 13 |  * |  * |  * |  * |
  361.  |----|    |    |    |    |    |    |   |----|    |    |    |    |    |    |
  362.  |  2 |  4 |  8 | 12 | 18 |  * |  * |   | 21 | 22 | 23 | 24 | 25 |  * |  * |
  363.  |---------|    |    |    |    |    |   |---------|    |    |    |    |    |
  364.  |  5    6 |  9 | 13 | 19 |  * | 31 |   | 31 | 32 | 33 | 34 | 35 |  * | 37 |
  365.  |--------------|    |    |    |    |   |--------------|    |    |    |    |
  366.  |  *   10   11 | 14 | 20 | 24 | 32 |   |  *   42   43 | 44 | 45 | 46 | 47 |
  367.  |-------------------|    |    |    |   |-------------------|    |    |    |
  368.  |  *   15   16   17 | 21 | 25 | 33 |   |  *   52   53   54 | 55 | 56 | 57 |
  369.  |------------------------|    |    |   |------------------------|    |    |
  370.  |  *    *    *   22   23 | 26 | 34 |   |  *    *    *   64   65 | 66 | 67 |
  371.  |-----------------------------|    |   |-----------------------------|    |
  372.  |  *    *   27   28   29   30 | 35 |   |  *    *   73   74   75   76 | 77 |
  373.  ------------------------------------   ------------------------------------
  374.  
  375. The elements are stored in the array V in a packed fashion following
  376. a fishbone design. The rows below the diagonal and the columns above
  377. the diagonal are stored with stride 1. Elements marked with * are
  378. not stored in V. If you look at V as a 2 dimensional array A, you
  379. see that A is symmetric in shape, but the corresponding symmetric
  380. values are not the same.
  381.  
  382.  
  383.  
  384.