Bug #10178: norm returned error message: 'division by zero' for some sparse matrices. 79/8779/2
Adeline CARNIS [Thu, 23 Aug 2012 07:38:02 +0000 (09:38 +0200)]
Change-Id: I40b8ec74a4dbdc3c678c53ee9cdd2fecad9ae242

scilab/CHANGES_5.4.X
scilab/modules/linear_algebra/tests/nonreg_tests/bug_10178.dia.ref [new file with mode: 0644]
scilab/modules/linear_algebra/tests/nonreg_tests/bug_10178.tst [new file with mode: 0644]
scilab/modules/overloading/macros/%sp_norm.sci

index 4c6e2bf..5fe6b4c 100644 (file)
@@ -78,6 +78,9 @@ Bug Fixes
                      give the possibility to use another norm than that imposed
                      (by default 2-norm).
 
+* Bug #10178 fixed - The norm function returned the message: "division by zero" for
+                     some sparse matrices.
+
 * Bug #11411 fixed - save function used unsigned char to store length of string.
                      Now it is an integer.
 
diff --git a/scilab/modules/linear_algebra/tests/nonreg_tests/bug_10178.dia.ref b/scilab/modules/linear_algebra/tests/nonreg_tests/bug_10178.dia.ref
new file mode 100644 (file)
index 0000000..f76f6ab
--- /dev/null
@@ -0,0 +1,40 @@
+// =============================================================================
+// Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
+// Copyright (C) 2012 - Scilab Enterprises - Adeline CARNIS
+//
+//  This file is distributed under the same license as the Scilab package.
+// =============================================================================
+// <-- CLI SHELL MODE -->
+// <-- Non-regression test for bug 10178 -->
+//
+// <-- Bugzilla URL -->
+// http://bugzilla.scilab.org/show_bug.cgi?id=10178
+//
+// <-- Short Description -->
+//    norm function falied for some sparse matrices.
+// =============================================================================
+A = [1 0 0 1 0;0 -1 -1 0 -1];
+AS = sparse(A);
+assert_checkalmostequal(norm(A), norm(AS));
+assert_checkalmostequal(norm(A'), norm(AS'));
+assert_checkequal(norm(A, 1), norm(AS, 1));
+assert_checkequal(norm(A, %inf), norm(AS, %inf));
+assert_checkequal(norm(A, 'fro'), norm(AS, 'fro'));
+A = [1 1 1 1 1;-1 -1 -1 -1 -1];
+AS = sparse(A);
+assert_checkequal(norm(A), norm(AS));
+assert_checkequal(norm(A)', norm(AS)');
+assert_checkequal(norm(A, 1), norm(AS, 1));
+assert_checkequal(norm(A, %inf), norm(AS, %inf));
+assert_checkequal(norm(A, 'fro'), norm(AS, 'fro'));
+A = [
+    0.    0.    0.    0.    0.  
+    0.    0.    0.    0.    0.  
+    0.    0.   -1.    0.    1.  
+    0.    0.    0.    0.    0.  
+    0.    0.    1.    0.   -1.  ];
+AS = sparse(A);
+assert_checkequal(norm(A), norm(AS));
+assert_checkequal(norm(A, 1), norm(AS, 1));
+assert_checkequal(norm(A, %inf), norm(AS, %inf));
+assert_checkequal(norm(A, 'fro'), norm(AS, 'fro'));
diff --git a/scilab/modules/linear_algebra/tests/nonreg_tests/bug_10178.tst b/scilab/modules/linear_algebra/tests/nonreg_tests/bug_10178.tst
new file mode 100644 (file)
index 0000000..606cc5e
--- /dev/null
@@ -0,0 +1,50 @@
+// =============================================================================
+// Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
+// Copyright (C) 2012 - Scilab Enterprises - Adeline CARNIS
+//
+//  This file is distributed under the same license as the Scilab package.
+// =============================================================================
+
+// <-- CLI SHELL MODE -->
+
+// <-- Non-regression test for bug 10178 -->
+//
+// <-- Bugzilla URL -->
+// http://bugzilla.scilab.org/show_bug.cgi?id=10178
+//
+// <-- Short Description -->
+//    norm function falied for some sparse matrices.
+// =============================================================================
+
+A = [1 0 0 1 0;0 -1 -1 0 -1];
+AS = sparse(A);
+
+assert_checkalmostequal(norm(A), norm(AS));
+assert_checkalmostequal(norm(A'), norm(AS'));
+assert_checkequal(norm(A, 1), norm(AS, 1));
+assert_checkequal(norm(A, %inf), norm(AS, %inf));
+assert_checkequal(norm(A, 'fro'), norm(AS, 'fro'));
+
+A = [1 1 1 1 1;-1 -1 -1 -1 -1];
+AS = sparse(A);
+
+assert_checkequal(norm(A), norm(AS));
+assert_checkequal(norm(A)', norm(AS)');
+assert_checkequal(norm(A, 1), norm(AS, 1));
+assert_checkequal(norm(A, %inf), norm(AS, %inf));
+assert_checkequal(norm(A, 'fro'), norm(AS, 'fro'));
+
+
+A = [
+    0.    0.    0.    0.    0.  
+    0.    0.    0.    0.    0.  
+    0.    0.   -1.    0.    1.  
+    0.    0.    0.    0.    0.  
+    0.    0.    1.    0.   -1.  ];
+
+AS = sparse(A);
+
+assert_checkequal(norm(A), norm(AS));
+assert_checkequal(norm(A, 1), norm(AS, 1));
+assert_checkequal(norm(A, %inf), norm(AS, %inf));
+assert_checkequal(norm(A, 'fro'), norm(AS, 'fro'));
index 8d4498e..06c4c67 100644 (file)
@@ -1,5 +1,6 @@
 // Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
 // Copyright (C) INRIA
