Addons/math/mt
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.
Install
Download
If you prefer jconsole or Qt Jconsole:
load 'pacman'
'update' jpkg ''
Updating server catalog...
Local JAL information has never been updated.
Installed addons are up to date, 1 addon are not yet installed.
The base library is up to date.
'install' jpkg 'math/mt'
Installing 1 package
Downloading math/mt...
Installing math/mt...
Done.
All available packages are installed and up to date.
The base library is up to date.
Or use Package Manager.
Pre-configure
NB. set paths to libs
LIB_mtbla_=: '/usr/lib64/openblas-pthreads/libopenblas.so.0' NB. system's BLAS
LIB_mtbli_=: (getenv 'HOME') , '/lib/libblis_threads=n.so' NB. user's BLIS
This step is necessary if an external library which performance you are going to compare with mt addon resides in non-standard place.
See also lapack2 pre-configuring if you are planning to use lapack2 addon, too.
Load
load 'math/mt'
Post-configure
NB. load pre-configured external library explicitly
load 'math/mt/external/blis/blis' NB. populates 'mtbli' locale
NB. enable multi-threading for BLIS library
thread_set_num_threads_cd_mtbli_ 8 NB. e.g. 8 threads
This step may be necessary for some external libraries involved in performance comparison.
Introspect
Environment
env_mt_ ''
Hardware
architecture: x86_64
cores: 8
maxthreads: 63
worker threads: 0
threads in pool# (#idle,#unfinished,#threads):
#0: 0 0 0
#1: 0 0 0
#2: 0 0 0
#3: 0 0 0
#4: 0 0 0
#5: 0 0 0
#6: 0 0 0
#7: 0 0 0
executing thread# (0 means master thread): 0
Thresholds to switch (+/ .*) to GEMM (built-in BLAS implementation)
(switches if threshold <= m*n*p, _1 = switch never)
for datatypes:
integer : _1
floating: _1
complex : _1
Interfaces to external libraries (optional)
LAPACK
library: liblapack.so.3
version: 3.12.0
BLAS
library: /usr/lib64/openblas-pthreads/libopenblas.so.0
version: OpenBLAS 0.3.27 DYNAMIC_ARCH NO_AFFINITY Haswell MAX_THREADS=64
BLIS
library: /home/user/lib/libblis_threads=n.so
version: 0.9.0
Some name
NB. get description of 'qnlen_mt_' monad
'd' info_mt_ 'qnlen'
Description:
Length
NB. get syntax of 'qnlen_mt_' monad
's' info_mt_ 'qnlen'
Syntax:
len=. qnlen q
where
q - a single quaternion or laminated quatertions
len - a scalar or (# q)-vector of non-negative numbers
NB. get both description and syntax of 'qnlen_mt_' monad
'ds' info_mt_ 'qnlen'
Description:
Length
Syntax:
len=. qnlen q
where
q - a single quaternion or laminated quatertions
len - a scalar or (# q)-vector of non-negative numbers
NB. get full information about 'qnlen_mt_' monad
info_mt_ 'qnlen'
---------------------------------------------------------
qnlen
Description:
Length
Syntax:
len=. qnlen q
where
q - a single quaternion or laminated quatertions
len - a scalar or (# q)-vector of non-negative numbers
Assertions:
0 *./@:<: qnlen q
Self-test
Testing is performed with random input data which are generated after the process is started.
Testing routines return the table of following columns:
- 'sentence' - 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
This table is also displayed in console while produced row-by-row.
To compare mt routines against alternatives, the following addons would be installed:
- math/lapack2
- math/misc
If your system have BLAS or BLIS libraries, their routines will be compared, too.
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 file.
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 for long execution time for matrices of big size (>500).
Currently more than 13,400 tests are supplied which cover functionality of all computing modules except scl.ijs.
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
A 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
A 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
Verify
Verifying is performed with pre-defined hard-coded input and output data which were selected thoroughly.
Verifying routine takes no arguments and returns result as a boolean scalar value:
verify_mt_ ''
1
0 verify_mt_ '' NB. run silently, the same as above
1
1 verify_mt_ '' NB. run verbosely, to see tests probed
<very long output...>
1
Verbose mode is useful when some test(s) is/are failed.
Currently more than 11,200 tests are supplied which cover functionality of 7 basic modules.
Benchmark
See Benchmarking sub-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
lapack2 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 lapack2 addon, mt is written on pure J and doesn't require LAPACK library 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 template.ijs file may be used as a seed to start with.
Directories structure
mt |
main staff |
Locales structure
mttmp |
temporal staff for testing: utilities, LAPACK cover verbs, other math addons |
Debug
Any verb can be switched into the debugging mode by applying a conjunction dbg (see dbg.ijs file). Under this mode, verb outputs in console some additional information. Verbosity of debug is defined by global noun DEBUG (see mt.ijs file). 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, locale name, input and output shapes and values
NB. debug (+/) with current verbosity level (2)
(+/ dbg_mt_ '+/')\ 100 20 3
┌──┬─────┬─────┬────┬─┬─┬───┐
│+/│MONAD│_ _ _│base│y│1│100│
└──┴─────┴─────┴────┴─┴─┴───┘
┌──┬───────┬────┬──────┬┬───┐
│+/│SUCCEED│base│result││100│
└──┴───────┴────┴──────┴┴───┘
┌──┬─────┬─────┬────┬─┬─┬──────┐
│+/│MONAD│_ _ _│base│y│2│100 20│
└──┴─────┴─────┴────┴─┴─┴──────┘
┌──┬───────┬────┬──────┬┬───┐
│+/│SUCCEED│base│result││120│
└──┴───────┴────┴──────┴┴───┘
┌──┬─────┬─────┬────┬─┬─┬────────┐
│+/│MONAD│_ _ _│base│y│3│100 20 3│
└──┴─────┴─────┴────┴─┴─┴────────┘
┌──┬───────┬────┬──────┬┬───┐
│+/│SUCCEED│base│result││123│
└──┴───────┴────┴──────┴┴───┘
100 120 123
NB. set verbosity level to 1, meaning output of rank and
NB. valency, locale name, input's and output's shapes without values
DEBUG_mt_=: 1
NB. debug (+/) with current verbosity level (1)
(+/ dbg_mt_ '+/')\ 100 20 3
┌──┬─────┬─────┬────┬─┬─┐
│+/│MONAD│_ _ _│base│y│1│
└──┴─────┴─────┴────┴─┴─┘
┌──┬───────┬────┬──────┬┐
│+/│SUCCEED│base│result││
└──┴───────┴────┴──────┴┘
┌──┬─────┬─────┬────┬─┬─┐
│+/│MONAD│_ _ _│base│y│2│
└──┴─────┴─────┴────┴─┴─┘
┌──┬───────┬────┬──────┬┐
│+/│SUCCEED│base│result││
└──┴───────┴────┴──────┴┘
┌──┬─────┬─────┬────┬─┬─┐
│+/│MONAD│_ _ _│base│y│3│
└──┴─────┴─────┴────┴─┴─┘
┌──┬───────┬────┬──────┬┐
│+/│SUCCEED│base│result││
└──┴───────┴────┴──────┴┘
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 icut struct.ijs bcut1 Roger Hui 2005-11-24 02:11 [JWiki] Essays/Block Matrix Inverse issym mm.ijs sym Roger Hui 2003-08-08 00:46:00 [Jgeneral] Testing whether array is symmetric see also [JWiki] Essays/Symmetric Array
upd struct.ijs Amend Oleg Kobchenko 2007-03-02 14: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
- [Jprogramming] Question about least squares with %. Igor Zhuravlov, 2021-05-08 03:58:58 UTC
- How to work with Matrix Market Exchange Formats