4 #include "gromacs/utility/real.h"
6 #include "../gmx_blas.h"
7 #include "../gmx_lapack.h"
10 F77_FUNC(dsytd2,DSYTD2)(const char * uplo,
17 int gmx_unused * info)
20 double taui,alpha,tmp;
22 const char ch=toupper(*uplo);
31 for(i=*n-1;i>=1;i--) {
34 F77_FUNC(dlarfg,DLARFG)(&i,&(a[i*(*lda)+(i-1)]),&(a[i*(*lda)+0]),&ti1,&taui);
35 e[i-1] = a[i*(*lda) + (i-1)];
36 if(fabs(taui)>GMX_DOUBLE_MIN) {
37 a[i*(*lda)+(i-1)] = 1.0;
40 F77_FUNC(dsymv,DSYMV)("U",&i,&taui,a,lda,&(a[i*(*lda)+0]),&ti1,&zero,tau,&ti1);
42 tmp = F77_FUNC(ddot,DDOT)(&i,tau,&ti1,&(a[i*(*lda)+0]),&ti1);
44 alpha = -0.5*taui*tmp;
46 F77_FUNC(daxpy,DAXPY)(&i,&alpha,&(a[i*(*lda)+0]),&ti1,tau,&ti1);
48 F77_FUNC(dsyr2,DSYR2)("U",&i,&minusone,&(a[i*(*lda)+0]),&ti1,tau,&ti1,a,lda);
50 a[i*(*lda)+(i-1)] = e[i-1];
64 ti2 = ( *n < i+2) ? *n : i+2;
66 F77_FUNC(dlarfg,DLARFG)(&ti1,&(a[(i-1)*(*lda)+(i)]),&(a[(i-1)*(*lda)+ti2-1]),&ti3,&taui);
68 e[i-1] = a[(i-1)*(*lda) + (i)];
70 if(fabs(taui)>GMX_DOUBLE_MIN) {
71 a[(i-1)*(*lda)+(i)] = 1.0;
75 F77_FUNC(dsymv,DSYMV)(uplo,&ti1,&taui,&(a[i*(*lda)+i]),lda,&(a[(i-1)*(*lda)+i]),
76 &ti2,&zero,&(tau[i-1]),&ti2);
78 tmp = F77_FUNC(ddot,DDOT)(&ti1,&(tau[i-1]),&ti2,&(a[(i-1)*(*lda)+i]),&ti2);
80 alpha = -0.5*taui*tmp;
82 F77_FUNC(daxpy,DAXPY)(&ti1,&alpha,&(a[(i-1)*(*lda)+i]),&ti2,&(tau[i-1]),&ti2);
84 F77_FUNC(dsyr2,DSYR2)(uplo,&ti1,&minusone,&(a[(i-1)*(*lda)+i]),&ti2,&(tau[i-1]),&ti2,
85 &(a[(i)*(*lda)+i]),lda);
87 a[(i-1)*(*lda)+(i)] = e[i-1];
90 d[i-1] = a[(i-1)*(*lda)+i-1];
93 d[*n-1] = a[(*n-1)*(*lda)+(*n-1)];