Complex left-divide and SVD give ireproducible results

View: New views
8 Messages — Rating Filter:   Alert me  

Complex left-divide and SVD give ireproducible results

by Dirk Laurie-2 :: Rate this Message:

Reply to Author | View Threaded | Show Only this Message

--------
Bug report for Octave 2.1.73 configured for i486-pc-linux-gnu    

Description:
-----------

The following operations give irreproducible results:
    z=A\b  when A is a nonsquare complex matrix
    [U,S,V]=svd(A) when A is a complex matrix

By 'irreproducible' is meant that the results are not identical
to the last bit when the operation is run more than once on
identical input, although the results in both cases are numerically
correct in the sense that b-A*z and A-U*S*V respectively are at
roundoff error level.  The two bugs are reported as one because
it is thought probable that they share a common cause.

A possible cause might be that an uninitialized value is used in
a comparison of which the result would normally not critically
affect the output, e.g. whether columns in a full-rank orthogonal
factorization are to be interchanged.


Repeat-By:
---------

The left-divide operation is so capricious that it suffices to
generate a random test case.  The following code exposes the bug:
 
A=rand(200,40)+1i*rand(200,40); b=rand(200,1); x0=A\b;
for iter=1:10, x=A\b; err(iter)=norm(x0-x,Inf); end
err

The otput should be all zeros.  If I save the above code in an M-file,
I get nonzero results most of the time, whether running the code in
an interactive shell or from the command line.

The SVD operation can also be exposed by a random test case,
albeit one in which the matrix is scaled down, as follows:

A=eps*(rand(40,8)+1i*rand(40,8)); [U0,S0,V0]=svd(A);
for iter=1:10,
    [U,S,V]=svd(A);
    err(iter)=max([norm(U-U0,Inf),norm(S-S0,Inf),norm(V-V0,Inf)]);
    end
err

The same comments as before apply to the output.

The same errors are also present, in fact manifesting more often, on
all more recent Octave releases that I have access to, including
Octave 3.0.0 as supplied with Ubuntu Hardy.


Configuration (please do not edit this section):
-----------------------------------------------

