nm10.htm

C++矩阵算法库

<TD width="25%">nl_ex.cpp</TD>
<TD width="50%">example of non-linear least squares</TD>
</TR>
<TR>
<TD width="25%">&nbsp;</TD>
<TD width="25%">nl_ex.txt</TD>
<TD width="50%">output from example</TD>
</TR>
<TR>
<TD width="25%">&nbsp;</TD>
<TD width="25%">garch.cpp</TD>
<TD width="50%">example of maximum-likelihood fit</TD>
</TR>
<TR>
<TD width="25%">&nbsp;</TD>
<TD width="25%">garch.dat</TD>
<TD width="50%">data file for garch.cpp</TD>
</TR>
<TR>
<TD width="25%">&nbsp;</TD>
<TD width="25%">garch.txt</TD>
<TD width="50%">output from example</TD>
</TR>
<TR>
<TD width="25%">&nbsp;</TD>
<TD width="25%">test_exc.cpp</TD>
<TD width="50%">demonstration exceptions</TD>
</TR>
<TR>
<TD width="25%">&nbsp;</TD>
<TD width="25%">test_exc.txt</TD>
<TD width="50%">output from test_exc.cpp</TD>
</TR>
<TR>
<TD width="25%"><B>Test files</B></TD>
<TD width="25%">tmt.h</TD>
<TD width="50%">header file for test files</TD>
</TR>
<TR>
<TD width="25%">&nbsp;</TD>
<TD width="25%">tmt*.cpp</TD>
<TD width="50%">test files (see the section on <A
HREF="#testing">testing</A>)</TD>
</TR>
<TR>
<TD width="25%">&nbsp;</TD>
<TD width="25%">tmt.txt</TD>
<TD width="50%">output from test files</TD>
</TR>
<TR>
<TD width="25%"><B>Make files</B></TD>
<TD width="25%">nm_gnu.mak</TD>
<TD width="50%">make file for Gnu G++</TD>
</TR>
<TR>
<TD width="25%">&nbsp;</TD>
<TD width="25%">nm_cc.mak</TD>
<TD width="50%">make file for AT&amp;T, Sun and HPUX</TD>
</TR>
<TR>
<TD width="25%">&nbsp;</TD>
<TD width="25%">nm_b55.mak</TD>
<TD width="50%">make file for Borland C++ 5.5</TD>
</TR>
<TR>
<TD width="25%">&nbsp;</TD>
<TD width="25%">newmat.lfl</TD>
<TD width="50%">library file list for use with <a href="genmake.htm">genmake</a></TD>
</TR>
<TR>
<TD width="25%">&nbsp;</TD>
<TD width="25%">nm_targ.txt</TD>
<TD width="50%">target file list for use with <a href="genmake.htm">genmake</a></TD>
</TR>
</TABLE>
<H2><A NAME="problem"></A>2.11 Problem report form
</H2>
<P CLASS="small"><A HREF="#refer">next</A> - <A HREF="#refer">skip</A> -
<A HREF="#starting">up</A> - <A HREF="#top">start</A></P>
<P>Copy and paste this to your editor; fill it out and email to
<B>robert (at) statsresearch.co.nz</B> </P>
<P>But first look in my web page <a href="http://www.robertnz.net/bugs">http://www.robertnz.net/bugs.htm</a> to see if the bug has
already been reported. </P>
<PRE> Version: ............... newmat10A (8 October 2002)
 Your email address: ....
 Today's date: ..........
 Your machine: ..........
 Operating system: ......
 Compiler &amp; version: ....
 Compiler options
   (eg GUI or console)...
 Describe the problem - attach examples if possible:





-----------------------------------------------------------
</PRE>

