Re: [eigen] Custom complex type

[ Thread Index | Date Index | More lists.tuxfamily.org/eigen Archives ]


Hi,

I don't really understand why you cannot specialize std::complex for
mpreal, e.g:

namespace std {
template<> class complex<mpreal> {
// copy paste here the content of your current mpcomplex
};
}

gael

On Tue, Jun 21, 2011 at 4:06 AM, Pavel Holoborodko
<pavel@xxxxxxxxxxxxxxx> wrote:
> I see, thank you for your prompt response.
> Yes, std::complex<mpreal> specialization will be much easier than
> modification of Eigen for custom complex type.
> I wanted mpcomplex to be based on low level C library, which uses POD types
> for real and imag parts (C pointers actually).
> It would be the most efficient to use them directly without wrapping to C++
> numeric class.
> Since it is impossible to specialize std::complex with POD - I had to write
> my own mpcomplex.
> Now, as Eigen heavily based on std::complex I have no other way than to
> make std::complex<mpreal> class specialization.
> I think it is acceptable and maybe even better solution in a long run (e.g.
> compatibility with other libraries).
>
> Thank you for your help.
> Pavel Holoborodko.
> On Tue, Jun 21, 2011 at 4:29 AM, Gael Guennebaud <gael.guennebaud@xxxxxxxxx>
> wrote:
>>
>> ok, got it. The pb is in gebp_traits
>> (src/Core/products/GeneralBlockPanelKernel.h) where I have
>> specializations on std::complex to identify complexes.
>>
>> Just a remark: cannot you simply specialize std::complex<mpreal> and
>> implement your optimization there ?
>>
>> gael
>>
>> On Mon, Jun 20, 2011 at 9:21 PM, Gael Guennebaud
>> <gael.guennebaud@xxxxxxxxx> wrote:
>> > On Mon, Jun 20, 2011 at 11:38 AM, Pavel Holoborodko
>> > <pavel@xxxxxxxxxxxxxxx> wrote:
>> >> In order to compile my Eigen-based code I had to add definitions for
>> >> various
>> >> kinds of conj_helper<>:
>> >> template<> struct conj_helper<mpc::mpcomplex, mpc::mpcomplex,...>
>> >> And
>> >> template<> struct scalar_product_traits<mpc::mpcomplex, mpc::mpcomplex>
>> >
>> > hm, I was afraid of that. This part of Eigen's internals requires some
>> > love. In theory, you should not have to do such specializations.
>> >
>> >> Code compiled well after that, but after my tests I noticed that in
>> >> expression (pseudo inverse calculation in SVD):
>> >> m_matrixV * invertedSingVals.asDiagonal() * m_matrixU.adjoint();
>> >> m_matrixU.adjoint() works just as plain transpose() despite that
>> >> m_matrixU
>> >> is complex.  invertedSingVals is of mpreal type.
>> >> If I use m_matrixU.adjoint().eval() instead everything is all right (as
>> >> well
>> >> as adjointInPlace()) .
>> >> I guess some mpcomplex-specialization is missing for proper matrix
>> >> expression template evaluation.
>> >> Could you advise what it could be?
>> >
>> > I really have to clean up all this mess with complexes!  Could you
>> > check whether m_matrixV * m_matrixU.adjoint() is equivalent to
>> > m_matrixV * m_matrixU.transpose()?
>> >
>> > gael
>> >
>> >> NumTraits<mpcomplex> is in attach. mpcomplex class as simple as:
>> >> class mpcomplex{
>> >> private:
>> >> mpreal re;
>> >> mpreal im;
>> >> // Standard API
>> >> .....
>> >> Thanks in advance.
>> >> Pavel Holoborodko.
>> >> P.S.
>> >> std::complex<mpreal> works nicely without any shaman dances.
>> >> On Wed, Jun 8, 2011 at 12:13 PM, Pavel Holoborodko
>> >> <pavel@xxxxxxxxxxxxxxx>
>> >> wrote:
>> >>>
>> >>> Hi Gael,
>> >>> Thanks for answers everyone.
>> >>> Gael, I'll try this way (custom mpcomplex + std comp. API) with hope
>> >>> that
>> >>> it will work with Eigen "in practice" too.
>> >>> I'll let you know about the results.
>> >>> On Wed, Jun 8, 2011 at 3:43 AM, Gael Guennebaud
>> >>> <gael.guennebaud@xxxxxxxxx> wrote:
>> >>>>
>> >>>> Hi Pavel,
>> >>>>
>> >>>> in theory you should be able to write your own mpcomplex class, and
>> >>>> define NumTraits<mpcomplex>::IsComplex to true as long as you
>> >>>> implement a STD compatible API.
>> >>>>
>> >>>> gael
>> >>>>
>> >>>> On Tue, Jun 7, 2011 at 8:04 AM, Pavel Holoborodko
>> >>>> <pavel@xxxxxxxxxxxxxxx>
>> >>>> wrote:
>> >>>> > Hello,
>> >>>> > Eigen's magicians, please advise on how to implement custom complex
>> >>>> > type for
>> >>>> > smooth integration with Eigen.
>> >>>> > I plan to create multi-precision C++ complex scalar class based on
>> >>>> > MPC: http://www.multiprecision.org/index.php?prog=mpc.
>> >>>> > MPC is cousin of MPFR from the same authors - http://mpfr.org.
>> >>>> > There are two possible ways to proceed:
>> >>>> > (1).  Add template specialization for std::complex<mpfr::mpreal>.
>> >>>> > (2).  Create own scalar class something like mpcomplex based on
>> >>>> > low-level
>> >>>> > MPC.
>> >>>> > Second option seems attractive from optimization point of view -
>> >>>> > minimum
>> >>>> > overhead, no conversion mpreal <-> MPC low level stuff is required.
>> >>>> > There is only one problem - will it be possible to integrate such
>> >>>> > complex
>> >>>> > type with Eigen?
>> >>>> > Documentation on NumTraits says:
>> >>>> > "An enum value IsComplex. It is equal to 1 if T is
>> >>>> > a std::complex type,
>> >>>> > and
>> >>>> > to 0 otherwise."
>> >>>> > Does it mean that only std::complex<T> custom complex types are
>> >>>> > supported in
>> >>>> > Eigen?
>> >>>> > If yes, do I need to provide something besides
>> >>>> > NumTraits<std::complex<mpreal> > specialization?
>> >>>> > Thanks in advance,
>> >>>> > Pavel.
>> >>>> >
>> >>>> >
>> >>>> >
>> >>>> >
>> >>>>
>> >>>>
>> >>>
>> >>
>> >>
>> >
>>
>>
>
>



Mail converted by MHonArc 2.6.19+ http://listengine.tuxfamily.org/