1 /* atou1.c CCMATH mathematics library source code.
2 *
3 * Copyright (C) 2000 Daniel A. Atkinson All rights reserved.
4 * This code may be redistributed under the terms of the GNU library
5 * public license (LGPL). ( See the lgpl.license file for details.)
6 * ------------------------------------------------------------------------
7 */
8 #include <stdlib.h>
atou1(double * a,int m,int n)9 void atou1(double *a,int m,int n)
10 { double *p0,*p,*q,*w;
11 int i,j,k,mm;
12 double s,h;
13 w=(double *)calloc(m,sizeof(double));
14 p0=a+n*n-1; i=n-1; mm=m-n;
15 if(mm==0){ *p0=1.; p0-=n+1; --i; ++mm;}
16 for(; i>=0 ;--i,++mm,p0-=n+1){
17 if(*p0!=0.){
18 for(j=0,p=p0+n; j<mm ;p+=n) w[j++]= *p;
19 h= *p0; *p0=1.-h;
20 for(j=0,p=p0+n; j<mm ;p+=n) *p= -h*w[j++];
21 for(k=i+1,q=p0+1; k<n ;++k){
22 for(j=0,p=q+n,s=0.; j<mm ;p+=n) s+=w[j++]* *p;
23 s*=h;
24 for(j=0,p=q+n; j<mm ;p+=n) *p-=s*w[j++];
25 *q++ = -s;
26 }
27 }
28 else{
29 *p0=1.;
30 for(j=0,p=p0+n,q=p0+1; j<mm ;++j,p+=n) *p= *q++ =0.;
31 }
32 }
33 free(w);
34 }
35