home *** CD-ROM | disk | FTP | other *** search
/ NetNews Usenet Archive 1992 #30 / NN_1992_30.iso / spool / comp / lang / cplus / 18050 < prev    next >
Encoding:
Text File  |  1992-12-16  |  7.1 KB  |  233 lines
  1. Xref: sparky comp.lang.c++:18050 sci.math.num-analysis:3602
  2. Newsgroups: comp.lang.c++,sci.math.num-analysis
  3. Path: sparky!uunet!munnari.oz.au!comp.vuw.ac.nz!kauri.vuw.ac.nz!robertd
  4. From: robertd@kauri.vuw.ac.nz (Robert Davies)
  5. Subject: Newmat06, a matrix package in C++
  6. Nntp-Posting-Host: kauri.vuw.ac.nz
  7. Message-ID: <BzCpGL.3M7@comp.vuw.ac.nz>
  8. Organization: Victoria University of Wellington
  9. Sender: news@comp.vuw.ac.nz (News Admin)
  10. Date: Wed, 16 Dec 1992 12:02:44 GMT
  11. Lines: 220
  12.  
  13. I have released a new version of my matrix package in C++. This version
  14. is called newmat06.
  15.  
  16. The package supports classes Matrix, UpperTriangularMatrix,
  17. LowerTriangularMatrix, DiagonalMatrix, SymmetricMatrix, BandMatrix,
  18. UpperBandMatrix, LowerBandMatrix, SymmetricBandMatrix, RowVector,
  19. ColumnVector.
  20.  
  21. Only one element type (float or double) is supported.
  22.  
  23. The package includes the operations *, +, -, inverse, transpose,
  24. conversion between types, submatrix, determinant, Cholesky
  25. decomposition, Householder triangularisation, singular value
  26. decomposition, symmetric eigenvalue analysis, fast Fourier transform,
  27. sort, print, an interface to "Numerical Recipes in C" programs, and an
  28. emulation of exceptions.
  29.  
  30. It is suitable for matrices in the size range 15 by 15 up to the maximum
  31. size that can conveniently be stored in a single block of memory.
  32.  
  33. The package works with AT&T C++ (versions 2.1 and 3), Gnu G++ (version
  34. 2.2), Borland C++ (version 3.1), Microsoft C++ (7.0) (with changes - see
  35. below), and mostly with Zortech (version 3.04).
  36.  
  37. It is available on Simtel20 (192.88.110.20) in the pd1:<msdos.cplusplus>
  38. directory in .zip format, and on the archives of comp.sources.misc,
  39. volume 34 issues 7 to 13 (shar format). Look for it on
  40. wuarchive.wustl.edu (128.252.135.4). It is also on Compuserve on the
  41. Borland bulletin board. And it is being distributed by the C users'
  42. group.
  43.  
  44. Here are some notes on the package.
  45.  
  46. The documentation says the package is good for matrices in the range 4 x
  47. 4 up to the maximum size that can be stored in a single block of memory.
  48. In fact the package becomes inefficient for matrices smaller than about
  49. 15 x 15 (the actual figure will depend on what you are doing and the
  50. speed of your floating point operations). The upper bound is set by the
  51. maximum size of an array and the maximum value of an "int".
  52.  
  53.  
  54. Gnu 2.3.x
  55. ---------
  56.  
  57. The latest version of Gnu picked up a couple of errors in the package.
  58. You should probably make these changes whatever compiler you are using.
  59.  
  60. In except.h change the destructor for class Janitor to
  61.  
  62.    virtual ~Janitor();
  63.  
  64. In newmat.h, class InvertedMatrix, delete the final const in the
  65. second definition of operator*
  66.  
  67. #ifndef TEMPS_DESTROYED_QUICKLY
  68.    SolvedMatrix operator*(const BaseMatrix&) const;  // inverse(A) * B
  69. #else
  70.    SolvedMatrix& operator*(const BaseMatrix&);       // inverse(A) * B
  71. #endif
  72.  
  73.  
  74. and make the corresponding changes in newmat6.cxx
  75.  
  76.  
  77. SolvedMatrix& InvertedMatrix::operator*(const BaseMatrix& bmx)
  78. {
  79.    REPORT
  80.    SolvedMatrix* x;
  81.    Try { x = new SolvedMatrix(bm, &bmx); MatrixErrorNoSpace(x); }
  82.    CatchAll { delete this; Throw(); }
  83.    delete this;                // since we are using bm rather than this
  84.    return *x;
  85. }
  86.  
  87.  
  88. I have not checked that G++ 2.3.x works as it is not yet loaded here.
  89. Please report any further problems to me. One correspondent reported a
  90. compiler error in tmte.cxx. Has anyone else had this error?
  91.  
  92.  
  93. Zortech 3.04
  94. ------------
  95.  
  96. 1. Replace controlw.h by the following
  97.  
  98.  
  99. //$$ controlw.h                Control word class
  100.  
  101. #ifndef CONTROL_WORD_LIB
  102. #define CONTROL_WORD_LIB 0
  103.  
  104. // for organising an int as a series of bits which indicate whether an
  105. // option is on or off.
  106.  
  107. class ControlWord
  108. {
  109. protected:
  110.    int cw;                                      // the control word
  111. public:
  112.    ControlWord() : cw(0) {}                     // do nothing
  113.    ControlWord(int i) : cw(i) {}                // load an integer
  114.  
  115.       // select specific bits (for testing at least one set)
  116.    ControlWord operator*(const ControlWord& i) const
  117.       { return ControlWord(cw & i.cw); }
  118.    void operator*=(const ControlWord& i)  { cw &= i.cw; }
  119.  
  120.       // set bits
  121.    ControlWord operator+(const ControlWord& i) const
  122.       { return ControlWord(cw | i.cw); }
  123.    void operator+=(const ControlWord& i)  { cw |= i.cw; }
  124.  
  125.       // reset bits
  126.    ControlWord operator-(const ControlWord& i) const
  127.       { return ControlWord(cw - (cw & i.cw)); }
  128.    void operator-=(const ControlWord& i) { cw -= (cw & i.cw); }
  129.  
  130.       // check if all of selected bits set or reset
  131.    Boolean operator>=(const ControlWord& i) const { return (cw & i.cw) == i.cw; }
  132.    Boolean operator<=(const ControlWord& i) const { return (cw & i.cw) == cw; }
  133.  
  134.       // flip selected bits
  135.    ControlWord operator^(const ControlWord& i) const
  136.       { return ControlWord(cw ^ i.cw); }
  137.    ControlWord operator~() const { return ControlWord(~cw); }
  138.  
  139.       // convert to integer
  140.    int operator+() const { return cw; }
  141.    int operator!() const { return cw==0; }
  142.    FREE_CHECK(ControlWord)
  143. };
  144.  
  145.  
  146. #endif
  147.  
  148.  
  149. 2. Do not #define TEMPS_DESTROYED_QUICKLY; do not #define Version21
  150.  
  151. 3. Do not use exceptions or the NRIC interface.
  152.  
  153. 4. Use the large model
  154.  
  155.  
  156. Microsoft C/C++, 7.0
  157. --------------------
  158.  
  159. Include the following block in include.h and delete the "const"s as for
  160. G++.
  161.  
  162. #define __MSC__
  163.  
  164.  
  165. #ifdef __MSC__                          // Microsoft
  166.    #include <stdlib.h>
  167.  
  168.    typedef int jmp_buf[9];
  169.    extern "C"
  170.    {
  171.       int __cdecl setjmp(jmp_buf);
  172.       void __cdecl longjmp(jmp_buf, int);
  173.    }
  174.  
  175.    #ifdef WANT_STREAM
  176.       #include <iostream.h>
  177.       #include <iomanip.h>
  178.    #endif
  179.    #ifdef WANT_MATH
  180.       #include <math.h>
  181.       #include <float.h>
  182.    #endif
  183. #endif
  184.  
  185.  
  186. In newmat.h you must include to following statement
  187.  
  188. void operator=(const Matrix& m) { operator=((const BaseMatrix&)m); }
  189.  
  190. in every class where operator=(const BaseMatrix&) is defined. "Matrix"
  191. is to be replaced by the class name of the class you are amending.
  192.  
  193. If you don't make this change Microsoft uses the default assignment
  194. operator. (What should it do, ARM doesn't make it clear?)
  195.  
  196. In jacobi.cxx, function void Jacobi(const SymmetricMatrix& X,
  197. DiagonalMatrix& D, SymmetricMatrix& A, Matrix& V, Boolean eivec), make
  198. the following changes:
  199.  
  200. Immediately after the initial { insert the line
  201.  
  202.    Real epsilon = FloatingPointPrecision::Epsilon();
  203.  
  204. change
  205.         if (i>4 && adp+g==adp && adq+g==adq) apq = 0.0;
  206. to
  207.         if (i>4 && g < epsilon*adp && g < epsilon*adq) apq = 0.0;
  208. and change
  209.            if (ah+g==ah) t = apq / h;
  210. to
  211.            if (g < epsilon*ah) t = apq / h;
  212.  
  213. Do #define TEMPS_DESTROYED_QUICKLY as the cleanup after an exception
  214. seems to work better.
  215.  
  216. Microsoft warns that setjmp and longjmp, which are used in the emulation
  217. of the exceptions, are not compatible with C++. So, be warned.
  218.  
  219.  
  220. AT&T (and others?)
  221. ------------------
  222.  
  223. If you are using this with interviews you may get a conflict with Catch.
  224. Either #undefine Catch or replace Catch with CATCH throughout my
  225. package.
  226.  
  227. I have assumed that math.h contains a function "abs". I understand I
  228. should be using iabs. So, if abs causes a problem, replace it with iabs.
  229. (In the HP unix workstation, the problem was the opposite, abs was
  230. defined in both math.h and stdlib.h. You need to delete it from one of
  231. these.)
  232.  
  233.