1 #ifdef HAVE_STDLIB_H 2 #include <stdlib.h> 3 #endif 4 #include "yagi.h" 5 mprove(double ** a,double ** alud,int n,int * indx,double * b,double * x)6void mprove(double **a,double **alud,int n,int *indx,double *b,double *x) 7 { 8 int j,i; 9 double sdp; 10 double *r; 11 12 r=dvector(1,n); 13 for (i=1;i<=n;i++) { 14 sdp = -b[i]; 15 for (j=1;j<=n;j++) sdp += a[i][j]*x[j]; 16 r[i]=sdp; 17 } 18 lubksb(alud,n,indx,r); 19 for (i=1;i<=n;i++) x[i] -= r[i]; 20 free_dvector(r,1,n); 21 } 22