This is gsl-ref.info, produced by makeinfo version 4.0 from gsl-ref.texi. INFO-DIR-SECTION Scientific software START-INFO-DIR-ENTRY * gsl-ref: (gsl-ref). GNU Scientific Library - Reference END-INFO-DIR-ENTRY This file documents the GNU Scientific Library. Copyright (C) 1996, 1997, 1998, 1999, 2000, 2001 The GSL Team. Permission is granted to copy, distribute and/or modify this document under the terms of the GNU Free Documentation License, Version 1.1 or any later version published by the Free Software Foundation; with no Invariant Sections, no Front-Cover Texts, and no Back-Cover Texts. A copy of the license is included in the section entitled "GNU Free Documentation License".  File: gsl-ref.info, Node: Vector views, Next: Copying vectors, Prev: Reading and writing vectors, Up: Vectors Vector views ------------ In addition to creating vectors from slices of blocks it is also possible to slice vectors and create vector views. For example, a subvector of another vector can be described with a view, or two views can be made which provide access to the even and odd elements of a vector. A vector view is a temporary object, stored on the stack, which can be used to operate on a subset of vector elements. Vector views can be defined for both constant and non-constant vectors, using separate types that preserve constness. A vector view has the type `gsl_vector_view' and a constant vector view has the type `gsl_vector_const_view'. In both cases the elements of the view can be accessed as a `gsl_vector' using the `vector' component of the view object. A pointer to a vector of type `gsl_vector *' or `const gsl_vector *' can be obtained by taking the address of this component with the `&' operator. - Function: gsl_vector_view gsl_vector_subvector (gsl_vector *V, size_t OFFSET, size_t N) - Function: gsl_vector_const_view gsl_vector_const_subvector (const gsl_vector * V, size_t OFFSET, size_t N) These functions return a vector view of a subvector of another vector V. The start of the new vector is offset by OFFSET elements from the start of the original vector. The new vector has N elements. Mathematically, the I-th element of the new vector V' is given by, v'(i) = v->data[(offset + i)*v->stride] where the index I runs from 0 to `n-1'. The `data' pointer of the returned vector struct is set to null if the combined parameters (OFFSET,N) overrun the end of the original vector. The new vector is only a view of the block underlying the original vector, V. The block containing the elements of V is not owned by the new vector. When the view goes out of scope the original vector V and its block will continue to exist. The original memory can only be deallocated by freeing the original vector. Of course, the original vector should not be deallocated while the view is still in use. The function `gsl_vector_const_subvector' is equivalent to `gsl_vector_subvector' but can be used for vectors which are declared `const'. - Function: gsl_vector gsl_vector_subvector_with_stride (gsl_vector *V, size_t OFFSET, size_t STRIDE, size_t N) - Function: gsl_vector_const_view gsl_vector_const_subvector_with_stride (const gsl_vector * V, size_t OFFSET, size_t STRIDE, size_t N) These functions return a vector view of a subvector of another vector V with an additional stride argument. The subvector is formed in the same way as for `gsl_vector_subvector' but the new vector has N elements with a step-size of STRIDE from one element to the next in the original vector. Mathematically, the I-th element of the new vector V' is given by, v'(i) = v->data[(offset + i*stride)*v->stride] where the index I runs from 0 to `n-1'. Note that subvector views give direct access to the underlying elements of the original vector. For example, the following code will zero the even elements of the vector `v' of length `n', while leaving the odd elements untouched, gsl_vector_view v_even = gsl_vector_subvector_with_stride (v, 0, 2, n/2); gsl_vector_set_zero (&v_even.vector); A vector view can be passed to any subroutine which takes a vector argument just as a directly allocated vector would be, using `&'VIEW`.vector'. For example, the following code computes the norm of odd elements of `v' using the BLAS routine DNRM2, gsl_vector_view v_odd = gsl_vector_subvector_with_stride (v, 1, 2, n/2); double r = gsl_blas_dnrm2 (&v_odd.vector); The function `gsl_vector_const_subvector_with_stride' is equivalent to `gsl_vector_subvector_with_stride' but can be used for vectors which are declared `const'. - Function: gsl_vector_view gsl_vector_complex_real (gsl_vector_complex *V) - Function: gsl_vector_const_view gsl_vector_complex_const_real (const gsl_vector_complex *V) These functions return a vector view of the real parts of the complex vector V. The function `gsl_vector_complex_const_real' is equivalent to `gsl_vector_complex_real' but can be used for vectors which are declared `const'. - Function: gsl_vector_view gsl_vector_complex_imag (gsl_vector_complex *V) - Function: gsl_vector_const_view gsl_vector_complex_const_imag (const gsl_vector_complex *V) These functions return a vector view of the imaginary parts of the complex vector V. The function `gsl_vector_complex_const_imag' is equivalent to `gsl_vector_complex_imag' but can be used for vectors which are declared `const'. - Function: gsl_vector_view gsl_vector_view_array (double *BASE, size_t N) - Function: gsl_vector_const_view gsl_vector_const_view_array (const double *BASE, size_t N) These functions return a vector view of an array. The start of the new vector is given by BASE and has N elements. Mathematically, the I-th element of the new vector V' is given by, v'(i) = base[i] where the index I runs from 0 to `n-1'. The array containing the elements of V is not owned by the new vector view. When the view goes out of scope the original array will continue to exist. The original memory can only be deallocated by freeing the original pointer BASE. Of course, the original array should not be deallocated while the view is still in use. The function `gsl_vector_const_view_array' is equivalent to `gsl_vector_view_array' but can be used for arrays which are declared `const'. - Function: gsl_vector_view gsl_vector_view_array_with_stride (double * BASE, size_t STRIDE, size_t N) - Function: gsl_vector_const_view gsl_vector_const_view_array_with_stride (const double * BASE, size_t STRIDE, size_t N) These functions return a vector view of an array BASE with an additional stride argument. The subvector is formed in the same way as for `gsl_vector_view_array' but the new vector has N elements with a step-size of STRIDE from one element to the next in the original array. Mathematically, the I-th element of the new vector V' is given by, v'(i) = base[i*stride] where the index I runs from 0 to `n-1'. Note that the view gives direct access to the underlying elements of the original array. A vector view can be passed to any subroutine which takes a vector argument just as a directly allocated vector would be, using `&'VIEW`.vector'. The function `gsl_vector_const_view_array_with_stride' is equivalent to `gsl_vector_view_array_with_stride' but can be used for arrays which are declared `const'.  File: gsl-ref.info, Node: Copying vectors, Next: Exchanging elements, Prev: Vector views, Up: Vectors Copying vectors --------------- Common operations on vectors such as addition and multiplication are available in the BLAS part of the library (*note BLAS Support::). However, it is useful to have a small number of utility functions which do not require the full BLAS code. The following functions fall into this category. - Function: int gsl_vector_memcpy (gsl_vector * DEST, const gsl_vector * SRC) This function copies the elements of the vector SRC into the vector DEST. The two vectors must have the same length. - Function: int gsl_vector_swap (gsl_vector * V, gsl_vector * W) This function exchanges the elements of the vectors V and W by copying. The two vectors must have the same length.  File: gsl-ref.info, Node: Exchanging elements, Next: Vector operations, Prev: Copying vectors, Up: Vectors Exchanging elements ------------------- The following function can be used to exchange, or permute, the elements of a vector. - Function: int gsl_vector_swap_elements (gsl_vector * V, size_t i, size_t j) This function exchanges the I-th and J-th elements of the vector V in-place. - Function: int gsl_vector_reverse (gsl_vector * V) This function reverses the order of the elements of the vector V.  File: gsl-ref.info, Node: Vector operations, Next: Finding maximum and minimum elements of vectors, Prev: Exchanging elements, Up: Vectors Vector operations ----------------- The following operations are only defined for real vectors. - Function: int gsl_vector_add (gsl_vector * A, const gsl_vector * B) This function adds the elements of vector B to the elements of vector A, a'_i = a_i + b_i. The two vectors must have the same length. - Function: int gsl_vector_sub (gsl_vector * A, const gsl_vector * B) This function subtracts the elements of vector B from the elements of vector A, a'_i = a_i - b_i. The two vectors must have the same length. - Function: int gsl_vector_mul (gsl_vector * A, const gsl_vector * B) This function multiplies the elements of vector A by the elements of vector B, a'_i = a_i * b_i. The two vectors must have the same length. - Function: int gsl_vector_div (gsl_vector * A, const gsl_vector * B) This function divides the elements of vector A by the elements of vector B, a'_i = a_i / b_i. The two vectors must have the same length. - Function: int gsl_vector_scale (gsl_vector * A, const double X) This function multiplies the elements of vector A by the constant factor X, a'_i = x a_i. - Function: int gsl_vector_add_constant (gsl_vector * A, const double X) This function adds the constant value X to the elements of the vector A, a'_i = a_i + x.  File: gsl-ref.info, Node: Finding maximum and minimum elements of vectors, Next: Vector properties, Prev: Vector operations, Up: Vectors Finding maximum and minimum elements of vectors ----------------------------------------------- - Function: double gsl_vector_max (const gsl_vector * V) This function returns the maximum value in the vector V. - Function: double gsl_vector_min (const gsl_vector * V) This function returns the minimum value in the vector V. - Function: void gsl_vector_minmax (const gsl_vector * V, double * MIN_OUT, double * MAX_OUT) This function returns the minimum and maximum values in the vector V, storing them in MIN_OUT and MAX_OUT. - Function: size_t gsl_vector_max_index (const gsl_vector * V) This function returns the index of the maximum value in the vector V. When there are several equal maximum elements then the lowest index is returned. - Function: size_t gsl_vector_min_index (const gsl_vector * V) This function returns the index of the minimum value in the vector V. When there are several equal minimum elements then the lowest index is returned. - Function: void gsl_vector_minmax_index (const gsl_vector * V, size_t * IMIN, size_t * IMAX) This function returns the indices of the minimum and maximum values in the vector V, storing them in IMIN and IMAX. When there are several equal minimum or maximum elements then the lowest indices are returned.  File: gsl-ref.info, Node: Vector properties, Next: Example programs for vectors, Prev: Finding maximum and minimum elements of vectors, Up: Vectors Vector properties ----------------- - Function: int gsl_vector_isnull (const gsl_vector * V) This function returns 1 if all the elements of the vector V are zero, and 0 otherwise.  File: gsl-ref.info, Node: Example programs for vectors, Prev: Vector properties, Up: Vectors Example programs for vectors ---------------------------- This program shows how to allocate, initialize and read from a vector using the functions `gsl_vector_alloc', `gsl_vector_set' and `gsl_vector_get'. #include #include int main (void) { int i; gsl_vector * v = gsl_vector_alloc (3); for (i = 0; i < 3; i++) { gsl_vector_set (v, i, 1.23 + i); } for (i = 0; i < 100; i++) { printf("v_%d = %g\n", i, gsl_vector_get (v, i)); } return 0; } Here is the output from the program. The final loop attempts to read outside the range of the vector `v', and the error is trapped by the range-checking code in `gsl_vector_get'. v_0 = 1.23 v_1 = 2.23 v_2 = 3.23 gsl: vector_source.c:12: ERROR: index out of range IOT trap/Abort (core dumped) The next program shows how to write a vector to a file. #include #include int main (void) { int i; gsl_vector * v = gsl_vector_alloc (100); for (i = 0; i < 100; i++) { gsl_vector_set (v, i, 1.23 + i); } { FILE * f = fopen("test.dat", "w"); gsl_vector_fprintf (f, v, "%.5g"); fclose (f); } return 0; } After running this program the file `test.dat' should contain the elements of `v', written using the format specifier `%.5g'. The vector could then be read back in using the function `gsl_vector_fscanf (f, v)'.  File: gsl-ref.info, Node: Matrices, Next: Vector and Matrix References and Further Reading, Prev: Vectors, Up: Vectors and Matrices Matrices ======== Matrices are defined by a `gsl_matrix' structure which describes a generalized slice of a block. Like a vector it represents a set of elements in an area of memory, but uses two indices instead of one. The `gsl_matrix' structure contains six components, the two dimensions of the matrix, a physical dimension, a pointer to the memory where the elements of the matrix are stored, DATA, a pointer to the block owned by the matrix BLOCK, if any, and an ownership flag, OWNER. The physical dimension determines the memory layout and can differ from the matrix dimension to allow the use of submatrices. The `gsl_matrix' structure is very simple and looks like this, typedef struct { size_t size1; size_t size2; size_t tda; double * data; gsl_block * block; int owner; } gsl_matrix; Matrices are stored in row-major order, meaning that each row of elements forms a contiguous block in memory. This is the standard "C-language ordering" of two-dimensional arrays. Note that FORTRAN stores arrays in column-major order. The number of rows is SIZE1. The range of valid row indices runs from 0 to `size1-1'. Similarly SIZE2 is the number of columns. The range of valid column indices runs from 0 to `size2-1'. The physical row dimension TDA, or "trailing dimension", specifies the size of a row of the matrix as laid out in memory. For example, in the following matrix SIZE1 is 3, SIZE2 is 4, and TDA is 8. The physical memory layout of the matrix begins in the top left hand-corner and proceeds from left to right along each row in turn. 00 01 02 03 XX XX XX XX 10 11 12 13 XX XX XX XX 20 21 22 23 XX XX XX XX Each unused memory location is represented by "`XX'". The pointer DATA gives the location of the first element of the matrix in memory. The pointer BLOCK stores the location of the memory block in which the elements of the matrix are located (if any). If the matrix owns this block then the OWNER field is set to one and the block will be deallocated when the matrix is freed. If the matrix is only a slice of a block owned by another object then the OWNER field is zero and any underlying block will not be freed. The functions for allocating and accessing matrices are defined in `gsl_matrix.h' * Menu: * Matrix allocation:: * Accessing matrix elements:: * Initializing matrix elements:: * Reading and writing matrices:: * Matrix views:: * Creating row and column views:: * Copying matrices:: * Copying rows and columns:: * Exchanging rows and columns:: * Matrix operations:: * Finding maximum and minimum elements of matrices:: * Matrix properties:: * Example programs for matrices::  File: gsl-ref.info, Node: Matrix allocation, Next: Accessing matrix elements, Up: Matrices Matrix allocation ----------------- The functions for allocating memory to a matrix follow the style of `malloc' and `free'. They also perform their own error checking. If there is insufficient memory available to allocate a vector then the functions call the GSL error handler (with an error number of `GSL_ENOMEM') in addition to returning a null pointer. Thus if you use the library error handler to abort your program then it isn't necessary to check every `alloc'. - Function: gsl_matrix * gsl_matrix_alloc (size_t N1, size_t N2) This function creates a matrix of size N1 rows by N2 columns, returning a pointer to a newly initialized matrix struct. A new block is allocated for the elements of the matrix, and stored in the BLOCK component of the matrix struct. The block is "owned" by the matrix, and will be deallocated when the matrix is deallocated. - Function: gsl_matrix * gsl_matrix_calloc (size_t N1, size_t N2) This function allocates memory for a matrix of size N1 rows by N2 columns and initializes all the elements of the matrix to zero. - Function: void gsl_matrix_free (gsl_matrix * M) This function frees a previously allocated matrix M. If the matrix was created using `gsl_matrix_alloc' then the block underlying the matrix will also be deallocated. If the matrix has been created from another object then the memory is still owned by that object and will not be deallocated.  File: gsl-ref.info, Node: Accessing matrix elements, Next: Initializing matrix elements, Prev: Matrix allocation, Up: Matrices Accessing matrix elements ------------------------- The functions for accessing the elements of a matrix use the same range checking system as vectors. You turn off range checking by recompiling your program with the preprocessor definition `GSL_RANGE_CHECK_OFF'. The elements of the matrix are stored in "C-order", where the second index moves continuously through memory. More precisely, the element accessed by the function `gsl_matrix_get(m,i,j)' and `gsl_matrix_set(m,i,j,x)' is m->data[i * m->tda + j] where TDA is the physical row-length of the matrix. - Function: double gsl_matrix_get (const gsl_matrix * M, size_t I, size_t J) This function returns the (I,J)th element of a matrix M. If I or J lie outside the allowed range of 0 to N1-1 and 0 to N2-1 then the error handler is invoked and 0 is returned. - Function: void gsl_matrix_set (gsl_matrix * M, size_t I, size_t J, double X) This function sets the value of the (I,J)th element of a matrix M to X. If I or J lies outside the allowed range of 0 to N1-1 and 0 to N2-1 then the error handler is invoked. - Function: double * gsl_matrix_ptr (gsl_matrix * M, size_t I, size_t J) - Function: const double * gsl_matrix_ptr (const gsl_matrix * M, size_t I, size_t J) These functions return a pointer to the (I,J)th element of a matrix M. If I or J lie outside the allowed range of 0 to N1-1 and 0 to N2-1 then the error handler is invoked and a null pointer is returned.  File: gsl-ref.info, Node: Initializing matrix elements, Next: Reading and writing matrices, Prev: Accessing matrix elements, Up: Matrices Initializing matrix elements ---------------------------- - Function: void gsl_matrix_set_all (gsl_matrix * M, double X) This function sets all the elements of the matrix M to the value X. - Function: void gsl_matrix_set_zero (gsl_matrix * M) This function sets all the elements of the matrix M to zero. - Function: void gsl_matrix_set_identity (gsl_matrix * M) This function sets the elements of the matrix M to the corresponding elements of the identity matrix, m(i,j) = \delta(i,j), i.e. a unit diagonal with all off-diagonal elements zero. This applies to both square and rectangular matrices.  File: gsl-ref.info, Node: Reading and writing matrices, Next: Matrix views, Prev: Initializing matrix elements, Up: Matrices Reading and writing matrices ---------------------------- The library provides functions for reading and writing matrices to a file as binary data or formatted text. - Function: int gsl_matrix_fwrite (FILE * STREAM, const gsl_matrix * M) This function writes the elements of the matrix M to the stream STREAM in binary format. The return value is 0 for success and `GSL_EFAILED' if there was a problem writing to the file. Since the data is written in the native binary format it may not be portable between different architectures. - Function: int gsl_matrix_fread (FILE * STREAM, gsl_matrix * M) This function reads into the matrix M from the open stream STREAM in binary format. The matrix M must be preallocated with the correct dimensions since the function uses the size of M to determine how many bytes to read. The return value is 0 for success and `GSL_EFAILED' if there was a problem reading from the file. The data is assumed to have been written in the native binary format on the same architecture. - Function: int gsl_matrix_fprintf (FILE * STREAM, const gsl_matrix * M, const char * FORMAT) This function writes the elements of the matrix M line-by-line to the stream STREAM using the format specifier FORMAT, which should be one of the `%g', `%e' or `%f' formats for floating point numbers and `%d' for integers. The function returns 0 for success and `GSL_EFAILED' if there was a problem writing to the file. - Function: int gsl_matrix_fscanf (FILE * STREAM, gsl_matrix * M) This function reads formatted data from the stream STREAM into the matrix M. The matrix M must be preallocated with the correct dimensions since the function uses the size of M to determine how many numbers to read. The function returns 0 for success and `GSL_EFAILED' if there was a problem reading from the file.  File: gsl-ref.info, Node: Matrix views, Next: Creating row and column views, Prev: Reading and writing matrices, Up: Matrices Matrix views ------------ A matrix view is a temporary object, stored on the stack, which can be used to operate on a subset of matrix elements. Matrix views can be defined for both constant and non-constant matrices using separate types that preserve constness. A matrix view has the type `gsl_matrix_view' and a constant matrix view has the type `gsl_matrix_const_view'. In both cases the elements of the view can by accessed using the `matrix' component of the view object. A pointer `gsl_matrix *' or `const gsl_matrix *' can be obtained by taking the address of the `matrix' component with the `&' operator. In addition to matrix views it is also possible to create vector views of a matrix, such as row or column views. - Function: gsl_matrix_view gsl_matrix_submatrix (gsl_matrix * M, size_t K1, size_t K2, size_t N1, size_t N2) - Function: gsl_matrix_const_view gsl_matrix_const_submatrix (const gsl_matrix * M, size_t K1, size_t K2, size_t N1, size_t N2) These functions return a matrix view of a submatrix of the matrix M. The upper-left element of the submatrix is the element (K1,K2) of the original matrix. The submatrix has N1 rows and N2 columns. The physical number of columns in memory given by TDA is unchanged. Mathematically, the (I,J)-th element of the new matrix is given by, m'(i,j) = m->data[(k1*m->tda + k1) + i*m->tda + j] where the index I runs from 0 to `n1-1' and the index J runs from 0 to `n2-1'. The `data' pointer of the returned matrix struct is set to null if the combined parameters (I,J,N1,N2,TDA) overrun the ends of the original matrix. The new matrix view is only a view of the block underlying the existing matrix, M. The block containing the elements of M is not owned by the new matrix view. When the view goes out of scope the original matrix M and its block will continue to exist. The original memory can only be deallocated by freeing the original matrix. Of course, the original matrix should not be deallocated while the view is still in use. The function `gsl_matrix_const_submatrix' is equivalent to `gsl_matrix_submatrix' but can be used for matrices which are declared `const'. - Function: gsl_matrix_view gsl_matrix_view_array (double * BASE, size_t N1, size_t N2) - Function: gsl_matrix_const_view gsl_matrix_const_view_array (const double * BASE, size_t N1, size_t N2) These functions return a matrix view of the array BASE. The matrix has N1 rows and N2 columns. The physical number of columns in memory is also given by N2. Mathematically, the (I,J)-th element of the new matrix is given by, m'(i,j) = base[i*n2 + j] where the index I runs from 0 to `n1-1' and the index J runs from 0 to `n2-1'. The new matrix is only a view of the array BASE. When the view goes out of scope the original array BASE will continue to exist. The original memory can only be deallocated by freeing the original array. Of course, the original array should not be deallocated while the view is still in use. The function `gsl_matrix_const_view_array' is equivalent to `gsl_matrix_view_array' but can be used for matrices which are declared `const'. - Function: gsl_matrix_view gsl_matrix_view_array_with_tda (double * BASE, size_t N1, size_t N2, size_t TDA) - Function: gsl_matrix_const_view gsl_matrix_const_view_array_with_tda (const double * BASE, size_t N1, size_t N2, size_t TDA) These functions return a matrix view of the array BASE with a physical number of columns TDA which may differ from corresponding the dimension of the matrix. The matrix has N1 rows and N2 columns, and the physical number of columns in memory is given by TDA. Mathematically, the (I,J)-th element of the new matrix is given by, m'(i,j) = base[i*tda + j] where the index I runs from 0 to `n1-1' and the index J runs from 0 to `n2-1'. The new matrix is only a view of the array BASE. When the view goes out of scope the original array BASE will continue to exist. The original memory can only be deallocated by freeing the original array. Of course, the original array should not be deallocated while the view is still in use. The function `gsl_matrix_const_view_array_with_tda' is equivalent to `gsl_matrix_view_array_with_tda' but can be used for matrices which are declared `const'. - Function: gsl_matrix_view gsl_matrix_view_vector (gsl_vector * V, size_t N1, size_t N2) - Function: gsl_matrix_const_view gsl_matrix_const_view_vector (const gsl_vector * V, size_t N1, size_t N2) These functions return a matrix view of the vector V. The matrix has N1 rows and N2 columns. The vector must have unit stride. The physical number of columns in memory is also given by N2. Mathematically, the (I,J)-th element of the new matrix is given by, m'(i,j) = v->data[i*n2 + j] where the index I runs from 0 to `n1-1' and the index J runs from 0 to `n2-1'. The new matrix is only a view of the vector V. When the view goes out of scope the original vector V will continue to exist. The original memory can only be deallocated by freeing the original vector. Of course, the original vector should not be deallocated while the view is still in use. The function `gsl_matrix_const_view_vector' is equivalent to `gsl_matrix_view_vector' but can be used for matrices which are declared `const'. - Function: gsl_matrix_view gsl_matrix_view_vector_with_tda (gsl_vector * V, size_t N1, size_t N2, size_t TDA) - Function: gsl_matrix_const_view gsl_matrix_const_view_vector_with_tda (const gsl_vector * V, size_t N1, size_t N2, size_t TDA) These functions return a matrix view of the vector V with a physical number of columns TDA which may differ from the corresponding matrix dimension. The vector must have unit stride. The matrix has N1 rows and N2 columns, and the physical number of columns in memory is given by TDA. Mathematically, the (I,J)-th element of the new matrix is given by, m'(i,j) = v->data[i*tda + j] where the index I runs from 0 to `n1-1' and the index J runs from 0 to `n2-1'. The new matrix is only a view of the vector V. When the view goes out of scope the original vector V will continue to exist. The original memory can only be deallocated by freeing the original vector. Of course, the original vector should not be deallocated while the view is still in use. The function `gsl_matrix_const_view_vector_with_tda' is equivalent to `gsl_matrix_view_vector_with_tda' but can be used for matrices which are declared `const'.  File: gsl-ref.info, Node: Creating row and column views, Next: Copying matrices, Prev: Matrix views, Up: Matrices Creating row and column views ----------------------------- In general there are two ways to access an object, by reference or by copying. The functions described in this section create vector views which allow access to a row or column of a matrix by reference. Modifying elements of the view is equivalent to modifying the matrix, since both the vector view and the matrix point to the same memory block. - Function: gsl_vector_view gsl_matrix_row (gsl_matrix * M, size_t I) - Function: gsl_vector_const_view gsl_matrix_const_row (const gsl_matrix * M, size_t I) These functions return a vector view of the I-th row of the matrix M. The `data' pointer of the new vector is set to null if I is out of range. The function `gsl_vector_const_row' is equivalent to `gsl_matrix_row' but can be used for matrices which are declared `const'. - Function: gsl_vector_view gsl_matrix_column (gsl_matrix * M, size_t J) - Function: gsl_vector_const_view gsl_matrix_const_column (const gsl_matrix * M, size_t J) These functions return a vector view of the J-th column of the matrix M. The `data' pointer of the new vector is set to null if J is out of range. The function `gsl_vector_const_column' equivalent to `gsl_matrix_column' but can be used for matrices which are declared `const'. - Function: gsl_vector_view gsl_matrix_diagonal (gsl_matrix * M) - Function: gsl_vector_const_view gsl_matrix_const_diagonal (const gsl_matrix * M) These functions returns a vector view of the diagonal of the matrix M. The matrix M is not required to be square. For a rectangular matrix the length of the diagonal is the same as the smaller dimension of the matrix. The function `gsl_matrix_const_diagonal' is equivalent to `gsl_matrix_diagonal' but can be used for matrices which are declared `const'. - Function: gsl_vector_view gsl_matrix_subdiagonal (gsl_matrix * M, size_t K) - Function: gsl_vector_const_view gsl_matrix_const_subdiagonal (const gsl_matrix * M, size_t K) These functions return a vector view of the K-th subdiagonal of the matrix M. The matrix M is not required to be square. The diagonal of the matrix corresponds to k = 0. The function `gsl_matrix_const_subdiagonal' is equivalent to `gsl_matrix_subdiagonal' but can be used for matrices which are declared `const'. - Function: gsl_vector_view gsl_matrix_superdiagonal (gsl_matrix * M, size_t K) - Function: gsl_vector_const_view gsl_matrix_const_superdiagonal (const gsl_matrix * M, size_t K) These functions return a vector view of the K-th superdiagonal of the matrix M. The matrix M is not required to be square. The diagonal of the matrix corresponds to k = 0. The function `gsl_matrix_const_superdiagonal' is equivalent to `gsl_matrix_superdiagonal' but can be used for matrices which are declared `const'.  File: gsl-ref.info, Node: Copying matrices, Next: Copying rows and columns, Prev: Creating row and column views, Up: Matrices Copying matrices ---------------- - Function: int gsl_matrix_memcpy (gsl_matrix * DEST, const gsl_matrix * SRC) This function copies the elements of the matrix SRC into the matrix DEST. The two matrices must have the same size. - Function: int gsl_matrix_swap (gsl_matrix * M1, const gsl_matrix * M2) This function exchanges the elements of the matrices M1 and M2 by copying. The two matrices must have the same size.  File: gsl-ref.info, Node: Copying rows and columns, Next: Exchanging rows and columns, Prev: Copying matrices, Up: Matrices Copying rows and columns ------------------------ The functions described in this section copy a row or column of a matrix into a vector. This allows the elements of the vector and the matrix to be modified independently. Note that if the matrix and the vector point to overlapping regions of memory then the result will be undefined. The same effect can be achieved with more generality using `gsl_vector_memcpy' with vector views of rows and columns. - Function: int gsl_matrix_get_row (gsl_vector * V, const gsl_matrix * M, size_t I) This function copies the elements of the I-th row of the matrix M into the vector V. The length of the vector must be the same as the length of the row. - Function: int gsl_matrix_get_col (gsl_vector * V, const gsl_matrix * M, size_t J) This function copies the elements of the I-th column of the matrix M into the vector V. The length of the vector must be the same as the length of the column. - Function: int gsl_matrix_set_row (gsl_matrix * M, size_t I, const gsl_vector * V) This function copies the elements of the vector V into the I-th row of the matrix M. The length of the vector must be the same as the length of the row. - Function: int gsl_matrix_set_col (gsl_matrix * M, size_t J, const gsl_vector * V) This function copies the elements of the vector V into the I-th column of the matrix M. The length of the vector must be the same as the length of the column.  File: gsl-ref.info, Node: Exchanging rows and columns, Next: Matrix operations, Prev: Copying rows and columns, Up: Matrices Exchanging rows and columns --------------------------- The following functions can be used to exchange the rows and columns of a matrix. - Function: int gsl_matrix_swap_rows (gsl_matrix * M, size_t I, size_t J) This function exchanges the I-th and J-th rows of the matrix M in-place. - Function: int gsl_matrix_swap_columns (gsl_matrix * M, size_t I, size_t J) This function exchanges the I-th and J-th columns of the matrix M in-place. - Function: int gsl_matrix_swap_rowcol (gsl_matrix * M, size_t I, size_t J) This function exchanges the I-th row and J-th column of the matrix M in-place. The matrix must be square for this operation to be possible. - Function: int gsl_matrix_transpose_memcpy (gsl_matrix * DEST, gsl_matrix * SRC) This function makes the matrix DEST the transpose of the matrix SRC by copying the elements of SRC into DEST. This function works for all matrices provided that the dimensions of the matrix DEST match the transposed dimensions of the matrix SRC. - Function: int gsl_matrix_transpose (gsl_matrix * M) This function replaces the matrix M by its transpose by copying the elements of the matrix in-place. The matrix must be square for this operation to be possible.  File: gsl-ref.info, Node: Matrix operations, Next: Finding maximum and minimum elements of matrices, Prev: Exchanging rows and columns, Up: Matrices Matrix operations ----------------- The following operations are only defined for real matrices. - Function: int gsl_matrix_add (gsl_matrix * A, const gsl_matrix * B) This function adds the elements of matrix B to the elements of matrix A, a'(i,j) = a(i,j) + b(i,j). The two matrices must have the same dimensions. - Function: int gsl_matrix_sub (gsl_matrix * A, const gsl_matrix * B) This function subtracts the elements of matrix B from the elements of matrix A, a'(i,j) = a(i,j) - b(i,j). The two matrices must have the same dimensions. - Function: int gsl_matrix_mul_elements (gsl_matrix * A, const gsl_matrix * B) This function multiplies the elements of matrix A by the elements of matrix B, a'(i,j) = a(i,j) * b(i,j). The two matrices must have the same dimensions. - Function: int gsl_matrix_div_elements (gsl_matrix * A, const gsl_matrix * B) This function divides the elements of matrix A by the elements of matrix B, a'(i,j) = a(i,j) / b(i,j). The two matrices must have the same dimensions. - Function: int gsl_matrix_scale (gsl_matrix * A, const double X) This function multiplies the elements of matrix A by the constant factor X, a'(i,j) = x a(i,j). - Function: int gsl_matrix_add_constant (gsl_matrix * A, const double X) This function adds the constant value X to the elements of the matrix A, a'(i,j) = a(i,j) + x.  File: gsl-ref.info, Node: Finding maximum and minimum elements of matrices, Next: Matrix properties, Prev: Matrix operations, Up: Matrices Finding maximum and minimum elements of matrices ------------------------------------------------ - Function: double gsl_matrix_max (const gsl_matrix * M) This function returns the maximum value in the matrix M. - Function: double gsl_matrix_min (const gsl_matrix * M) This function returns the minimum value in the matrix M. - Function: void gsl_matrix_minmax (const gsl_matrix * M, double * MIN_OUT, double * MAX_OUT) This function returns the minimum and maximum values in the matrix M, storing them in MIN_OUT and MAX_OUT. - Function: void gsl_matrix_max_index (const gsl_matrix * M, size_t * IMAX, size_t * JMAX) This function returns the indices of the maximum value in the matrix M, storing them in IMAX and JMAX. When there are several equal maximum elements then the first element found is returned. - Function: void gsl_matrix_min_index (const gsl_matrix * M, size_t * IMAX, size_t * JMAX) This function returns the indices of the minimum value in the matrix M, storing them in IMAX and JMAX. When there are several equal minimum elements then the first element found is returned. - Function: void gsl_matrix_minmax_index (const gsl_matrix * M, size_t * IMIN, size_t * IMAX) This function returns the indices of the minimum and maximum values in the matrix M, storing them in (IMIN,JMIN) and (IMAX,JMAX). When there are several equal minimum or maximum elements then the first elements found are returned.  File: gsl-ref.info, Node: Matrix properties, Next: Example programs for matrices, Prev: Finding maximum and minimum elements of matrices, Up: Matrices Matrix properties ----------------- - Function: int gsl_matrix_isnull (const gsl_matrix * M) This function returns 1 if all the elements of the matrix M are zero, and 0 otherwise.  File: gsl-ref.info, Node: Example programs for matrices, Prev: Matrix properties, Up: Matrices Example programs for matrices ----------------------------- The program below shows how to allocate, initialize and read from a matrix using the functions `gsl_matrix_alloc', `gsl_matrix_set' and `gsl_matrix_get'. #include #include int main (void) { int i, j; gsl_matrix * m = gsl_matrix_alloc (10, 3); for (i = 0; i < 10; i++) for (j = 0; j < 3; j++) gsl_matrix_set (m, i, j, 0.23 + 100*i + j); for (i = 0; i < 100; i++) for (j = 0; j < 3; j++) printf("m(%d,%d) = %g\n", i, j, gsl_matrix_get (m, i, j)); return 0; } Here is the output from the program. The final loop attempts to read outside the range of the matrix `m', and the error is trapped by the range-checking code in `gsl_matrix_get'. m(0,0) = 0.23 m(0,1) = 1.23 m(0,2) = 2.23 m(1,0) = 100.23 m(1,1) = 101.23 m(1,2) = 102.23 ... m(9,2) = 902.23 gsl: matrix_source.c:13: ERROR: first index out of range IOT trap/Abort (core dumped) The next program shows how to write a matrix to a file. #include #include int main (void) { int i, j, k = 0; gsl_matrix * m = gsl_matrix_alloc (100, 100); gsl_matrix * a = gsl_matrix_alloc (100, 100); for (i = 0; i < 100; i++) for (j = 0; j < 100; j++) gsl_matrix_set (m, i, j, 0.23 + i + j); { FILE * f = fopen("test.dat", "w"); gsl_matrix_fwrite (f, m); fclose (f); } { FILE * f = fopen("test.dat", "r"); gsl_matrix_fread (f, a); fclose (f); } for (i = 0; i < 100; i++) for (j = 0; j < 100; j++) { double mij = gsl_matrix_get(m, i, j); double aij = gsl_matrix_get(a, i, j); if (mij != aij) k++; } printf("differences = %d (should be zero)\n", k); return (k > 0); } After running this program the file `test.dat' should contain the elements of `m', written in binary format. The matrix which is read back in using the function `gsl_matrix_fread' should be exactly equal to the original matrix. The following program demonstrates the use of vector views. The program computes the column-norms of a matrix. #include #include #include #include int main (void) { size_t i,j; gsl_matrix *m = gsl_matrix_alloc (10, 10); for (i = 0; i < 10; i++) for (j = 0; j < 10; j++) gsl_matrix_set (m, i, j, sin (i) + cos (j)); for (j = 0; j < 10; j++) { gsl_vector_view column = gsl_matrix_column (m, j); double d; d = gsl_blas_dnrm2 (&column.vector); printf ("matrix column %d, norm = %g\n", j, d); } gsl_matrix_free (m); } Here is the output of the program, which can be confirmed using GNU OCTAVE, $ ./a.out matrix column 0, norm = 4.31461 matrix column 1, norm = 3.1205 matrix column 2, norm = 2.19316 matrix column 3, norm = 3.26114 matrix column 4, norm = 2.53416 matrix column 5, norm = 2.57281 matrix column 6, norm = 4.20469 matrix column 7, norm = 3.65202 matrix column 8, norm = 2.08524 matrix column 9, norm = 3.07313 octave> m = sin(0:9)' * ones(1,10) + ones(10,1) * cos(0:9); octave> sqrt(sum(m.^2)) ans = 4.3146 3.1205 2.1932 3.2611 2.5342 2.5728 4.2047 3.6520 2.0852 3.0731  File: gsl-ref.info, Node: Vector and Matrix References and Further Reading, Prev: Matrices, Up: Vectors and Matrices References and Further Reading ============================== The block, vector and matrix objects in GSL follow the `valarray' model of C++. A description of this model can be found in the following reference, B. Stroustrup, `The C++ Programming Language' (3rd Ed), Section 22.4 Vector Arithmetic. Addison-Wesley 1997, ISBN 0-201-88954-4.  File: gsl-ref.info, Node: Permutations, Next: Combinations, Prev: Vectors and Matrices, Up: Top Permutations ************ This chapter describes functions for creating and manipulating permutations. A permutation p is represented by an array of n integers in the range 0 .. n-1, where each value p_i occurs once and only once. The application of a permutation p to a vector v yields a new vector v' where v'_i = v_{p_i}. For example, the array (0,1,3,2) represents a permutation which exchanges the last two elements of a four element vector. The corresponding identity permutation is (0,1,2,3). Note that the permutations produced by the linear algebra routines correspond to the exchange of matrix columns, and so should be considered as applying to row-vectors in the form v' = v P rather than column-vectors, when permuting the elements of a vector. The functions described in this chapter are defined in the header file `gsl_permutation.h'. * Menu: * The Permutation struct:: * Permutation allocation:: * Accessing permutation elements:: * Permutation properties:: * Permutation functions:: * Applying Permutations:: * Reading and writing permutations:: * Permutation Examples:: * Permutation References and Further Reading::  File: gsl-ref.info, Node: The Permutation struct, Next: Permutation allocation, Up: Permutations The Permutation struct ====================== A permutation is stored by a structure containing two components, the size of the permutation and a pointer to the permutation array. The elements of the permutation array are all of type `size_t'. The `gsl_permutation' structure looks like this, typedef struct { size_t size; size_t * data; } gsl_permutation;  File: gsl-ref.info, Node: Permutation allocation, Next: Accessing permutation elements, Prev: The Permutation struct, Up: Permutations Permutation allocation ====================== - Function: gsl_permutation * gsl_permutation_alloc (size_t N) This function allocates memory for a new permutation of size N. The permutation is not initialized and its elements are undefined. Use the function `gsl_permutation_calloc' if you want to create a permutation which is initialized to the identity. A null pointer is returned if insufficient memory is available to create the permutation. - Function: gsl_permutation * gsl_permutation_calloc (size_t N) This function allocates memory for a new permutation of size N and initializes it to the identity. A null pointer is returned if insufficient memory is available to create the permutation. - Function: void gsl_permutation_init (gsl_permutation * P) This function initializes the permutation P to the identity, i.e. (0,1,2,...,n-1). - Function: void gsl_permutation_free (gsl_permutation * P) This function frees all the memory used by the permutation P. - Function: int gsl_permutation_memcpy (gsl_permutation * DEST, const gsl_permutation * SRC) This function copies the elements of the permutation SRC into the permutation DEST. The two permutations must have the same size.  File: gsl-ref.info, Node: Accessing permutation elements, Next: Permutation properties, Prev: Permutation allocation, Up: Permutations Accessing permutation elements ============================== The following functions can be used to access and manipulate permutations. - Function: size_t gsl_permutation_get (const gsl_permutation * P, const size_t I) This function returns the value of the I-th element of the permutation P. If I lies outside the allowed range of 0 to N-1 then the error handler is invoked and 0 is returned. - Function: int gsl_permutation_swap (gsl_permutation * P, const size_t I, const size_t J) This function exchanges the I-th and J-th elements of the permutation P.  File: gsl-ref.info, Node: Permutation properties, Next: Permutation functions, Prev: Accessing permutation elements, Up: Permutations Permutation properties ====================== - Function: size_t gsl_permutation_size (const gsl_permutation * P) This function returns the size of the permutation P. - Function: size_t * gsl_permutation_data (const gsl_permutation * P) This function returns a pointer to the array of elements in the permutation P. - Function: int gsl_permutation_valid (gsl_permutation * P) This function checks that the permutation P is valid. The N elements should contain each of the numbers 0 .. N-1 once and only once.