uname output:     Linux gawie 2.6.22-14-generic #1 SMP Tue Feb 12 07:42:25 UTC 2008 i686 GNU/Linux
configure opts:   '--prefix=/usr' '--datadir=/usr/share' '--libdir=/usr/lib' '--libexecdir=/usr/lib' '--infodir=/usr/share/info' '--mandir=/usr/share/man' '--with-blas=-lblas-3' '--with-lapack=-llapack-3' '--with-hdf5' '--with-fftw' '--with-f77=/usr/bin/gfortran' '--enable-shared' '--enable-rpath' '--disable-static' '--build' 'i486-linux-gnu' 'build_alias=i486-linux-gnu' 'CC=/usr/bin/gcc' 'CXX=/usr/bin/g++' 'F77=/usr/bin/gfortran'
Fortran compiler: /usr/bin/gfortran
FFLAGS:           -O2
F2C:              
F2CFLAGS:        
FLIBS:            -L/usr/lib/gcc/i486-linux-gnu/4.1.2 -L/usr/lib/gcc/i486-linux-gnu/4.1.2/../../../../lib -L/lib/../lib -L/usr/lib/../lib -lhdf5 -lz -lgfortranbegin -lgfortran -lm
CPPFLAGS:        
INCFLAGS:         -I. -I. -I./liboctave -I./src -I./libcruft/misc
C compiler:       /usr/bin/gcc, version 4.1.2 20061103 (prerelease) (Ubuntu 4.1.1-18ubuntu2)
CFLAGS:           -O2
CPICFLAG:         -fPIC
C++ compiler:     /usr/bin/g++, version 4.1.2
CXXFLAGS:         -O2
CXXPICFLAG:       -fPIC
LD_CXX:           /usr/bin/g++
LDFLAGS:          -s
LIBFLAGS:         -L.
RLD_FLAG:         -Wl,-rpath -Wl,/usr/lib/octave-2.1.73
BLAS_LIBS:        -llapack-3 -lblas-3
FFTW_LIBS:        -lfftw3
LIBS:             -lreadline  -lncurses -ldl -lhdf5 -lz -lm
LEXLIB:          
LIBPLPLOT:        
LIBDLFCN:        
LIBGLOB:          %OCTAVE_CONF_LIBGLOB%
SED:              /bin/sed
DEFS:

  -DPACKAGE_NAME="" -DPACKAGE_TARNAME="" -DPACKAGE_VERSION=""
  -DPACKAGE_STRING="" -DPACKAGE_BUGREPORT="" -DOCTAVE_SOURCE=1
  -D_GNU_SOURCE=1 -DSTDC_HEADERS=1 -DHAVE_SYS_TYPES_H=1
  -DHAVE_SYS_STAT_H=1 -DHAVE_STDLIB_H=1 -DHAVE_STRING_H=1
  -DHAVE_MEMORY_H=1 -DHAVE_STRINGS_H=1 -DHAVE_INTTYPES_H=1 -DHAVE_STDINT_H=1
  -DHAVE_UNISTD_H=1 -DSEPCHAR=':' -DSEPCHAR_STR=":" -D__NO_MATH_INLINES=1
  -DCXX_NEW_FRIEND_TEMPLATE_DECL=1 -DCXX_ISO_COMPLIANT_LIBRARY=1
  -DCXX_ABI=gnu_v3 -DHAVE_LIBM=1 -DHAVE_HDF5_H=1 -DHAVE_HDF5=1
  -DHAVE_H5GGET_NUM_OBJS=1 -DHAVE_FFTW3=1 -DHAVE_IEEE754_DATA_FORMAT=1
  -DF77_FUNC(name,NAME)=name ## _ -DF77_FUNC_(name,NAME)=name ## _
  -DHAVE_BLAS=1 -DHAVE_GETHOSTNAME=1 -DHAVE_GETPWNAM=1 -DHAVE_DEV_T=1
  -DHAVE_INO_T=1 -DHAVE_NLINK_T=1 -DHAVE_NLINK_T=1 -DHAVE_LONG_LONG_INT=1
  -DHAVE_UNSIGNED_LONG_LONG_INT=1 -DHAVE_SIGSET_T=1 -DHAVE_SIG_ATOMIC_T=1
  -DSIZEOF_SHORT=2 -DSIZEOF_INT=4 -DSIZEOF_LONG=4 -DSIZEOF_LONG_LONG=8
  -DHAVE_ALLOCA_H=1 -DHAVE_ALLOCA=1 -DNPOS=std::string::npos
  -DHAVE_PLACEMENT_DELETE=1 -DHAVE_DYNAMIC_AUTO_ARRAYS=1 -DSTDC_HEADERS=1
  -DHAVE_DIRENT_H=1 -DTIME_WITH_SYS_TIME=1 -DHAVE_SYS_WAIT_H=1
  -DHAVE_ASSERT_H=1 -DHAVE_CURSES_H=1 -DHAVE_DLFCN_H=1 -DHAVE_FCNTL_H=1
  -DHAVE_FLOAT_H=1 -DHAVE_GRP_H=1 -DHAVE_LIMITS_H=1 -DHAVE_MEMORY_H=1
  -DHAVE_NCURSES_H=1 -DHAVE_POLL_H=1 -DHAVE_PWD_H=1 -DHAVE_STDLIB_H=1
  -DHAVE_STRING_H=1 -DHAVE_SYS_IOCTL_H=1 -DHAVE_SYS_PARAM_H=1
  -DHAVE_SYS_POLL_H=1 -DHAVE_SYS_RESOURCE_H=1 -DHAVE_SYS_SELECT_H=1
  -DHAVE_SYS_STAT_H=1 -DHAVE_SYS_TIME_H=1 -DHAVE_SYS_TIMES_H=1
  -DHAVE_SYS_TYPES_H=1 -DHAVE_SYS_UTSNAME_H=1 -DHAVE_TERMCAP_H=1
  -DHAVE_UNISTD_H=1 -DHAVE_SSTREAM=1 -DHAVE_TERMIO_H=1
  -DHAVE_SGTTY_H=1 -DHAVE_ATEXIT=1 -DHAVE_BASENAME=1 -DHAVE_BCOPY=1
  -DHAVE_BZERO=1 -DHAVE_CANONICALIZE_FILE_NAME=1 -DHAVE_DUP2=1
  -DHAVE_ENDGRENT=1 -DHAVE_ENDPWENT=1 -DHAVE_EXECVP=1 -DHAVE_FCNTL=1
  -DHAVE_FORK=1 -DHAVE_GETCWD=1 -DHAVE_GETEGID=1 -DHAVE_GETEUID=1
  -DHAVE_GETGID=1 -DHAVE_GETGRENT=1 -DHAVE_GETGRGID=1 -DHAVE_GETGRNAM=1
  -DHAVE_GETPGRP=1 -DHAVE_GETPID=1 -DHAVE_GETPPID=1 -DHAVE_GETPWENT=1
  -DHAVE_GETPWUID=1 -DHAVE_GETTIMEOFDAY=1 -DHAVE_GETUID=1 -DHAVE_GETWD=1
  -DHAVE_KILL=1 -DHAVE_LINK=1 -DHAVE_LOCALTIME_R=1 -DHAVE_LSTAT=1
  -DHAVE_MEMMOVE=1 -DHAVE_MKDIR=1 -DHAVE_MKFIFO=1 -DHAVE_MKSTEMP=1
  -DHAVE_ON_EXIT=1 -DHAVE_PIPE=1 -DHAVE_POLL=1 -DHAVE_PUTENV=1
  -DHAVE_RAISE=1 -DHAVE_READLINK=1 -DHAVE_RENAME=1 -DHAVE_RINDEX=1
  -DHAVE_RMDIR=1 -DHAVE_ROUND=1 -DHAVE_SELECT=1 -DHAVE_SETGRENT=1
  -DHAVE_SETPWENT=1 -DHAVE_SETVBUF=1 -DHAVE_SIGACTION=1 -DHAVE_SIGLONGJMP=1
  -DHAVE_SIGPENDING=1 -DHAVE_SIGPROCMASK=1 -DHAVE_SIGSUSPEND=1 -DHAVE_STAT=1
  -DHAVE_STRCASECMP=1 -DHAVE_STRDUP=1 -DHAVE_STRERROR=1 -DHAVE_STRFTIME=1
  -DHAVE_STRNCASECMP=1 -DHAVE_STRPTIME=1 -DHAVE_SYMLINK=1 -DHAVE_TEMPNAM=1
  -DHAVE_UMASK=1 -DHAVE_UNLINK=1 -DHAVE_USLEEP=1 -DHAVE_VFPRINTF=1
  -DHAVE_VSPRINTF=1 -DHAVE_VSNPRINTF=1 -DHAVE_WAITPID=1 -DHAVE_LIBDL=1
  -DHAVE_DLOPEN=1 -DHAVE_DLSYM=1 -DHAVE_DLERROR=1 -DHAVE_DLCLOSE=1
  -DHAVE_DLOPEN_API=1 -DENABLE_DYNAMIC_LINKING=1 -DHAVE_TIMEVAL=1
  -DHAVE_FINITE=1 -DHAVE_ISNAN=1 -DHAVE_ISINF=1 -DHAVE_COPYSIGN=1
  -DHAVE_DECL_SIGNBIT=1 -DHAVE_ACOSH=1 -DHAVE_ASINH=1 -DHAVE_ATANH=1
  -DHAVE_ERF=1 -DHAVE_ERFC=1 -DHAVE_STRUCT_STAT_ST_BLKSIZE=1
  -DHAVE_STRUCT_STAT_ST_BLOCKS=1 -DHAVE_STRUCT_STAT_ST_RDEV=1
  -DHAVE_STRUCT_TM_TM_ZONE=1 -DHAVE_TM_ZONE=1 -DUSE_READLINE=1
  -DEXCEPTION_IN_MATH=1 -DRETSIGTYPE=void -DSYS_SIGLIST_DECLARED=1
  -DHAVE_SYS_SIGLIST=1 -DHAVE_POSIX_SIGNALS=1 -DHAVE_GETRUSAGE=1
  -DHAVE_TIMES=1 -DGNUPLOT_BINARY="gnuplot" -DGNUPLOT_HAS_FRAMES=

