Linear Algebra Routines

[Abel]

Back when the 8087 coprocessor made its debut, I had occasion to put together asm coding for some fundamental linear algebraic routines to solve problems that otherwise required access to a central computer with a library of Fortran linpack/eispack/lapack routines.  Fortran source code libraries are available online at http://www.netlib.org and offer a reliable starting point for algorithm translations to other languages.

To get in a little programming exercise, I tried my hand at seeing how difficult it was to adapt the Pentium to solving such problems, which so far involve the manipulation of 2-dimensional matrices of real, double-precision floating point numbers.  Asm programming appears to have two clear advantages over other languages.  The addressing modes, e.g. qword ptr [-8+8*ecx+esi], let one access A[i-1,j] with ecx=i and [esi]=A[0,j] so that only register variables are required for indexing and i and j values never occupy memory/stack addresses.  Variable column addressing is less convenient but, fortunately, coding can be arranged so that most looping is over common column indices.  Similarly, the floating point stack may be exploited to hold intermediate values without need for external memory.

Presently, I've worked out codings to factor a square matrix into the product of a lower and upper triangular matrix (LU factorization) and to find the unitary transform diagonalizing a symmetric, square matrix (Ad = ~Q*A*Q).  The former routine provides a basis for finding determinants, inverses, and solving sets of linear equations.  The latter is basically the eigenvalue/vector problem.  My humble 120MHz/48M box solves 300 simultaneous equations for 300 unknowns in 5.0 sec and takes 5.5 sec to find all the eigenvectors and eigenvalues of a 100x100 matrix (0.3 sec for values only).  Unfortunately, the algorithms scale as N cubed for large N and N~1000 seems a likely upper limit wrt both CPU time and memory consumption, if one isn't interested in all-night calculations.

To me, these seem rather essential tools for data processing, yet I've come up with zilch searching the net for equivalent, open source coding.  Any links would be appreciated.  Otherwise, if there's anyone who follows what the h##! I'm talking about and interested in testing and refining a wip, I'll be happy to zip commented asm code for the current DLL and a couple BASIC text codings I've used for tests, with the GNU/GPL restrictions to open source application only.  No experience here with optimization and have no idea what can be accomplished for 30 pages of floating point opcode other than the default allocation alignment on 8-byte boundaries.

Abel

[Biterider]

Hello Abel
ObjAsm32 has 2 objects called Array and Matrix that lets you perform complex mathematical computations on nxn matices. NaN has written 3 examples on how to use them (LU factorization, complex matrix handling, matrix inversion, etc).
Take a look at this work, it can be an interesting starting point or perhaps we can share our code.

Regards,

Biterider

Our link: http://objasm32.tripod.com/index.htm

[Abel]

Hi Byterider,
Many thanks for the link.  So many ways to skin the poor cat!  As an amateur programmer, what I'm looking for is open source code I'm free to fold, spindle, mutilate and use as I will.  In return, I'm obliged to clearly reference my sources and offer what bug fixes and refinements of significance I can.  As I understand your EULA, our interests aren't in full alignment.

For any amateurs of like-mind, I've zipped the pack below.  The matrix routines are all in one DLL, but individual procedures can be extracted as separate library functions.  The easiest way to test algebraic software is to plug the answer back in the original equations and check that the equalities are satisfied.  To simplify matters, there's also a routine to fill matrices with random, floating-point integers.  Checks should then regurgitate screens of 0's to 15 or so decimal places if all's hunky-dory.  For optimization, remember that 1 cycle saved in the innermost loop of a 100x100 matrix is worth a million cycles outermost (N^3).  There's also a local version of random.inc from my include folder.  It's Pentium only (.586) and outside Hutch's ,486 spec, but the rdtsc opcode is just too convenient to avoid when a 'random' seed is called for.

Abel

http://prove66.home.att.net/files/asm/Matpack2.zip

[Biterider]

Hi Abel!
Thanks for your replay. I don't think that we are missaligned about the EULA. Can you specify the points you mean?
In the next days I'll check the content of your link.

