\File{gaussj.c},{12:26},{May 15 1997}
\L{\LB{\C{}\/* Slightly modified version of the stock numerical recipes}}
\L{\LB{}\Tab{8}{Gauss Jordan maxtrix solver *\/\CE{}}}
\L{\LB{}}
\L{\LB{\K{\#include} \<math.h\>}}
\L{\LB{\K{\#define} NRANSI}}
\L{\LB{\K{\#include} \S{}\"nrutil.h\"\SE{}}}
\L{\LB{\K{\#define} SWAP(a,b) \{temp=(a);(a)=(b);(b)=temp;\}}}
\L{\LB{\K{\#include} \<stdio.h\>}}
\L{\LB{}}
\L{\LB{ \K{void} gaussj(\K{float} a[6][6], \K{int} n, \K{float} b[6][2], \K{int} m) }}
\L{\LB{\{}}
\L{\LB{}\Tab{8}{\K{int} *indxc,*indxr,*ipiv;}}
\L{\LB{}\Tab{8}{\K{int} i,icol,irow,j,k,l,ll;}}
\L{\LB{}\Tab{8}{\K{float} big,dum,pivinv,temp;}}
\L{\LB{}}
\L{\LB{}\Tab{8}{indxc=ivector(1,n);}}
\L{\LB{}\Tab{8}{indxr=ivector(1,n);}}
\L{\LB{}\Tab{8}{ipiv=ivector(1,n);}}
\L{\LB{}\Tab{8}{\K{for} (j=1;j\<=n;j++) ipiv[j]=0;}}
\L{\LB{}\Tab{8}{\K{for} (i=1;i\<=n;i++) \{}}
\L{\LB{}\Tab{16}{big=0.0;}}
\L{\LB{}\Tab{16}{\K{for} (j=1;j\<=n;j++)}}
\L{\LB{}\Tab{24}{\K{if} (ipiv[j] != 1)}}
\L{\LB{}\Tab{32}{\K{for} (k=1;k\<=n;k++) \{}}
\L{\LB{}\Tab{40}{\K{if} (ipiv[k] == 0) \{}}
\L{\LB{}\Tab{48}{\K{if} (fabs(a[j][k]) \>= big) \{}}
\L{\LB{}\Tab{56}{big=fabs(a[j][k]);}}
\L{\LB{}\Tab{56}{irow=j;}}
\L{\LB{}\Tab{56}{icol=k;}}
\L{\LB{}\Tab{48}{\}}}
\L{\LB{}\Tab{40}{\} \K{else} \K{if} (ipiv[k] \> 1) nrerror(\S{}\"gaussj: Singular Matrix\-1\"\SE{});}}
\L{\LB{}\Tab{32}{\}}}
\L{\LB{}\Tab{16}{++(ipiv[icol]);}}
\L{\LB{}\Tab{16}{\K{if} (irow != icol) \{}}
\L{\LB{}\Tab{24}{\K{for} (l=1;l\<=n;l++) SWAP(a[irow][l],a[icol][l])}}
\L{\LB{}\Tab{24}{\K{for} (l=1;l\<=m;l++) SWAP(b[irow][l],b[icol][l])}}
\L{\LB{}\Tab{16}{\}}}
\L{\LB{}\Tab{16}{indxr[i]=irow;}}
\L{\LB{}\Tab{16}{indxc[i]=icol;}}
\L{\LB{}\Tab{16}{\K{if} (a[icol][icol] == 0.0) nrerror(\S{}\"gaussj: Singular Matrix\-2\"\SE{});}}
\L{\LB{}\Tab{16}{pivinv=1.0\/a[icol][icol];}}
\L{\LB{}\Tab{16}{a[icol][icol]=1.0;}}
\L{\LB{}\Tab{16}{\K{for} (l=1;l\<=n;l++) a[icol][l] *= pivinv;}}
\L{\LB{}\Tab{16}{\K{for} (l=1;l\<=m;l++) b[icol][l] *= pivinv;}}
\L{\LB{}\Tab{16}{\K{for} (ll=1;ll\<=n;ll++)}}
\L{\LB{}\Tab{24}{\K{if} (ll != icol) \{}}
\L{\LB{}\Tab{32}{dum=a[ll][icol];}}
\L{\LB{}\Tab{32}{a[ll][icol]=0.0;}}
\L{\LB{}\Tab{32}{\K{for} (l=1;l\<=n;l++) a[ll][l] \-= a[icol][l]*dum;}}
\L{\LB{}\Tab{32}{\K{for} (l=1;l\<=m;l++) b[ll][l] \-= b[icol][l]*dum;}}
\L{\LB{}\Tab{24}{\}}}
\L{\LB{}\Tab{8}{\}}}
\L{\LB{}\Tab{8}{\K{for} (l=n;l\>=1;l\-\-) \{}}
\L{\LB{}\Tab{16}{\K{if} (indxr[l] != indxc[l])}}
\L{\LB{}\Tab{24}{\K{for} (k=1;k\<=n;k++)}}
\L{\LB{}\Tab{32}{SWAP(a[k][indxr[l]],a[k][indxc[l]]);}}
\L{\LB{}\Tab{8}{\}}}
\L{\LB{}\Tab{8}{free\_ivector(ipiv,1,n);}}
\L{\LB{}\Tab{8}{free\_ivector(indxr,1,n);}}
\L{\LB{}\Tab{8}{free\_ivector(indxc,1,n);}}
\L{\LB{\}}}
\L{\LB{\K{\#undef} SWAP}}
\L{\LB{\K{\#undef} NRANSI}}