User-preferences (please do not edit this section):
--------------------------------------------------

  DEFAULT_EXEC_PATH = "/usr/lib/octave/2.1.73/site/exec/i486-pc-linux-gnu:/usr/lib/octave/site/exec/i486-pc-linux-gnu:/usr/lib/octave/2.1.73/exec/i486-pc-linux-gnu:/usr/bin"
  DEFAULT_LOADPATH = ".:/usr/lib/octave/2.1.73/site/oct/i486-pc-linux-gnu//:/usr/lib/octave/site/oct/api-v13/i486-pc-linux-gnu//:/usr/lib/octave/site/oct/i486-pc-linux-gnu//:/usr/share/octave/2.1.73/site/m//:/usr/share/octave/site/api-v13/m//:/usr/share/octave/site/m//:/usr/lib/octave/2.1.73/oct/i486-pc-linux-gnu//:/usr/share/octave/2.1.73/m//"
  EDITOR = "emacs"
  EXEC_PATH = ":/home/dirk/bin:/usr/local/sbin:/usr/local/bin:/usr/sbin:/usr/bin:/sbin:/bin:/usr/bin/X11:/usr/games"
  IMAGEPATH = ".:/usr/share/octave/2.1.73/imagelib//"
  INFO_FILE = "/usr/share/info/octave2.1.info"
  INFO_PROGRAM = "info"
  LOADPATH = "/home/dirk/octave/lib//::/usr/local/share/octave/site-m//:"
  PAGER = "pager"
  PS1 = "> "
  PS2 =
  PS4 = "+ "
  automatic_replot = 1
  beep_on_error = 1
  completion_append_char = " "
  crash_dumps_octave_core = 1
  default_save_format = "ascii"
  echo_executing_commands = 0
  fixed_point_format = 0
  gnuplot_binary = "gnuplot"
  gnuplot_command_end = "\n"
  gnuplot_command_plot = "pl"
  gnuplot_command_replot = "rep"
  gnuplot_command_splot = "sp"
  gnuplot_command_title = "t"
  gnuplot_command_using = "u"
  gnuplot_command_with = "w"
  gnuplot_has_frames = 1
  history_file = "/home/dirk/.octave_hist"
  history_size = 1024
  ignore_function_time_stamp = "system"
  max_recursion_depth = 256
  output_max_field_width = 20
  output_precision = 15
  page_output_immediately = 0
  page_screen_output = 1
  print_answer_id_name = 1
  print_empty_dimensions = 1
  print_rhs_assign_val = 0
