相关代码 |
|
dojohnavg(M,Mout,n)/* John's averaging process for matrices */ register double *M,*Mout; register int n; { register int i,j,iM,jM; register double term,diag; double tr; double trace(); tr=trace(M,n)*2/n; for(i=iM=0;i { diag=M[iM+i]; for(j=jM=0;j { if(i==j)continue; term=M[iM+j]*tr/(diag+M[jM+j]); Mout[iM+j]=term; } } term=tr/2; for(i=0;i *Mout=term; return; }