Addons/math/mt

From J Wiki
Jump to: navigation, search
User Guide | Installation | Development | Categories | Git | Build Log

math/mt - Matrix Toolbox.

mt is a set of routines for solving some problems in matrix algebra: transforming, decomposing, reducing to condensed form, factorization, equation solving, function applying, condition number estimating. It's based mostly on LAPACK algorithms.

Browse history, source and examples using GitHub.

Installation

Use JAL/Package Manager.

Load math/mt addon with the following line:

   load 'math/mt'

Addon self-testing

Testing routines output the table of following columns:

  • 'algorithm' - command used to test, either verb or bonded verb
  • 'rcond' - reciprocal of condition number of a square matrix, or NaN
  • 'rel fwd err' - relative forward error if possible, or NaN
  • 'rel bwd err' - relative backward error if possible, or NaN
  • 'time, sec.' - execution time, in seconds
  • 'space' - space required to execute command, in bytes

While testing, output is directed into:

  • console (always)
  • global noun TESTLOG (always)
  • log file (optional, depends on global noun TESTLOGFILE value)

To compare mt routines against alternatives, the following addons would be installed:

  • math/lapack of latest version
  • math/misc

First of all, let's define appropriate random array generator:

   NB. Monad to generate an array of shape y
   NB. with elements distributed as U(0,1)
   rag=. ?@$&0

For various generator examples see rand.ijs .

Self-test the entire add-on

   NB. use rag to test addon by random matrices of size 150x150,
   NB. see mt.ijs
   rag test_mt_ 150 150

Be ready of long execution time for matrices of big size (>500).

Self-test the group of algorithms

   NB. use rag to test low-level algorithms by random matrices of size 150x150,
   NB. see mt.ijs
   rag testlow_mt_ 150 150

   NB. use rag to test mid-level algorithms by random matrices of size 150x150,
   NB. see mt.ijs
   rag testmid_mt_ 150 150

   NB. use rag to test high-level algorithms by random matrices of size 150x150,
   NB. see mt.ijs
   rag testhigh_mt_ 150 150

Self-test particular file

   NB. use rag to test xxtrfxxxxx verbs (triangular factorization)
   NB. by random matrices of size 150x150, see trf.ijs
   rag testtrf_mt_ 150 150

Self-test particular set of verbs by pre-defined matrix

   NB. generate random general matrix of size 150x150
   A=. rag 150 150

   NB. use A to test getrfxxxx (triangular factorization of
   NB. a general matrix) by pre-defined matrix, see trf.ijs
   testgetrf_mt_ A

Addon benchmarking

See Benchmarking page.

Example usage

Let's make some arrays to play with:

   A=. rag 5 5
   X=. rag 5 3
   B=. A mp_mt_ X

Now perform some calculations:

   NB. do triangular factorization: P*L1*U=A, see trf.ijs
   'ip L1U'=. getrfpl1u_mt_ A
   L1=. trl1_mt_ L1U
   U=. tru_mt_ L1U
   P=. ip2P_mt_ ip
   A -: P mp_mt_ L1 mp_mt_ U
1

   NB. solve equation: A*X=B, see sv.ijs
   Xapprox=. A gesvax_mt_ B
   X -: Xapprox
1

   NB. do QR factorization: Q*R=A, see qf.ijs
   QfR=. geqrf_mt_ A
   Q=. ungqr_mt_ QfR
   R=. tru_mt_ QfR
   A -: Q mp_mt_ R