<H2><A NAME="refer"></A>3. Reference manual</H2>
<P CLASS="small"><A HREF="#constr">next</A> - <A HREF="#error">skip</A> -
<A HREF="#top">up</A> - <A HREF="#top">start</A></P>
<TABLE WIDTH="100%">
<TR>
<TD VALIGN="TOP" ALIGN="LEFT" WIDTH="50%"> <A HREF="#constr">3.1
Constructors </A><BR>
<A HREF="#elements">3.2 Accessing elements </A><BR>
<A HREF="#copy">3.3 Assignment and copying </A><BR>
<A HREF="#entering">3.4 Entering values </A><BR>
<A HREF="#unary">3.5 Unary operations </A><BR>
<A HREF="#binary">3.6 Binary operations </A><BR>
<A HREF="#matscal">3.7 Matrix and scalar ops </A><BR>
<A HREF="#scalar1">3.8 Scalar functions - size &amp; shape </A><BR>
<A HREF="#scalar2">3.9 Scalar functions - maximum &amp; minimum </A><BR>
<A HREF="#scalar3">3.10 Scalar functions - numerical </A><BR>
<A HREF="#submat">3.11 Submatrices </A><BR>
<A HREF="#dimen">3.12 Change dimension </A><BR>
<A HREF="#ch_type">3.13 Change type </A><BR>
<A HREF="#solve">3.14 Multiple matrix solve </A><BR>
<A HREF="#memory">3.15 Memory management </A><BR>
<A HREF="#efficien">3.16 Efficiency </A></TD>
<TD VALIGN="TOP" ALIGN="LEFT" WIDTH="50%">  <A
HREF="#output">3.17 Output </A><BR>
<A HREF="#unspec">3.18 Accessing unspecified type </A><BR>
<A HREF="#cholesky">3.19 Cholesky decomposition </A><BR>
<A HREF="#qr">3.20 QR decomposition </A><BR>
<A HREF="#svd">3.21 Singular value decomposition </A><BR>
<A HREF="#evalues">3.22 Eigenvalue decomposition </A><BR>
<A HREF="#sorting">3.23 Sorting </A><BR>
<A HREF="#fft">3.24 Fast Fourier transform </A><BR>
<A HREF="#trigtran">3.25 Fast trigonometric transforms </A><BR>
<A HREF="#nric">3.26 Numerical recipes in C </A><BR>
<A HREF="#except">3.27 Exceptions </A><BR>
<A HREF="#cleanup">3.28 Cleanup following exception </A><BR>
<A HREF="#nonlin">3.29 Non-linear applications </A><BR>
<A HREF="#stl">3.30 Standard template library </A><BR>
<A HREF="#namesp">3.31 Namespace </A></TD>
</TR>
</TABLE>
<H2><A NAME="constr"></A>3.1 Constructors</H2>
<P CLASS="small"><A HREF="#elements">next</A> - <A HREF="#elements">skip</A> -
<A HREF="#refer">up</A> - <A HREF="#top">start</A></P>
<P>To construct an <I>m</I> x <I>n</I> matrix, <TT>A</TT>, (<I>m</I> and
<I>n</I> are integers) use </P>
<PRE>    Matrix A(m,n);
</PRE>

<P>The UpperTriangularMatrix, LowerTriangularMatrix, SymmetricMatrix and
DiagonalMatrix types are square. To construct an <I>n</I> x <I>n</I> matrix
use, for example </P>
<PRE>    UpperTriangularMatrix UT(n);
    LowerTriangularMatrix LT(n);
    SymmetricMatrix S(n);
    DiagonalMatrix D(n);
</PRE>

<P>Band matrices need to include bandwidth information in their constructors. 
</P>
<PRE>    BandMatrix BM(n, lower, upper);
    UpperBandMatrix UB(n, upper);
    LowerBandMatrix LB(n, lower);
    SymmetricBandMatrix SB(n, lower);
</PRE>

<P>The integers <I>upper</I> and <I>lower</I> are the number of non-zero
diagonals above and below the diagonal (<I>excluding</I> the diagonal)
respectively. </P>
<P>The RowVector and ColumnVector types take just one argument in their
constructors: </P>
<PRE>    RowVector RV(n);
    ColumnVector CV(n);
</PRE>

<P><b>These constructors do <EM>not</EM> initialise the elements of the matrices.
</b>
To set all the elements to zero use, for example, </P>
<PRE>    Matrix A(m, n); A = 0.0;
</PRE>

<P>The IdentityMatrix takes one argument in its constructor specifying its 
dimension.</P>
<pre>    IdentityMatrix I(n);</pre>
<p>The value of the diagonal elements <b>is</b> set to 1 by default, but you can 
change this value as with other matrix types. </p>

<P>You can also construct vectors and matrices without specifying the
dimension. For example </P>
<PRE>    Matrix A;
</PRE>

<P>In this case the dimension must be set by an <A HREF="#copy">assignment
statement</A> or a <A HREF="#dimen">re-size statement</A>. </P>
<P>You can also use a constructor to set a matrix equal to another matrix or
matrix expression. </P>
<PRE>    Matrix A = UT;
    Matrix A = UT * LT;
</PRE>

