wkb@cbnews.ATT.COM (Wm. Keith Brummett) (04/21/89)
I have need of a small, fast routine in C language that will invert
matrices of order <= 4. (Actually, I'll take a large, slow one that
does order 2 or 3). Does anyone have something like this that they
would care to send me? If it matters, it will be used in a small
graphics project.
BTW, can anyone tell me why it is that every language except C seems
to have standard subroutines to do this?
Thanks in advance!
--
| Wm. Keith Brummett (614) 860-3187 AT&T, Room 2B-252 |
| att!cbnews!wkb Cornet 353 6200 East Broad St. |
| or, wkb@cbnews.ATT.COM Columbus, OH 43213-1550 |
`----------------------------------------------------------------------'k1wolcot@pikes.Colorado.EDU (Kenneth A. Wolcott) (04/21/89)
In article <5785@cbnews.ATT.COM>, wkb@cbnews.ATT.COM (Wm. Keith Brummett) writes: [ ... lines deleted ... ] > BTW, can anyone tell me why it is that every language except C seems > to have standard subroutines to do this? [ ... lines deleted ... ] I did not see any matrix inversion routines "BUILT-IN" since I started in Pascal in 1981 -- let's see -- I've looked at, used or read about as many as twenty different PC-based or mini-computer Pascal compilers -- none had matrix inversion routines! Ada doesn't have one (and people accuse Ada of being the "kitchen-sink" language) -- I must have missed it when I looked at Prime Pl/1 -- I don't think any of our 7 FORTRAN compilers have it... I think that you are after a LIBRARY function (not BUILT-IN) for matrix inversion. If so, that is a installation-dependent characteristic! Someone at site #1 has a FORTRAN compiler on a PC with that library -- someone at site #2 has the same compiler and the same hardware, but doesn't have that library. Same thing with "C". Same thing with Ada, Pascal, Modula-2, PL/1, COBOL -- you name the language -- you'll find different libraries installed. ANSI never (I don't think so) did standardize the libraries for each language! How could it be done? It might be real nice to be able to go to six different sites and get identical library support! Kenneth A. Wolcott University of Colorado at Denver Computing Services Advising & Operations "k1wolcot@pikes.colorado.edu" or "kwolcott@pikes.colorado.edu" or "kwolcott@copper.colorado.edu" or "kwolcott@orphan.colorado.edu" or "kawolcott@cudenver.bitnet" ==========================================================================
gwyn@smoke.BRL.MIL (Doug Gwyn) (04/21/89)
In article <5785@cbnews.ATT.COM> wkb@cbnews.ATT.COM (Wm. Keith Brummett) writes: > I have need of a small, fast routine in C language that will invert > matrices of order <= 4. /* The following do not attempt to handle singularities or ill-conditioning. */ void Inv1( double a[1][1], double b[1][1] ) { b[0][0] = 1.0 / a[0][0]; } void Inv2( double a[2][2], double b[2][2] ) { register double d = a[0][0] * a[1][1] - a[0][1] * a[1][0]; b[0][0] = a[1][1] / d; b[0][1] = -a[1][0] / d; b[1][0] = -a[0][1] / d; b[1][1] = a[0][0] / d; } void Inv3( double a[3][3], double b[3][3] ) { double m00 = a[1][1] * a[2][2] - a[2][1] * a[1][2]; double m01 = a[1][2] * a[2][0] - a[2][2] * a[1][0]; double m02 = a[1][0] * a[2][1] - a[2][0] * a[1][1]; register double d = a[0][0] * m00 + a[0][1] * m01 + a[0][2] * m02; b[0][0] = m00 / d; b[0][1] = m01 / d; b[0][2] = m01 / d; b[1][0] = (a[2][1] * a[0][2] - a[0][1] * a[2][2]) / d; b[1][1] = (a[2][2] * a[0][0] - a[0][2] * a[2][0]) / d; b[1][2] = (a[2][0] * a[0][1] - a[0][0] * a[2][1]) / d; b[2][0] = (a[0][1] * a[1][2] - a[1][1] * a[0][2]) / d; b[2][1] = (a[0][2] * a[1][0] - a[1][2] * a[0][0]) / d; b[2][2] = (a[0][0] * a[1][1] - a[1][0] * a[0][1]) / d; } void Inv4( double a[4][4], double b[4][4] ) { /* XXX -- you provide this yourself, I'm getting tired */ } > BTW, can anyone tell me why it is that every language except C seems > to have standard subroutines to do this? Very few languages have standard matrix inversion functions. APL is about the only one I know of. Certainly not Fortran or Pascal.
vevea@paideia.uchicago.edu (Jack L. Vevea) (04/22/89)
In article <10087@smoke.BRL.MIL> gwyn@brl.arpa (Doug Gwyn) writes: >In article <5785@cbnews.ATT.COM> wkb@cbnews.ATT.COM (Wm. Keith Brummett) writes: >> I have need of a small, fast routine in C language that will invert >> matrices of order <= 4. While we're on this subject, I have need of a long, slow matrix inversion routine that is at least callable from C. I say "long, slow" because I need something I can apply to matrices with dimension in the hundreds, so it must be sensitive to the possibility of rounding error, and pick the optimal pivots. It would be nice if, in addition, it could track its own precision in some way. If anyone could at least point me in the direction of a good algorithm, I would be eternally grateful. (Apologies if this doesn't belong in comp.lang.c)
faustus@dogwood.Berkeley.EDU (Wayne A. Christopher) (04/24/89)
In article <2846@tank.uchicago.edu>, vevea@paideia.uchicago.edu (Jack L. Vevea) writes: > While we're on this subject, I have need of a long, slow matrix > inversion routine that is at least callable from C. I say > "long, slow" because I need something I can apply to matrices with dimension > in the hundreds, so it must be sensitive to the possibility of rounding > error, and pick the optimal pivots. There are a lot of algorithms used to do this sort of thing with sparse matrices in circuit simulation. In particular, the Spice 3 circuit simulator available from Berkeley has a sparse matrix manipulation package that you could probably use. Wayne
kurtk@tekcae.CAX.TEK.COM (Kurt Krueger) (04/24/89)
There are two 'obvious' places to look for code or algorithms: The Collected Alogrithmns of the ACM (CACM). Any technical/college library should have this. The articles usually contain code (but stuff that is well known is likely to be in older articles in Fortran or Pascal) and references to the algorithmn used. Any book on numerical methods. Note: "C" is not well suited for writing code for dealing with an arbitrary 2 (or higher) dimension array. Brush up on your pointers. You get to do the row/column location by hand. A good source would be old Fortran code (Fortran 2 or IBM 1130 Fortran IV) that handles 2-d arrays as a strung out 1-d array. Another "trick" is to NOT pass a 2-d array, but pass a vector of row pointers. You can then use matrix[row][column] in your function as if you have a "real" 2-d array. But don't expect blinding execution speeds.
georg@cbnewsh.ATT.COM (georg.k.karawas) (04/28/89)
The book 'Numerical Recipes in C' should contain a routine for matrix inversion. As it was already pointed in other articles, it is a hassle to work with pointers. A better solution might be to take a 'canned' FORTRAN routine from a standard package, like LINPACK, and call it from C. Be careful though that C stores arrays row-by-row. UNIX f77 compilers can produce object code which can be linked with a C calling program.