/* fd.f -- translated by f2c (version 20000817).
   You must link the resulting object file with the libraries:
	-lf2c -lm   (in that order)
*/

#include "f2c.h"
#include "blaswrap.h"

/* Table of constant values */

static integer c__1 = 1;

/* Copyright (C) 2000, International Business Machines */
/* Corporation and others.  All Rights Reserved. */
/* Subroutine */ int fd_(lafla, pntint, values, poly, n, q, neqcon, lptint, 
	lpoly)
doublereal *lafla, *pntint, *values, *poly;
integer *n, *q, *neqcon, *lptint, *lpoly;
{
    /* System generated locals */
    integer i__1, i__2, i__3, i__4;

    /* Local variables */
    static integer i__, k, l, dd;
    extern /* Subroutine */ int evalnp_();
    static integer np1, np2, qp1;
    static doublereal val;


/*  *********************************************************************** */
/*  THIS SUBROUTINE COMPUTED FINITE DIFFERENCES FOR MULTIVARIATE */
/*  INTERPOLATION, USING 'Q' INTERPOLATION POINTS AND 'Q' NEWTON */
/*  FUNDAMENTAL POLYNOMIALS. */

/*  PARAMETERS: */

/*   N      (INPUT)  DIMENTION OF THE PROBLEM */

/*   Q      (INPUT)  NUMBER OF INTERPOLATION POINTS */

/*   PNTINT (INPUT)  LIST OF  'NIND' DATA POINTS. THE I-TH POINT OCCUPIES */
/*                   POSITIONS ( I - 1 ) * N + 1 TO I * N. */
/*   VALUES (INPUT)  VALUES OF THE OBJECTIVE FUNCTION AT THE 'NIND' */
/*                   DATA POINTS CONTAINED IN POINTS. */
/*   POLY   (INPUT)  THE ARRAY CONTAINING COEFFICIENTS OF NEWTON FUNDAMENTAL */
/*                   POLYNOMIALS (AS COMPUTED BY 'NBUILD') */
/*   NEQCON (INPUT)  NUMBER OF LINEARLY INDEP. EQUALITY CONSTRAINTS */

/*   LAFLA  (OUTPUT) THE ARRAY OF THE  FINITE DIFFERENCES */

/*  ************************************************************************** */


/*  LOCAL VARIABLES */

    /* Parameter adjustments */
    --values;
    --lafla;
    --pntint;
    --poly;

    /* Function Body */
    np1 = *n + 1;
    np2 = *n + 2;
    dd = np1 * np2 / (float)2.;
    qp1 = *q + 1;

/*  INITIALISE THE LAFLA, THE FINITE DIFFERENCE OPERATORS ON F */

    i__1 = *q;
    for (i__ = 1; i__ <= i__1; ++i__) {
	lafla[i__] = values[i__];
/* L10: */
    }

/*  START THE MAJOR ITERATION */

    i__1 = *q;
    for (i__ = 1; i__ <= i__1; ++i__) {

/*  FOR ALL MULTIINDICES ALFA (WHOSE CARDINALITY IS GREATER THAN OR EQUAL */
/*  TO I), UPDATE THE ASSOCIATED FUNCTION L_ALFA USING */
/*  THE NEWTON FUNDAMENTAL POLYNOMIALS OF DEGREE EQUAL TO THE */
/*  CARDINALITY OF ALFA_{I-1} */

	if (i__ == 2) {

/*  ALL FINITE DIFFERENCES CORRESPONDING TO DEGREES >= 1 ARE UPDATED. */
/*  THE NEWTON FUNDAMENTAL POLYNOMIALS OF DEGREE CARDINALITY OF ALFA_I = 0 */
/*  ARE USED IN THE UPDATE. THE UPDATES ARE SKIPPED FOR DUMMY POLYNOMIALS */
/*  WHOSE INDEX IS ALWAYS HIGHER THAN NP1-NEQCON */

/* Computing MIN */
	    i__3 = *q, i__4 = np1 - *neqcon;
	    i__2 = min(i__3,i__4);
	    for (l = 2; l <= i__2; ++l) {
		evalnp_(&val, &pntint[1], &l, &poly[1], &c__1, n, lptint, 
			lpoly);
		lafla[l] -= lafla[1] * val;
/* L110: */
	    }
	    i__2 = *q;
	    for (l = np2; l <= i__2; ++l) {
		evalnp_(&val, &pntint[1], &l, &poly[1], &c__1, n, lptint, 
			lpoly);
		lafla[l] -= lafla[1] * val;
/* L115: */
	    }
	} else {
	    if (i__ == np2) {

/*  ALL FINITE DIFFERENCES CORRESPONDING TO DEGREES = 2 ARE UPDATED. */
/*  THE NEWTON FUNDAMENTAL POLYNOMIALS OF DEGREE CARDINALITY OF ALFA_I = 1 */
/*  ARE USED IN THE UPDATE */

		i__2 = *q;
		for (l = np2; l <= i__2; ++l) {
		    i__3 = np1 - *neqcon;
		    for (k = 2; k <= i__3; ++k) {
			evalnp_(&val, &pntint[1], &l, &poly[1], &k, n, lptint,
				 lpoly);
			lafla[l] -= lafla[k] * val;
/* L130: */
		    }
/* L120: */
		}
	    }
	}
/* L100: */
    }
    return 0;
} /* fd_ */