Re: [eigen] mixed-precision Cholesky?

[Gael Guennebaud <gael.guennebaud@xxxxxxxxx>]

if you want to give it a try, you can checkout
http://bitbucket.org/ggael/eigen-complex. Here LLT and LDLT solving
with a real matrix and a complex rhs (vector only atm) works fine.

gael

On Tue, Jul 13, 2010 at 2:15 PM, Gael Guennebaud
<gael.guennebaud@xxxxxxxxx> wrote:
> ah ok, I got confused because now Dynamic==-1, instead of 3331...
> Anyway, as explained in the other thread the changes I were talking
> about are happening in the eigen-complex fork. I'd prefer to have
> something relatively stable before merging to the default branch.
>
> gael
>
> On Tue, Jul 13, 2010 at 4:06 AM, Manoj Rajagopalan <rmanoj@xxxxxxxxx> wrote:
>>
>> I am using the default branch. I just pulled the entire eigen repo into a
>> clean directory, performed 'hg branch default' followed by a make install -
>> same error list indicating an invalid attempt to access a pointer at this
>> line, for a CwiseUnaryOp<ei_scalar_cast_op<...> > expression.
>>
>> Shall I still try the experiment you suggested?
>>
>> -- Manoj
>>
>>
>> On Monday 12 July 2010 10:50:09 am Gael Guennebaud wrote:
>>> oh, I see you are using the 2.0 branch... indeed in this very special
>>> case calling the high level routines would probably work though it
>>> will skip the efficient matrix-vector implementation.
>>>
>>> For instance, try to replace the line containing the
>>> ei_cache_friendly_product_colmajor_times_vector by what is commented
>>> below this line, something like:
>>>
>>> //                            if (IsLowerTriangular)
>>> //           other.col(c).end(size-endBlock) += (lhs.block(endBlock,
>>> startBlock, size-endBlock, endBlock-startBlock)
>>> //                                           *
>>> other.col(c).block(startBlock,endBlock-startBlock)).lazy();
>>> //                            else
>>> //           other.col(c).end(size-endBlock) += (lhs.block(endBlock,
>>> startBlock, size-endBlock, endBlock-startBlock)
>>> //                                           *
>>> other.col(c).block(startBlock,endBlock-startBlock)).lazy();
>>>
>>>
>>> Beside I'm currently working on a true support for mixing real and
>>> complex with vectorization and optimal performance. See the other
>>> threads on that topic. So this should work out of the box soon.
>>>
>>> gael
>>>
>>> On Mon, Jul 12, 2010 at 3:56 PM, Manoj Rajagopalan <rmanoj@xxxxxxxxx> wrote:
>>> > On Thursday 10 June 2010 05:26:30 am Gael Guennebaud wrote:
>>> >> On Thu, Jun 10, 2010 at 1:37 AM, Manoj Rajagopalan <rmanoj@xxxxxxxxx>
>> wrote:
>>> >> > Hi,
>>> >> >
>>> >> >   I have situation where, due to memory pressure, I need to store a
>>> >> > symmetric
>>> >> > matrix as real-values but it must solve a complex RHS for a complex
>>> >> > solution
>>> >> > vector. Currently I am storing the matrix as complex-type so that
>>> >> > llt() and ldlt() can be called without compilation errors but is the
>>> >> > above situation realizable in Eigen?
>>> >> >
>>> >> >  When invoking
>>> >> > real_matrix.selfadjointView<Lower>().ldlt().solveInPlace(complex_vecto
>>> >> >r)
>>> >> >
>>> >> >  I get a gcc 4.2.4 error saying:
>>> >> >
>>> >> > no matching function for call
>>> >> > to
>>> >> > ‘ei_cache_friendly_product_colmajor_times_vector(Eigen::DenseIndex&,
>>> >> > double*, Eigen::DenseIndex,
>>> >> > Eigen::VectorBlock<Eigen::Matrix<std::complex<double>, 33331, 1, 0,
>>> >> > 33331, 1>, 33331>, std::complex<double>*, std::complex<double>)’
>>> >> >
>>> >> >
>>> >> > no matching function for call
>>> >> > to ‘ei_cache_friendly_product_rowmajor_times_vector(double*,
>>> >> > Eigen::DenseIndex, std::complex<double>*, Eigen::DenseIndex&,
>>> >> > Eigen::VectorBlock<Eigen::Matrix<std::complex<double>, 33331, 1, 0,
>>> >> > 33331, 1>, 33331>&, std::complex<double>)’
>>> >> >
>>> >> >
>>> >> >
>>> >> >   Would it make sense to introduce a cast<> on LLT and LDLT so that
>>> >> > such a mixed-precision operation becomes possible?
>>> >>
>>> >> The idea is rather to allow mixing types where it makes sense, as for
>>> >> products, triangular solves and the like, but I don't know when this
>>> >> will be available.
>>> >>
>>> >> gael
>>> >
>>> > Will replacing the call to ei_cache_friendly_product-*() with a
>>> > high-level Eigen GEMM routine be a simple fix for the problem? When I
>>> > call
>>> > solveInPlace() on
>>> >   real_matrix.cast<complex_type>().selfadjointView<Lower>().ldlt(),
>>> > this is the point where I get a compile error since
>>> > ei_cache_friendly_product*() expect each operand to have a valid
>>> > coeffRef() which doesn't make sense for a type-cast matrix.
>>> >
>>> > I wanted to try this myself but I got lost in the template branches that
>>> > specialize and optimize :-(
>>> >
>>> > thanks,
>>> > Manoj
>>
>>
>>
>>
>