1.1 --- /dev/null Thu Jan 01 00:00:00 1970 +0000
1.2 +++ b/src/glphbm.c Mon Dec 06 13:09:21 2010 +0100
1.3 @@ -0,0 +1,519 @@
1.4 +/* glphbm.c */
1.5 +
1.6 +/***********************************************************************
1.7 +* This code is part of GLPK (GNU Linear Programming Kit).
1.8 +*
1.9 +* Copyright (C) 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2008,
1.10 +* 2009, 2010 Andrew Makhorin, Department for Applied Informatics,
1.11 +* Moscow Aviation Institute, Moscow, Russia. All rights reserved.
1.12 +* E-mail: <mao@gnu.org>.
1.13 +*
1.14 +* GLPK is free software: you can redistribute it and/or modify it
1.15 +* under the terms of the GNU General Public License as published by
1.16 +* the Free Software Foundation, either version 3 of the License, or
1.17 +* (at your option) any later version.
1.18 +*
1.19 +* GLPK is distributed in the hope that it will be useful, but WITHOUT
1.20 +* ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
1.21 +* or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public
1.22 +* License for more details.
1.23 +*
1.24 +* You should have received a copy of the GNU General Public License
1.25 +* along with GLPK. If not, see <http://www.gnu.org/licenses/>.
1.26 +***********************************************************************/
1.27 +
1.28 +#define _GLPSTD_ERRNO
1.29 +#define _GLPSTD_STDIO
1.30 +#include "glphbm.h"
1.31 +#include "glpenv.h"
1.32 +
1.33 +/***********************************************************************
1.34 +* NAME
1.35 +*
1.36 +* hbm_read_mat - read sparse matrix in Harwell-Boeing format
1.37 +*
1.38 +* SYNOPSIS
1.39 +*
1.40 +* #include "glphbm.h"
1.41 +* HBM *hbm_read_mat(const char *fname);
1.42 +*
1.43 +* DESCRIPTION
1.44 +*
1.45 +* The routine hbm_read_mat reads a sparse matrix in the Harwell-Boeing
1.46 +* format from a text file whose name is the character string fname.
1.47 +*
1.48 +* Detailed description of the Harwell-Boeing format recognised by this
1.49 +* routine is given in the following report:
1.50 +*
1.51 +* I.S.Duff, R.G.Grimes, J.G.Lewis. User's Guide for the Harwell-Boeing
1.52 +* Sparse Matrix Collection (Release I), TR/PA/92/86, October 1992.
1.53 +*
1.54 +* RETURNS
1.55 +*
1.56 +* If no error occured, the routine hbm_read_mat returns a pointer to
1.57 +* a data structure containing the matrix. In case of error the routine
1.58 +* prints an appropriate error message and returns NULL. */
1.59 +
1.60 +struct dsa
1.61 +{ /* working area used by routine hbm_read_mat */
1.62 + const char *fname;
1.63 + /* name of input text file */
1.64 + FILE *fp;
1.65 + /* stream assigned to input text file */
1.66 + int seqn;
1.67 + /* card sequential number */
1.68 + char card[80+1];
1.69 + /* card image buffer */
1.70 + int fmt_p;
1.71 + /* scale factor */
1.72 + int fmt_k;
1.73 + /* iterator */
1.74 + int fmt_f;
1.75 + /* format code */
1.76 + int fmt_w;
1.77 + /* field width */
1.78 + int fmt_d;
1.79 + /* number of decimal places after point */
1.80 +};
1.81 +
1.82 +/***********************************************************************
1.83 +* read_card - read next data card
1.84 +*
1.85 +* This routine reads the next 80-column card from the input text file
1.86 +* and stores its image into the character string card. If the card was
1.87 +* read successfully, the routine returns zero, otherwise non-zero. */
1.88 +
1.89 +static int read_card(struct dsa *dsa)
1.90 +{ int k, c;
1.91 + dsa->seqn++;
1.92 + memset(dsa->card, ' ', 80), dsa->card[80] = '\0';
1.93 + k = 0;
1.94 + for (;;)
1.95 + { c = fgetc(dsa->fp);
1.96 + if (ferror(dsa->fp))
1.97 + { xprintf("%s:%d: read error - %s\n", dsa->fname, dsa->seqn,
1.98 + strerror(errno));
1.99 + return 1;
1.100 + }
1.101 + if (feof(dsa->fp))
1.102 + { if (k == 0)
1.103 + xprintf("%s:%d: unexpected EOF\n", dsa->fname,
1.104 + dsa->seqn);
1.105 + else
1.106 + xprintf("%s:%d: missing final LF\n", dsa->fname,
1.107 + dsa->seqn);
1.108 + return 1;
1.109 + }
1.110 + if (c == '\r') continue;
1.111 + if (c == '\n') break;
1.112 + if (iscntrl(c))
1.113 + { xprintf("%s:%d: invalid control character 0x%02X\n",
1.114 + dsa->fname, dsa->seqn, c);
1.115 + return 1;
1.116 + }
1.117 + if (k == 80)
1.118 + { xprintf("%s:%d: card image too long\n", dsa->fname,
1.119 + dsa->seqn);
1.120 + return 1;
1.121 + }
1.122 + dsa->card[k++] = (char)c;
1.123 + }
1.124 + return 0;
1.125 +}
1.126 +
1.127 +/***********************************************************************
1.128 +* scan_int - scan integer value from the current card
1.129 +*
1.130 +* This routine scans an integer value from the current card, where fld
1.131 +* is the name of the field, pos is the position of the field, width is
1.132 +* the width of the field, val points to a location to which the scanned
1.133 +* value should be stored. If the value was scanned successfully, the
1.134 +* routine returns zero, otherwise non-zero. */
1.135 +
1.136 +static int scan_int(struct dsa *dsa, char *fld, int pos, int width,
1.137 + int *val)
1.138 +{ char str[80+1];
1.139 + xassert(1 <= width && width <= 80);
1.140 + memcpy(str, dsa->card + pos, width), str[width] = '\0';
1.141 + if (str2int(strspx(str), val))
1.142 + { xprintf("%s:%d: field `%s' contains invalid value `%s'\n",
1.143 + dsa->fname, dsa->seqn, fld, str);
1.144 + return 1;
1.145 + }
1.146 + return 0;
1.147 +}
1.148 +
1.149 +/***********************************************************************
1.150 +* parse_fmt - parse Fortran format specification
1.151 +*
1.152 +* This routine parses the Fortran format specification represented as
1.153 +* character string which fmt points to and stores format elements into
1.154 +* appropriate static locations. Should note that not all valid Fortran
1.155 +* format specifications may be recognised. If the format specification
1.156 +* was recognised, the routine returns zero, otherwise non-zero. */
1.157 +
1.158 +static int parse_fmt(struct dsa *dsa, char *fmt)
1.159 +{ int k, s, val;
1.160 + char str[80+1];
1.161 + /* first character should be left parenthesis */
1.162 + if (fmt[0] != '(')
1.163 +fail: { xprintf("hbm_read_mat: format `%s' not recognised\n", fmt);
1.164 + return 1;
1.165 + }
1.166 + k = 1;
1.167 + /* optional scale factor */
1.168 + dsa->fmt_p = 0;
1.169 + if (isdigit((unsigned char)fmt[k]))
1.170 + { s = 0;
1.171 + while (isdigit((unsigned char)fmt[k]))
1.172 + { if (s == 80) goto fail;
1.173 + str[s++] = fmt[k++];
1.174 + }
1.175 + str[s] = '\0';
1.176 + if (str2int(str, &val)) goto fail;
1.177 + if (toupper((unsigned char)fmt[k]) != 'P') goto iter;
1.178 + dsa->fmt_p = val, k++;
1.179 + if (!(0 <= dsa->fmt_p && dsa->fmt_p <= 255)) goto fail;
1.180 + /* optional comma may follow scale factor */
1.181 + if (fmt[k] == ',') k++;
1.182 + }
1.183 + /* optional iterator */
1.184 + dsa->fmt_k = 1;
1.185 + if (isdigit((unsigned char)fmt[k]))
1.186 + { s = 0;
1.187 + while (isdigit((unsigned char)fmt[k]))
1.188 + { if (s == 80) goto fail;
1.189 + str[s++] = fmt[k++];
1.190 + }
1.191 + str[s] = '\0';
1.192 + if (str2int(str, &val)) goto fail;
1.193 +iter: dsa->fmt_k = val;
1.194 + if (!(1 <= dsa->fmt_k && dsa->fmt_k <= 255)) goto fail;
1.195 + }
1.196 + /* format code */
1.197 + dsa->fmt_f = toupper((unsigned char)fmt[k++]);
1.198 + if (!(dsa->fmt_f == 'D' || dsa->fmt_f == 'E' ||
1.199 + dsa->fmt_f == 'F' || dsa->fmt_f == 'G' ||
1.200 + dsa->fmt_f == 'I')) goto fail;
1.201 + /* field width */
1.202 + if (!isdigit((unsigned char)fmt[k])) goto fail;
1.203 + s = 0;
1.204 + while (isdigit((unsigned char)fmt[k]))
1.205 + { if (s == 80) goto fail;
1.206 + str[s++] = fmt[k++];
1.207 + }
1.208 + str[s] = '\0';
1.209 + if (str2int(str, &dsa->fmt_w)) goto fail;
1.210 + if (!(1 <= dsa->fmt_w && dsa->fmt_w <= 255)) goto fail;
1.211 + /* optional number of decimal places after point */
1.212 + dsa->fmt_d = 0;
1.213 + if (fmt[k] == '.')
1.214 + { k++;
1.215 + if (!isdigit((unsigned char)fmt[k])) goto fail;
1.216 + s = 0;
1.217 + while (isdigit((unsigned char)fmt[k]))
1.218 + { if (s == 80) goto fail;
1.219 + str[s++] = fmt[k++];
1.220 + }
1.221 + str[s] = '\0';
1.222 + if (str2int(str, &dsa->fmt_d)) goto fail;
1.223 + if (!(0 <= dsa->fmt_d && dsa->fmt_d <= 255)) goto fail;
1.224 + }
1.225 + /* last character should be right parenthesis */
1.226 + if (!(fmt[k] == ')' && fmt[k+1] == '\0')) goto fail;
1.227 + return 0;
1.228 +}
1.229 +
1.230 +/***********************************************************************
1.231 +* read_int_array - read array of integer type
1.232 +*
1.233 +* This routine reads an integer array from the input text file, where
1.234 +* name is array name, fmt is Fortran format specification that controls
1.235 +* reading, n is number of array elements, val is array of integer type.
1.236 +* If the array was read successful, the routine returns zero, otherwise
1.237 +* non-zero. */
1.238 +
1.239 +static int read_int_array(struct dsa *dsa, char *name, char *fmt,
1.240 + int n, int val[])
1.241 +{ int k, pos;
1.242 + char str[80+1];
1.243 + if (parse_fmt(dsa, fmt)) return 1;
1.244 + if (!(dsa->fmt_f == 'I' && dsa->fmt_w <= 80 &&
1.245 + dsa->fmt_k * dsa->fmt_w <= 80))
1.246 + { xprintf(
1.247 + "%s:%d: can't read array `%s' - invalid format `%s'\n",
1.248 + dsa->fname, dsa->seqn, name, fmt);
1.249 + return 1;
1.250 + }
1.251 + for (k = 1, pos = INT_MAX; k <= n; k++, pos++)
1.252 + { if (pos >= dsa->fmt_k)
1.253 + { if (read_card(dsa)) return 1;
1.254 + pos = 0;
1.255 + }
1.256 + memcpy(str, dsa->card + dsa->fmt_w * pos, dsa->fmt_w);
1.257 + str[dsa->fmt_w] = '\0';
1.258 + strspx(str);
1.259 + if (str2int(str, &val[k]))
1.260 + { xprintf(
1.261 + "%s:%d: can't read array `%s' - invalid value `%s'\n",
1.262 + dsa->fname, dsa->seqn, name, str);
1.263 + return 1;
1.264 + }
1.265 + }
1.266 + return 0;
1.267 +}
1.268 +
1.269 +/***********************************************************************
1.270 +* read_real_array - read array of real type
1.271 +*
1.272 +* This routine reads a real array from the input text file, where name
1.273 +* is array name, fmt is Fortran format specification that controls
1.274 +* reading, n is number of array elements, val is array of real type.
1.275 +* If the array was read successful, the routine returns zero, otherwise
1.276 +* non-zero. */
1.277 +
1.278 +static int read_real_array(struct dsa *dsa, char *name, char *fmt,
1.279 + int n, double val[])
1.280 +{ int k, pos;
1.281 + char str[80+1], *ptr;
1.282 + if (parse_fmt(dsa, fmt)) return 1;
1.283 + if (!(dsa->fmt_f != 'I' && dsa->fmt_w <= 80 &&
1.284 + dsa->fmt_k * dsa->fmt_w <= 80))
1.285 + { xprintf(
1.286 + "%s:%d: can't read array `%s' - invalid format `%s'\n",
1.287 + dsa->fname, dsa->seqn, name, fmt);
1.288 + return 1;
1.289 + }
1.290 + for (k = 1, pos = INT_MAX; k <= n; k++, pos++)
1.291 + { if (pos >= dsa->fmt_k)
1.292 + { if (read_card(dsa)) return 1;
1.293 + pos = 0;
1.294 + }
1.295 + memcpy(str, dsa->card + dsa->fmt_w * pos, dsa->fmt_w);
1.296 + str[dsa->fmt_w] = '\0';
1.297 + strspx(str);
1.298 + if (strchr(str, '.') == NULL && strcmp(str, "0"))
1.299 + { xprintf("%s(%d): can't read array `%s' - value `%s' has no "
1.300 + "decimal point\n", dsa->fname, dsa->seqn, name, str);
1.301 + return 1;
1.302 + }
1.303 + /* sometimes lower case letters appear */
1.304 + for (ptr = str; *ptr; ptr++)
1.305 + *ptr = (char)toupper((unsigned char)*ptr);
1.306 + ptr = strchr(str, 'D');
1.307 + if (ptr != NULL) *ptr = 'E';
1.308 + /* value may appear with decimal exponent but without letters
1.309 + E or D (for example, -123.456-012), so missing letter should
1.310 + be inserted */
1.311 + ptr = strchr(str+1, '+');
1.312 + if (ptr == NULL) ptr = strchr(str+1, '-');
1.313 + if (ptr != NULL && *(ptr-1) != 'E')
1.314 + { xassert(strlen(str) < 80);
1.315 + memmove(ptr+1, ptr, strlen(ptr)+1);
1.316 + *ptr = 'E';
1.317 + }
1.318 + if (str2num(str, &val[k]))
1.319 + { xprintf(
1.320 + "%s:%d: can't read array `%s' - invalid value `%s'\n",
1.321 + dsa->fname, dsa->seqn, name, str);
1.322 + return 1;
1.323 + }
1.324 + }
1.325 + return 0;
1.326 +}
1.327 +
1.328 +HBM *hbm_read_mat(const char *fname)
1.329 +{ struct dsa _dsa, *dsa = &_dsa;
1.330 + HBM *hbm = NULL;
1.331 + dsa->fname = fname;
1.332 + xprintf("hbm_read_mat: reading matrix from `%s'...\n",
1.333 + dsa->fname);
1.334 + dsa->fp = fopen(dsa->fname, "r");
1.335 + if (dsa->fp == NULL)
1.336 + { xprintf("hbm_read_mat: unable to open `%s' - %s\n",
1.337 + dsa->fname, strerror(errno));
1.338 + goto fail;
1.339 + }
1.340 + dsa->seqn = 0;
1.341 + hbm = xmalloc(sizeof(HBM));
1.342 + memset(hbm, 0, sizeof(HBM));
1.343 + /* read the first heading card */
1.344 + if (read_card(dsa)) goto fail;
1.345 + memcpy(hbm->title, dsa->card, 72), hbm->title[72] = '\0';
1.346 + strtrim(hbm->title);
1.347 + xprintf("%s\n", hbm->title);
1.348 + memcpy(hbm->key, dsa->card+72, 8), hbm->key[8] = '\0';
1.349 + strspx(hbm->key);
1.350 + xprintf("key = %s\n", hbm->key);
1.351 + /* read the second heading card */
1.352 + if (read_card(dsa)) goto fail;
1.353 + if (scan_int(dsa, "totcrd", 0, 14, &hbm->totcrd)) goto fail;
1.354 + if (scan_int(dsa, "ptrcrd", 14, 14, &hbm->ptrcrd)) goto fail;
1.355 + if (scan_int(dsa, "indcrd", 28, 14, &hbm->indcrd)) goto fail;
1.356 + if (scan_int(dsa, "valcrd", 42, 14, &hbm->valcrd)) goto fail;
1.357 + if (scan_int(dsa, "rhscrd", 56, 14, &hbm->rhscrd)) goto fail;
1.358 + xprintf("totcrd = %d; ptrcrd = %d; indcrd = %d; valcrd = %d; rhsc"
1.359 + "rd = %d\n", hbm->totcrd, hbm->ptrcrd, hbm->indcrd,
1.360 + hbm->valcrd, hbm->rhscrd);
1.361 + /* read the third heading card */
1.362 + if (read_card(dsa)) goto fail;
1.363 + memcpy(hbm->mxtype, dsa->card, 3), hbm->mxtype[3] = '\0';
1.364 + if (strchr("RCP", hbm->mxtype[0]) == NULL ||
1.365 + strchr("SUHZR", hbm->mxtype[1]) == NULL ||
1.366 + strchr("AE", hbm->mxtype[2]) == NULL)
1.367 + { xprintf("%s:%d: matrix type `%s' not recognised\n",
1.368 + dsa->fname, dsa->seqn, hbm->mxtype);
1.369 + goto fail;
1.370 + }
1.371 + if (scan_int(dsa, "nrow", 14, 14, &hbm->nrow)) goto fail;
1.372 + if (scan_int(dsa, "ncol", 28, 14, &hbm->ncol)) goto fail;
1.373 + if (scan_int(dsa, "nnzero", 42, 14, &hbm->nnzero)) goto fail;
1.374 + if (scan_int(dsa, "neltvl", 56, 14, &hbm->neltvl)) goto fail;
1.375 + xprintf("mxtype = %s; nrow = %d; ncol = %d; nnzero = %d; neltvl ="
1.376 + " %d\n", hbm->mxtype, hbm->nrow, hbm->ncol, hbm->nnzero,
1.377 + hbm->neltvl);
1.378 + /* read the fourth heading card */
1.379 + if (read_card(dsa)) goto fail;
1.380 + memcpy(hbm->ptrfmt, dsa->card, 16), hbm->ptrfmt[16] = '\0';
1.381 + strspx(hbm->ptrfmt);
1.382 + memcpy(hbm->indfmt, dsa->card+16, 16), hbm->indfmt[16] = '\0';
1.383 + strspx(hbm->indfmt);
1.384 + memcpy(hbm->valfmt, dsa->card+32, 20), hbm->valfmt[20] = '\0';
1.385 + strspx(hbm->valfmt);
1.386 + memcpy(hbm->rhsfmt, dsa->card+52, 20), hbm->rhsfmt[20] = '\0';
1.387 + strspx(hbm->rhsfmt);
1.388 + xprintf("ptrfmt = %s; indfmt = %s; valfmt = %s; rhsfmt = %s\n",
1.389 + hbm->ptrfmt, hbm->indfmt, hbm->valfmt, hbm->rhsfmt);
1.390 + /* read the fifth heading card (optional) */
1.391 + if (hbm->rhscrd <= 0)
1.392 + { strcpy(hbm->rhstyp, "???");
1.393 + hbm->nrhs = 0;
1.394 + hbm->nrhsix = 0;
1.395 + }
1.396 + else
1.397 + { if (read_card(dsa)) goto fail;
1.398 + memcpy(hbm->rhstyp, dsa->card, 3), hbm->rhstyp[3] = '\0';
1.399 + if (scan_int(dsa, "nrhs", 14, 14, &hbm->nrhs)) goto fail;
1.400 + if (scan_int(dsa, "nrhsix", 28, 14, &hbm->nrhsix)) goto fail;
1.401 + xprintf("rhstyp = `%s'; nrhs = %d; nrhsix = %d\n",
1.402 + hbm->rhstyp, hbm->nrhs, hbm->nrhsix);
1.403 + }
1.404 + /* read matrix structure */
1.405 + hbm->colptr = xcalloc(1+hbm->ncol+1, sizeof(int));
1.406 + if (read_int_array(dsa, "colptr", hbm->ptrfmt, hbm->ncol+1,
1.407 + hbm->colptr)) goto fail;
1.408 + hbm->rowind = xcalloc(1+hbm->nnzero, sizeof(int));
1.409 + if (read_int_array(dsa, "rowind", hbm->indfmt, hbm->nnzero,
1.410 + hbm->rowind)) goto fail;
1.411 + /* read matrix values */
1.412 + if (hbm->valcrd <= 0) goto done;
1.413 + if (hbm->mxtype[2] == 'A')
1.414 + { /* assembled matrix */
1.415 + hbm->values = xcalloc(1+hbm->nnzero, sizeof(double));
1.416 + if (read_real_array(dsa, "values", hbm->valfmt, hbm->nnzero,
1.417 + hbm->values)) goto fail;
1.418 + }
1.419 + else
1.420 + { /* elemental (unassembled) matrix */
1.421 + hbm->values = xcalloc(1+hbm->neltvl, sizeof(double));
1.422 + if (read_real_array(dsa, "values", hbm->valfmt, hbm->neltvl,
1.423 + hbm->values)) goto fail;
1.424 + }
1.425 + /* read right-hand sides */
1.426 + if (hbm->nrhs <= 0) goto done;
1.427 + if (hbm->rhstyp[0] == 'F')
1.428 + { /* dense format */
1.429 + hbm->nrhsvl = hbm->nrow * hbm->nrhs;
1.430 + hbm->rhsval = xcalloc(1+hbm->nrhsvl, sizeof(double));
1.431 + if (read_real_array(dsa, "rhsval", hbm->rhsfmt, hbm->nrhsvl,
1.432 + hbm->rhsval)) goto fail;
1.433 + }
1.434 + else if (hbm->rhstyp[0] == 'M' && hbm->mxtype[2] == 'A')
1.435 + { /* sparse format */
1.436 + /* read pointers */
1.437 + hbm->rhsptr = xcalloc(1+hbm->nrhs+1, sizeof(int));
1.438 + if (read_int_array(dsa, "rhsptr", hbm->ptrfmt, hbm->nrhs+1,
1.439 + hbm->rhsptr)) goto fail;
1.440 + /* read sparsity pattern */
1.441 + hbm->rhsind = xcalloc(1+hbm->nrhsix, sizeof(int));
1.442 + if (read_int_array(dsa, "rhsind", hbm->indfmt, hbm->nrhsix,
1.443 + hbm->rhsind)) goto fail;
1.444 + /* read values */
1.445 + hbm->rhsval = xcalloc(1+hbm->nrhsix, sizeof(double));
1.446 + if (read_real_array(dsa, "rhsval", hbm->rhsfmt, hbm->nrhsix,
1.447 + hbm->rhsval)) goto fail;
1.448 + }
1.449 + else if (hbm->rhstyp[0] == 'M' && hbm->mxtype[2] == 'E')
1.450 + { /* elemental format */
1.451 + hbm->rhsval = xcalloc(1+hbm->nrhsvl, sizeof(double));
1.452 + if (read_real_array(dsa, "rhsval", hbm->rhsfmt, hbm->nrhsvl,
1.453 + hbm->rhsval)) goto fail;
1.454 + }
1.455 + else
1.456 + { xprintf("%s:%d: right-hand side type `%c' not recognised\n",
1.457 + dsa->fname, dsa->seqn, hbm->rhstyp[0]);
1.458 + goto fail;
1.459 + }
1.460 + /* read starting guesses */
1.461 + if (hbm->rhstyp[1] == 'G')
1.462 + { hbm->nguess = hbm->nrow * hbm->nrhs;
1.463 + hbm->sguess = xcalloc(1+hbm->nguess, sizeof(double));
1.464 + if (read_real_array(dsa, "sguess", hbm->rhsfmt, hbm->nguess,
1.465 + hbm->sguess)) goto fail;
1.466 + }
1.467 + /* read solution vectors */
1.468 + if (hbm->rhstyp[2] == 'X')
1.469 + { hbm->nexact = hbm->nrow * hbm->nrhs;
1.470 + hbm->xexact = xcalloc(1+hbm->nexact, sizeof(double));
1.471 + if (read_real_array(dsa, "xexact", hbm->rhsfmt, hbm->nexact,
1.472 + hbm->xexact)) goto fail;
1.473 + }
1.474 +done: /* reading has been completed */
1.475 + xprintf("hbm_read_mat: %d cards were read\n", dsa->seqn);
1.476 + fclose(dsa->fp);
1.477 + return hbm;
1.478 +fail: /* something wrong in Danish kingdom */
1.479 + if (hbm != NULL)
1.480 + { if (hbm->colptr != NULL) xfree(hbm->colptr);
1.481 + if (hbm->rowind != NULL) xfree(hbm->rowind);
1.482 + if (hbm->rhsptr != NULL) xfree(hbm->rhsptr);
1.483 + if (hbm->rhsind != NULL) xfree(hbm->rhsind);
1.484 + if (hbm->values != NULL) xfree(hbm->values);
1.485 + if (hbm->rhsval != NULL) xfree(hbm->rhsval);
1.486 + if (hbm->sguess != NULL) xfree(hbm->sguess);
1.487 + if (hbm->xexact != NULL) xfree(hbm->xexact);
1.488 + xfree(hbm);
1.489 + }
1.490 + if (dsa->fp != NULL) fclose(dsa->fp);
1.491 + return NULL;
1.492 +}
1.493 +
1.494 +/***********************************************************************
1.495 +* NAME
1.496 +*
1.497 +* hbm_free_mat - free sparse matrix in Harwell-Boeing format
1.498 +*
1.499 +* SYNOPSIS
1.500 +*
1.501 +* #include "glphbm.h"
1.502 +* void hbm_free_mat(HBM *hbm);
1.503 +*
1.504 +* DESCRIPTION
1.505 +*
1.506 +* The hbm_free_mat routine frees all the memory allocated to the data
1.507 +* structure containing a sparse matrix in the Harwell-Boeing format. */
1.508 +
1.509 +void hbm_free_mat(HBM *hbm)
1.510 +{ if (hbm->colptr != NULL) xfree(hbm->colptr);
1.511 + if (hbm->rowind != NULL) xfree(hbm->rowind);
1.512 + if (hbm->rhsptr != NULL) xfree(hbm->rhsptr);
1.513 + if (hbm->rhsind != NULL) xfree(hbm->rhsind);
1.514 + if (hbm->values != NULL) xfree(hbm->values);
1.515 + if (hbm->rhsval != NULL) xfree(hbm->rhsval);
1.516 + if (hbm->sguess != NULL) xfree(hbm->sguess);
1.517 + if (hbm->xexact != NULL) xfree(hbm->xexact);
1.518 + xfree(hbm);
1.519 + return;
1.520 +}
1.521 +
1.522 +/* eof */