# return_last_computed_value = <no value or error in displaying it>
  save_precision = 15
  saving_history = 1
  sighup_dumps_octave_core = 1
  sigterm_dumps_octave_core = 1
  silent_functions = 0
  split_long_rows = 1
  string_fill_char = " "
  struct_levels_to_print = 2
  suppress_verbose_help_message = 1
  warn_assign_as_truth_value = 1
  warn_divide_by_zero = 1
  warn_empty_list_elements = 0
  warn_fortran_indexing = 0
  warn_function_name_clash = 0
  warn_future_time_stamp = 1
  warn_imag_to_real = 0
  warn_missing_semicolon = 0
  warn_neg_dim_as_zero = 0
  warn_num_to_str = 1
  warn_resize_on_range_error = 0
  warn_separator_insert = 0
  warn_single_quote_string = 0
  warn_str_to_num = 0
  warn_undefined_return_values = 1
  warn_variable_switch_label = 0
_______________________________________________
Bug-octave mailing list
Bug-octave@...
https://www.cae.wisc.edu/mailman/listinfo/bug-octave

Re: Complex left-divide and SVD give ireproducible results

by David Bateman :: Rate this Message:

Reply to Author | View Threaded | Show Only this Message

Dirk Laurie wrote:

> --------
> Bug report for Octave 2.1.73 configured for i486-pc-linux-gnu    
>
> Description:
> -----------
>
> The following operations give irreproducible results:
>     z=A\b  when A is a nonsquare complex matrix
>     [U,S,V]=svd(A) when A is a complex matrix
>
> By 'irreproducible' is meant that the results are not identical
> to the last bit when the operation is run more than once on
> identical input, although the results in both cases are numerically
> correct in the sense that b-A*z and A-U*S*V respectively are at
> roundoff error level.  The two bugs are reported as one because
> it is thought probable that they share a common cause.
>
> A possible cause might be that an uninitialized value is used in
> a comparison of which the result would normally not critically
> affect the output, e.g. whether columns in a full-rank orthogonal
> factorization are to be interchanged.
>  

