octave-maintainers
[Top][All Lists]
Advanced

[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

Re: Seg fault with test_sparse.m


From: David Bateman
Subject: Re: Seg fault with test_sparse.m
Date: Fri, 18 Aug 2006 20:50:16 +0200
User-agent: Mozilla Thunderbird 1.0.6-7.6.20060mdk (X11/20050322)

Joe Koski wrote:
> on 8/17/06 3:30 PM, David Bateman at address@hidden wrote:
> 
> 
>>Joe Koski wrote:
>>
>>>on 8/16/06 2:47 PM, David Bateman at address@hidden wrote:
>>>
>>>
>>>
>>>>Joe Koski wrote:
>>>>
>>>>
>>>>>Hi all,
>>>>>
>>>>>After getting a successful build of octave-2.9.7 on my Mac G5 with OS X
>>>>>10.4.7, Xcode tools 2.3 (gcc-4.0.1) and gfortran, I installed dejaGNU and
>>>>>ran make check. I got one IO error, but things were going well until I hit
>>>>>test_sparse. There I got a seg fault. I have attached the terminal output
>>>>>and fntest.log output in case you want to see them.
>>>>>
>>>>>I then moved to the test directory and ran octave. Here's what I got:
>>>>>
>>>>>octave:1> test_sparse
>>>>>writing test output to sptest.log
>>>>> 
>>>>>
>>>>>
>>>>>>>>>>/Tools/octave-2.9.7/test/test_sparse.m
>>>>>>>>>>         
>>>>>
>>>>>warning: Calculating Sparse Matrix Type
>>>>>panic: Bus error -- stopping myself...
>>>>>Bus error
>>>>>Joe-Koskis-Computer:/Tools/octave-2.9.7/test jakoski$
>>>>>
>>>>>I have also attached the sptest.log for your information.
>>>>>
>>>>>To my untrained eye, the script seems to be running through quite a few
>>>>>tests before the failure. Things are almost working.
>>>>>
>>>>>Anyway, I'll be available to try things, but maybe it's time for someone
>>>>>with sparse matrix knowledge (definitely not me) to try this on a Mac. I'll
>>>>>be available to assist anyone who wants to try building octave-2.9.7 or
>>>>>later on a Mac under OS X 10.4. I can cut many hours from the build process
>>>>>with my experience from the last few days.
>>>>>
>>>>>The other problem that I have is with the legend.m routine in octave-forge.
>>>>>The problem is persistent with any script that uses legend. As I said, I
>>>>>can
>>>>>file a bug report, if it is an unknown problem.
>>>>>
>>>>>Meanwhile, I'll go back to octave-2.1.73 for my routine stuff. Let me know
>>>>>when or if I can help or try something.
>>>>>
>>>>>As an aside, I've been around complicated software and complicated
>>>>>engineering systems for over 40 years, and I have seldom seen things that
>>>>>went cleanly from start to finish. There's always something, and the
>>>>>problem
>>>>>never seems to be the thing that you were worrying about.
>>>>>
>>>>>Joe
>>>>>
>>>>>
>>>>>
>>>>>
>>>>>
>>>>> 
>>>>
>>>>I tried looking at this in June with Marius Schamschula, but couldn't
>>>>identify the exact test causing the issues. Firstly, is UFSparse
>>>>correctly installed? If so which version? In any case, I'd promised
>>>>myself I'd steal my wife's powerbook for a few hours to test this out
>>>>one of these days, though I'd be happy if some one else could identify
>>>>exactly which test is failing. From your logs it seems that something like
>>>>
>>>>sm = sprandn(8,8,0.2)
>>>>fm = rand(8,8);
>>>>am = sm .^ fm;
>>>>am2 = am(find(am != 0)
>>>>
>>>>might fail for your. If it does I have a chance of finding the problem
>>>>without bothering my wife :-)
>>>>
>>>>D.
>>>>
>>>
>>>David,
>>>
>>>You are so correct. Your test fails on my system:
>>>
>>>Joe-Koskis-Computer:~/Codes/octave_test_files jakoski$ octave
>>>GNU Octave, version 2.9.7 (powerpc-apple-darwin8.7.0).
>>>Copyright (C) 2006 John W. Eaton.
>>>This is free software; see the source code for copying conditions.
>>>There is ABSOLUTELY NO WARRANTY; not even for MERCHANTIBILITY or
>>>FITNESS FOR A PARTICULAR PURPOSE.  For details, type `warranty'.
>>>
>>>Additional information about Octave is available at http://www.octave.org.
>>>
>>>Please contribute if you find this software useful.
>>>For more information, visit http://www.octave.org/help-wanted.html
>>>
>>>Report bugs to <address@hidden> (but first, please read
>>>http://www.octave.org/bugs.html to learn how to write a helpful report).
>>>
>>>octave:1> sm = sprandn(8,8,0.2)
>>>sm =
>>>
>>>Compressed Column Sparse (rows = 8, cols = 8, nnz = 12)
>>>
>>>  (1, 1) -> -1.2817
>>>  (6, 1) ->  1.0076
>>>  (7, 1) -> -0.65782
>>>  (5, 2) -> -1.2952
>>>  (7, 3) -> -0.66335
>>>  (1, 4) -> -0.64771
>>>  (4, 4) ->  0.63989
>>>  (7, 4) -> -0.40588
>>>  (8, 6) -> -0.82241
>>>  (5, 8) ->  1.1687
>>>  (6, 8) ->  1.4763
>>>  (7, 8) ->  0.24860
>>>
>>>octave:2> fm = rand(8,8);
>>>octave:3> am = sm .^ fm;
>>>panic: Bus error -- stopping myself...
>>>
>>>My UFSparse came directly from the Florida website:
>>>
>>>  http://www.cise.ufl.edu/research/sparse/umfpack/
>>>
>>>Which shows in the README as version 5.0. Of course there are seven other
>>>libraries that must also be built before all the externals are satisfied
>>>(according to the octave configure script, anyway): AMD-2.0, COLAMD-2.5,
>>>CCOLAMD-2.5, CHOLMOD-1.1.1, CXSparse-2.0.1, metis-4.0, and CAMD-2.1. Since
>>>all these are in separate directories, my fear is that a header file in one
>>>library doesn't match the same header in another library. Just a thought.
>>>Unfortunately, none of the libraries are easy to identify without opening
>>>each individual README file. Of course metis-4.0 comes from
>>>
>>>  http://glaros.dtc.umn.edu/gkhome/views/metis
>>>
>>>If you want me to rebuild and retest, I'm available. Just let me know. I
>>>just installed the Apple Xcode-2.4 developer tools which were released about
>>>a week ago. I haven't tried building anything with the new tools yet,
>>>though. I'm trying to keep the latest version of my OS and Xcode tools on
>>>the premise that it's the target configuration for octave.
>>>
>>>Joe
>>>
>>
>>Joe,
>>
>>I have no idea why the existing code in sparse-xpow.cc is causing issues
>>on big-endian Mac's. However, looking at the code I think it needs a
>>major overhaul in any case as the existing code is much slower than it
>>needs to be. In any case ".^" is a pretty stupid operator to use on a
>>sparse matrix, so I don't see that it should be optimized too much. In
>>any case, can you try the attached patch on your Mac and see if the
>>short test I sent you now works. Can you also then run "make check"
>>again and find the next place (if any) where there are issues :-)
>>
>>D.
>>
> 
> 
> David,
> 
> I applied the patch and did a complete rebuild (27 minutes with make -j2).
> Things are much better. I did the rebuild because I wasn't sure if my new
> Xcode-2.4 was totally compatible with the previous Xcode-2.3, and to make
> sure I could still build octave with 2.4.
> 
> Now, no seg fault, instead "make check" gives
> 
> test_sparse.m .......................................... PASS 1156/1176 FAIL
> 20
> 
> Then, as before, I moved to /test and ran test_sparse directly from octave.
> I'm forwarding the sptest.log to you separately for your examination. I'm
> not familiar enough with the tests to understand the output.

Ok, will look at this soon to see where the 20 failing tests are..


> 
> On the simple test that you sent, it runs through as follows.
<snip>
> I guess my only question is whether the NaNs should be there.

I don't think any element should be NaN, I wonder what

sm(find(am == NaN))
fm(find(am == NaN))

gives. One possible cause for this I can see is that the full matrix
versions of these operators used the std::pow function rather than "pow"
internally. Perhaps the behaviour is different in the two cases for a
Mac. Can you try the attached version of the patch instead and see if
this helps.

> 
> I agree that a sparse matrix with each element to a different non-sparse
> power (at least that's my interpretation of the test that caused the
> original problem) is not often done, especially by me.

Yeah, its basically a pretty dumb operator for sparse matrices, but it
should work correctly and consistently with full matrices.


> 
> If we wait a couple of years, most Macs will be little-endian, and this
> problem will be history, but there are always a few PowerPCs, Irixes, etc.
> out there.

There'll always be some BE machine somewhere that some one will want to
run octave on, so better to address this on a MaC now as that is one
very common BE platform.

> 
> Let me know what to try next. I'm available for testing.

Apart from the attached patch I'll propose something after I read the
log file.

D.


> 
> Joe
> 
> 
> 

Index: src/sparse-xpow.cc
===================================================================
RCS file: /cvs/octave/src/sparse-xpow.cc,v
retrieving revision 1.8
diff -c -r1.8 sparse-xpow.cc
*** src/sparse-xpow.cc  3 May 2006 19:32:47 -0000       1.8
--- src/sparse-xpow.cc  18 Aug 2006 18:45:37 -0000
***************
*** 249,255 ****
          for (octave_idx_type i = 0; i < nr; i++)
            {
              OCTAVE_QUIT;
!             result (i, j) = pow (atmp, b(i,j));
            }
        }
  
--- 249,255 ----
          for (octave_idx_type i = 0; i < nr; i++)
            {
              OCTAVE_QUIT;
!             result (i, j) = std::pow (atmp, b(i,j));
            }
        }
  
***************
*** 264,270 ****
          for (octave_idx_type i = 0; i < nr; i++)
            {
              OCTAVE_QUIT;
!             result (i, j) = pow (a, b(i,j));
            }
        }
  
--- 264,270 ----
          for (octave_idx_type i = 0; i < nr; i++)
            {
              OCTAVE_QUIT;
!             result (i, j) = std::pow (a, b(i,j));
            }
        }
  
***************
*** 289,295 ****
        for (octave_idx_type i = 0; i < nr; i++)
        {
          OCTAVE_QUIT;
!         result (i, j) = pow (atmp, b(i,j));
        }
      }
  
--- 289,295 ----
        for (octave_idx_type i = 0; i < nr; i++)
        {
          OCTAVE_QUIT;
!         result (i, j) = std::pow (atmp, b(i,j));
        }
      }
  
***************
*** 315,321 ****
  
        if (static_cast<int> (b) != b && a.any_element_is_negative ())
        {
!         ComplexMatrix result (nr, nc, Complex (pow (0.0, b)));
  
          // FIXME -- avoid apparent GNU libm bug by
          // converting A and B to complex instead of just A.
--- 315,321 ----
  
        if (static_cast<int> (b) != b && a.any_element_is_negative ())
        {
!         ComplexMatrix result (nr, nc, Complex (std::pow (0.0, b)));
  
          // FIXME -- avoid apparent GNU libm bug by
          // converting A and B to complex instead of just A.
***************
*** 328,347 ****
              
                Complex atmp (a.data (i));
                
!               result (a.ridx(i), j) = pow (atmp, btmp);
              }
  
          retval = octave_value (result);
        }
        else
        {
!         Matrix result (nr, nc, (pow (0.0, b)));
  
          for (octave_idx_type j = 0; j < nc; j++)
            for (octave_idx_type i = a.cidx(j); i < a.cidx(j+1); i++)
              {
                OCTAVE_QUIT;
!               result (a.ridx(i), j) = pow (a.data (i), b);
              }
  
          retval = octave_value (result);
--- 328,347 ----
              
                Complex atmp (a.data (i));
                
!               result (a.ridx(i), j) = std::pow (atmp, btmp);
              }
  
          retval = octave_value (result);
        }
        else
        {
!         Matrix result (nr, nc, (std::pow (0.0, b)));
  
          for (octave_idx_type j = 0; j < nc; j++)
            for (octave_idx_type i = a.cidx(j); i < a.cidx(j+1); i++)
              {
                OCTAVE_QUIT;
!               result (a.ridx(i), j) = std::pow (a.data (i), b);
              }
  
          retval = octave_value (result);
***************
*** 361,367 ****
          Complex atmp (a.data (i));
          Complex btmp (b);
  
!         result.data (i) = pow (atmp, btmp);
        }
  
        result.maybe_compress (true);
--- 361,367 ----
          Complex atmp (a.data (i));
          Complex btmp (b);
  
!         result.data (i) = std::pow (atmp, btmp);
        }
  
        result.maybe_compress (true);
***************
*** 375,381 ****
        for (octave_idx_type i = 0; i < nz; i++)
        {
          OCTAVE_QUIT;
!         result.data (i) = pow (a.data (i), b);
        }
  
        result.maybe_compress (true);
--- 375,381 ----
        for (octave_idx_type i = 0; i < nz; i++)
        {
          OCTAVE_QUIT;
!         result.data (i) = std::pow (a.data (i), b);
        }
  
        result.maybe_compress (true);
***************
*** 406,480 ****
  
    int convert_to_complex = 0;
    for (octave_idx_type j = 0; j < nc; j++)
!     for (octave_idx_type i = 0; i < nr; i++)
        {
!       OCTAVE_QUIT;
!       double atmp = a (i, j);
!       double btmp = b (i, j);
!       if (atmp < 0.0 && static_cast<int> (btmp) != btmp)
          {
!           convert_to_complex = 1;
!           goto done;
          }
        }
  
  done:
  
!   octave_idx_type nel = 0;
!   for (octave_idx_type j = 0; j < nc; j++) 
!     for (octave_idx_type i = 0; i < nr; i++)
!       if (!(a.elem (i, j) == 0. && b.elem (i, j) != 0.))
!       nel++;
  
    if (convert_to_complex)
      {
!       SparseComplexMatrix complex_result (nr, nc, nel);
  
-       octave_idx_type ii = 0;
-       complex_result.cidx(0) = 0;
        for (octave_idx_type j = 0; j < nc; j++)
        {
!         for (octave_idx_type i = 0; i < nr; i++)
            {
              OCTAVE_QUIT;
!             Complex atmp (a (i, j));
!             Complex btmp (b (i, j));
!             Complex tmp =  pow (atmp, btmp);
!             if (tmp != 0.)
!               {
!                 complex_result.data (ii) = tmp;
!                 complex_result.ridx (ii++) = i;
!               }
            }
-         complex_result.cidx (j+1) = ii;
        }
!       complex_result.maybe_compress ();
! 
        retval = complex_result;
      }
    else
      {
!       SparseMatrix result (nr, nc, nel);
!       octave_idx_type ii = 0;
  
-       result.cidx (0) = 0;
        for (octave_idx_type j = 0; j < nc; j++)
        {
!         for (octave_idx_type i = 0; i < nr; i++)
            {
              OCTAVE_QUIT;
!             double tmp = pow (a (i, j), b (i, j));
!             if (tmp != 0.)
!               {
!                 result.data (ii) = tmp;
!                 result.ridx (ii++) = i;
!               }
            }
-         result.cidx (j+1) = ii;
        }
! 
!       result.maybe_compress ();
! 
        retval = result;
      }
  
--- 406,461 ----
  
    int convert_to_complex = 0;
    for (octave_idx_type j = 0; j < nc; j++)
!     for (octave_idx_type i = a.cidx(j); i < a.cidx(j+1); i++)
        {
!       if (a.data(i) < 0.0)
          {
!           double btmp = b (a.ridx(i), j);
!           if (static_cast<int> (btmp) != btmp)
!             {
!               convert_to_complex = 1;
!               goto done;
!             }
          }
        }
  
  done:
  
!   // This is a dumb operator for sparse matrices anyway, and there is
!   // no sensible way to handle the 0.^0 versus the 0.^x cases. Therefore
!   // allocate a full matrix filled for the 0.^0 case and shrink it later
!   // as needed
  
    if (convert_to_complex)
      {
!       SparseComplexMatrix complex_result (nr, nc, Complex(1.0, 0.0));
  
        for (octave_idx_type j = 0; j < nc; j++)
        {
!         for (octave_idx_type i = a.cidx(j); i < a.cidx(j+1); i++)
            {
              OCTAVE_QUIT;
!             complex_result.xelem(a.ridx(i), j) = std::pow (a.data(i), 
!                                                       b (a.ridx(i), j));
            }
        }
!       complex_result.maybe_compress (true);
        retval = complex_result;
      }
    else
      {
!       SparseMatrix result (nr, nc, 1.0);
  
        for (octave_idx_type j = 0; j < nc; j++)
        {
!         for (octave_idx_type i = a.cidx(j); i < a.cidx(j+1); i++)
            {
              OCTAVE_QUIT;
!             result.xelem(a.ridx(i), j) = std::pow (a.data(i), 
!                                                    b (a.ridx(i), j));
            }
        }
!       result.maybe_compress (true);
        retval = result;
      }
  
***************
*** 498,504 ****
        for (octave_idx_type i = 0; i < nz; i++)
        {
          OCTAVE_QUIT;
!         result.data (i) = pow (Complex (a.data (i)), b);
        }
    
        result.maybe_compress (true);
--- 479,485 ----
        for (octave_idx_type i = 0; i < nz; i++)
        {
          OCTAVE_QUIT;
!         result.data (i) = std::pow (Complex (a.data (i)), b);
        }
    
        result.maybe_compress (true);
***************
*** 525,556 ****
        return octave_value ();
      }
  
!   octave_idx_type nel = 0;
!   for (octave_idx_type j = 0; j < nc; j++) 
!     for (octave_idx_type i = 0; i < nr; i++)
!       if (!(a.elem (i, j) == 0. && b.elem (i, j) != 0.))
!       nel++;
! 
!   SparseComplexMatrix result (nr, nc, nel);
!   octave_idx_type ii = 0;
! 
!   result.cidx(0) = 0;
    for (octave_idx_type j = 0; j < nc; j++)
      {
!       for (octave_idx_type i = 0; i < nr; i++)
        {
          OCTAVE_QUIT;
!         Complex tmp = pow (Complex (a (i, j)), b (i, j));
!         if (tmp != 0.)
!           {
!             result.data (ii) = tmp; 
!             result.ridx (ii++) = i; 
!           }
        }
-       result.cidx (j+1) = ii;
      }
  
!   result.maybe_compress ();
  
    return result;
  }
--- 506,522 ----
        return octave_value ();
      }
  
!   SparseComplexMatrix result (nr, nc, Complex(1.0, 0.0));
    for (octave_idx_type j = 0; j < nc; j++)
      {
!       for (octave_idx_type i = a.cidx(j); i < a.cidx(j+1); i++)
        {
          OCTAVE_QUIT;
!         result.xelem(a.ridx(i), j) = std::pow (a.data(i), b (a.ridx(i), j));
        }
      }
  
!   result.maybe_compress (true);
  
    return result;
  }
***************
*** 571,579 ****
          OCTAVE_QUIT;
          double btmp = b (i, j);
          if (xisint (btmp))
!           result (i, j) = pow (a, static_cast<int> (btmp));
          else
!           result (i, j) = pow (a, btmp);
        }
      }
  
--- 537,545 ----
          OCTAVE_QUIT;
          double btmp = b (i, j);
          if (xisint (btmp))
!           result (i, j) = std::pow (a, static_cast<int> (btmp));
          else
!           result (i, j) = std::pow (a, btmp);
        }
      }
  
***************
*** 592,598 ****
      for (octave_idx_type i = 0; i < nr; i++)
        {
        OCTAVE_QUIT;
!       result (i, j) = pow (a, b (i, j));
        }
  
    return result;
--- 558,564 ----
      for (octave_idx_type i = 0; i < nr; i++)
        {
        OCTAVE_QUIT;
!       result (i, j) = std::pow (a, b (i, j));
        }
  
    return result;
***************
*** 609,615 ****
        octave_idx_type nr = a.rows ();
        octave_idx_type nc = a.cols ();
  
!       ComplexMatrix result (nr, nc, Complex (pow (0.0, b)));
  
        if (xisint (b))
        {
--- 575,581 ----
        octave_idx_type nr = a.rows ();
        octave_idx_type nc = a.cols ();
  
!       ComplexMatrix result (nr, nc, Complex (std::pow (0.0, b)));
  
        if (xisint (b))
        {
***************
*** 618,624 ****
              {
                OCTAVE_QUIT;
                result (a.ridx(i), j) = 
!                 pow (a.data (i), static_cast<int> (b));
              }
        }
        else
--- 584,590 ----
              {
                OCTAVE_QUIT;
                result (a.ridx(i), j) = 
!                 std::pow (a.data (i), static_cast<int> (b));
              }
        }
        else
***************
*** 627,633 ****
            for (octave_idx_type i = a.cidx(j); i < a.cidx(j+1); i++)
              {
                OCTAVE_QUIT;
!               result (a.ridx(i), j) = pow (a.data (i), b);
              }
        }  
  
--- 593,599 ----
            for (octave_idx_type i = a.cidx(j); i < a.cidx(j+1); i++)
              {
                OCTAVE_QUIT;
!               result (a.ridx(i), j) = std::pow (a.data (i), b);
              }
        }  
  
***************
*** 644,650 ****
          for (octave_idx_type i = 0; i < nz; i++)
            {
              OCTAVE_QUIT;
!             result.data (i) = pow (a.data (i), static_cast<int> (b));
            }
        }
        else
--- 610,616 ----
          for (octave_idx_type i = 0; i < nz; i++)
            {
              OCTAVE_QUIT;
!             result.data (i) = std::pow (a.data (i), static_cast<int> (b));
            }
        }
        else
***************
*** 652,658 ****
          for (octave_idx_type i = 0; i < nz; i++)
            {
              OCTAVE_QUIT;
!             result.data (i) = pow (a.data (i), b);
            }
        }  
  
--- 618,624 ----
          for (octave_idx_type i = 0; i < nz; i++)
            {
              OCTAVE_QUIT;
!             result.data (i) = std::pow (a.data (i), b);
            }
        }  
  
***************
*** 680,717 ****
        return octave_value ();
      }
  
!   octave_idx_type nel = 0;
!   for (octave_idx_type j = 0; j < nc; j++) 
!     for (octave_idx_type i = 0; i < nr; i++)
!       if (!(a.elem (i, j) == 0. && b.elem (i, j) != 0.))
!       nel++;
! 
!   SparseComplexMatrix result (nr, nc, nel);
!   octave_idx_type ii = 0;
! 
!   result.cidx (0) = 0;
    for (octave_idx_type j = 0; j < nc; j++)
      {
!       for (octave_idx_type i = 0; i < nr; i++)
        {
          OCTAVE_QUIT;
!         double btmp = b (i, j);
          Complex tmp;
  
          if (xisint (btmp))
!           tmp = pow (a (i, j), static_cast<int> (btmp));
          else
!           tmp = pow (a (i, j), btmp);
!         if (tmp != 0.)
!           {
!             result.data (ii) = tmp; 
!             result.ridx (ii++) = i; 
!           }
        }
-       result.cidx (j+1) = ii;
      }
  
!   result.maybe_compress ();
  
    return result;
  }
--- 646,669 ----
        return octave_value ();
      }
  
!   SparseComplexMatrix result (nr, nc, Complex(1.0, 0.0));
    for (octave_idx_type j = 0; j < nc; j++)
      {
!       for (octave_idx_type i = a.cidx(j); i < a.cidx(j+1); i++)
        {
          OCTAVE_QUIT;
!         double btmp = b (a.ridx(i), j);
          Complex tmp;
  
          if (xisint (btmp))
!           result.xelem(a.ridx(i), j) = std::pow (a.data (i), 
!                                             static_cast<int> (btmp));
          else
!           result.xelem(a.ridx(i), j) = std::pow (a.data (i), btmp);
        }
      }
  
!   result.maybe_compress (true);
  
    return result;
  }
***************
*** 735,741 ****
        for (octave_idx_type i = 0; i < nz; i++)
        {
          OCTAVE_QUIT;
!         result.data (i) = pow (a.data (i), b);
        }
  
        result.maybe_compress (true);
--- 687,693 ----
        for (octave_idx_type i = 0; i < nz; i++)
        {
          OCTAVE_QUIT;
!         result.data (i) = std::pow (a.data (i), b);
        }
  
        result.maybe_compress (true);
***************
*** 762,790 ****
        return octave_value ();
      }
  
!   octave_idx_type nel = 0;
!   for (octave_idx_type j = 0; j < nc; j++) 
!     for (octave_idx_type i = 0; i < nr; i++)
!       if (!(a.elem (i, j) == 0. && b.elem (i, j) != 0.))
!       nel++;
! 
!   SparseComplexMatrix result (nr, nc, nel);
!   octave_idx_type ii = 0;
! 
!   result.cidx (0) = 0;
    for (octave_idx_type j = 0; j < nc; j++) 
      {
!       for (octave_idx_type i = 0; i < nr; i++)
        {
          OCTAVE_QUIT;
!         Complex tmp = pow (a (i, j), b (i, j));
!         if (tmp != 0.)
!           {
!             result.data (ii) = tmp;
!             result.ridx (ii++) = i;
!           }
        }
-       result.cidx (j+1) = ii;
      }
    result.maybe_compress (true);
  
--- 714,727 ----
        return octave_value ();
      }
  
!   SparseComplexMatrix result (nr, nc, Complex(1.0, 0.0));
    for (octave_idx_type j = 0; j < nc; j++) 
      {
!       for (octave_idx_type i = a.cidx(j); i < a.cidx(j+1); i++)
        {
          OCTAVE_QUIT;
!         result.xelem(a.ridx(i), j) = std::pow (a.data (i), b (a.ridx(i), j));
        }
      }
    result.maybe_compress (true);
  

reply via email to

[Prev in Thread] Current Thread [Next in Thread]