3
<!-- This HTML file has been created by texi2html 1.54+ (gsl)
4
from ../gsl-ref.texi -->
6
<TITLE>GNU Scientific Library -- Reference Manual - Vectors and Matrices</TITLE>
7
<!-- <LINK rel="stylesheet" title="Default Style Sheet" href="/css/texinfo.css" type="text/css"> -->
8
<link href="gsl-ref_9.html" rel=Next>
9
<link href="gsl-ref_7.html" rel=Previous>
10
<link href="gsl-ref_toc.html" rel=ToC>
14
<p>Go to the <A HREF="gsl-ref_1.html">first</A>, <A HREF="gsl-ref_7.html">previous</A>, <A HREF="gsl-ref_9.html">next</A>, <A HREF="gsl-ref_50.html">last</A> section, <A HREF="gsl-ref_toc.html">table of contents</A>.
18
<H1><A NAME="SEC152" HREF="gsl-ref_toc.html#TOC152">Vectors and Matrices</A></H1>
26
The functions described in this chapter provide a simple vector and
27
matrix interface to ordinary C arrays. The memory management of these
28
arrays is implemented using a single underlying type, known as a
29
block. By writing your functions in terms of vectors and matrices you
30
can pass a single structure containing both data and dimensions as an
31
argument without needing additional function parameters. The structures
32
are compatible with the vector and matrix formats used by BLAS
39
<H2><A NAME="SEC153" HREF="gsl-ref_toc.html#TOC153">Data types</A></H2>
42
All the functions are available for each of the standard data-types.
43
The versions for <CODE>double</CODE> have the prefix <CODE>gsl_block</CODE>,
44
<CODE>gsl_vector</CODE> and <CODE>gsl_matrix</CODE>. Similarly the versions for
45
single-precision <CODE>float</CODE> arrays have the prefix
46
<CODE>gsl_block_float</CODE>, <CODE>gsl_vector_float</CODE> and
47
<CODE>gsl_matrix_float</CODE>. The full list of available types is given
55
gsl_block_long_double long double
57
gsl_block_uint unsigned int
59
gsl_block_ulong unsigned long
61
gsl_block_ushort unsigned short
63
gsl_block_uchar unsigned char
64
gsl_block_complex complex double
65
gsl_block_complex_float complex float
66
gsl_block_complex_long_double complex long double
70
Corresponding types exist for the <CODE>gsl_vector</CODE> and
71
<CODE>gsl_matrix</CODE> functions.
77
<H2><A NAME="SEC154" HREF="gsl-ref_toc.html#TOC154">Blocks</A></H2>
80
For consistency all memory is allocated through a <CODE>gsl_block</CODE>
81
structure. The structure contains two components, the size of an area of
82
memory and a pointer to the memory. The <CODE>gsl_block</CODE> structure looks
96
Vectors and matrices are made by <I>slicing</I> an underlying block. A
97
slice is a set of elements formed from an initial offset and a
98
combination of indices and step-sizes. In the case of a matrix the
99
step-size for the column index represents the row-length. The step-size
100
for a vector is known as the <I>stride</I>.
104
The functions for allocating and deallocating blocks are defined in
105
<TT>'gsl_block.h'</TT>
111
<H3><A NAME="SEC155" HREF="gsl-ref_toc.html#TOC155">Block allocation</A></H3>
114
The functions for allocating memory to a block follow the style of
115
<CODE>malloc</CODE> and <CODE>free</CODE>. In addition they also perform their own
116
error checking. If there is insufficient memory available to allocate a
117
block then the functions call the GSL error handler (with an error
118
number of <CODE>GSL_ENOMEM</CODE>) in addition to returning a null
119
pointer. Thus if you use the library error handler to abort your program
120
then it isn't necessary to check every <CODE>alloc</CODE>.
125
<DT><U>Function:</U> gsl_block * <B>gsl_block_alloc</B> <I>(size_t <VAR>n</VAR>)</I>
126
<DD><A NAME="IDX808"></A>
127
This function allocates memory for a block of <VAR>n</VAR> double-precision
128
elements, returning a pointer to the block struct. The block is not
129
initialized and so the values of its elements are undefined. Use the
130
function <CODE>gsl_block_calloc</CODE> if you want to ensure that all the
131
elements are initialized to zero.
135
A null pointer is returned if insufficient memory is available to create
142
<DT><U>Function:</U> gsl_block * <B>gsl_block_calloc</B> <I>(size_t <VAR>n</VAR>)</I>
143
<DD><A NAME="IDX809"></A>
144
This function allocates memory for a block and initializes all the
145
elements of the block to zero.
151
<DT><U>Function:</U> void <B>gsl_block_free</B> <I>(gsl_block * <VAR>b</VAR>)</I>
152
<DD><A NAME="IDX810"></A>
153
This function frees the memory used by a block <VAR>b</VAR> previously
154
allocated with <CODE>gsl_block_alloc</CODE> or <CODE>gsl_block_calloc</CODE>.
160
<H3><A NAME="SEC156" HREF="gsl-ref_toc.html#TOC156">Reading and writing blocks</A></H3>
163
The library provides functions for reading and writing blocks to a file
164
as binary data or formatted text.
169
<DT><U>Function:</U> int <B>gsl_block_fwrite</B> <I>(FILE * <VAR>stream</VAR>, const gsl_block * <VAR>b</VAR>)</I>
170
<DD><A NAME="IDX811"></A>
171
This function writes the elements of the block <VAR>b</VAR> to the stream
172
<VAR>stream</VAR> in binary format. The return value is 0 for success and
173
<CODE>GSL_EFAILED</CODE> if there was a problem writing to the file. Since the
174
data is written in the native binary format it may not be portable
175
between different architectures.
181
<DT><U>Function:</U> int <B>gsl_block_fread</B> <I>(FILE * <VAR>stream</VAR>, gsl_block * <VAR>b</VAR>)</I>
182
<DD><A NAME="IDX812"></A>
183
This function reads into the block <VAR>b</VAR> from the open stream
184
<VAR>stream</VAR> in binary format. The block <VAR>b</VAR> must be preallocated
185
with the correct length since the function uses the size of <VAR>b</VAR> to
186
determine how many bytes to read. The return value is 0 for success and
187
<CODE>GSL_EFAILED</CODE> if there was a problem reading from the file. The
188
data is assumed to have been written in the native binary format on the
195
<DT><U>Function:</U> int <B>gsl_block_fprintf</B> <I>(FILE * <VAR>stream</VAR>, const gsl_block * <VAR>b</VAR>, const char * <VAR>format</VAR>)</I>
196
<DD><A NAME="IDX813"></A>
197
This function writes the elements of the block <VAR>b</VAR> line-by-line to
198
the stream <VAR>stream</VAR> using the format specifier <VAR>format</VAR>, which
199
should be one of the <CODE>%g</CODE>, <CODE>%e</CODE> or <CODE>%f</CODE> formats for
200
floating point numbers and <CODE>%d</CODE> for integers. The function returns
201
0 for success and <CODE>GSL_EFAILED</CODE> if there was a problem writing to
208
<DT><U>Function:</U> int <B>gsl_block_fscanf</B> <I>(FILE * <VAR>stream</VAR>, gsl_block * <VAR>b</VAR>)</I>
209
<DD><A NAME="IDX814"></A>
210
This function reads formatted data from the stream <VAR>stream</VAR> into the
211
block <VAR>b</VAR>. The block <VAR>b</VAR> must be preallocated with the correct
212
length since the function uses the size of <VAR>b</VAR> to determine how many
213
numbers to read. The function returns 0 for success and
214
<CODE>GSL_EFAILED</CODE> if there was a problem reading from the file.
220
<H3><A NAME="SEC157" HREF="gsl-ref_toc.html#TOC157">Example programs for blocks</A></H3>
223
The following program shows how to allocate a block,
227
<PRE class="example">
228
#include <stdio.h>
229
#include <gsl/gsl_block.h>
234
gsl_block * b = gsl_block_alloc (100);
236
printf ("length of block = %u\n", b->size);
237
printf ("block data address = %#x\n", b->data);
245
Here is the output from the program,
249
<PRE class="example">
250
length of block = 100
251
block data address = 0x804b0d8
256
<H2><A NAME="SEC158" HREF="gsl-ref_toc.html#TOC158">Vectors</A></H2>
258
<A NAME="IDX815"></A>
259
<A NAME="IDX816"></A>
263
Vectors are defined by a <CODE>gsl_vector</CODE> structure which describes a
264
slice of a block. Different vectors can be created which point to the
265
same block. A vector slice is a set of equally-spaced elements of an
270
The <CODE>gsl_vector</CODE> structure contains five components, the
271
<I>size</I>, the <I>stride</I>, a pointer to the memory where the elements
272
are stored, <VAR>data</VAR>, a pointer to the block owned by the vector,
273
<VAR>block</VAR>, if any, and an ownership flag, <VAR>owner</VAR>. The structure
274
is very simple and looks like this,
278
<PRE class="example">
290
The <VAR>size</VAR> is simply the number of vector elements. The range of
291
valid indices runs from 0 to <CODE>size-1</CODE>. The <VAR>stride</VAR> is the
292
step-size from one element to the next in physical memory, measured in
293
units of the appropriate datatype. The pointer <VAR>data</VAR> gives the
294
location of the first element of the vector in memory. The pointer
295
<VAR>block</VAR> stores the location of the memory block in which the vector
296
elements are located (if any). If the vector owns this block then the
297
<VAR>owner</VAR> field is set to one and the block will be deallocated when the
298
vector is freed. If the vector points to a block owned by another
299
object then the <VAR>owner</VAR> field is zero and any underlying block will not be
304
The functions for allocating and accessing vectors are defined in
305
<TT>'gsl_vector.h'</TT>
311
<H3><A NAME="SEC159" HREF="gsl-ref_toc.html#TOC159">Vector allocation</A></H3>
314
The functions for allocating memory to a vector follow the style of
315
<CODE>malloc</CODE> and <CODE>free</CODE>. In addition they also perform their own
316
error checking. If there is insufficient memory available to allocate a
317
vector then the functions call the GSL error handler (with an error
318
number of <CODE>GSL_ENOMEM</CODE>) in addition to returning a null
319
pointer. Thus if you use the library error handler to abort your program
320
then it isn't necessary to check every <CODE>alloc</CODE>.
325
<DT><U>Function:</U> gsl_vector * <B>gsl_vector_alloc</B> <I>(size_t <VAR>n</VAR>)</I>
326
<DD><A NAME="IDX817"></A>
327
This function creates a vector of length <VAR>n</VAR>, returning a pointer to
328
a newly initialized vector struct. A new block is allocated for the
329
elements of the vector, and stored in the <VAR>block</VAR> component of the
330
vector struct. The block is "owned" by the vector, and will be
331
deallocated when the vector is deallocated.
337
<DT><U>Function:</U> gsl_vector * <B>gsl_vector_calloc</B> <I>(size_t <VAR>n</VAR>)</I>
338
<DD><A NAME="IDX818"></A>
339
This function allocates memory for a vector of length <VAR>n</VAR> and
340
initializes all the elements of the vector to zero.
346
<DT><U>Function:</U> void <B>gsl_vector_free</B> <I>(gsl_vector * <VAR>v</VAR>)</I>
347
<DD><A NAME="IDX819"></A>
348
This function frees a previously allocated vector <VAR>v</VAR>. If the
349
vector was created using <CODE>gsl_vector_alloc</CODE> then the block
350
underlying the vector will also be deallocated. If the vector has been
351
created from another object then the memory is still owned by that
352
object and will not be deallocated.
358
<H3><A NAME="SEC160" HREF="gsl-ref_toc.html#TOC160">Accessing vector elements</A></H3>
360
<A NAME="IDX820"></A>
361
<A NAME="IDX821"></A>
362
<A NAME="IDX822"></A>
363
<A NAME="IDX823"></A>
364
<A NAME="IDX824"></A>
368
Unlike FORTRAN compilers, C compilers do not usually provide
369
support for range checking of vectors and matrices. Range checking is
370
available in the GNU C Compiler extension <CODE>checkergcc</CODE> but it is
371
not available on every platform. The functions <CODE>gsl_vector_get</CODE>
372
and <CODE>gsl_vector_set</CODE> can perform portable range checking for you
373
and report an error if you attempt to access elements outside the
378
The functions for accessing the elements of a vector or matrix are
379
defined in <TT>'gsl_vector.h'</TT> and declared <CODE>extern inline</CODE> to
380
eliminate function-call overhead. You must compile your program with
381
the macro <CODE>HAVE_INLINE</CODE> defined to use these functions. If
382
necessary you can turn off range checking completely without modifying
383
any source files by recompiling your program with the preprocessor
384
definition <CODE>GSL_RANGE_CHECK_OFF</CODE>. Provided your compiler supports
385
inline functions the effect of turning off range checking is to replace
386
calls to <CODE>gsl_vector_get(v,i)</CODE> by <CODE>v->data[i*v->stride]</CODE> and
387
calls to <CODE>gsl_vector_set(v,i,x)</CODE> by <CODE>v->data[i*v->stride]=x</CODE>.
388
Thus there should be no performance penalty for using the range checking
389
functions when range checking is turned off.
394
<DT><U>Function:</U> double <B>gsl_vector_get</B> <I>(const gsl_vector * <VAR>v</VAR>, size_t <VAR>i</VAR>)</I>
395
<DD><A NAME="IDX825"></A>
396
This function returns the <VAR>i</VAR>-th element of a vector <VAR>v</VAR>. If
397
<VAR>i</VAR> lies outside the allowed range of 0 to <VAR>n-1</VAR> then the error
398
handler is invoked and 0 is returned.
404
<DT><U>Function:</U> void <B>gsl_vector_set</B> <I>(gsl_vector * <VAR>v</VAR>, size_t <VAR>i</VAR>, double <VAR>x</VAR>)</I>
405
<DD><A NAME="IDX826"></A>
406
This function sets the value of the <VAR>i</VAR>-th element of a vector
407
<VAR>v</VAR> to <VAR>x</VAR>. If <VAR>i</VAR> lies outside the allowed range of 0 to
408
<VAR>n-1</VAR> then the error handler is invoked.
414
<DT><U>Function:</U> double * <B>gsl_vector_ptr</B> <I>(gsl_vector * <VAR>v</VAR>, size_t <VAR>i</VAR>)</I>
415
<DD><A NAME="IDX827"></A>
416
<DT><U>Function:</U> const double * <B>gsl_vector_const_ptr</B> <I>(const gsl_vector * <VAR>v</VAR>, size_t <VAR>i</VAR>)</I>
417
<DD><A NAME="IDX828"></A>
418
These functions return a pointer to the <VAR>i</VAR>-th element of a vector
419
<VAR>v</VAR>. If <VAR>i</VAR> lies outside the allowed range of 0 to <VAR>n-1</VAR>
420
then the error handler is invoked and a null pointer is returned.
426
<H3><A NAME="SEC161" HREF="gsl-ref_toc.html#TOC161">Initializing vector elements</A></H3>
428
<A NAME="IDX829"></A>
429
<A NAME="IDX830"></A>
434
<DT><U>Function:</U> void <B>gsl_vector_set_all</B> <I>(gsl_vector * <VAR>v</VAR>, double <VAR>x</VAR>)</I>
435
<DD><A NAME="IDX831"></A>
436
This function sets all the elements of the vector <VAR>v</VAR> to the value
443
<DT><U>Function:</U> void <B>gsl_vector_set_zero</B> <I>(gsl_vector * <VAR>v</VAR>)</I>
444
<DD><A NAME="IDX832"></A>
445
This function sets all the elements of the vector <VAR>v</VAR> to zero.
451
<DT><U>Function:</U> int <B>gsl_vector_set_basis</B> <I>(gsl_vector * <VAR>v</VAR>, size_t <VAR>i</VAR>)</I>
452
<DD><A NAME="IDX833"></A>
453
This function makes a basis vector by setting all the elements of the
454
vector <VAR>v</VAR> to zero except for the <VAR>i</VAR>-th element which is set to
461
<H3><A NAME="SEC162" HREF="gsl-ref_toc.html#TOC162">Reading and writing vectors</A></H3>
464
The library provides functions for reading and writing vectors to a file
465
as binary data or formatted text.
470
<DT><U>Function:</U> int <B>gsl_vector_fwrite</B> <I>(FILE * <VAR>stream</VAR>, const gsl_vector * <VAR>v</VAR>)</I>
471
<DD><A NAME="IDX834"></A>
472
This function writes the elements of the vector <VAR>v</VAR> to the stream
473
<VAR>stream</VAR> in binary format. The return value is 0 for success and
474
<CODE>GSL_EFAILED</CODE> if there was a problem writing to the file. Since the
475
data is written in the native binary format it may not be portable
476
between different architectures.
482
<DT><U>Function:</U> int <B>gsl_vector_fread</B> <I>(FILE * <VAR>stream</VAR>, gsl_vector * <VAR>v</VAR>)</I>
483
<DD><A NAME="IDX835"></A>
484
This function reads into the vector <VAR>v</VAR> from the open stream
485
<VAR>stream</VAR> in binary format. The vector <VAR>v</VAR> must be preallocated
486
with the correct length since the function uses the size of <VAR>v</VAR> to
487
determine how many bytes to read. The return value is 0 for success and
488
<CODE>GSL_EFAILED</CODE> if there was a problem reading from the file. The
489
data is assumed to have been written in the native binary format on the
496
<DT><U>Function:</U> int <B>gsl_vector_fprintf</B> <I>(FILE * <VAR>stream</VAR>, const gsl_vector * <VAR>v</VAR>, const char * <VAR>format</VAR>)</I>
497
<DD><A NAME="IDX836"></A>
498
This function writes the elements of the vector <VAR>v</VAR> line-by-line to
499
the stream <VAR>stream</VAR> using the format specifier <VAR>format</VAR>, which
500
should be one of the <CODE>%g</CODE>, <CODE>%e</CODE> or <CODE>%f</CODE> formats for
501
floating point numbers and <CODE>%d</CODE> for integers. The function returns
502
0 for success and <CODE>GSL_EFAILED</CODE> if there was a problem writing to
509
<DT><U>Function:</U> int <B>gsl_vector_fscanf</B> <I>(FILE * <VAR>stream</VAR>, gsl_vector * <VAR>v</VAR>)</I>
510
<DD><A NAME="IDX837"></A>
511
This function reads formatted data from the stream <VAR>stream</VAR> into the
512
vector <VAR>v</VAR>. The vector <VAR>v</VAR> must be preallocated with the correct
513
length since the function uses the size of <VAR>v</VAR> to determine how many
514
numbers to read. The function returns 0 for success and
515
<CODE>GSL_EFAILED</CODE> if there was a problem reading from the file.
521
<H3><A NAME="SEC163" HREF="gsl-ref_toc.html#TOC163">Vector views</A></H3>
524
In addition to creating vectors from slices of blocks it is also
525
possible to slice vectors and create vector views. For example, a
526
subvector of another vector can be described with a view, or two views
527
can be made which provide access to the even and odd elements of a
532
A vector view is a temporary object, stored on the stack, which can be
533
used to operate on a subset of vector elements. Vector views can be
534
defined for both constant and non-constant vectors, using separate types
535
that preserve constness. A vector view has the type
536
<CODE>gsl_vector_view</CODE> and a constant vector view has the type
537
<CODE>gsl_vector_const_view</CODE>. In both cases the elements of the view
538
can be accessed as a <CODE>gsl_vector</CODE> using the <CODE>vector</CODE> component
539
of the view object. A pointer to a vector of type <CODE>gsl_vector *</CODE>
540
or <CODE>const gsl_vector *</CODE> can be obtained by taking the address of
541
this component with the <CODE>&</CODE> operator.
545
When using this pointer, it is important to ensure that the view itself
546
remains in scope--the simplest way to do so is by always writing the
547
pointer as <CODE>&</CODE><VAR>view</VAR><CODE>.vector</CODE>, and not storing this value
548
in another pointer variable.
553
<DT><U>Function:</U> gsl_vector_view <B>gsl_vector_subvector</B> <I>(gsl_vector *<VAR>v</VAR>, size_t <VAR>offset</VAR>, size_t <VAR>n</VAR>)</I>
554
<DD><A NAME="IDX838"></A>
555
<DT><U>Function:</U> gsl_vector_const_view <B>gsl_vector_const_subvector</B> <I>(const gsl_vector * <VAR>v</VAR>, size_t <VAR>offset</VAR>, size_t <VAR>n</VAR>)</I>
556
<DD><A NAME="IDX839"></A>
557
These functions return a vector view of a subvector of another vector
558
<VAR>v</VAR>. The start of the new vector is offset by <VAR>offset</VAR> elements
559
from the start of the original vector. The new vector has <VAR>n</VAR>
560
elements. Mathematically, the <VAR>i</VAR>-th element of the new vector
561
<VAR>v'</VAR> is given by,
565
<PRE class="example">
566
v'(i) = v->data[(offset + i)*v->stride]
570
where the index <VAR>i</VAR> runs from 0 to <CODE>n-1</CODE>.
574
The <CODE>data</CODE> pointer of the returned vector struct is set to null if
575
the combined parameters (<VAR>offset</VAR>,<VAR>n</VAR>) overrun the end of the
580
The new vector is only a view of the block underlying the original
581
vector, <VAR>v</VAR>. The block containing the elements of <VAR>v</VAR> is not
582
owned by the new vector. When the view goes out of scope the original
583
vector <VAR>v</VAR> and its block will continue to exist. The original
584
memory can only be deallocated by freeing the original vector. Of
585
course, the original vector should not be deallocated while the view is
590
The function <CODE>gsl_vector_const_subvector</CODE> is equivalent to
591
<CODE>gsl_vector_subvector</CODE> but can be used for vectors which are
592
declared <CODE>const</CODE>.
599
<DT><U>Function:</U> gsl_vector_view <B>gsl_vector_subvector_with_stride</B> <I>(gsl_vector *<VAR>v</VAR>, size_t <VAR>offset</VAR>, size_t <VAR>stride</VAR>, size_t <VAR>n</VAR>)</I>
600
<DD><A NAME="IDX840"></A>
601
<DT><U>Function:</U> gsl_vector_const_view <B>gsl_vector_const_subvector_with_stride</B> <I>(const gsl_vector * <VAR>v</VAR>, size_t <VAR>offset</VAR>, size_t <VAR>stride</VAR>, size_t <VAR>n</VAR>)</I>
602
<DD><A NAME="IDX841"></A>
603
These functions return a vector view of a subvector of another vector
604
<VAR>v</VAR> with an additional stride argument. The subvector is formed in
605
the same way as for <CODE>gsl_vector_subvector</CODE> but the new vector has
606
<VAR>n</VAR> elements with a step-size of <VAR>stride</VAR> from one element to
607
the next in the original vector. Mathematically, the <VAR>i</VAR>-th element
608
of the new vector <VAR>v'</VAR> is given by,
612
<PRE class="example">
613
v'(i) = v->data[(offset + i*stride)*v->stride]
617
where the index <VAR>i</VAR> runs from 0 to <CODE>n-1</CODE>.
621
Note that subvector views give direct access to the underlying elements
622
of the original vector. For example, the following code will zero the
623
even elements of the vector <CODE>v</CODE> of length <CODE>n</CODE>, while leaving the
624
odd elements untouched,
628
<PRE class="example">
629
gsl_vector_view v_even
630
= gsl_vector_subvector_with_stride (v, 0, 2, n/2);
631
gsl_vector_set_zero (&v_even.vector);
635
A vector view can be passed to any subroutine which takes a vector
636
argument just as a directly allocated vector would be, using
637
<CODE>&</CODE><VAR>view</VAR><CODE>.vector</CODE>. For example, the following code
638
computes the norm of odd elements of <CODE>v</CODE> using the BLAS
643
<PRE class="example">
644
gsl_vector_view v_odd
645
= gsl_vector_subvector_with_stride (v, 1, 2, n/2);
646
double r = gsl_blas_dnrm2 (&v_odd.vector);
650
The function <CODE>gsl_vector_const_subvector_with_stride</CODE> is equivalent
651
to <CODE>gsl_vector_subvector_with_stride</CODE> but can be used for vectors
652
which are declared <CODE>const</CODE>.
658
<DT><U>Function:</U> gsl_vector_view <B>gsl_vector_complex_real</B> <I>(gsl_vector_complex *<VAR>v</VAR>)</I>
659
<DD><A NAME="IDX842"></A>
660
<DT><U>Function:</U> gsl_vector_const_view <B>gsl_vector_complex_const_real</B> <I>(const gsl_vector_complex *<VAR>v</VAR>)</I>
661
<DD><A NAME="IDX843"></A>
662
These functions return a vector view of the real parts of the complex
667
The function <CODE>gsl_vector_complex_const_real</CODE> is equivalent to
668
<CODE>gsl_vector_complex_real</CODE> but can be used for vectors which are
669
declared <CODE>const</CODE>.
675
<DT><U>Function:</U> gsl_vector_view <B>gsl_vector_complex_imag</B> <I>(gsl_vector_complex *<VAR>v</VAR>)</I>
676
<DD><A NAME="IDX844"></A>
677
<DT><U>Function:</U> gsl_vector_const_view <B>gsl_vector_complex_const_imag</B> <I>(const gsl_vector_complex *<VAR>v</VAR>)</I>
678
<DD><A NAME="IDX845"></A>
679
These functions return a vector view of the imaginary parts of the
680
complex vector <VAR>v</VAR>.
684
The function <CODE>gsl_vector_complex_const_imag</CODE> is equivalent to
685
<CODE>gsl_vector_complex_imag</CODE> but can be used for vectors which are
686
declared <CODE>const</CODE>.
692
<DT><U>Function:</U> gsl_vector_view <B>gsl_vector_view_array</B> <I>(double *<VAR>base</VAR>, size_t <VAR>n</VAR>)</I>
693
<DD><A NAME="IDX846"></A>
694
<DT><U>Function:</U> gsl_vector_const_view <B>gsl_vector_const_view_array</B> <I>(const double *<VAR>base</VAR>, size_t <VAR>n</VAR>)</I>
695
<DD><A NAME="IDX847"></A>
696
These functions return a vector view of an array. The start of the new
697
vector is given by <VAR>base</VAR> and has <VAR>n</VAR> elements. Mathematically,
698
the <VAR>i</VAR>-th element of the new vector <VAR>v'</VAR> is given by,
702
<PRE class="example">
707
where the index <VAR>i</VAR> runs from 0 to <CODE>n-1</CODE>.
711
The array containing the elements of <VAR>v</VAR> is not owned by the new
712
vector view. When the view goes out of scope the original array will
713
continue to exist. The original memory can only be deallocated by
714
freeing the original pointer <VAR>base</VAR>. Of course, the original array
715
should not be deallocated while the view is still in use.
719
The function <CODE>gsl_vector_const_view_array</CODE> is equivalent to
720
<CODE>gsl_vector_view_array</CODE> but can be used for arrays which are
721
declared <CODE>const</CODE>.
727
<DT><U>Function:</U> gsl_vector_view <B>gsl_vector_view_array_with_stride</B> <I>(double * <VAR>base</VAR>, size_t <VAR>stride</VAR>, size_t <VAR>n</VAR>)</I>
728
<DD><A NAME="IDX848"></A>
729
<DT><U>Function:</U> gsl_vector_const_view <B>gsl_vector_const_view_array_with_stride</B> <I>(const double * <VAR>base</VAR>, size_t <VAR>stride</VAR>, size_t <VAR>n</VAR>)</I>
730
<DD><A NAME="IDX849"></A>
731
These functions return a vector view of an array <VAR>base</VAR> with an
732
additional stride argument. The subvector is formed in the same way as
733
for <CODE>gsl_vector_view_array</CODE> but the new vector has <VAR>n</VAR> elements
734
with a step-size of <VAR>stride</VAR> from one element to the next in the
735
original array. Mathematically, the <VAR>i</VAR>-th element of the new
736
vector <VAR>v'</VAR> is given by,
740
<PRE class="example">
741
v'(i) = base[i*stride]
745
where the index <VAR>i</VAR> runs from 0 to <CODE>n-1</CODE>.
749
Note that the view gives direct access to the underlying elements of the
750
original array. A vector view can be passed to any subroutine which
751
takes a vector argument just as a directly allocated vector would be,
752
using <CODE>&</CODE><VAR>view</VAR><CODE>.vector</CODE>.
756
The function <CODE>gsl_vector_const_view_array_with_stride</CODE> is
757
equivalent to <CODE>gsl_vector_view_array_with_stride</CODE> but can be used
758
for arrays which are declared <CODE>const</CODE>.
765
<H3><A NAME="SEC164" HREF="gsl-ref_toc.html#TOC164">Copying vectors</A></H3>
768
Common operations on vectors such as addition and multiplication are
769
available in the BLAS part of the library (see section <A HREF="gsl-ref_12.html#SEC212">BLAS Support</A>). However, it is useful to have a small number of utility
770
functions which do not require the full BLAS code. The following
771
functions fall into this category.
776
<DT><U>Function:</U> int <B>gsl_vector_memcpy</B> <I>(gsl_vector * <VAR>dest</VAR>, const gsl_vector * <VAR>src</VAR>)</I>
777
<DD><A NAME="IDX850"></A>
778
This function copies the elements of the vector <VAR>src</VAR> into the
779
vector <VAR>dest</VAR>. The two vectors must have the same length.
785
<DT><U>Function:</U> int <B>gsl_vector_swap</B> <I>(gsl_vector * <VAR>v</VAR>, gsl_vector * <VAR>w</VAR>)</I>
786
<DD><A NAME="IDX851"></A>
787
This function exchanges the elements of the vectors <VAR>v</VAR> and <VAR>w</VAR>
788
by copying. The two vectors must have the same length.
795
<H3><A NAME="SEC165" HREF="gsl-ref_toc.html#TOC165">Exchanging elements</A></H3>
798
The following function can be used to exchange, or permute, the elements
804
<DT><U>Function:</U> int <B>gsl_vector_swap_elements</B> <I>(gsl_vector * <VAR>v</VAR>, size_t <VAR>i</VAR>, size_t <VAR>j</VAR>)</I>
805
<DD><A NAME="IDX852"></A>
806
This function exchanges the <VAR>i</VAR>-th and <VAR>j</VAR>-th elements of the
807
vector <VAR>v</VAR> in-place.
813
<DT><U>Function:</U> int <B>gsl_vector_reverse</B> <I>(gsl_vector * <VAR>v</VAR>)</I>
814
<DD><A NAME="IDX853"></A>
815
This function reverses the order of the elements of the vector <VAR>v</VAR>.
822
<H3><A NAME="SEC166" HREF="gsl-ref_toc.html#TOC166">Vector operations</A></H3>
825
The following operations are only defined for real vectors.
830
<DT><U>Function:</U> int <B>gsl_vector_add</B> <I>(gsl_vector * <VAR>a</VAR>, const gsl_vector * <VAR>b</VAR>)</I>
831
<DD><A NAME="IDX854"></A>
832
This function adds the elements of vector <VAR>b</VAR> to the elements of
833
vector <VAR>a</VAR>, a'_i = a_i + b_i. The two vectors must have the
840
<DT><U>Function:</U> int <B>gsl_vector_sub</B> <I>(gsl_vector * <VAR>a</VAR>, const gsl_vector * <VAR>b</VAR>)</I>
841
<DD><A NAME="IDX855"></A>
842
This function subtracts the elements of vector <VAR>b</VAR> from the elements of
843
vector <VAR>a</VAR>, a'_i = a_i - b_i. The two vectors must have the
850
<DT><U>Function:</U> int <B>gsl_vector_mul</B> <I>(gsl_vector * <VAR>a</VAR>, const gsl_vector * <VAR>b</VAR>)</I>
851
<DD><A NAME="IDX856"></A>
852
This function multiplies the elements of vector <VAR>a</VAR> by the elements of
853
vector <VAR>b</VAR>, a'_i = a_i * b_i. The two vectors must have the
860
<DT><U>Function:</U> int <B>gsl_vector_div</B> <I>(gsl_vector * <VAR>a</VAR>, const gsl_vector * <VAR>b</VAR>)</I>
861
<DD><A NAME="IDX857"></A>
862
This function divides the elements of vector <VAR>a</VAR> by the elements of
863
vector <VAR>b</VAR>, a'_i = a_i / b_i. The two vectors must have the
870
<DT><U>Function:</U> int <B>gsl_vector_scale</B> <I>(gsl_vector * <VAR>a</VAR>, const double <VAR>x</VAR>)</I>
871
<DD><A NAME="IDX858"></A>
872
This function multiplies the elements of vector <VAR>a</VAR> by the constant
873
factor <VAR>x</VAR>, a'_i = x a_i.
879
<DT><U>Function:</U> int <B>gsl_vector_add_constant</B> <I>(gsl_vector * <VAR>a</VAR>, const double <VAR>x</VAR>)</I>
880
<DD><A NAME="IDX859"></A>
881
This function adds the constant value <VAR>x</VAR> to the elements of the
882
vector <VAR>a</VAR>, a'_i = a_i + x.
888
<H3><A NAME="SEC167" HREF="gsl-ref_toc.html#TOC167">Finding maximum and minimum elements of vectors</A></H3>
892
<DT><U>Function:</U> double <B>gsl_vector_max</B> <I>(const gsl_vector * <VAR>v</VAR>)</I>
893
<DD><A NAME="IDX860"></A>
894
This function returns the maximum value in the vector <VAR>v</VAR>.
900
<DT><U>Function:</U> double <B>gsl_vector_min</B> <I>(const gsl_vector * <VAR>v</VAR>)</I>
901
<DD><A NAME="IDX861"></A>
902
This function returns the minimum value in the vector <VAR>v</VAR>.
908
<DT><U>Function:</U> void <B>gsl_vector_minmax</B> <I>(const gsl_vector * <VAR>v</VAR>, double * <VAR>min_out</VAR>, double * <VAR>max_out</VAR>)</I>
909
<DD><A NAME="IDX862"></A>
910
This function returns the minimum and maximum values in the vector
911
<VAR>v</VAR>, storing them in <VAR>min_out</VAR> and <VAR>max_out</VAR>.
917
<DT><U>Function:</U> size_t <B>gsl_vector_max_index</B> <I>(const gsl_vector * <VAR>v</VAR>)</I>
918
<DD><A NAME="IDX863"></A>
919
This function returns the index of the maximum value in the vector <VAR>v</VAR>.
920
When there are several equal maximum elements then the lowest index is
927
<DT><U>Function:</U> size_t <B>gsl_vector_min_index</B> <I>(const gsl_vector * <VAR>v</VAR>)</I>
928
<DD><A NAME="IDX864"></A>
929
This function returns the index of the minimum value in the vector <VAR>v</VAR>.
930
When there are several equal minimum elements then the lowest index is
937
<DT><U>Function:</U> void <B>gsl_vector_minmax_index</B> <I>(const gsl_vector * <VAR>v</VAR>, size_t * <VAR>imin</VAR>, size_t * <VAR>imax</VAR>)</I>
938
<DD><A NAME="IDX865"></A>
939
This function returns the indices of the minimum and maximum values in
940
the vector <VAR>v</VAR>, storing them in <VAR>imin</VAR> and <VAR>imax</VAR>. When
941
there are several equal minimum or maximum elements then the lowest
942
indices are returned.
948
<H3><A NAME="SEC168" HREF="gsl-ref_toc.html#TOC168">Vector properties</A></H3>
952
<DT><U>Function:</U> int <B>gsl_vector_isnull</B> <I>(const gsl_vector * <VAR>v</VAR>)</I>
953
<DD><A NAME="IDX866"></A>
954
This function returns 1 if all the elements of the vector <VAR>v</VAR> are
955
zero, and 0 otherwise.
961
<H3><A NAME="SEC169" HREF="gsl-ref_toc.html#TOC169">Example programs for vectors</A></H3>
964
This program shows how to allocate, initialize and read from a vector
965
using the functions <CODE>gsl_vector_alloc</CODE>, <CODE>gsl_vector_set</CODE> and
966
<CODE>gsl_vector_get</CODE>.
970
<PRE class="example">
971
#include <stdio.h>
972
#include <gsl/gsl_vector.h>
978
gsl_vector * v = gsl_vector_alloc (3);
980
for (i = 0; i < 3; i++)
982
gsl_vector_set (v, i, 1.23 + i);
985
for (i = 0; i < 100; i++)
987
printf ("v_%d = %g\n", i, gsl_vector_get (v, i));
995
Here is the output from the program. The final loop attempts to read
996
outside the range of the vector <CODE>v</CODE>, and the error is trapped by
997
the range-checking code in <CODE>gsl_vector_get</CODE>.
1001
<PRE class="example">
1005
gsl: vector_source.c:12: ERROR: index out of range
1006
IOT trap/Abort (core dumped)
1010
The next program shows how to write a vector to a file.
1014
<PRE class="example">
1015
#include <stdio.h>
1016
#include <gsl/gsl_vector.h>
1022
gsl_vector * v = gsl_vector_alloc (100);
1024
for (i = 0; i < 100; i++)
1026
gsl_vector_set (v, i, 1.23 + i);
1030
FILE * f = fopen ("test.dat", "w");
1031
gsl_vector_fprintf (f, v, "%.5g");
1039
After running this program the file <TT>'test.dat'</TT> should contain the
1040
elements of <CODE>v</CODE>, written using the format specifier
1041
<CODE>%.5g</CODE>. The vector could then be read back in using the function
1042
<CODE>gsl_vector_fscanf (f, v)</CODE> as follows:
1046
<PRE class="example">
1047
#include <stdio.h>
1048
#include <gsl/gsl_vector.h>
1054
gsl_vector * v = gsl_vector_alloc (10);
1057
FILE * f = fopen ("test.dat", "r");
1058
gsl_vector_fscanf (f, v);
1062
for (i = 0; i < 10; i++)
1064
printf ("%g\n", gsl_vector_get(v, i));
1073
<H2><A NAME="SEC170" HREF="gsl-ref_toc.html#TOC170">Matrices</A></H2>
1075
<A NAME="IDX867"></A>
1076
<A NAME="IDX868"></A>
1077
<A NAME="IDX869"></A>
1078
<A NAME="IDX870"></A>
1082
Matrices are defined by a <CODE>gsl_matrix</CODE> structure which describes a
1083
generalized slice of a block. Like a vector it represents a set of
1084
elements in an area of memory, but uses two indices instead of one.
1088
The <CODE>gsl_matrix</CODE> structure contains six components, the two
1089
dimensions of the matrix, a physical dimension, a pointer to the memory
1090
where the elements of the matrix are stored, <VAR>data</VAR>, a pointer to
1091
the block owned by the matrix <VAR>block</VAR>, if any, and an ownership
1092
flag, <VAR>owner</VAR>. The physical dimension determines the memory layout
1093
and can differ from the matrix dimension to allow the use of
1094
submatrices. The <CODE>gsl_matrix</CODE> structure is very simple and looks
1099
<PRE class="example">
1112
Matrices are stored in row-major order, meaning that each row of
1113
elements forms a contiguous block in memory. This is the standard
1114
"C-language ordering" of two-dimensional arrays. Note that FORTRAN
1115
stores arrays in column-major order. The number of rows is <VAR>size1</VAR>.
1116
The range of valid row indices runs from 0 to <CODE>size1-1</CODE>. Similarly
1117
<VAR>size2</VAR> is the number of columns. The range of valid column indices
1118
runs from 0 to <CODE>size2-1</CODE>. The physical row dimension <VAR>tda</VAR>, or
1119
<I>trailing dimension</I>, specifies the size of a row of the matrix as
1124
For example, in the following matrix <VAR>size1</VAR> is 3, <VAR>size2</VAR> is 4,
1125
and <VAR>tda</VAR> is 8. The physical memory layout of the matrix begins in
1126
the top left hand-corner and proceeds from left to right along each row
1131
<PRE class="example">
1132
00 01 02 03 XX XX XX XX
1133
10 11 12 13 XX XX XX XX
1134
20 21 22 23 XX XX XX XX
1138
Each unused memory location is represented by "<CODE>XX</CODE>". The
1139
pointer <VAR>data</VAR> gives the location of the first element of the matrix
1140
in memory. The pointer <VAR>block</VAR> stores the location of the memory
1141
block in which the elements of the matrix are located (if any). If the
1142
matrix owns this block then the <VAR>owner</VAR> field is set to one and the
1143
block will be deallocated when the matrix is freed. If the matrix is
1144
only a slice of a block owned by another object then the <VAR>owner</VAR> field is
1145
zero and any underlying block will not be freed.
1149
The functions for allocating and accessing matrices are defined in
1150
<TT>'gsl_matrix.h'</TT>
1156
<H3><A NAME="SEC171" HREF="gsl-ref_toc.html#TOC171">Matrix allocation</A></H3>
1159
The functions for allocating memory to a matrix follow the style of
1160
<CODE>malloc</CODE> and <CODE>free</CODE>. They also perform their own error
1161
checking. If there is insufficient memory available to allocate a vector
1162
then the functions call the GSL error handler (with an error number of
1163
<CODE>GSL_ENOMEM</CODE>) in addition to returning a null pointer. Thus if you
1164
use the library error handler to abort your program then it isn't
1165
necessary to check every <CODE>alloc</CODE>.
1170
<DT><U>Function:</U> gsl_matrix * <B>gsl_matrix_alloc</B> <I>(size_t <VAR>n1</VAR>, size_t <VAR>n2</VAR>)</I>
1171
<DD><A NAME="IDX871"></A>
1172
This function creates a matrix of size <VAR>n1</VAR> rows by <VAR>n2</VAR>
1173
columns, returning a pointer to a newly initialized matrix struct. A new
1174
block is allocated for the elements of the matrix, and stored in the
1175
<VAR>block</VAR> component of the matrix struct. The block is "owned" by the
1176
matrix, and will be deallocated when the matrix is deallocated.
1182
<DT><U>Function:</U> gsl_matrix * <B>gsl_matrix_calloc</B> <I>(size_t <VAR>n1</VAR>, size_t <VAR>n2</VAR>)</I>
1183
<DD><A NAME="IDX872"></A>
1184
This function allocates memory for a matrix of size <VAR>n1</VAR> rows by
1185
<VAR>n2</VAR> columns and initializes all the elements of the matrix to zero.
1191
<DT><U>Function:</U> void <B>gsl_matrix_free</B> <I>(gsl_matrix * <VAR>m</VAR>)</I>
1192
<DD><A NAME="IDX873"></A>
1193
This function frees a previously allocated matrix <VAR>m</VAR>. If the
1194
matrix was created using <CODE>gsl_matrix_alloc</CODE> then the block
1195
underlying the matrix will also be deallocated. If the matrix has been
1196
created from another object then the memory is still owned by that
1197
object and will not be deallocated.
1203
<H3><A NAME="SEC172" HREF="gsl-ref_toc.html#TOC172">Accessing matrix elements</A></H3>
1205
<A NAME="IDX874"></A>
1206
<A NAME="IDX875"></A>
1210
The functions for accessing the elements of a matrix use the same range
1211
checking system as vectors. You can turn off range checking by recompiling
1212
your program with the preprocessor definition
1213
<CODE>GSL_RANGE_CHECK_OFF</CODE>.
1217
The elements of the matrix are stored in "C-order", where the second
1218
index moves continuously through memory. More precisely, the element
1219
accessed by the function <CODE>gsl_matrix_get(m,i,j)</CODE> and
1220
<CODE>gsl_matrix_set(m,i,j,x)</CODE> is
1224
<PRE class="example">
1225
m->data[i * m->tda + j]
1229
where <VAR>tda</VAR> is the physical row-length of the matrix.
1234
<DT><U>Function:</U> double <B>gsl_matrix_get</B> <I>(const gsl_matrix * <VAR>m</VAR>, size_t <VAR>i</VAR>, size_t <VAR>j</VAR>)</I>
1235
<DD><A NAME="IDX876"></A>
1236
This function returns the (i,j)-th element of a matrix
1237
<VAR>m</VAR>. If <VAR>i</VAR> or <VAR>j</VAR> lie outside the allowed range of 0 to
1238
<VAR>n1-1</VAR> and 0 to <VAR>n2-1</VAR> then the error handler is invoked and 0
1245
<DT><U>Function:</U> void <B>gsl_matrix_set</B> <I>(gsl_matrix * <VAR>m</VAR>, size_t <VAR>i</VAR>, size_t <VAR>j</VAR>, double <VAR>x</VAR>)</I>
1246
<DD><A NAME="IDX877"></A>
1247
This function sets the value of the (i,j)-th element of a
1248
matrix <VAR>m</VAR> to <VAR>x</VAR>. If <VAR>i</VAR> or <VAR>j</VAR> lies outside the
1249
allowed range of 0 to <VAR>n1-1</VAR> and 0 to <VAR>n2-1</VAR> then the error
1256
<DT><U>Function:</U> double * <B>gsl_matrix_ptr</B> <I>(gsl_matrix * <VAR>m</VAR>, size_t <VAR>i</VAR>, size_t <VAR>j</VAR>)</I>
1257
<DD><A NAME="IDX878"></A>
1258
<DT><U>Function:</U> const double * <B>gsl_matrix_const_ptr</B> <I>(const gsl_matrix * <VAR>m</VAR>, size_t <VAR>i</VAR>, size_t <VAR>j</VAR>)</I>
1259
<DD><A NAME="IDX879"></A>
1260
These functions return a pointer to the (i,j)-th element of a
1261
matrix <VAR>m</VAR>. If <VAR>i</VAR> or <VAR>j</VAR> lie outside the allowed range of
1262
0 to <VAR>n1-1</VAR> and 0 to <VAR>n2-1</VAR> then the error handler is invoked
1263
and a null pointer is returned.
1269
<H3><A NAME="SEC173" HREF="gsl-ref_toc.html#TOC173">Initializing matrix elements</A></H3>
1271
<A NAME="IDX880"></A>
1272
<A NAME="IDX881"></A>
1273
<A NAME="IDX882"></A>
1274
<A NAME="IDX883"></A>
1275
<A NAME="IDX884"></A>
1276
<A NAME="IDX885"></A>
1277
<A NAME="IDX886"></A>
1278
<A NAME="IDX887"></A>
1283
<DT><U>Function:</U> void <B>gsl_matrix_set_all</B> <I>(gsl_matrix * <VAR>m</VAR>, double <VAR>x</VAR>)</I>
1284
<DD><A NAME="IDX888"></A>
1285
This function sets all the elements of the matrix <VAR>m</VAR> to the value
1292
<DT><U>Function:</U> void <B>gsl_matrix_set_zero</B> <I>(gsl_matrix * <VAR>m</VAR>)</I>
1293
<DD><A NAME="IDX889"></A>
1294
This function sets all the elements of the matrix <VAR>m</VAR> to zero.
1300
<DT><U>Function:</U> void <B>gsl_matrix_set_identity</B> <I>(gsl_matrix * <VAR>m</VAR>)</I>
1301
<DD><A NAME="IDX890"></A>
1302
This function sets the elements of the matrix <VAR>m</VAR> to the
1303
corresponding elements of the identity matrix, m(i,j) =
1304
\delta(i,j), i.e. a unit diagonal with all off-diagonal elements zero.
1305
This applies to both square and rectangular matrices.
1311
<H3><A NAME="SEC174" HREF="gsl-ref_toc.html#TOC174">Reading and writing matrices</A></H3>
1314
The library provides functions for reading and writing matrices to a file
1315
as binary data or formatted text.
1320
<DT><U>Function:</U> int <B>gsl_matrix_fwrite</B> <I>(FILE * <VAR>stream</VAR>, const gsl_matrix * <VAR>m</VAR>)</I>
1321
<DD><A NAME="IDX891"></A>
1322
This function writes the elements of the matrix <VAR>m</VAR> to the stream
1323
<VAR>stream</VAR> in binary format. The return value is 0 for success and
1324
<CODE>GSL_EFAILED</CODE> if there was a problem writing to the file. Since the
1325
data is written in the native binary format it may not be portable
1326
between different architectures.
1332
<DT><U>Function:</U> int <B>gsl_matrix_fread</B> <I>(FILE * <VAR>stream</VAR>, gsl_matrix * <VAR>m</VAR>)</I>
1333
<DD><A NAME="IDX892"></A>
1334
This function reads into the matrix <VAR>m</VAR> from the open stream
1335
<VAR>stream</VAR> in binary format. The matrix <VAR>m</VAR> must be preallocated
1336
with the correct dimensions since the function uses the size of <VAR>m</VAR> to
1337
determine how many bytes to read. The return value is 0 for success and
1338
<CODE>GSL_EFAILED</CODE> if there was a problem reading from the file. The
1339
data is assumed to have been written in the native binary format on the
1346
<DT><U>Function:</U> int <B>gsl_matrix_fprintf</B> <I>(FILE * <VAR>stream</VAR>, const gsl_matrix * <VAR>m</VAR>, const char * <VAR>format</VAR>)</I>
1347
<DD><A NAME="IDX893"></A>
1348
This function writes the elements of the matrix <VAR>m</VAR> line-by-line to
1349
the stream <VAR>stream</VAR> using the format specifier <VAR>format</VAR>, which
1350
should be one of the <CODE>%g</CODE>, <CODE>%e</CODE> or <CODE>%f</CODE> formats for
1351
floating point numbers and <CODE>%d</CODE> for integers. The function returns
1352
0 for success and <CODE>GSL_EFAILED</CODE> if there was a problem writing to
1359
<DT><U>Function:</U> int <B>gsl_matrix_fscanf</B> <I>(FILE * <VAR>stream</VAR>, gsl_matrix * <VAR>m</VAR>)</I>
1360
<DD><A NAME="IDX894"></A>
1361
This function reads formatted data from the stream <VAR>stream</VAR> into the
1362
matrix <VAR>m</VAR>. The matrix <VAR>m</VAR> must be preallocated with the correct
1363
dimensions since the function uses the size of <VAR>m</VAR> to determine how many
1364
numbers to read. The function returns 0 for success and
1365
<CODE>GSL_EFAILED</CODE> if there was a problem reading from the file.
1371
<H3><A NAME="SEC175" HREF="gsl-ref_toc.html#TOC175">Matrix views</A></H3>
1374
A matrix view is a temporary object, stored on the stack, which can be
1375
used to operate on a subset of matrix elements. Matrix views can be
1376
defined for both constant and non-constant matrices using separate types
1377
that preserve constness. A matrix view has the type
1378
<CODE>gsl_matrix_view</CODE> and a constant matrix view has the type
1379
<CODE>gsl_matrix_const_view</CODE>. In both cases the elements of the view
1380
can by accessed using the <CODE>matrix</CODE> component of the view object. A
1381
pointer <CODE>gsl_matrix *</CODE> or <CODE>const gsl_matrix *</CODE> can be obtained
1382
by taking the address of the <CODE>matrix</CODE> component with the <CODE>&</CODE>
1383
operator. In addition to matrix views it is also possible to create
1384
vector views of a matrix, such as row or column views.
1389
<DT><U>Function:</U> gsl_matrix_view <B>gsl_matrix_submatrix</B> <I>(gsl_matrix * <VAR>m</VAR>, size_t <VAR>k1</VAR>, size_t <VAR>k2</VAR>, size_t <VAR>n1</VAR>, size_t <VAR>n2</VAR>)</I>
1390
<DD><A NAME="IDX895"></A>
1391
<DT><U>Function:</U> gsl_matrix_const_view <B>gsl_matrix_const_submatrix</B> <I>(const gsl_matrix * <VAR>m</VAR>, size_t <VAR>k1</VAR>, size_t <VAR>k2</VAR>, size_t <VAR>n1</VAR>, size_t <VAR>n2</VAR>)</I>
1392
<DD><A NAME="IDX896"></A>
1393
These functions return a matrix view of a submatrix of the matrix
1394
<VAR>m</VAR>. The upper-left element of the submatrix is the element
1395
(<VAR>k1</VAR>,<VAR>k2</VAR>) of the original matrix. The submatrix has <VAR>n1</VAR>
1396
rows and <VAR>n2</VAR> columns. The physical number of columns in memory
1397
given by <VAR>tda</VAR> is unchanged. Mathematically, the
1398
(i,j)-th element of the new matrix is given by,
1402
<PRE class="example">
1403
m'(i,j) = m->data[(k1*m->tda + k1) + i*m->tda + j]
1407
where the index <VAR>i</VAR> runs from 0 to <CODE>n1-1</CODE> and the index <VAR>j</VAR>
1408
runs from 0 to <CODE>n2-1</CODE>.
1412
The <CODE>data</CODE> pointer of the returned matrix struct is set to null if
1413
the combined parameters (<VAR>i</VAR>,<VAR>j</VAR>,<VAR>n1</VAR>,<VAR>n2</VAR>,<VAR>tda</VAR>)
1414
overrun the ends of the original matrix.
1418
The new matrix view is only a view of the block underlying the existing
1419
matrix, <VAR>m</VAR>. The block containing the elements of <VAR>m</VAR> is not
1420
owned by the new matrix view. When the view goes out of scope the
1421
original matrix <VAR>m</VAR> and its block will continue to exist. The
1422
original memory can only be deallocated by freeing the original matrix.
1423
Of course, the original matrix should not be deallocated while the view
1428
The function <CODE>gsl_matrix_const_submatrix</CODE> is equivalent to
1429
<CODE>gsl_matrix_submatrix</CODE> but can be used for matrices which are
1430
declared <CODE>const</CODE>.
1437
<DT><U>Function:</U> gsl_matrix_view <B>gsl_matrix_view_array</B> <I>(double * <VAR>base</VAR>, size_t <VAR>n1</VAR>, size_t <VAR>n2</VAR>)</I>
1438
<DD><A NAME="IDX897"></A>
1439
<DT><U>Function:</U> gsl_matrix_const_view <B>gsl_matrix_const_view_array</B> <I>(const double * <VAR>base</VAR>, size_t <VAR>n1</VAR>, size_t <VAR>n2</VAR>)</I>
1440
<DD><A NAME="IDX898"></A>
1441
These functions return a matrix view of the array <VAR>base</VAR>. The
1442
matrix has <VAR>n1</VAR> rows and <VAR>n2</VAR> columns. The physical number of
1443
columns in memory is also given by <VAR>n2</VAR>. Mathematically, the
1444
(i,j)-th element of the new matrix is given by,
1448
<PRE class="example">
1449
m'(i,j) = base[i*n2 + j]
1453
where the index <VAR>i</VAR> runs from 0 to <CODE>n1-1</CODE> and the index <VAR>j</VAR>
1454
runs from 0 to <CODE>n2-1</CODE>.
1458
The new matrix is only a view of the array <VAR>base</VAR>. When the view
1459
goes out of scope the original array <VAR>base</VAR> will continue to exist.
1460
The original memory can only be deallocated by freeing the original
1461
array. Of course, the original array should not be deallocated while
1462
the view is still in use.
1466
The function <CODE>gsl_matrix_const_view_array</CODE> is equivalent to
1467
<CODE>gsl_matrix_view_array</CODE> but can be used for matrices which are
1468
declared <CODE>const</CODE>.
1475
<DT><U>Function:</U> gsl_matrix_view <B>gsl_matrix_view_array_with_tda</B> <I>(double * <VAR>base</VAR>, size_t <VAR>n1</VAR>, size_t <VAR>n2</VAR>, size_t <VAR>tda</VAR>)</I>
1476
<DD><A NAME="IDX899"></A>
1477
<DT><U>Function:</U> gsl_matrix_const_view <B>gsl_matrix_const_view_array_with_tda</B> <I>(const double * <VAR>base</VAR>, size_t <VAR>n1</VAR>, size_t <VAR>n2</VAR>, size_t <VAR>tda</VAR>)</I>
1478
<DD><A NAME="IDX900"></A>
1479
These functions return a matrix view of the array <VAR>base</VAR> with a
1480
physical number of columns <VAR>tda</VAR> which may differ from the corresponding
1481
dimension of the matrix. The matrix has <VAR>n1</VAR> rows and <VAR>n2</VAR>
1482
columns, and the physical number of columns in memory is given by
1483
<VAR>tda</VAR>. Mathematically, the (i,j)-th element of the new
1488
<PRE class="example">
1489
m'(i,j) = base[i*tda + j]
1493
where the index <VAR>i</VAR> runs from 0 to <CODE>n1-1</CODE> and the index <VAR>j</VAR>
1494
runs from 0 to <CODE>n2-1</CODE>.
1498
The new matrix is only a view of the array <VAR>base</VAR>. When the view
1499
goes out of scope the original array <VAR>base</VAR> will continue to exist.
1500
The original memory can only be deallocated by freeing the original
1501
array. Of course, the original array should not be deallocated while
1502
the view is still in use.
1506
The function <CODE>gsl_matrix_const_view_array_with_tda</CODE> is equivalent
1507
to <CODE>gsl_matrix_view_array_with_tda</CODE> but can be used for matrices
1508
which are declared <CODE>const</CODE>.
1514
<DT><U>Function:</U> gsl_matrix_view <B>gsl_matrix_view_vector</B> <I>(gsl_vector * <VAR>v</VAR>, size_t <VAR>n1</VAR>, size_t <VAR>n2</VAR>)</I>
1515
<DD><A NAME="IDX901"></A>
1516
<DT><U>Function:</U> gsl_matrix_const_view <B>gsl_matrix_const_view_vector</B> <I>(const gsl_vector * <VAR>v</VAR>, size_t <VAR>n1</VAR>, size_t <VAR>n2</VAR>)</I>
1517
<DD><A NAME="IDX902"></A>
1518
These functions return a matrix view of the vector <VAR>v</VAR>. The matrix
1519
has <VAR>n1</VAR> rows and <VAR>n2</VAR> columns. The vector must have unit
1520
stride. The physical number of columns in memory is also given by
1521
<VAR>n2</VAR>. Mathematically, the (i,j)-th element of the new
1526
<PRE class="example">
1527
m'(i,j) = v->data[i*n2 + j]
1531
where the index <VAR>i</VAR> runs from 0 to <CODE>n1-1</CODE> and the index <VAR>j</VAR>
1532
runs from 0 to <CODE>n2-1</CODE>.
1536
The new matrix is only a view of the vector <VAR>v</VAR>. When the view
1537
goes out of scope the original vector <VAR>v</VAR> will continue to exist.
1538
The original memory can only be deallocated by freeing the original
1539
vector. Of course, the original vector should not be deallocated while
1540
the view is still in use.
1544
The function <CODE>gsl_matrix_const_view_vector</CODE> is equivalent to
1545
<CODE>gsl_matrix_view_vector</CODE> but can be used for matrices which are
1546
declared <CODE>const</CODE>.
1553
<DT><U>Function:</U> gsl_matrix_view <B>gsl_matrix_view_vector_with_tda</B> <I>(gsl_vector * <VAR>v</VAR>, size_t <VAR>n1</VAR>, size_t <VAR>n2</VAR>, size_t <VAR>tda</VAR>)</I>
1554
<DD><A NAME="IDX903"></A>
1555
<DT><U>Function:</U> gsl_matrix_const_view <B>gsl_matrix_const_view_vector_with_tda</B> <I>(const gsl_vector * <VAR>v</VAR>, size_t <VAR>n1</VAR>, size_t <VAR>n2</VAR>, size_t <VAR>tda</VAR>)</I>
1556
<DD><A NAME="IDX904"></A>
1557
These functions return a matrix view of the vector <VAR>v</VAR> with a
1558
physical number of columns <VAR>tda</VAR> which may differ from the
1559
corresponding matrix dimension. The vector must have unit stride. The
1560
matrix has <VAR>n1</VAR> rows and <VAR>n2</VAR> columns, and the physical number
1561
of columns in memory is given by <VAR>tda</VAR>. Mathematically, the
1562
(i,j)-th element of the new matrix is given by,
1566
<PRE class="example">
1567
m'(i,j) = v->data[i*tda + j]
1571
where the index <VAR>i</VAR> runs from 0 to <CODE>n1-1</CODE> and the index <VAR>j</VAR>
1572
runs from 0 to <CODE>n2-1</CODE>.
1576
The new matrix is only a view of the vector <VAR>v</VAR>. When the view
1577
goes out of scope the original vector <VAR>v</VAR> will continue to exist.
1578
The original memory can only be deallocated by freeing the original
1579
vector. Of course, the original vector should not be deallocated while
1580
the view is still in use.
1584
The function <CODE>gsl_matrix_const_view_vector_with_tda</CODE> is equivalent
1585
to <CODE>gsl_matrix_view_vector_with_tda</CODE> but can be used for matrices
1586
which are declared <CODE>const</CODE>.
1593
<H3><A NAME="SEC176" HREF="gsl-ref_toc.html#TOC176">Creating row and column views</A></H3>
1596
In general there are two ways to access an object, by reference or by
1597
copying. The functions described in this section create vector views
1598
which allow access to a row or column of a matrix by reference.
1599
Modifying elements of the view is equivalent to modifying the matrix,
1600
since both the vector view and the matrix point to the same memory
1606
<DT><U>Function:</U> gsl_vector_view <B>gsl_matrix_row</B> <I>(gsl_matrix * <VAR>m</VAR>, size_t <VAR>i</VAR>)</I>
1607
<DD><A NAME="IDX905"></A>
1608
<DT><U>Function:</U> gsl_vector_const_view <B>gsl_matrix_const_row</B> <I>(const gsl_matrix * <VAR>m</VAR>, size_t <VAR>i</VAR>)</I>
1609
<DD><A NAME="IDX906"></A>
1610
These functions return a vector view of the <VAR>i</VAR>-th row of the matrix
1611
<VAR>m</VAR>. The <CODE>data</CODE> pointer of the new vector is set to null if
1612
<VAR>i</VAR> is out of range.
1616
The function <CODE>gsl_vector_const_row</CODE> is equivalent to
1617
<CODE>gsl_matrix_row</CODE> but can be used for matrices which are declared
1624
<DT><U>Function:</U> gsl_vector_view <B>gsl_matrix_column</B> <I>(gsl_matrix * <VAR>m</VAR>, size_t <VAR>j</VAR>)</I>
1625
<DD><A NAME="IDX907"></A>
1626
<DT><U>Function:</U> gsl_vector_const_view <B>gsl_matrix_const_column</B> <I>(const gsl_matrix * <VAR>m</VAR>, size_t <VAR>j</VAR>)</I>
1627
<DD><A NAME="IDX908"></A>
1628
These functions return a vector view of the <VAR>j</VAR>-th column of the
1629
matrix <VAR>m</VAR>. The <CODE>data</CODE> pointer of the new vector is set to
1630
null if <VAR>j</VAR> is out of range.
1634
The function <CODE>gsl_vector_const_column</CODE> equivalent to
1635
<CODE>gsl_matrix_column</CODE> but can be used for matrices which are declared
1641
<A NAME="IDX909"></A>
1642
<A NAME="IDX910"></A>
1644
<DT><U>Function:</U> gsl_vector_view <B>gsl_matrix_diagonal</B> <I>(gsl_matrix * <VAR>m</VAR>)</I>
1645
<DD><A NAME="IDX911"></A>
1646
<DT><U>Function:</U> gsl_vector_const_view <B>gsl_matrix_const_diagonal</B> <I>(const gsl_matrix * <VAR>m</VAR>)</I>
1647
<DD><A NAME="IDX912"></A>
1648
These functions returns a vector view of the diagonal of the matrix
1649
<VAR>m</VAR>. The matrix <VAR>m</VAR> is not required to be square. For a
1650
rectangular matrix the length of the diagonal is the same as the smaller
1651
dimension of the matrix.
1655
The function <CODE>gsl_matrix_const_diagonal</CODE> is equivalent to
1656
<CODE>gsl_matrix_diagonal</CODE> but can be used for matrices which are
1657
declared <CODE>const</CODE>.
1662
<A NAME="IDX913"></A>
1663
<A NAME="IDX914"></A>
1665
<DT><U>Function:</U> gsl_vector_view <B>gsl_matrix_subdiagonal</B> <I>(gsl_matrix * <VAR>m</VAR>, size_t <VAR>k</VAR>)</I>
1666
<DD><A NAME="IDX915"></A>
1667
<DT><U>Function:</U> gsl_vector_const_view <B>gsl_matrix_const_subdiagonal</B> <I>(const gsl_matrix * <VAR>m</VAR>, size_t <VAR>k</VAR>)</I>
1668
<DD><A NAME="IDX916"></A>
1669
These functions return a vector view of the <VAR>k</VAR>-th subdiagonal of
1670
the matrix <VAR>m</VAR>. The matrix <VAR>m</VAR> is not required to be square.
1671
The diagonal of the matrix corresponds to k = 0.
1675
The function <CODE>gsl_matrix_const_subdiagonal</CODE> is equivalent to
1676
<CODE>gsl_matrix_subdiagonal</CODE> but can be used for matrices which are
1677
declared <CODE>const</CODE>.
1682
<A NAME="IDX917"></A>
1683
<A NAME="IDX918"></A>
1685
<DT><U>Function:</U> gsl_vector_view <B>gsl_matrix_superdiagonal</B> <I>(gsl_matrix * <VAR>m</VAR>, size_t <VAR>k</VAR>)</I>
1686
<DD><A NAME="IDX919"></A>
1687
<DT><U>Function:</U> gsl_vector_const_view <B>gsl_matrix_const_superdiagonal</B> <I>(const gsl_matrix * <VAR>m</VAR>, size_t <VAR>k</VAR>)</I>
1688
<DD><A NAME="IDX920"></A>
1689
These functions return a vector view of the <VAR>k</VAR>-th superdiagonal of
1690
the matrix <VAR>m</VAR>. The matrix <VAR>m</VAR> is not required to be square. The
1691
diagonal of the matrix corresponds to k = 0.
1695
The function <CODE>gsl_matrix_const_superdiagonal</CODE> is equivalent to
1696
<CODE>gsl_matrix_superdiagonal</CODE> but can be used for matrices which are
1697
declared <CODE>const</CODE>.
1704
<H3><A NAME="SEC177" HREF="gsl-ref_toc.html#TOC177">Copying matrices</A></H3>
1708
<DT><U>Function:</U> int <B>gsl_matrix_memcpy</B> <I>(gsl_matrix * <VAR>dest</VAR>, const gsl_matrix * <VAR>src</VAR>)</I>
1709
<DD><A NAME="IDX921"></A>
1710
This function copies the elements of the matrix <VAR>src</VAR> into the
1711
matrix <VAR>dest</VAR>. The two matrices must have the same size.
1717
<DT><U>Function:</U> int <B>gsl_matrix_swap</B> <I>(gsl_matrix * <VAR>m1</VAR>, gsl_matrix * <VAR>m2</VAR>)</I>
1718
<DD><A NAME="IDX922"></A>
1719
This function exchanges the elements of the matrices <VAR>m1</VAR> and
1720
<VAR>m2</VAR> by copying. The two matrices must have the same size.
1726
<H3><A NAME="SEC178" HREF="gsl-ref_toc.html#TOC178">Copying rows and columns</A></H3>
1729
The functions described in this section copy a row or column of a matrix
1730
into a vector. This allows the elements of the vector and the matrix to
1731
be modified independently. Note that if the matrix and the vector point
1732
to overlapping regions of memory then the result will be undefined. The
1733
same effect can be achieved with more generality using
1734
<CODE>gsl_vector_memcpy</CODE> with vector views of rows and columns.
1739
<DT><U>Function:</U> int <B>gsl_matrix_get_row</B> <I>(gsl_vector * <VAR>v</VAR>, const gsl_matrix * <VAR>m</VAR>, size_t <VAR>i</VAR>)</I>
1740
<DD><A NAME="IDX923"></A>
1741
This function copies the elements of the <VAR>i</VAR>-th row of the matrix
1742
<VAR>m</VAR> into the vector <VAR>v</VAR>. The length of the vector must be the
1743
same as the length of the row.
1749
<DT><U>Function:</U> int <B>gsl_matrix_get_col</B> <I>(gsl_vector * <VAR>v</VAR>, const gsl_matrix * <VAR>m</VAR>, size_t <VAR>j</VAR>)</I>
1750
<DD><A NAME="IDX924"></A>
1751
This function copies the elements of the <VAR>j</VAR>-th column of the matrix
1752
<VAR>m</VAR> into the vector <VAR>v</VAR>. The length of the vector must be the
1753
same as the length of the column.
1759
<DT><U>Function:</U> int <B>gsl_matrix_set_row</B> <I>(gsl_matrix * <VAR>m</VAR>, size_t <VAR>i</VAR>, const gsl_vector * <VAR>v</VAR>)</I>
1760
<DD><A NAME="IDX925"></A>
1761
This function copies the elements of the vector <VAR>v</VAR> into the
1762
<VAR>i</VAR>-th row of the matrix <VAR>m</VAR>. The length of the vector must be
1763
the same as the length of the row.
1769
<DT><U>Function:</U> int <B>gsl_matrix_set_col</B> <I>(gsl_matrix * <VAR>m</VAR>, size_t <VAR>j</VAR>, const gsl_vector * <VAR>v</VAR>)</I>
1770
<DD><A NAME="IDX926"></A>
1771
This function copies the elements of the vector <VAR>v</VAR> into the
1772
<VAR>j</VAR>-th column of the matrix <VAR>m</VAR>. The length of the vector must be
1773
the same as the length of the column.
1779
<H3><A NAME="SEC179" HREF="gsl-ref_toc.html#TOC179">Exchanging rows and columns</A></H3>
1782
The following functions can be used to exchange the rows and columns of
1788
<DT><U>Function:</U> int <B>gsl_matrix_swap_rows</B> <I>(gsl_matrix * <VAR>m</VAR>, size_t <VAR>i</VAR>, size_t <VAR>j</VAR>)</I>
1789
<DD><A NAME="IDX927"></A>
1790
This function exchanges the <VAR>i</VAR>-th and <VAR>j</VAR>-th rows of the matrix
1791
<VAR>m</VAR> in-place.
1797
<DT><U>Function:</U> int <B>gsl_matrix_swap_columns</B> <I>(gsl_matrix * <VAR>m</VAR>, size_t <VAR>i</VAR>, size_t <VAR>j</VAR>)</I>
1798
<DD><A NAME="IDX928"></A>
1799
This function exchanges the <VAR>i</VAR>-th and <VAR>j</VAR>-th columns of the
1800
matrix <VAR>m</VAR> in-place.
1806
<DT><U>Function:</U> int <B>gsl_matrix_swap_rowcol</B> <I>(gsl_matrix * <VAR>m</VAR>, size_t <VAR>i</VAR>, size_t <VAR>j</VAR>)</I>
1807
<DD><A NAME="IDX929"></A>
1808
This function exchanges the <VAR>i</VAR>-th row and <VAR>j</VAR>-th column of the
1809
matrix <VAR>m</VAR> in-place. The matrix must be square for this operation to
1816
<DT><U>Function:</U> int <B>gsl_matrix_transpose_memcpy</B> <I>(gsl_matrix * <VAR>dest</VAR>, const gsl_matrix * <VAR>src</VAR>)</I>
1817
<DD><A NAME="IDX930"></A>
1818
This function makes the matrix <VAR>dest</VAR> the transpose of the matrix
1819
<VAR>src</VAR> by copying the elements of <VAR>src</VAR> into <VAR>dest</VAR>. This
1820
function works for all matrices provided that the dimensions of the matrix
1821
<VAR>dest</VAR> match the transposed dimensions of the matrix <VAR>src</VAR>.
1827
<DT><U>Function:</U> int <B>gsl_matrix_transpose</B> <I>(gsl_matrix * <VAR>m</VAR>)</I>
1828
<DD><A NAME="IDX931"></A>
1829
This function replaces the matrix <VAR>m</VAR> by its transpose by copying
1830
the elements of the matrix in-place. The matrix must be square for this
1831
operation to be possible.
1837
<H3><A NAME="SEC180" HREF="gsl-ref_toc.html#TOC180">Matrix operations</A></H3>
1840
The following operations are defined for real and complex matrices.
1845
<DT><U>Function:</U> int <B>gsl_matrix_add</B> <I>(gsl_matrix * <VAR>a</VAR>, const gsl_matrix * <VAR>b</VAR>)</I>
1846
<DD><A NAME="IDX932"></A>
1847
This function adds the elements of matrix <VAR>b</VAR> to the elements of
1848
matrix <VAR>a</VAR>, a'(i,j) = a(i,j) + b(i,j). The two matrices must have the
1855
<DT><U>Function:</U> int <B>gsl_matrix_sub</B> <I>(gsl_matrix * <VAR>a</VAR>, const gsl_matrix * <VAR>b</VAR>)</I>
1856
<DD><A NAME="IDX933"></A>
1857
This function subtracts the elements of matrix <VAR>b</VAR> from the elements of
1858
matrix <VAR>a</VAR>, a'(i,j) = a(i,j) - b(i,j). The two matrices must have the
1865
<DT><U>Function:</U> int <B>gsl_matrix_mul_elements</B> <I>(gsl_matrix * <VAR>a</VAR>, const gsl_matrix * <VAR>b</VAR>)</I>
1866
<DD><A NAME="IDX934"></A>
1867
This function multiplies the elements of matrix <VAR>a</VAR> by the elements of
1868
matrix <VAR>b</VAR>, a'(i,j) = a(i,j) * b(i,j). The two matrices must have the
1875
<DT><U>Function:</U> int <B>gsl_matrix_div_elements</B> <I>(gsl_matrix * <VAR>a</VAR>, const gsl_matrix * <VAR>b</VAR>)</I>
1876
<DD><A NAME="IDX935"></A>
1877
This function divides the elements of matrix <VAR>a</VAR> by the elements of
1878
matrix <VAR>b</VAR>, a'(i,j) = a(i,j) / b(i,j). The two matrices must have the
1885
<DT><U>Function:</U> int <B>gsl_matrix_scale</B> <I>(gsl_matrix * <VAR>a</VAR>, const double <VAR>x</VAR>)</I>
1886
<DD><A NAME="IDX936"></A>
1887
This function multiplies the elements of matrix <VAR>a</VAR> by the constant
1888
factor <VAR>x</VAR>, a'(i,j) = x a(i,j).
1894
<DT><U>Function:</U> int <B>gsl_matrix_add_constant</B> <I>(gsl_matrix * <VAR>a</VAR>, const double <VAR>x</VAR>)</I>
1895
<DD><A NAME="IDX937"></A>
1896
This function adds the constant value <VAR>x</VAR> to the elements of the
1897
matrix <VAR>a</VAR>, a'(i,j) = a(i,j) + x.
1903
<H3><A NAME="SEC181" HREF="gsl-ref_toc.html#TOC181">Finding maximum and minimum elements of matrices</A></H3>
1906
The following operations are only defined for real matrices.
1911
<DT><U>Function:</U> double <B>gsl_matrix_max</B> <I>(const gsl_matrix * <VAR>m</VAR>)</I>
1912
<DD><A NAME="IDX938"></A>
1913
This function returns the maximum value in the matrix <VAR>m</VAR>.
1919
<DT><U>Function:</U> double <B>gsl_matrix_min</B> <I>(const gsl_matrix * <VAR>m</VAR>)</I>
1920
<DD><A NAME="IDX939"></A>
1921
This function returns the minimum value in the matrix <VAR>m</VAR>.
1927
<DT><U>Function:</U> void <B>gsl_matrix_minmax</B> <I>(const gsl_matrix * <VAR>m</VAR>, double * <VAR>min_out</VAR>, double * <VAR>max_out</VAR>)</I>
1928
<DD><A NAME="IDX940"></A>
1929
This function returns the minimum and maximum values in the matrix
1930
<VAR>m</VAR>, storing them in <VAR>min_out</VAR> and <VAR>max_out</VAR>.
1936
<DT><U>Function:</U> void <B>gsl_matrix_max_index</B> <I>(const gsl_matrix * <VAR>m</VAR>, size_t * <VAR>imax</VAR>, size_t * <VAR>jmax</VAR>)</I>
1937
<DD><A NAME="IDX941"></A>
1938
This function returns the indices of the maximum value in the matrix
1939
<VAR>m</VAR>, storing them in <VAR>imax</VAR> and <VAR>jmax</VAR>. When there are
1940
several equal maximum elements then the first element found is returned,
1941
searching in row-major order.
1947
<DT><U>Function:</U> void <B>gsl_matrix_min_index</B> <I>(const gsl_matrix * <VAR>m</VAR>, size_t * <VAR>imax</VAR>, size_t * <VAR>jmax</VAR>)</I>
1948
<DD><A NAME="IDX942"></A>
1949
This function returns the indices of the minimum value in the matrix
1950
<VAR>m</VAR>, storing them in <VAR>imax</VAR> and <VAR>jmax</VAR>. When there are
1951
several equal minimum elements then the first element found is returned,
1952
searching in row-major order.
1958
<DT><U>Function:</U> void <B>gsl_matrix_minmax_index</B> <I>(const gsl_matrix * <VAR>m</VAR>, size_t * <VAR>imin</VAR>, size_t * <VAR>imax</VAR>)</I>
1959
<DD><A NAME="IDX943"></A>
1960
This function returns the indices of the minimum and maximum values in
1961
the matrix <VAR>m</VAR>, storing them in (<VAR>imin</VAR>,<VAR>jmin</VAR>) and
1962
(<VAR>imax</VAR>,<VAR>jmax</VAR>). When there are several equal minimum or maximum
1963
elements then the first elements found are returned, searching in
1970
<H3><A NAME="SEC182" HREF="gsl-ref_toc.html#TOC182">Matrix properties</A></H3>
1974
<DT><U>Function:</U> int <B>gsl_matrix_isnull</B> <I>(const gsl_matrix * <VAR>m</VAR>)</I>
1975
<DD><A NAME="IDX944"></A>
1976
This function returns 1 if all the elements of the matrix <VAR>m</VAR> are
1977
zero, and 0 otherwise.
1983
<H3><A NAME="SEC183" HREF="gsl-ref_toc.html#TOC183">Example programs for matrices</A></H3>
1986
The program below shows how to allocate, initialize and read from a matrix
1987
using the functions <CODE>gsl_matrix_alloc</CODE>, <CODE>gsl_matrix_set</CODE> and
1988
<CODE>gsl_matrix_get</CODE>.
1992
<PRE class="example">
1993
#include <stdio.h>
1994
#include <gsl/gsl_matrix.h>
2000
gsl_matrix * m = gsl_matrix_alloc (10, 3);
2002
for (i = 0; i < 10; i++)
2003
for (j = 0; j < 3; j++)
2004
gsl_matrix_set (m, i, j, 0.23 + 100*i + j);
2006
for (i = 0; i < 10; i++)
2007
for (j = 0; j < 3; j++)
2008
printf ("m(%d,%d) = %g\n", i, j,
2009
gsl_matrix_get (m, i, j));
2016
Here is the output from the program. The final loop attempts to read
2017
outside the range of the matrix <CODE>m</CODE>, and the error is trapped by
2018
the range-checking code in <CODE>gsl_matrix_get</CODE>.
2022
<PRE class="example">
2031
gsl: matrix_source.c:13: ERROR: first index out of range
2032
IOT trap/Abort (core dumped)
2036
The next program shows how to write a matrix to a file.
2040
<PRE class="example">
2041
#include <stdio.h>
2042
#include <gsl/gsl_matrix.h>
2048
gsl_matrix * m = gsl_matrix_alloc (100, 100);
2049
gsl_matrix * a = gsl_matrix_alloc (100, 100);
2051
for (i = 0; i < 100; i++)
2052
for (j = 0; j < 100; j++)
2053
gsl_matrix_set (m, i, j, 0.23 + i + j);
2056
FILE * f = fopen ("test.dat", "wb");
2057
gsl_matrix_fwrite (f, m);
2062
FILE * f = fopen ("test.dat", "rb");
2063
gsl_matrix_fread (f, a);
2067
for (i = 0; i < 100; i++)
2068
for (j = 0; j < 100; j++)
2070
double mij = gsl_matrix_get (m, i, j);
2071
double aij = gsl_matrix_get (a, i, j);
2072
if (mij != aij) k++;
2075
printf ("differences = %d (should be zero)\n", k);
2081
After running this program the file <TT>'test.dat'</TT> should contain the
2082
elements of <CODE>m</CODE>, written in binary format. The matrix which is read
2083
back in using the function <CODE>gsl_matrix_fread</CODE> should be exactly
2084
equal to the original matrix.
2088
The following program demonstrates the use of vector views. The program
2089
computes the column-norms of a matrix.
2093
<PRE class="example">
2094
#include <math.h>
2095
#include <stdio.h>
2096
#include <gsl/gsl_matrix.h>
2097
#include <gsl/gsl_blas.h>
2104
gsl_matrix *m = gsl_matrix_alloc (10, 10);
2106
for (i = 0; i < 10; i++)
2107
for (j = 0; j < 10; j++)
2108
gsl_matrix_set (m, i, j, sin (i) + cos (j));
2110
for (j = 0; j < 10; j++)
2112
gsl_vector_view column = gsl_matrix_column (m, j);
2115
d = gsl_blas_dnrm2 (&column.vector);
2117
printf ("matrix column %d, norm = %g\n", j, d);
2120
gsl_matrix_free (m);
2127
Here is the output of the program, which can be confirmed using GNU
2132
<PRE class="example">
2134
matrix column 0, norm = 4.31461
2135
matrix column 1, norm = 3.1205
2136
matrix column 2, norm = 2.19316
2137
matrix column 3, norm = 3.26114
2138
matrix column 4, norm = 2.53416
2139
matrix column 5, norm = 2.57281
2140
matrix column 6, norm = 4.20469
2141
matrix column 7, norm = 3.65202
2142
matrix column 8, norm = 2.08524
2143
matrix column 9, norm = 3.07313
2147
<PRE class="example">
2148
octave> m = sin(0:9)' * ones(1,10)
2149
+ ones(10,1) * cos(0:9);
2150
octave> sqrt(sum(m.^2))
2153
4.3146 3.1205 2.1932 3.2611 2.5342 2.5728
2154
4.2047 3.6520 2.0852 3.0731
2159
<H2><A NAME="SEC184" HREF="gsl-ref_toc.html#TOC184">References and Further Reading</A></H2>
2162
The block, vector and matrix objects in GSL follow the <CODE>valarray</CODE>
2163
model of C++. A description of this model can be found in the following
2168
<UL class="itemize">
2172
<CITE>The C++ Programming Language</CITE> (3rd Ed),
2173
Section 22.4 Vector Arithmetic.
2174
Addison-Wesley 1997, ISBN 0-201-88954-4.
2178
<p>Go to the <A HREF="gsl-ref_1.html">first</A>, <A HREF="gsl-ref_7.html">previous</A>, <A HREF="gsl-ref_9.html">next</A>, <A HREF="gsl-ref_50.html">last</A> section, <A HREF="gsl-ref_toc.html">table of contents</A>.