[94] in mathematical software users group
Re: Multidimensional Matrix Operations
daemon@ATHENA.MIT.EDU (Sinan Karahan)
Thu May 20 20:48:37 1993
Date: Thu, 20 May 93 17:50:48 PDT
From: sinan@isi.com (Sinan Karahan)
To: msug@Athena.MIT.EDU, sjw@bongo.cc.utexas.edu
> From sjw@bongo.cc.utexas.edu Thu May 20 09:17:04 1993
> Return-Path: <sjw@bongo.cc.utexas.edu>
> Received: from isi.com (HOPPER) by hardrock (4.1/inc/isi-1.6s)
> id AA08450; Thu, 20 May 93 09:17:02 PDT
> Received: by isi.com (5.57/Ultrix3.0-C)
> id AA05558; Thu, 20 May 93 09:05:34 -0700
> Received: from bongo.cc.utexas.edu by Athena.MIT.EDU with SMTP
> id AA04578; Thu, 20 May 93 12:02:34 EDT
> Received: by bongo.cc.utexas.edu (5.57/Ultrix3.0/cc-ix.mc-1.5)
> id AA14353; Thu, 20 May 93 11:02:16 -0500
> Date: Thu, 20 May 93 11:02:16 -0500
> From: sjw@bongo.cc.utexas.edu (Sue J. Worden)
> Message-Id: <9305201602.AA14353@bongo.cc.utexas.edu>
> To: msug@Athena.MIT.EDU
> Subject: Multidimensional Matrix Operations
> Status: R
>
> > To: sjw@bongo.cc.utexas.edu (Sue J. Worden)
> > Cc: msug@Athena.MIT.EDU
> > Subject: Re: 3D FFT
> > Date: Thu, 20 May 93 11:01:10 EDT
> > From: "Reid M. Pinchback" <reidmp@Athena.MIT.EDU>
>
> > What kinds of mathematical operations do people perform on
> > multidimensional matrices? [ ... ]
>
> Well, I'm not sure, to tell you the truth. The user who
> inquired here was doing some kind of biomechanic modeling
> using tensors, and he said he needed to multiply two 3D
> matrices. I'm not sure what multiplication *means* in a
> case like that. For one thing, it seems ambiguous. I've
> been hoping that if I saw some code, it would all become
> crystal clear... :-)
>
> Sue Worden
> Mathematical Services
> Computation Center, Mail Code 12700
> University of Texas at Austin
> Austin, Texas USA 78712-1110
> voice: (512) 471-3359
> fax: (512) 471-1582
> email: sjw@bongo.cc.utexas.edu
>
I did some investigation on this topic some time ago, I thought I
would share it with you.
CURRENT STATE OF THE ART:
Currently, there is no commercial mathematics/controls/signal
processing software with multi-dimensional objects. Multi-dimensional
matrices can be handled to a certain extent using Kronecker Algebra
techniques but this requires the introduction of some tedious and
complicated notation.
For a summary of these techniques used in Control Systems, see
[Brewer1978]. In matrix-based numerical software packages like Xmath,
Matrixx (Integrated Systems) and Matlab (Mathworks), the Kronecker
multiplication function is useful in these applications, though in a
very limited way. Whenever there is a need to handle data that
requires higher dimensional arrays, one has to resort to inefficient
methods such as concatenating the sub- matrices in an array, and
referencing to the data by index- offsetting.
There are many problems in applied mathematics, control engineering,
signal processing, fluid mechanics, solid mechanics, physics and other
areas of applied mathematics which can be cast in multi-indexed
notation. Some of these applications have traditionally used tensors.
SOME EXAMPLES:
CONTROL SYSTEMS
1) Lyapunov Equation
This is a matrix equation of the form LX + XN = Y. It can be solved
by column-stacking the unknown matrix X and the right-hand side matrix
Y. When this is done, after applying some Kronecker-algebraic
identities one obtains:
(I*L + N'*I) cs(X) = cs(Y)
where * is the Kronecker multiplication operator, ' stands for
transpose, and cs(.) is the column-stack operator for matrices. This
gives rise to n^2 linear equations.
2) Partial Fraction Expansion of Rational Matrix Transfer Functions
The collection of the residue matrices can be easily identified with
three indices, corresponding to inputs, outputs, and eigenvalues.
Treated this way, the computational procedure becomes a less
cumbersome task because of the more efficient notation. See [Brewer,
Tait, Vetter] for the above two applications and some other useful
examples.
3) Sensitivity analysis of systems with varying or uncertain parameters
For these applications, one needs to take derivatives of matrices with
respect to vectors or other matrices. Organization of this data is a
very tedious task without multi-dimensional objects. For a good
example of some earlier attempts in this, along with the treatment of
multi-variable Taylor expansions, see [Vetter]. This is a good example
of the kind of applications that led to the development of polyx
theory (see appendix).
4) Taylor expansions of nonlinear functions that represent system,
input, and output nonlinearities in terms of higher degree
polynomials, and polynomial fitting of these data,
5) Computation of "directional derivatives" of vector fields to
calculate linearizing or stabilizing control laws. For even the
simplest cases (second- and third-order expansions) the amount of data
that needs to be manipulated is of order n^3 and n^4, respectively,
where n is the dimension of the system.
MULTI-DIMENSIONAL SIGNAL PROCESSING
In many signal processing applications, there is a need to express
patterned matrices in a compact and efficient way. There have been
some attempts to use Kronecker products, direct sums, and unitary
transformations to handle these cases. The drawbacks of these are:
a) The notation is awkward and difficult to use,
b) The expression, and ultimately, the storage of the data
require very large amounts of memory, most of which is
sparse and unnecessary.
Multi-dimensional raster-scan mappings of data, and filter bank
applications are good examples. See, for instance, [Regalia and
Mitra, 1988].
OTHER APPLICATIONS
Multi-variate interpolation.
Numerical Optimization (calculation of gradients, hessians, etc.)
Elasticity, Fluid Mechanics, Physics and all other engineering/
scientific disciplines that use tensors.
Neural Networks, where large amounts of data need to be
treated in a structured way.
APPENDIX
POLYXA: A NEW PARADIGM FOR MULTI-INDEXED OBJECTS
A number of Control Scinece and Applied Mathematics
researchers (J. W. Brewer, G. R. Tait, W. J. Vetter) have recently
developed a new method for dealing with multi-indexed quantites.
Vetter has coined the term polyx for POLY-indeXed quantity. Tensors
and matrices are sub-classes of polyxa (plural form of polyx) but not
exhaustive examples. They have given these types of structures a
systematic notation (which is traditionally one of the biggest
obstacles in dealing with multi-indexed objects).
In this appendix, I will try to present a very short summary
of their work, which, I believe, can serve as a very useful guide for
subsequent development of Xmath features using multi-dimensional
matrices. The material presented here is yet to be published, and it
is mainly quoted from [Brewer, Tait, Vetter]. If you want further
information, I can provide some references.
The notation for a polyx is a name with two strings of (free-running)
indices. For example,
iknm
P
ab
The mapping of these polyxa to matrices is important, because only in
this way, one is able to employ the many matrix functions and
operations available to us. The polyx rule chosen here is consistent
with the tensor literature (covariant tensors, i.e. subscripts are
mapped to row vectors and contravariant tensors, i.e. superscripts are
mapped into columns). In order to illustrate this rule, we notice
that the elements of a matrix constitute a polyx and the mapping rule
requires that the simple polyx element
i
M
k
maps to the matrix element in the ith row and kth column of the matrix
M, for fixed i and k.
Polyx equality is used to indicate that two polyxa have the same
information content though the structures might be different. For
instance,
i ik
P = Q
k
Structural Equality means equality in the matrix sense: information
content as well as structural location.
To illustrate with examples, lets consider how two identity polyxa
maps to matrices:
2 | 1 0 |
I = | |
2 | 0 1 |
is a 2x2 matrix, and
22
I = [1 0 0 1]
is a 1x4 vector. These two objects aare Polyx-equivalent, but not
structural-equivalent.
A polyx multiplication is an operation that forms a larger polyx from
two smaller ones:
abcd ab cd
P = R H
efgh ef gh
This operation has a resemblance to kronecker multiplication of
matricesin the matrix domain.
A polyx contraction of two polyxa is the sum over a common index of
the product of polyx elements. As in tensor notation, a common index
among two polyxa automatically implies a sum over the range of the
indices. Example:
ia ---- i ra
P = \ R H
kb / kr b
----
r
The above operation has some correspondence to the usual matrix
product in the matrix domain.
There are many more operations, identities and equivalence relations
defined for polyxa and their mappings to matrix domain. The
systematic treatment of multi-dimensional data types in this framework
makes the task of creating polyx applications much easier. To my
knowledge, this is the only study of its kind on multi-indexed
quantities.
If anyone is interested in receiving copies of the unpublished work
mentioned below, I can help them get in touch with the authors
(though I cannot vouch for the availibility of the manuscripts).
REFERENCES
[Brewer1978] Brewer, J.W., "Kronecker Products and Matrix Calculus in
System Theory," IEEE Trans. on Circuits and Systems, Vol. CAS-25, No.
9, pp. 772-81Sept. 1978.
[Brewer, Tait, Vetter] Brewer, J.W., G.R. Tait, and W.J. Vetter,
"Polyxa: Applications to System and Control Theory," unpublished work.
[Regalia, Mitra] Regalia, A. P., and Sanjit K. Mitra, "Kronecker
Products, Unitary Matrices, and Signal Processing Applications," SIAM
Review, Vol.31, No. 4, pp. 586-613, Dec. 1989.
[Tait, Brewer, Vetter] Tait, G.R., Brewer, J.W., and W.J. Vetter,
"Polyx Methods: Theoretical Foundations and Examples," unpublished
work.
[Vetter] Vetter, W.J., "Matrix Calculus Operations and Taylor
Expansions," SIAM Review, Vol. 15, No. 2, pp. 352-369, Apr. 1973.
===============================================================================
Dr. Sinan Karahan sinan@isi.com
Integrated Systems Inc (408) 980-1500 x340
3260 Jay Street, Santa Clara, CA 95054 (408) 980-0400 (fax)
===============================================================================