+// Copyright (C) Scilab Enterprises - Adeline CARNIS
 // 
 // This file must be used under the terms of the CeCILL.
 // This source file is licensed as described in the file COPYING, which
 
 function res=%sp_norm(S,flag)
 
-[lhs,rhs]=argn(0)
-if rhs==1 then flag=2;end //norm(S)=norm(S,2)
-[m,n]=size(S)
+    [lhs,rhs]=argn(0)
+    if rhs==1 then flag=2;end //norm(S)=norm(S,2)
+    [m,n]=size(S)
 
-if m==1|n==1 then //vector norm
-  [ij,v]=spget(S);
-  res=norm(v,flag);
-  return
-end
+    if m==1|n==1 then //vector norm
+        [ij,v]=spget(S);
+        res=norm(v,flag);
+        return
+    end
 
-select flag
-case 1 then
-  res=max(ones(1,m)*abs(S))
-case 2 then
-  if m<n then
-    S=S*S'
-  elseif m>n then
-    S=S'*S
-  end
-  tol=1.d-5;
-  x = sum(abs(S),1)';
-  res = norm(x);
-  if res==0 then return,end
-  x = x/res;res0 = 0;
-  while abs(res-res0) > tol*res
-    res0 = res;   Sx = S*x; res = norm(Sx);
-    x = S'*Sx;
-    x = x/norm(x);
-  end
-  if m<>n then res=sqrt(res),end
-case %inf then
-  res=max(abs(S)*ones(n,1))
-case 'inf' then
-  res=max(abs(S)*ones(n,1))
-case 'fro' then
-  [ij,v]=spget(S);
-  res=sqrt(sum(abs(v.*v)))
-else
-  [ij,v]=spget(S);
-  res=norm(v,flag);
-end
+    select flag
+    case 1 then
+        res=max(ones(1,m)*abs(S))
+    case 2 then
+        if m<n then
+            S1=S*S'
+        elseif m>n then
+            S1=S'*S
+        else
+            S1 = S;
+        end
+
+        tol=%eps;
+        x = sum(abs(S1),1)';
+        res = norm(x);
+        if res==0 then return,end
+        x = x/res;res0 = 0;
+        while abs(res-res0) > tol*res
+            res0 = res;   Sx = S1*x;
+            
+            // Bug #10178: norm failed for some sparse matrix
+            // If Sx = 0, we had "division by zero" with x/norm(x)
+            // So, use to sum(abs(S).^2).^(1/2)
+            if Sx == 0 then
+                res = sum(abs(S).^2).^(1/2);
+                return
+            end
+            // End Bug #10178
+            
+            res = norm(Sx);
+            x = S1'*Sx;
+            x = x/norm(x);
+        end
+        if m<>n then res=sqrt(res),end
+    case %inf then
+        res=max(abs(S)*ones(n,1))
+    case 'inf' then
+        res=max(abs(S)*ones(n,1))
+    case 'fro' then
+        [ij,v]=spget(S);
+        res=sqrt(sum(abs(v.*v)))
+    else
+        [ij,v]=spget(S);
+        res=norm(v,flag);
+    end
 endfunction