Note that 2.1.73 is an older release and you should probably be using
3.0.1, but that doesn't change this discussion. In fact due to the
change in the use of the lapack xGELSS function to use xGELSY instead in
recent versions of Octave, with your test case I can no longer exposes
the issue. There has been many e-mails about issues similar to this to
the Octave lists., and its related to many many different operators and
functions in Octave of floating point values.

Basically, I think you just have to live with it for the latest
compilers and if you use an optimized blas. My take on the reason for
such issues is that the x86 processors have 80 bit internal registers
and optimized blas implementations and the 4.0.x release of the gnu tool
chain and later make use of these additional bits of precision. That is
between operations on the x86 registers, the floating-point values don't
have to be stored back to the 64-bit value that contains it saving the
cost of a memory copy from and back to the registers. This can give
significant improvement in the performance and a marginal (and
unpredictable) improvement in the accuracy.

However, during a context switch on the x86 processors, all of the
registers are swapped back to memory and those 64-bit values and because
you can't control when the context switches might happen when the matrix
sizes start to get large you can get unpredictable differences in the
numerical solutions.

So what can you do if you absolutely have to always have the same
numerical values

* Don't use optimized version of the blas/lapack
* Use a gnu toolchain earlier than 4.0.x or use the -ffloat-store option
to g++, gcc and gfortran.

Though the cost in terms of speed will be horrendous. Also, if you
algorithm is so sensitivity to those last couple of bits of accuracy you
probably need to reconsider your algorithm. In any case I don't consider
this to be a bug in Octave as such.

Cheers
David


--
David Bateman                                David.Bateman@...
Motorola Labs - Paris                        +33 1 69 35 48 04 (Ph)
Parc Les Algorithmes, Commune de St Aubin    +33 6 72 01 06 33 (Mob)
91193 Gif-Sur-Yvette FRANCE                  +33 1 69 35 77 01 (Fax)

The information contained in this communication has been classified as:

[x] General Business Information
[ ] Motorola Internal Use Only
[ ] Motorola Confidential Proprietary

_______________________________________________
Bug-octave mailing list
Bug-octave@...
https://www.cae.wisc.edu/mailman/listinfo/bug-octave

Complex left-divide and SVD give ireproducible results

by John W. Eaton :: Rate this Message:

Reply to Author | View Threaded | Show Only this Message

On 28-Jul-2008, Dirk Laurie wrote:

| --------
| Bug report for Octave 2.1.73 configured for i486-pc-linux-gnu    
|
| Description:
| -----------
|
| The following operations give irreproducible results:
|     z=A\b  when A is a nonsquare complex matrix
|     [U,S,V]=svd(A) when A is a complex matrix
|
| By 'irreproducible' is meant that the results are not identical
| to the last bit when the operation is run more than once on
| identical input, although the results in both cases are numerically
| correct in the sense that b-A*z and A-U*S*V respectively are at
| roundoff error level.  The two bugs are reported as one because
| it is thought probable that they share a common cause.
|
| A possible cause might be that an uninitialized value is used in
| a comparison of which the result would normally not critically
| affect the output, e.g. whether columns in a full-rank orthogonal
| factorization are to be interchanged.
|
|
| Repeat-By:
| ---------
|
| The left-divide operation is so capricious that it suffices to
| generate a random test case.  The following code exposes the bug:
|  
| A=rand(200,40)+1i*rand(200,40); b=rand(200,1); x0=A\b;
| for iter=1:10, x=A\b; err(iter)=norm(x0-x,Inf); end
| err
|
| The otput should be all zeros.  If I save the above code in an M-file,
| I get nonzero results most of the time, whether running the code in
| an interactive shell or from the command line.
|
| The SVD operation can also be exposed by a random test case,
| albeit one in which the matrix is scaled down, as follows:
|
| A=eps*(rand(40,8)+1i*rand(40,8)); [U0,S0,V0]=svd(A);
| for iter=1:10,
|     [U,S,V]=svd(A);
|     err(iter)=max([norm(U-U0,Inf),norm(S-S0,Inf),norm(V-V0,Inf)]);
|     end
| err
|
| The same comments as before apply to the output.
|
| The same errors are also present, in fact manifesting more often, on
| all more recent Octave releases that I have access to, including
| Octave 3.0.0 as supplied with Ubuntu Hardy.

