Mark Starnes wrote:
>  Hi everyone,
> 
>  I've been trying to get this simple example to work for a couple of days 
>  now, to no avail.  I noticed at 
>  http://projects.scipy.org/scipy/scipy/ticket/329 that there seem to be 
>  come problems with csc and complex values but I'm hoping they're 
>  unrelated and someone can help me sort this problem!  I couldn't get the 
>  csc matrix to work at all here - the csr matrix was having a go and 
>  getting the right answer but now fails too!  Anyway, here's the code.  
>  I'm using Python2.4.3 on winxp with Scipy0.5.2, Numpy 1.0.1.
> 
>  Any help will be appreciated.
> 
>  Best regards,
> 
>  Mark Starnes.
> 
>  ***
> 
>  import numpy as n
>  import scipy as s
>  from scipy.linsolve import spsolve
>  from scipy import linalg as la
>  from numpy import zeros,complex64
> 
>  # solve kx=f for x, where
>  #
>  #/                    \  /    \    /       > #| 1+1i  2+2i  3+3i   | |   2  |   | -5+3i  |
>  #|                    | |      |   |        |
>  #| 2     3     4+4i   | | 3+2i | = | 1-6i   |
>  #|                    | |      |   |        |
>  #| 4+9i 5+10i  2+11i  | |  -3  |   | -3+25i |
>  #\                    /  \    /    \       /
>  #
> 
> 
>  k0=n.array([[1+1j,2+2j,3+3j],[2,3,4+4j],[4+9j,5+10j,2+11j]])
>  f0=[-5+3j,1-6j,-3+25j]
> 
>  # Test numpy routine.
> 
>  k=n.mat(k0)
>  f=n.mat(f0).T
>  fsol=s.dot(la.inv(k),f)
>  print fsol
>  print k*fsol
>  # result should = [2, 3+2j, -3].  Numpy passes.
> 
>  # Test Scipy sparse routines.
>  k=s.sparse.csc_matrix((3,3),dtype='F')
>  k[0,0]=1+1j ; k[0,1]=2+2j  ; k[0,2]=3+3j ;
>  k[1,0]=2    ; k[1,1]=3     ; k[1,2]=4+4j ;
>  k[2,0]=4+9j ; k[2,1]=5+10j ; k[2,2]=2+11j;
> 
>  f=zeros((3),complex64)
>  f[0]=-5+3j
>  f[1]=1-6j
>  f[2]=-3+25j
> 
>  #fsol2=spsolve(k,f)   # this crashes python with typeD, AND with type F
>  #print fsol2
>  #print k.dot(fsol2)
> 
>  k=s.sparse.csr_matrix((3,3),dtype='F')
>  k[0,0]=1+1j ; k[0,1]=2+2j  ; k[0,2]=3+3j ;
>  k[1,0]=2    ; k[1,1]=3     ; k[1,2]=4+4j ;
>  k[2,0]=4+9j ; k[2,1]=5+10j ; k[2,2]=2+11j;
> 
>  fsol3=spsolve(k,f)   # this crashes python with typeD, originally ok 
>  with type F
>  print fsol3         
>  print k.dot(fsol3)  
> 
>  _______________________________________________
>  SciPy-user mailing list
>  SciPy-user@[...].org
>  http://projects.scipy.org/mailman/listinfo/scipy-user
>    
Hi Mark,

This seems to be a an UMFPACK problem. It works fine for me with SuperLU.
I am using the latest svn versions and UMFPACK v4.4

Nils

With UMFPACK

python -i starnes.py
x [[ 2. +1.11022302e-16j]
 [ 3. +2.00000000e+00j]
 [-3. -7.28583860e-17j]]
f [[-5. +3.j]
 [ 1. -6.j]
 [-3.+25.j]]
Traceback (most recent call last):
  File "starnes.py", line 47, in ?
    fsol2=spsolve(k,f)   #
  File "/usr/lib64/python2.4/site-packages/scipy/linsolve/linsolve.py",
line 74, in spsolve
    return umf.linsolve( umfpack.UMFPACK_A, mat, b, autoTranspose = True )
  File
"/usr/lib64/python2.4/site-packages/scipy/linsolve/umfpack/umfpack.py",
line 561, in linsolve
    self.numeric( mtx )
  File
"/usr/lib64/python2.4/site-packages/scipy/linsolve/umfpack/umfpack.py",
line 374, in numeric
    self.symbolic( mtx )
  File
"/usr/lib64/python2.4/site-packages/scipy/linsolve/umfpack/umfpack.py",
line 356, in symbolic
    raise RuntimeError, '%s failed with %s' % (self.funs.symbolic,
RuntimeError: <function umfpack_zi_symbolic at 0x2aaaadff7500>  failed
with UMFPACK_ERROR_invalid_matrix

Without UMFPACK

python -i starnes.py
x [[ 2. +1.11022302e-16j]
 [ 3. +2.00000000e+00j]
 [-3. -7.28583860e-17j]]
f [[-5. +3.j]
 [ 1. -6.j]
 [-3.+25.j]]
Use minimum degree ordering on A'+A.
x [ 2. -7.58729969e-17j  3. +2.00000000e+00j -3. -2.77108197e-16j]
f [-5. +3.j  1. -6.j -3.+25.j]
Use minimum degree ordering on A'+A.
x [ 2. -3.44633288e-16j  3. +2.00000000e+00j -3. -2.91704595e-16j]
f [-5. +3.j  1. -6.j -3.+25.j]


Nils