<P>Only conversions that don't lose information are supported - eg you cannot
convert an upper triangular matrix into a diagonal matrix using =. </P>
<H2><A NAME="elements"></A>3.2 Accessing elements 
</H2>
<P CLASS="small"><A HREF="#copy">next</A> - <A HREF="#copy">skip</A> -
<A HREF="#refer">up</A> - <A HREF="#top">start</A></P>
<P>Elements are accessed by expressions of the form <TT>A(i,j)</TT> where
<I>i</I> and <I>j</I> run from 1 to the appropriate dimension. Access elements
of vectors with just one argument. Diagonal matrices can accept one or two
subscripts. </P>
<P>This is different from the earliest version of the package in which the
subscripts ran from 0 to one less than the appropriate dimension. Use
<TT>A.element(i,j)</TT> if you want this earlier convention. </P>
<P><TT>A(i,j)</TT> and <TT>A.element(i,j)</TT> can appear on either side of an
= sign. </P>
<P>If you activate the <TT>#define SETUP_C_SUBSCRIPTS</TT> in
<TT>include.h</TT> you can also access elements using the traditional C style
notation. That is <TT>A[i][j]</TT> for matrices (except diagonal) and
<TT>V[i]</TT> for vectors and diagonal matrices. The subscripts start at zero
(i.e. like element) and there is <I>no</I> range checking. Because of the
possibility of confusing <TT>V(i)</TT> and <TT>V[i]</TT>, I suggest you do
<I>not</I> activate this option unless you really want to use it. </P>
<P>Symmetric matrices are stored as lower triangular matrices. It is important
to remember this if you are using the <TT>A[i][j]</TT> method of accessing
elements. Make sure the first subscript is greater than or equal to the second
subscript. However, if you are using the <TT>A(i,j)</TT> method the program
will swap <TT>i</TT> and <TT>j</TT> if necessary; so it doesn't matter if you
think of the storage as being in the upper triangle (but it <I>does</I> matter
in some other situations such as when <A HREF="#entering">entering</A> data).</P>
<P>The IdentityMatrix type does not support element access. 
</P>
<H2><A NAME="copy"></A>3.3 Assignment and copying</H2>
<P CLASS="small"><A HREF="#entering">next</A> - <A HREF="#entering">skip</A> -
<A HREF="#refer">up</A> - <A HREF="#top">start</A></P>
<P>The operator <TT>=</TT> is used for copying matrices, converting matrices,
or evaluating expressions. For example </P>
<PRE>    A = B;  A = L;  A = L * U;
</PRE>

<P>Only conversions that don't lose information are supported. The dimensions
of the matrix on the left hand side are adjusted to those of the matrix or
expression on the right hand side. Elements on the right hand side which are
not present on the left hand side are set to zero. </P>
<P>The operator <TT>&lt;&lt;</TT> can be used in place of <TT>=</TT> where it
is permissible for information to be lost. </P>
<P>For example </P>
<PRE>    SymmetricMatrix S; Matrix A;
    ......
    S &lt;&lt; A.t() * A;
</PRE>

<P>is acceptable whereas </P>
<PRE>    S = A.t() * A;                            // error
</PRE>

<P>will cause a runtime error since the package does not (yet?) recognise
<TT>A.t()*A</TT> as symmetric. </P>
<P>Note that you can <I>not</I> use <TT>&lt;&lt;</TT> with constructors. For
example </P>
<PRE>    SymmetricMatrix S &lt;&lt; A.t() * A;           // error
</PRE>

<P>does <I>not</I> work. </P>
<P>Also note that <TT>&lt;&lt;</TT> cannot be used to load values from a full
matrix into a band matrix, since it will be unable to determine the bandwidth
of the band matrix. </P>
<P>A third copy routine is used in a similar role to <TT>=</TT>. Use </P>
<PRE>    A.Inject(D);
</PRE>

<P>to copy the elements of <TT>D</TT> to the corresponding elements of
<TT>A</TT> but leave the elements of <TT>A</TT> unchanged if there is no
corresponding element of <TT>D</TT> (the <TT>=</TT> operator would set them to
0). This is useful, for example, for setting the diagonal elements of a matrix
without disturbing the rest of the matrix. Unlike <TT>=</TT> and
<TT>&lt;&lt;</TT>, Inject does not reset the dimensions of <TT>A</TT>, which
must match those of <TT>D</TT>. Inject does <I>not</I> test for no loss of
information. </P>
<P>You cannot replace <TT>D</TT> by a matrix expressio