Of the systems I have easy access to, I'm only able to reproduce this
on one with Octave 2.1.73.  The failure happens on a 32-bit AMD system
running Debian.  I'm unable to reproduce the problem with Octave 3.0.0
on the same system.  The difference appears to be that Octave 3.0 is
linked with

  liblapack.so.3gf => /usr/lib/liblapack.so.3gf (0xa633d000)
  libblas.so.3gf => /usr/lib/libblas.so.3gf (0xa627c000)

which are provided by

  ii  libblas3gf                 1.2-1.5
  ii  liblapack3gf               3.1.1-0.4

and Octave 2.1.73 is linked with

  liblapack.so.3 => /usr/lib/atlas/liblapack.so.3 (0xa6b2e000)
  libblas.so.3 => /usr/lib/atlas/libblas.so.3 (0xa67cf000)

which are provided by

  ii  blas                       1.1-14
  ii  lapack3                    3.0.20000531a-6.1

I'd guess that the problem is in the linear algebra package, not
Octave.

jwe
_______________________________________________
Bug-octave mailing list
Bug-octave@...
https://www.cae.wisc.edu/mailman/listinfo/bug-octave

Re: Complex left-divide and SVD give ireproducible results

by Sergei Steshenko-2 :: Rate this Message:

Reply to Author | View Threaded | Show Only this Message




--- On Mon, 7/28/08, David Bateman <David.Bateman@...> wrote:

> From: David Bateman <David.Bateman@...>
> Subject: Re: Complex left-divide and SVD give ireproducible results
> To: "Dirk Laurie" <dpl@...>
> Cc: bug@...
> Date: Monday, July 28, 2008, 8:24 AM
> Dirk Laurie wrote:
> > --------
> > Bug report for Octave 2.1.73 configured for
> i486-pc-linux-gnu    
...

> Basically, I think you just have to live with it for the
> latest
> compilers and if you use an optimized blas. My take on the
> reason for
> such issues is that the x86 processors have 80 bit internal
> registers
> and optimized blas implementations and the 4.0.x release of
> the gnu tool
> chain and later make use of these additional bits of
> precision. That is
> between operations on the x86 registers, the floating-point
> values don't
> have to be stored back to the 64-bit value that contains it
> saving the
> cost of a memory copy from and back to the registers. This
> can give
> significant improvement in the performance and a marginal
> (and
> unpredictable) improvement in the accuracy.
>
> However, during a context switch on the x86 processors, all
> of the
> registers are swapped back to memory and those 64-bit
> values and because
> you can't control when the context switches might
> happen when the matrix
> sizes start to get large you can get unpredictable
> differences in the
> numerical solutions.
>
> So what can you do if you absolutely have to always have
> the same
> numerical values
...
> Cheers
> David
>
>
> --
> David Bateman                              

Won't using SSE2 arithmetic fix the problem ?

I think, though I'm not 100% sure, that SSE2 FP registers aew 64 bit,
not 80 bit.

--Sergei.


     
_______________________________________________
Bug-octave mailing list
Bug-octave@...
https://www.cae.wisc.edu/mailman/listinfo/bug-octave

Re: Complex left-divide and SVD give ireproducible results

by Miroslaw Kwasniak :: Rate this Message:

Reply to Author | View Threaded | Show Only this Message

