🚀 go-pugleaf

g++ version 1.36.0-, sun 4/110, SunOS 4.0
The following file causes the abort() at line 1201 of cplus-decl.c to
be called from function finish_file() in gcc-cc1plus, resulting in
compiler death from signal 6.

--------------- begin lin_solve.c ----------------
/*
routine to solve the d-dimensional linear system cx = b using Gaussian
elimination with scaled column pivoting. The algorithm used here
follows that given in Burden and Faires, pp. 326-328.
*/

/*
Should only include qflt.h and remes.h when variable-sized arrays work
*/
#include "poly.h"
#include "remes.h"

void lin_solve(qflt c[][max_deg + 2], qflt x[], qflt b[], int d)

{
/* Static to avoid construction every time */
static	qflt	a[max_deg + 2][max_deg + 3],        /* augmented matrix */
		s[max_deg + 2],			    /* scale factors */
		absa, mx;
int	i, j, k, p,
	rp[d];              /* row pointers */

/*
Copy the entire mess into an augmented matrix for two reasons:
   1) We don't want to sully the originals.
   2) We want nice uniformity in the algorithm.
At the same time, we compute the maximal row entries.
*/

for (i = 0; i < d; i++)
  {
    s[i] = 0.0;
    for (j = 0; j < d; j++)
      {
	a[i][j] = c[i][j];
	if ((absa = fabs(a[i][j])) > s[i]) s[i] = absa;
      }
    if (s[i] == qflt_zero) printf("lin_solve: singular!\n");
    a[i][d] = b[i];
    rp[i] = i;
  }

for (i = 0; i < d - 1; i++)
  {
    /*
    Find pivot row.
    */
    mx = qflt_zero;
    for (j = i; j < d; j++)
      if (mx < (absa = fabs(a[rp[j]][i] / s[rp[i]])))
	{
	  p = j;
	  mx = absa;
	}
    if ((mx == qflt_zero) || (a[rp[p]][i] == qflt_zero))
      printf("lin_solve: singular!\n");
    if (rp[i] != rp[p])
      {
	k = rp[i];
	rp[i] = rp[p];
	rp[p] = k;
      }
    for (j = i + 1; j < d; j++)
      {
	absa = a[rp[j]][i] / a[rp[i]][i];
	for (k = 0; k <= d; k++)
	  a[rp[j]][k] -= absa * a[rp[i]][k];
      }
  }
if (a[rp[d-1]][d-1] == qflt_zero) printf("lin_solve: singular!\n");
x[d-1] = a[rp[d-1]][d] / a[rp[d-1]][d-1];
for (i = d - 2; i >= 0; i--)
  {
    absa = qflt_zero;
    for (j = i + 1; j < d; j++)
      absa += a[rp[i]][j] * x[j];
    x[i] = (a[rp[i]][d] - absa) / a[rp[i]][i];
  }
}
-------------- end lin_solve.c -------------
-------------- begin poly.h    -------------
/*
poly.h

class declaration for polynomial functions
*/

#include "../qflt.h"

const int max_deg = 3;

class p_func;

class poly
{
  int d;
  qflt coeff[max_deg + 1];

 public:
  inline poly(int);
  qflt& operator [] (int);
  qflt operator () (qflt&);
  poly operator * (poly&);
  poly operator + (poly&);
  poly operator - ();
  poly operator - (poly&);
  friend ostream& operator << (ostream&, poly&);
  friend istream& operator >> (istream&, poly&);
  friend poly poly_deriv(poly&);
  friend qflt& p_func::operator [] (int);
  friend p_func::p_func(int, int);
  friend qflt p_func::coeff_ptl(qflt&, int);
  friend p_func p_deriv(p_func&);
  friend void bld_nguess(pqf, pqf, p_func&, int, qflt[], qflt&);
  friend int main(int, char *[]);
  friend void bld_r_iguess(pqf, pqf, p_func&, int, qflt[], qflt&, qflt&);
};

/*
Construct a polynomial of given degree.
*/
inline poly::poly(int deg)
{
d = deg;
for (int i = 0; i <= d; i++)
  coeff[i] = qflt_zero;
}
----------------- end poly.h ------------------
----------------- begin qflt.h ----------------
overload fabs;

#include <math.h>
double quiet_NaN(long);   /* Should be in math.h but isn't */
#include <stream.h>

inline int min(int a, int b)
{
  return a < b ? a : b;
}

inline int max(int a, int b)
{
  return a > b ? a : b;
}

// Max with absolute value
inline double maxabs(double& a, double& b)
{
  if (fabs(a) > fabs(b)) return a;
  else return b;
}

// Min with absolute value
inline double minabs(double& a, double& b)
{
  if (fabs(a) < fabs(b)) return a;
  else return b;
}

class qflt
{
  double high, low;
  void split(double, double&, double&);
  qflt exdblmul(double&, double&);
 public:
  qflt(double, double = 0.0);
  inline qflt();
  qflt operator + (qflt&);
  qflt operator - (qflt&);
  qflt operator - ();
  qflt operator * (qflt&);
  qflt operator / (qflt&);
  inline void operator *= (qflt&);
  inline int operator == (qflt&);
  inline int operator != (qflt&);
  inline int operator < (qflt&);
  inline int operator > (qflt&);
  inline int operator <= (qflt&);
  inline int operator >= (qflt&);
  inline double operator double();
  friend ostream& operator << (ostream&, qflt&);
  friend istream& operator >> (istream&, qflt&);
};

extern qflt qflt_zero;
extern qflt qflt_one;

inline qflt::qflt()
{
}

inline double qflt::operator double()
{
  return(high + low);
}

inline void qflt::operator *= (qflt& a)
{
  *this = *this * a;
}

inline int qflt::operator == (qflt& a)
{
  return (high == a.high) && (low == a.low);
}

inline int qflt::operator != (qflt& a)
{
  return (high != a.high) || (low != a.low);
}

inline int qflt::operator < (qflt& a)
{
  return (high < a.high) || ((high == a.high) && (low < a.low));
}

inline int qflt::operator > (qflt& a)
{
  return (high > a.high) || ((high == a.high) && (low > a.low));
}

inline int qflt::operator <= (qflt& a)
{
  return (high < a.high) || ((high == a.high) && (low <= a.low));
}

inline int qflt::operator >= (qflt& a)
{
  return (high > a.high) || ((high == a.high) && (low >= a.low));
}

inline qflt fabs(qflt& x)
{
  return x < qflt_zero ? -x : x;
}

/* Declaration for pointer to qflt function */
typedef qflt (*pqf)(qflt&);
------------------- end qflt.h ---------------------
------------------- begin remes.h ------------------
const qflt macheps = qflt(1.0e-25, 0.0);
const int  max_dstep = 5;
const int  max_iter  = 50;
------------------- end remes.h --------------------

I have not tried to figure out why this is happening...

	-- Robert