#include <stdio.h>#include <stdlib.h>#include <math.h>#include "math_util.h"#define MAX_ITT_CM 5000void cauchy( fun, dfun, x, n, tol )long double (*fun)( long double * );void (*dfun)( long double *, long double * );long double *x, tol;int n;{   long double dx, ndx, *grad, *h, *x0, *x1, f0, f1, alpha, err1;   int i, j;   long double line_search( long double(*)(), long double, long double *,                             long double *, int, long double );      if( tol <= 0.0 ) tol = 1.0e-4;   if( n <= 0 ) return;   grad = a1d_allo_dbl( n );   h = a1d_allo_dbl( n );   x0 = a1d_allo_dbl( n );      vec_copy( x0, x, n );      f0 = (*fun)( x0 );   for( i = 0; i < MAX_ITT_CM; i++ ) {      printf(" itt = %d, f = %Le\n",i,f0);      (*dfun)( grad, x0 );      err1 = enorm( grad, n );      if( err1 <= tol ) {         vec_copy( x, x0, n );         a1d_free_dbl( x0 );         a1d_free_dbl( h );         a1d_free_dbl( grad );         printf("***  Norm(grad f) <= eps\n");         return;      }      for( j = 0; j < n; j++ )         h[j] = -grad[j];               alpha = line_search( fun, f0, x0, h, n, tol );            ndx = 0.0;            for( j = 0; j < n; j++ ) {         dx = alpha*h[j];         x0[j] += dx;         ndx += dx*dx;      }      ndx = sqrt( ndx );               f1 = (*fun)( x0 );      /*      if( ndx <= tol/10 ) {         vec_copy( x, x0, n );         a1d_free_dbl( x0 );         a1d_free_dbl( h );         a1d_free_dbl( grad );         printf("***  |x1-x0| <= 0.1*eps\n");         return;      }*/         f0 = f1;      }   printf("Error: Too many Cauchy iterations\n");   vec_copy( x, x0, n );   a1d_free_dbl( x0 );   a1d_free_dbl( h );   a1d_free_dbl( grad );   return;   }