On Mon, Jul 28, 2008 at 11:38:26AM -0400, John W. Eaton wrote:

> Of the systems I have easy access to, I'm only able to reproduce this
> on one with Octave 2.1.73.  The failure happens on a 32-bit AMD system
> running Debian.  I'm unable to reproduce the problem with Octave 3.0.0
> on the same system.  The difference appears to be that Octave 3.0 is
> linked with
>
>   liblapack.so.3gf => /usr/lib/liblapack.so.3gf (0xa633d000)
>   libblas.so.3gf => /usr/lib/libblas.so.3gf (0xa627c000)
>
> which are provided by
>
>   ii  libblas3gf                 1.2-1.5
>   ii  liblapack3gf               3.1.1-0.4

Reproducible with 3.0.1 on debian lenny amd32

-  octave3.0 1:3.0.1-4
        liblapack.so.3gf => /usr/lib/atlas/liblapack.so.3gf (0xb638f000)
        libblas.so.3gf => /usr/lib/atlas/libblas.so.3gf (0xb6013000)
- libatlas3gf-base 3.6.0-21.5

This error disappears when I run with libraries provided by:
  libblas3gf     1.2-1.6  
  liblapack3gf   3.1.1-0.5  

LD_PRELOAD=/usr/lib/libblas.so.3gf.0:/usr/lib/liblapack.so.3gf.0 \
/usr/bin/octave-3.0.1

I found an anothor issue - inefficiency of debian's atlas on my Athlon XP
1600+:

A=rand(200,40)+1i*rand(200,40); b=rand(200,1);
t0=cputime;
x0=A\b;for iter=1:1000, x=A\b; err(iter)=norm(x0-x,Inf); end;
sum(err~=0),cputime-t0

atlas ~10.5 seconds
blas+lapack  ~9.2 seconds

> I'd guess that the problem is in the linear algebra package, not
> Octave.

Yes but I can't remove atlas on my machine because harminv forced it.

Mirek
_______________________________________________
Bug-octave mailing list
Bug-octave@...
https://www.cae.wisc.edu/mailman/listinfo/bug-octave

Re: Complex left-divide and SVD give ireproducible results

by dbateman :: Rate this Message:

Reply to Author | View Threaded | Show Only this Message


Sergei Steshenko-2 wrote:
Won't using SSE2 arithmetic fix the problem ?

I think, though I'm not 100% sure, that SSE2 FP registers aew 64 bit,
not 80 bit.
Hey I didn't know that. In that case I suppose yes, but are you really sure that all floating point operations use SSE2/SSE3 ? Even if a few don't and a context switch happens in between then the same situation might arise, but just less frequently.

Cheers
D.

Re: Complex left-divide and SVD give ireproducible results

by John W. Eaton :: Rate this Message:

Reply to Author | View Threaded | Show Only this Message

On 28-Jul-2008, Miroslaw Kwasniak wrote:

| Yes but I can't remove atlas on my machine because harminv forced it.

That sounds like a packaging problem, so this is probably the wrong
place to report it.

jwe
_______________________________________________
Bug-octave mailing list
Bug-octave@...
https://www.cae.wisc.edu/mailman/listinfo/bug-octave

Re: Complex left-divide and SVD give ireproducible results

by Thomas Weber-8 :: Rate this Message:

Reply to Author | View Threaded | Show Only this Message

On 28/07/08 21:50 +0200, Miroslaw Kwasniak wrote:
> On Mon, Jul 28, 2008 at 11:38:26AM -0400, John W. Eaton wrote:
>
> Yes but I can't remove atlas on my machine because harminv forced it.

That's actually pretty unlikely[1]. The package probably prefers atlas, but
it should be possible to deinstall it and install blas/lapack.

If not, file a bug against harminv.

[1] Atlas is not available on all architectures Debian supports, so most
maintainers use OR alternatives for provdiding the dependencies, meaning
atlas will be preferred, but you can install blas and deinstall atlas
easily.

        Thomas
_______________________________________________
Bug-octave mailing list
Bug-octave@...
https://www.cae.wisc.edu/mailman/listinfo/bug-octave