Best regards,

Biterider

[NaN]

Ya.. you got me curious over this point as well?  What is bothering you about our verison of the EULA?  In short as I remember it (didnt reread it for this post) is that your free to do what ever, just dont sell the source as your own.  The exe's are your own and free to do what you will with.

The basic reason is to protect our work (phyiscally in macro code and objects) while having the intent of our work remainll free for people to use and take advantage of (build exe's and working projects).  We just dont want someone writing books on our work, listing source and making a 'source cd' from OA32 with out special permissions. 

So as i understand your needs, your should be free and clear to build any Matrix work you want...
Regards,
:NaN:

[Abel]

Biterider & NaN,

My reservations are with regard to the EULA statement:

"2. You are completely free to distribute your own developed source, object or executable code for commercial purposes, but NOT the source code of the SOFTWARE PRODUCT."

That leaves ambiguous whether I may distribute all, a portion, or a modified portion, of the SOFTWARE PRODUCT for noncommercial purposes but the NOT tells me, if in doubt no.

Most forum members are source-code oriented.  My interpretation of the EULA precludes me from incorporating in posted source code a procedure derived from your library.  As a hypothetical, suppose the library contains a routine for plotting a series of points in a window using the PolyBezier API function to create a smoothed appearance.  This involves assigning values to two intermediate points per interval based on arbitrary constraints - continuity, slope, curvature, ...  The constraints I want to use differ from the library function's built-in defaults.  Am I permitted to post/distribute a zip containing modified source code with due acknowledgement of its derivative origin?

Libraries serve both as blackbox repositories and starting points from which one  builds on the work of others, e.g. www.netlib.org for FORTRAN users at universities et al.  Whether I'm in Cambridge MA or Cambridge UK, I can access a local library and find the same code.  Proprietary libraries, on the other hand, are single source sites generally tied to a specific software product be it commercial or otherwise - and they have every right to be as generous or restrictive in the use of their content as they deem fit. 

My background lies in the hard sciences of an age when free dissemination of ideas and information was the 'sine qua non' and problem-solving a collegial endeavour.  That's the qualification I place on my own efforts and any like-minded individual or library is welcome to whatever they find of value with the understanding that anyone is free to modify, build-upon, and contribute to the code as they see fit but make no proprietary claim to derivative content.

Regards,
Abel

[NaN]

That is in short our Basic intent.  I see your point, and confusion.  But I can assure you we only want to protect our invested time with the project.  Not take advantage of other's time.

I can't speak for BiteRider, as its more his product than my own.  But I can say there *should* be no reason why you can't modify and distribute any source.  The *core* routines where the bulk of the scripting resides is the focus of the "SOFTWARE PRODUCT" as I understand it from our discussions.  This is the definitions of macros like "Object" and "Method" etc.  The Nutz and bolts of the system is what we were aiming to protect from commercial redistribution (How objects are oganized in memory and instanciated is entirely his property as a result of his hard work in research and design).  Objects and Object sources are free to post and discuss and modify etc to your hearts content.  They are all useless with out the core scripts that know what to do when terms such as "Object" and "Method" is found.   Your can go ahead and make a book on objects you've written, or modified from our library of existing objects.  Who ever buys the book would have to download Biteriders package to make any use of them... as your NOT alowed to publish the core routines that defined the macros used in your Object Sources.

BiteRider is off seing the world at the moment... Last I heard he is in south america (far away from home), but when he is back I think this point should be revisited.   Im no lawyer and to me our EULA fit the bill, but your interpretation is definitely a persective I hadnt considered.  I think the best approach would be to list the core components in the package that we would deem protected from comercial redistribution. 

I certaining dont like the idea of an EULA that essentially leaves the impression that all source your create is banned from distribution with out written concent.  That would go 100% against our intent in producing the package in the first place...

Thank you for being specific on this point.  Its a bit of an eye opener ;)
Regards,
:NaN: