Re: LU decomposition of complex dense matrices

David Lidsky wrote:
> /*
>   Sample Output 
>   [(1,0),(0,-2),]
> 	Largest element in the vector z is at location 1
>  */
> // How can we make it give the desired answer:
> // "Largest element in the vector z is at location 2"

Your code works fine and '1' is the correct result. The first value has 
index 0, the second has index 1.

Btw. you should apply the appended patch to mtl.h which fixes the 'bug' 
in max_index (which now returns the index of the maximum). That means 
max_index does not work for complex any more (complex is not 
LessThanComparable). You should use max_abs_index, which returns the 
index of the element with max. magnitude. If your STL does not contain 
an function "double& abs(complex<double>)" just add this to your code:

namespace std {
double abs(const complex<double>& z)
{ return (z.real()*z.real()+z.imag()*z.imag()); }
}; // namespace std

In my tests LU did not work for fixed size matrices, maybe some one 
should look over this. For dense matrices its always best to use CLAPACK.

regards
Gunter Winkler

*** mtl.h	Mon Sep 30 14:48:36 2002
--- mtl.h.old	Mon Sep 30 14:52:58 2002
***************
*** 666,672 ****
  //!component: function
  //!definition: mtl.h
  //!complexity: O(n)
! // The location (index) of the element with the maximum value.
  //!example: max_index.cc
  //!typereqs: <tt>Vec::value_type</tt> must be LessThanComparible.
  template <class Vec>
--- 666,672 ----
  //!component: function
  //!definition: mtl.h
  //!complexity: O(n)
! // The location (index) of the element with the maximum absolute value.
  //!example: max_index.cc
  //!typereqs: <tt>Vec::value_type</tt> must be LessThanComparible.
  template <class Vec>
***************
*** 674,680 ****
  max_index(const Vec& x)
  {
    typename Vec::const_iterator maxi =
!     mtl_algo::max_element(x.begin(), x.end());
    return maxi.index();
  }
  
--- 674,680 ----
  max_index(const Vec& x)
  {
    typename Vec::const_iterator maxi =
!     mtl_algo::max_element(x.begin(), x.end(), abs_cmp());
    return maxi.index();
  }