alpar@9: /* glpmat.h (linear algebra routines) */ alpar@9: alpar@9: /*********************************************************************** alpar@9: * This code is part of GLPK (GNU Linear Programming Kit). alpar@9: * alpar@9: * Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008, alpar@9: * 2009, 2010, 2011 Andrew Makhorin, Department for Applied Informatics, alpar@9: * Moscow Aviation Institute, Moscow, Russia. All rights reserved. alpar@9: * E-mail: . alpar@9: * alpar@9: * GLPK is free software: you can redistribute it and/or modify it alpar@9: * under the terms of the GNU General Public License as published by alpar@9: * the Free Software Foundation, either version 3 of the License, or alpar@9: * (at your option) any later version. alpar@9: * alpar@9: * GLPK is distributed in the hope that it will be useful, but WITHOUT alpar@9: * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY alpar@9: * or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public alpar@9: * License for more details. alpar@9: * alpar@9: * You should have received a copy of the GNU General Public License alpar@9: * along with GLPK. If not, see . alpar@9: ***********************************************************************/ alpar@9: alpar@9: #ifndef GLPMAT_H alpar@9: #define GLPMAT_H alpar@9: alpar@9: /*********************************************************************** alpar@9: * FULL-VECTOR STORAGE alpar@9: * alpar@9: * For a sparse vector x having n elements, ne of which are non-zero, alpar@9: * the full-vector storage format uses two arrays x_ind and x_vec, which alpar@9: * are set up as follows: alpar@9: * alpar@9: * x_ind is an integer array of length [1+ne]. Location x_ind[0] is alpar@9: * not used, and locations x_ind[1], ..., x_ind[ne] contain indices of alpar@9: * non-zero elements in vector x. alpar@9: * alpar@9: * x_vec is a floating-point array of length [1+n]. Location x_vec[0] alpar@9: * is not used, and locations x_vec[1], ..., x_vec[n] contain numeric alpar@9: * values of ALL elements in vector x, including its zero elements. alpar@9: * alpar@9: * Let, for example, the following sparse vector x be given: alpar@9: * alpar@9: * (0, 1, 0, 0, 2, 3, 0, 4) alpar@9: * alpar@9: * Then the arrays are: alpar@9: * alpar@9: * x_ind = { X; 2, 5, 6, 8 } alpar@9: * alpar@9: * x_vec = { X; 0, 1, 0, 0, 2, 3, 0, 4 } alpar@9: * alpar@9: * COMPRESSED-VECTOR STORAGE alpar@9: * alpar@9: * For a sparse vector x having n elements, ne of which are non-zero, alpar@9: * the compressed-vector storage format uses two arrays x_ind and x_vec, alpar@9: * which are set up as follows: alpar@9: * alpar@9: * x_ind is an integer array of length [1+ne]. Location x_ind[0] is alpar@9: * not used, and locations x_ind[1], ..., x_ind[ne] contain indices of alpar@9: * non-zero elements in vector x. alpar@9: * alpar@9: * x_vec is a floating-point array of length [1+ne]. Location x_vec[0] alpar@9: * is not used, and locations x_vec[1], ..., x_vec[ne] contain numeric alpar@9: * values of corresponding non-zero elements in vector x. alpar@9: * alpar@9: * Let, for example, the following sparse vector x be given: alpar@9: * alpar@9: * (0, 1, 0, 0, 2, 3, 0, 4) alpar@9: * alpar@9: * Then the arrays are: alpar@9: * alpar@9: * x_ind = { X; 2, 5, 6, 8 } alpar@9: * alpar@9: * x_vec = { X; 1, 2, 3, 4 } alpar@9: * alpar@9: * STORAGE-BY-ROWS alpar@9: * alpar@9: * For a sparse matrix A, which has m rows, n columns, and ne non-zero alpar@9: * elements the storage-by-rows format uses three arrays A_ptr, A_ind, alpar@9: * and A_val, which are set up as follows: alpar@9: * alpar@9: * A_ptr is an integer array of length [1+m+1] also called "row pointer alpar@9: * array". It contains the relative starting positions of each row of A alpar@9: * in the arrays A_ind and A_val, i.e. element A_ptr[i], 1 <= i <= m, alpar@9: * indicates where row i begins in the arrays A_ind and A_val. If all alpar@9: * elements in row i are zero, then A_ptr[i] = A_ptr[i+1]. Location alpar@9: * A_ptr[0] is not used, location A_ptr[1] must contain 1, and location alpar@9: * A_ptr[m+1] must contain ne+1 that indicates the position after the alpar@9: * last element in the arrays A_ind and A_val. alpar@9: * alpar@9: * A_ind is an integer array of length [1+ne]. Location A_ind[0] is not alpar@9: * used, and locations A_ind[1], ..., A_ind[ne] contain column indices alpar@9: * of (non-zero) elements in matrix A. alpar@9: * alpar@9: * A_val is a floating-point array of length [1+ne]. Location A_val[0] alpar@9: * is not used, and locations A_val[1], ..., A_val[ne] contain numeric alpar@9: * values of non-zero elements in matrix A. alpar@9: * alpar@9: * Non-zero elements of matrix A are stored contiguously, and the rows alpar@9: * of matrix A are stored consecutively from 1 to m in the arrays A_ind alpar@9: * and A_val. The elements in each row of A may be stored in any order alpar@9: * in A_ind and A_val. Note that elements with duplicate column indices alpar@9: * are not allowed. alpar@9: * alpar@9: * Let, for example, the following sparse matrix A be given: alpar@9: * alpar@9: * | 11 . 13 . . . | alpar@9: * | 21 22 . 24 . . | alpar@9: * | . 32 33 . . . | alpar@9: * | . . 43 44 . 46 | alpar@9: * | . . . . . . | alpar@9: * | 61 62 . . . 66 | alpar@9: * alpar@9: * Then the arrays are: alpar@9: * alpar@9: * A_ptr = { X; 1, 3, 6, 8, 11, 11; 14 } alpar@9: * alpar@9: * A_ind = { X; 1, 3; 4, 2, 1; 2, 3; 4, 3, 6; 1, 2, 6 } alpar@9: * alpar@9: * A_val = { X; 11, 13; 24, 22, 21; 32, 33; 44, 43, 46; 61, 62, 66 } alpar@9: * alpar@9: * PERMUTATION MATRICES alpar@9: * alpar@9: * Let P be a permutation matrix of the order n. It is represented as alpar@9: * an integer array P_per of length [1+n+n] as follows: if p[i,j] = 1, alpar@9: * then P_per[i] = j and P_per[n+j] = i. Location P_per[0] is not used. alpar@9: * alpar@9: * Let A' = P*A. If i-th row of A corresponds to i'-th row of A', then alpar@9: * P_per[i'] = i and P_per[n+i] = i'. alpar@9: * alpar@9: * References: alpar@9: * alpar@9: * 1. Gustavson F.G. Some basic techniques for solving sparse systems of alpar@9: * linear equations. In Rose and Willoughby (1972), pp. 41-52. alpar@9: * alpar@9: * 2. Basic Linear Algebra Subprograms Technical (BLAST) Forum Standard. alpar@9: * University of Tennessee (2001). */ alpar@9: alpar@9: #define check_fvs _glp_mat_check_fvs alpar@9: int check_fvs(int n, int nnz, int ind[], double vec[]); alpar@9: /* check sparse vector in full-vector storage format */ alpar@9: alpar@9: #define check_pattern _glp_mat_check_pattern alpar@9: int check_pattern(int m, int n, int A_ptr[], int A_ind[]); alpar@9: /* check pattern of sparse matrix */ alpar@9: alpar@9: #define transpose _glp_mat_transpose alpar@9: void transpose(int m, int n, int A_ptr[], int A_ind[], double A_val[], alpar@9: int AT_ptr[], int AT_ind[], double AT_val[]); alpar@9: /* transpose sparse matrix */ alpar@9: alpar@9: #define adat_symbolic _glp_mat_adat_symbolic alpar@9: int *adat_symbolic(int m, int n, int P_per[], int A_ptr[], int A_ind[], alpar@9: int S_ptr[]); alpar@9: /* compute S = P*A*D*A'*P' (symbolic phase) */ alpar@9: alpar@9: #define adat_numeric _glp_mat_adat_numeric alpar@9: void adat_numeric(int m, int n, int P_per[], alpar@9: int A_ptr[], int A_ind[], double A_val[], double D_diag[], alpar@9: int S_ptr[], int S_ind[], double S_val[], double S_diag[]); alpar@9: /* compute S = P*A*D*A'*P' (numeric phase) */ alpar@9: alpar@9: #define min_degree _glp_mat_min_degree alpar@9: void min_degree(int n, int A_ptr[], int A_ind[], int P_per[]); alpar@9: /* minimum degree ordering */ alpar@9: alpar@9: #define amd_order1 _glp_mat_amd_order1 alpar@9: void amd_order1(int n, int A_ptr[], int A_ind[], int P_per[]); alpar@9: /* approximate minimum degree ordering (AMD) */ alpar@9: alpar@9: #define symamd_ord _glp_mat_symamd_ord alpar@9: void symamd_ord(int n, int A_ptr[], int A_ind[], int P_per[]); alpar@9: /* approximate minimum degree ordering (SYMAMD) */ alpar@9: alpar@9: #define chol_symbolic _glp_mat_chol_symbolic alpar@9: int *chol_symbolic(int n, int A_ptr[], int A_ind[], int U_ptr[]); alpar@9: /* compute Cholesky factorization (symbolic phase) */ alpar@9: alpar@9: #define chol_numeric _glp_mat_chol_numeric alpar@9: int chol_numeric(int n, alpar@9: int A_ptr[], int A_ind[], double A_val[], double A_diag[], alpar@9: int U_ptr[], int U_ind[], double U_val[], double U_diag[]); alpar@9: /* compute Cholesky factorization (numeric phase) */ alpar@9: alpar@9: #define u_solve _glp_mat_u_solve alpar@9: void u_solve(int n, int U_ptr[], int U_ind[], double U_val[], alpar@9: double U_diag[], double x[]); alpar@9: /* solve upper triangular system U*x = b */ alpar@9: alpar@9: #define ut_solve _glp_mat_ut_solve alpar@9: void ut_solve(int n, int U_ptr[], int U_ind[], double U_val[], alpar@9: double U_diag[], double x[]); alpar@9: /* solve lower triangular system U'*x = b */ alpar@9: alpar@9: #endif alpar@9: alpar@9: /* eof */