1 /* derivatives - for Bayes code in XLISP-STAT and S */
2 /* XLISP-STAT 2.1 Copyright (c) 1990, by Luke Tierney */
3 /* Additions to Xlisp 2.1, Copyright (c) 1989 by David Michael Betz */
4 /* You may give out copies of this software; for conditions see the */
5 /* file COPYING included with this distribution. */
8 extern double macheps();
15 numergrad(int n
, double *x
, double *grad
, double *fsum
,
16 int ffun(double *, double *, double *, double **),
17 double h
, double *typx
)
21 double old_xi
, f1
, f2
, hi
;
23 for (i
= 0; i
< n
; i
++) {
25 hi
= (typx
!= nil
) ? typx
[i
] * h
: h
;
27 (*ffun
)(x
, &f1
, nil
, nil
);
29 (*ffun
)(x
, &f2
, nil
, nil
);
31 grad
[i
] = (f1
- f2
) / (2.0 * hi
);
37 numerhess(int n
, double *x
, double **hess
, double f
, double *fsum
,
38 int *ffun(double *, double *, double *, double *),
39 double h
, double *typx
)
42 double old_xi
, old_xj
, f1
, f2
, hi
, hj
;
44 for (i
= 0; i
< n
; i
++) {
45 hi
= (typx
!= nil
) ? typx
[i
] * h
: h
;
46 hess
[i
][i
] = (fsum
[i
] - 2 * f
) / (hi
* hi
);
47 for (j
= i
+ 1; j
< n
; j
++) {
48 hj
= (typx
!= nil
) ? typx
[j
] * h
: h
;
53 (*ffun
)(x
, &f1
, nil
, nil
);
56 (*ffun
)(x
, &f2
, nil
, nil
);
59 hess
[i
][j
] = (2 * f
+ f1
+ f2
- fsum
[i
] - fsum
[j
]) / (2.0 * hi
* hj
);
60 hess
[j
][i
] = hess
[i
][j
];