test_binv

(tmp_for_tar/misc/test_binv.m)

Help text

 Check whether bonv works

 Last modified: April 2001

Cross-Reference Information

binv             tmp_for_tar/misc/binv.m

Listing of file tmp_for_tar/misc/test_binv.m

##
## Check whether bonv works


## Author:        Etienne Grossmann  <etienne@isr.ist.utl.pt>
## Last modified: April 2001

Q = floor(30*rand(1))+30 ;
P = Q+floor(30*rand(1))+30 ;

a = randn(P) ;
invA = inv(a(1:Q,1:Q)) ;

b1 = inv(a) ;
b2 = binv(a,invA) ;

if any( abs(b1(:)-b2(:)) > sqrt(eps) ),
  printf("not ok : binv fails\n") ;
else
  printf("ok : binv works\n") ;
end


Q = 400 ;
P = 2*Q ; 
printf("making a  %i x %i  random matrix, please wait ...\n",P,P) ;
a = rand(P) ;

printf("inverting it a  %i x %i  subblock ...\n",Q,Q) ;

invA = inv(a(1:Q,1:Q)) ;

for R=[1,2,5,10,20,50,100,200],
  printf("size A is %-4i, size D is %-4i ... ",Q,R) ;
  mytic();
  b1 = inv(a(1:Q+R,1:Q+R)) ;
  printf("inv() : %-8.3g .... ",mytic()) ;
  mytic() ;
  b2 = binv(a(1:Q+R,1:Q+R),invA) ;
  printf("binv() : %-8.3g.   ",mytic()) ;
  printf("max err : %-8.3g\n",max(abs(b1(:)-b2(:)))) ;
end