This is the mail archive of the
gsl-discuss@sources.redhat.com
mailing list for the GSL project.
Re: Matrix matters.
- From: Brian Gough <bjg at network-theory dot co dot uk>
- To: Toby White <tow at theor dot ch dot cam dot ac dot uk>
- Cc: gsl-discuss at sources dot redhat dot com
- Date: Mon, 18 Feb 2002 20:55:42 +0000 (GMT)
- Subject: Re: Matrix matters.
- References: <20020214194610.GA498@theor.ch.cam.ac.uk>
Toby White writes:
> I had cause to be playing around with the gsl_matrix sources, and
> the following things arose:
> 1) the gsl_matrix_*_swap function is prototyped as
> int gsl_matrix_swap(gsl_matrix * m1, const gsl_matrix * m2);
> I assume this should be
> int gsl_matrix_swap(gsl_matrix * m1, gsl_matrix * m2);
> (ie no const on the second matrix, since it will be altered.
Thanks, I've added this correction.
> 2) the arithmetic operations _scale, _add_constant, and
> _add_diagonal all take the constant as a double. Surely the
> constant should be in the same form as the matrix; ie, instead of
> int gsl_matrix_long_double_scale (gsl_matrix_long_double * a,
> const double x);
> it should be
> int gsl_matrix_long_double_scale (gsl_matrix_long_double * a,
> const long double x);
> and so on. This is a trivial change; diff at
> http://www-theor.ch.cam.ac.uk/people/tow/gsl/mat.constant.diff
Since the API is frozen I will have to leave it. I agree that it's
not ideal.
The motivation was that integer multiplication and division is not
very useful for scientific computing. Typically one wants to scale
integer arrays by a real factor, to adjust the dynamic range.
So I used double for that argument throughout, because the template
headers can only handle a fixed type for that.
> 3) the arithmetic operations were not implemented for complex
> matrices. This is relatively easy to do, and I've done it for my
> own use, at least. If you're interested, the code is at the same
> place as above, files mat.complex.diff and oper_complex_source.c
Thanks. I'd suggest in mul and div,
-- use temporary variables for the real/imag parts, it will make the
code easier to read (c.f. gsl_complex_divide), and move constants
outside the loop.
-- base the functions on the versions in complex/math.c, I think there
is a trick to prevent overflow in the division.
-- scale seems to have bug compared with mul (overwriting real part)
If you can do something like that I'll add them in.
Brian