~logan/ubuntu/trusty/suitesparse/4.2.1-3ubuntu1

« back to all changes in this revision

Viewing changes to CXSparse/MATLAB/Test/test_qr.m

  • Committer: Bazaar Package Importer
  • Author(s): Christophe Prud'homme
  • Date: 2007-05-29 09:36:29 UTC
  • mfrom: (1.1.1 upstream)
  • Revision ID: james.westby@ubuntu.com-20070529093629-zowquo0b7slkk6nc
Tags: 3.0.0-2
* suitesparse builds properly twice in a row
* Bug fix: "suitesparse - FTBFS: Broken build depens: libgfortran1-dev",
  thanks to Bastian Blank (Closes: #426349).
* Bug fix: "suitesparse_3.0.0-1: FTBFS: build-depends on
  libgfortran1-dev", thanks to Steve Langasek (Closes: #426354).

Show diffs side-by-side

added added

removed removed

Lines of Context:
 
1
function test_qr
 
2
%TEST_QR test various QR factorization methods
 
3
%
 
4
% Example:
 
5
%   test_qr
 
6
% See also: testall
 
7
 
 
8
%   Copyright 2006-2007, Timothy A. Davis.
 
9
%   http://www.cise.ufl.edu/research/sparse
 
10
 
 
11
index = UFget ;
 
12
[ignore f] = sort (max (index.nrows,index.ncols)) ;
 
13
 
 
14
% f = 276 
 
15
% f = 706
 
16
f = f (1:100) ;
 
17
 
 
18
for i = f 
 
19
 
 
20
    % Prob = UFget (i,index)
 
21
    Prob = UFget (i) ;
 
22
    disp (Prob) ;
 
23
    A = Prob.A ;
 
24
    [m n] = size (A) ;
 
25
    if (m < n)
 
26
        A = A' ;
 
27
    end
 
28
    [m n] = size (A) ;
 
29
    if (sprank (A) < n | ~isreal (A))                                       %#ok
 
30
        continue ;
 
31
    end
 
32
 
 
33
    [V,beta,p,R1,q] = cs_qr(A) ;
 
34
    A = A (p,q) ;
 
35
    parent = etree (A, 'col') ;                                             %#ok
 
36
 
 
37
    R0 = qr (A) ;
 
38
    R2 = qr_givens (full (A)) ;
 
39
    R3 = qr_givens_full (full (A)) ;
 
40
 
 
41
    subplot (2,2,1) ; cspy (R0) ; title ('matlab') ;
 
42
    subplot (2,2,2) ; cspy (R3) ; title ('qr-full') ;
 
43
    subplot (2,2,3) ; cspy (R2) ; title ('qr-givens') ;
 
44
    subplot (2,2,4) ; cspy (R1) ; title ('cs-qr') ;
 
45
 
 
46
    e0 = norm (A'*A-R0'*R0,1) / norm (A,1) ;
 
47
    e1 = norm (A'*A-R1'*R1,1) / norm (A,1) ;
 
48
    e2 = norm (A'*A-R2'*R2,1) / norm (A,1) ;
 
49
    e3 = norm (A'*A-R3'*R3,1) / norm (A,1) ;
 
50
    fprintf ('error %6.2e %6.2e %6.2e %6.2e\n', e0, e1, e2, e3) ;
 
51
    drawnow
 
52
    if (e1 > e0*1e3 | e2 > e0*1e3)                                          %#ok
 
53
        pause
 
54
    end
 
55
end