1

   NB. raise A to integer powers: A^i.30, see pow.ijs
   Apows=. (i. 30) gepow_mt_ A                            NB. via repeated squaring, O(n^3 * log(30)) FLOPs
   ApowsSF=. (idmat_mt_ # A) 0} |.!.0 mp_mt_/\ 30 # ,: A  NB. straightforward, O(n^3 * 30) FLOPs
   Apows -: ApowsSF
1

   NB. find matrix exponential: E:=exp(A), see exp.ijs
   E=. geexp_mt_ A                                        NB. by scaling and squaring method
   Eapprox=. +/ Apows % ! i. 30                           NB. by finite summation
   E -: Eapprox
1

mt for users

mt follows LAPACK's name conventions whenever possible, e.g. getrf verb of mt corresponds to both DGETRF and ZGETRF (xGETRF) subroutines of LAPACK. Floating point datatype input usually forces J-actor (verb/adverb/conjunction) to behave like a Dxxxxx counterpart in LAPACK. Complex datatype input usually forces J-actor to behave like a Zxxxxx counterpart.

Most J-actors in mt has particular specialization. E.g. larftbc corresponds to LAPACK subroutine xLARFT('B','C',...), larftfr corresponds to xLARFT('F','R',...), and so on.

When referenced, names may be grouped by wildcard 'x', e.g.:

  • xxexp means: geexp, diexp, heexp
  • getrfxxxx means: getrflu1p, getrfpl1u, getrfpu1l, getrful1p

Similarity of J-actor in mt to its correspondent counterpart in LAPACK has some levels:

  • implementation: input and output storage layouts are the same, algorighm is essentially the same; e.g. gesvax implements xGESV
  • model: input or output storage layouts may differ, or algorighm is differrent; e.g. getrf models xGETRF
  • simulation: input or output storage layouts are differrent, algorighm is essentially differrent; e.g. gecon1 simulates xGECON('1')
  • similar: provides similar functionality, but is completely different; e.g. hetrfpl is similar to xHETRF('L')
  • new: provides new functionality; e.g. geexp

Currently only the limited set of routines is implemented, see the MATRIX.

mt for developers

LAPACK J addon uses original LAPACK library as backend. So, it is rather simple to add new interface for any LAPACK subroutine. But, there is excessive data pre- and post-processing, and addon user is limited to data types and storage layouts provided by Fortran.

As opposed to LAPACK J addon, mt is written on pure J and doesn't require LAPACK for underlying calculations. So, any J data types and storage layouts can be used. But, adding new stuff requires programming of algorithms from scratch.

The file template.ijs may be used as a seed to start with.

Debug

Any verb can be switched into the debugging mode by applying a conjunction dbg (file dbg.ijs). Under this mode, verb outputs in console some additional information. Verbosity of debug is defined by global noun DEBUG (file mt.ijs). E.g.:

   NB. statement to debug
   +/\ 100 20 3
100 120 123

   NB. default verbosity level is 2, meaning output of rank and
   NB. valency, input's and output's shapes and values

   NB. debug (+/) with current verbosity level (2)
   (+/ dbg_mt_ '+/')\ 100 20 3
┌───┬────────────────┬─┬─┬───┐
│dbg│+/ [MONAD] _ _ _│y│1│100│
└───┴────────────────┴─┴─┴───┘
┌───┬──────────┬┬───┐
│dbg│+/ SUCCESS││100│
└───┴──────────┴┴───┘
┌───┬────────────────┬─┬─┬──────┐
│dbg│+/ [MONAD] _ _ _│y│2│100 20│
└───┴────────────────┴─┴─┴──────┘
┌───┬──────────┬┬───┐
│dbg│+/ SUCCESS││120│
└───┴──────────┴┴───┘
┌───┬────────────────┬─┬─┬────────┐
│dbg│+/ [MONAD] _ _ _│y│3│100 20 3│
└───┴────────────────┴─┴─┴────────┘
┌───┬──────────┬┬───┐
│dbg│+/ SUCCESS││123│
└───┴──────────┴┴───┘
100 120 123

   NB. set verbosity level to 1, meaning output of rank and
   NB. valency, input's and output's shapes without values
   DEBUG_mt_=: 1

   NB. debug (+/) with current verbosity level (1)
   (+/ dbg_mt_ '+/')\ 100 20 3
┌───┬────────────────┬─┬─┐
│dbg│+/ [MONAD] _ _ _│y│1│
└───┴────────────────┴─┴─┘
┌───┬──────────┬┐
│dbg│+/ SUCCESS││
└───┴──────────┴┘
┌───┬────────────────┬─┬─┐
│dbg│+/ [MONAD] _ _ _│y│2│
└───┴────────────────┴─┴─┘
┌───┬──────────┬┐
│dbg│+/ SUCCESS││
└───┴──────────┴┘
┌───┬────────────────┬─┬─┐
│dbg│+/ [MONAD] _ _ _│y│3│
└───┴────────────────┴─┴─┘
┌───┬──────────┬┐
│dbg│+/ SUCCESS││
└───┴──────────┴┘
100 120 123

   NB. set verbosity level to 0, meaning
   NB. execute original verb unmangled
   DEBUG_mt_=: 0

   NB. debug (+/) with current verbosity level (0)
   (+/ dbg_mt_ '+/')\ 100 20 3
100 120 123

Authors

Authorship of mt addon code belongs to various authors. All the authors are listed below:

Name File Orignial name Author Published Link
ag util.ijs appl Henry Rich 2005-10-22 06:37:12 [Jforum gerund apply]
icut struct.ijs bcut1 Roger Hui 2005-11-24 03:53:19 [JWiki Essays/Block Matrix Inverse]
upd struct.ijs Amend Oleg Kobchenko 2007-03-02 22:15:54 [Jprogramming Transform to Amend]

All the rest code was written by Igor Zhuravlov.

See Also

  • [Jgeneral GPL as addon's license], Igor Zhuravlov, 2010-05-27 13:31:28 HKT
  • [Jprogramming New addon: math/mt (Matrix toolbox)], Igor Zhuravlov, 2010-06-05 19:21:44 HKT
  • [Jprogramming Eigen values in J - example] Igor Zhuravlov, 2013-04-09 07:51:53 UTC