Studio/LAPACK

From J Wiki
Jump to navigation Jump to search

Lab: LAPACK

LAPACK (Linear Algebra Package) is a set of routines for solving systems of simultaneous linear equations, least-squares solutions of linear systems of equations, eigenvalue problems, and singular value problems.

The associated matrix factorizations, LU, Cholesky, QR, SVD, Schur, generalized Schur, are also provided, as are related computations such as reordering of the Schur factorizations and estimating condition numbers.

The LAPACK routines are accessed through the dll call mechanism as documented in User Manual chapter Dlls and Memory Management.

All definitions for LAPACK are in locale jlapack.

The scripts and other files for LAPACK addon are in directory:

   jpath '~addons/math/lapack'
D:\Math\j602\addons\math\lapack

The first step in using LAPACK is to load the lapack.ijs script to define the jlapack locale.

   load 'math/lapack'

There are hundreds of LAPACK routines but you only load interfaces for the ones you need. The interface is defined in a script with the name of the routine.

geev is a LAPACK routine for computing eigenvalues of a real or complex matrix.

If you want to use geev you need to load the geev script to make the necessary definitions in the jlapack locale.

   load 'math/lapack/geev'

   ]m=: ?.4 4$100
46 55 79 52
54 39 60 57
60 94 46 78
13 18 51 92
   ,.geev_jlapack_ m  NB. compute eigenvalues and eigenvectors
+----------------------------------------------------------------+
|_0.359503   _0.276552j0.152898   _0.276552j_0.152898    0.368957|
|_0.437463             0.840512              0.840512    0.269732|
|_0.505807  _0.408889j_0.157841    _0.408889j0.157841 _0.00397621|
|_0.650802 _0.0117517j0.0379341 _0.0117517j_0.0379341   _0.889438|
+----------------------------------------------------------------+
|219.662 _25.077j6.83493 _25.077j_6.83493 53.4921                |
+----------------------------------------------------------------+
|_0.531971  _0.630475j_0.191732  _0.630475j0.191732  0.509658    |
|_0.473031   0.0850361j0.190114 0.0850361j_0.190114  0.311404    |
|_0.602175              0.68181             0.68181  0.316593    |
|_0.361432 _0.238794j_0.0218802 _0.238794j0.0218802 _0.736916    |
+----------------------------------------------------------------+

   >1{geev_jlapack_ m  NB. just the eigenvalues
219.662 _25.077j6.83493 _25.077j_6.83493 53.4921

If you are going to be using geev a lot, you could define it in the base locale so that you do not have to use the full name.

   geev=: geev_jlapack_
   >1{geev m
219.662 _25.077j6.83493 _25.077j_6.83493 53.4921

The LAPACK routines are documented in the doc subdirectory with a file with the same name, but a "lap" suffix.

This documentation is directly from the LAPACK library and is for both the mathematician and the programmer. Unfortunately the programmer is assumed to be a Fortran programmer. This means you will have to put in some extra effort to understand the programming interface details.

docs_jlapack_ shows a window with a list of the lap files.

   docs_jlapack_ ''  NB. display list of lap files
Docs.png

JSvnAddons:math/lapack/doc

   load 'math/lapack/geev'
   >1{geev_jlapack_ 1j3+m        NB. eigenvalues of complex array
223.781j12.2838 53.362j_0.313279 _24.9688j_6.86053 _25.174j6.89006

gesvd computes the singular value decomposition (SVD) and left and right singular vectors of a real or complex matrix,

   load 'math/lapack/gesvd'
   ,.gesvd_jlapack_ m           NB. SVD
+------------------------------------------+
|_0.504658   _0.106691  0.625094  _0.585828|
|_0.455858 _0.00284226  0.383999   0.802951|
|_0.605257   _0.497406 _0.619549 _0.0490927|
|_0.413736    0.860928 _0.279215 _0.0983119|
+------------------------------------------+
|229.101     0      0       0              |
|      0 57.51      0       0              |
|      0     0 39.178       0              |
|      0     0      0 15.3781              |
+------------------------------------------+
|_0.390765 _0.412338  0.221745  0.792533   |
|_0.479597  _0.64751 _0.354979 _0.474034   |
|_0.507034  0.216093   0.75765 _0.349554   |
|_0.600172  0.603335 _0.500785  0.158098   |
+------------------------------------------+

The LAPACK routines are quite fast. For example:

   mat=: ?.100 100$1000

   time=: 6!:2

   NB. time in seconds to compute SVD:
   time 'gesvd_jlapack_ mat'
0.0380163


   load 'math/lapack/potrf'
   mp=: +/ . *
   potrf_jlapack_ m mp |: m   NB. Cholesky decomposition
118.684       0        0      0
103.914 22.0871        0      0
131.609 19.7358  54.0874      0
87.6357 27.2287 _1.42499 55.987


   load 'math/lapack/geqrf'
   ,.geqrf_jlapack_ m           NB. QR factorization
+----------------------------------------+
|_0.490334 0.00198915  0.515248 _0.702914|
|_0.575609  _0.749862 _0.231658  0.229598|
|_0.639566   0.657569 _0.348576  0.192493|
|_0.138573  0.0728396  0.747896  0.645091|
+----------------------------------------+
|_93.8136 _112.031  _109.76 _120.942     |
|       0  33.9874 _10.8716  15.3529     |
|       0        0  48.9134   55.206     |
|       0        0        0  50.8984     |
+----------------------------------------+


   load 'math/lapack/gesv'
   m1=: ?.4 3$30
   gesv_jlapack_ m;m1         NB. solves m * x = m1
 0.627233  0.413047 _0.437923
_0.277309 _0.229798  0.293523
_0.273575 _0.080829  0.562991
 0.269455   0.09662 _0.220684


   load 'math/lapack/getrf'
   getrf_jlapack_ m           NB. LU decomposition
+-----------------------------+------------------------+-------+
|       1         0        0 0|60    94     46       78|3 2 4 4|
|     0.9         1        0 0| 0 _45.6   18.6    _13.2|       |
|0.216667 0.0519006        1 0| 0     0 40.068  75.7851|       |
|0.766667  0.374269 0.917738 1| 0     0      0 _72.4105|       |
+-----------------------------+------------------------+-------+

routines_jlapack_ lists all the scripts in the lapack directory.

   routines_jlapack_ ''
Routines.png

JSvnAddons:math/lapack

As you can see, interfaces are defined for only a few LAPACK routines.

Compare the list of available routines in docs_jlapack_ with the list of implemented interfaces in routines_jlapack_ .

With a little effort you should be able to add an interface for the LAPACK routine you require. You need to study the lap documentation for interfaces already provided, and by following the pattern in their interfaces, write a script for the interface you want.

If you write LAPACK interfaces, please send them to us so they can be included in the package.

The calculation of eigenvalues using LAPACK is an integral part of a new demo of eigenpictures. These show the action of a matrix on unit vectors, and highlight the eigenvalues, which appear as a scaling factor times the eigenvectors.

Access this demo from the menu Studio|Demos|eigenpictures, or by loading it directly:

load '~system\examples\graphics\isigraph\eigenpic.ijs'

For more information on LAPACK: http://www.netlib.org