Loading ...
Sorry, an error occurred while loading the content.
 

Re: [ublas-dev] inner_prod() for std::complex working?

Expand Messages
  • Kresimir Fresl
    ... [...] ... Inner product is in uBLAS defined as inner_prod (v1, v2) = sum (v1 [i] * v2 [i]), that is, as v1^T v2 both for real and complex vectors -- it is
    Message 1 of 5 , Feb 1, 2003
      pobdrah wrote:

      > I've just joined this group. I am testing uBLAS and doing some speed
      > benchmark. In doing that, I found that inner_prod() works ok for
      > <float> but doesn't compile for std::complex<float>.
      [...]
      > // testing inner_prod() for <std::complex<float> >
      > typedef vector<std::complex<float> >::value_type value_type;
      > typedef vector<std::complex<float> >::size_type size_type;
      > typedef type_traits<value_type>::real_type real_type;
      > vector<std::complex<float> > v1(1), v2(1), v3(1);
      > real_type f = inner_prod(v1, v2); // type conversion error!
      >
      > Am I missing anything here, or it is a bug?

      Inner product is in uBLAS defined as

      inner_prod (v1, v2) = sum (v1 [i] * v2 [i]),

      that is, as

      v1^T v2

      both for real and complex vectors -- it is
      the equivalent of `[cz]dotu' from BLAS.

      The result type of the inner_prod() for two complex
      vectors is a complex number, not real; e.g.:

      v1 = { (1,-1), (1,-1), (1,-1), (1,-1), (1,-1), (1,-1) }
      v2 = { (0,1), (0,1), (0,1), (0,1), (0,1), (0,1) }

      inner_prod (v1, v2) = (-6, -6)
      inner_prod (v1, v1) = (0, -12)
      inner_prod (v2, v2) = (-6, 0)

      As the last example shows, imaginary part can sometimes
      be zero, but in general it is not.

      If you need

      v1^H v2

      (i.e. `usual' definition of inner product for complex vectors,
      i.e. equivalent of `[cz]dotc' from BLAS), you can write

      inner_prod (conj (v1), v2)

      e.g.:

      inner_prod (conj (v1), v2) = (6,-6)
      inner_prod (conj (v1), v1) = (12,0)
      inner_prod (conj (v2), v2) = (6,0)

      but the result is again complex number.

      If you calculate v^H v, imaginary part is always zero
      (second and third example), so you can say that the result
      is real -- but this is special case.

      Vector norms are real numbers both for real and complex
      vectors; e.g. ||v||_2, that is, ublas::norm_2 (v), is defined as

      sqrt (sum (v [i] * conj (v [i])))

      i.e.

      sqrt (inner_prod (v, conj (v))

      Sincerely,

      fres


      PS. Joerg, there is a small and completely unimportant typo
      in `vector_expression.hpp', line 1855: last `)' is missing in the
      comment:

      // inner_prod (v1, v2) = sum (v1 [i] * v2 [i]
    • jhr.walter@t-online.de
      ... [snip the excellent explanation] ... Fixed. How did COMO detect this one? ;-) Thanks, Joerg P.S.:I ve just uploaded my current version containing mostly
      Message 2 of 5 , Feb 1, 2003
        > pobdrah wrote:
        >
        > > I've just joined this group. I am testing uBLAS and doing some speed
        > > benchmark. In doing that, I found that inner_prod() works ok for
        > > <float> but doesn't compile for std::complex<float>.
        > [...]
        > > // testing inner_prod() for <std::complex<float> >
        > > typedef vector<std::complex<float> >::value_type value_type;
        > > typedef vector<std::complex<float> >::size_type size_type;
        > > typedef type_traits<value_type>::real_type real_type;
        > > vector<std::complex<float> > v1(1), v2(1), v3(1);
        > > real_type f = inner_prod(v1, v2); // type conversion error!
        > >
        > > Am I missing anything here, or it is a bug?

        [snip the excellent explanation]

        > PS. Joerg, there is a small and completely unimportant typo
        > in `vector_expression.hpp', line 1855: last `)' is missing in the
        > comment:
        >
        > // inner_prod (v1, v2) = sum (v1 [i] * v2 [i]

        Fixed. How did COMO detect this one? ;-)

        Thanks,

        Joerg

        P.S.:I've just uploaded my current version containing mostly documentation
        updates. Due to some bug reports regarding the examples in the documentation
        I've now compiled every example and hopefully fixed the corresponding bugs
        (BTW should I add the compilable examples to the distribution and to boost
        CVS?). I've also changed the tool for editing the HTML documentation: I'm
        using Mozilla Composer instead of MS Frontpage now.
      • Kresimir Fresl
        Hi Joerg, ... Compiler magick ;o) ... with a little help from a human friend. ... I think yes -- `examples in `test? and `bench? directories are not very
        Message 3 of 5 , Feb 2, 2003
          Hi Joerg,

          >>PS. Joerg, there is a small and completely unimportant typo
          >>in `vector_expression.hpp', line 1855: last `)' is missing in the
          >>comment:
          >>
          >>// inner_prod (v1, v2) = sum (v1 [i] * v2 [i]

          > Fixed. How did COMO detect this one? ;-)

          Compiler magick ;o) ... with a little help from a human friend.

          > P.S.:I've just uploaded my current version containing mostly documentation
          > updates. Due to some bug reports regarding the examples in the documentation
          > I've now compiled every example and hopefully fixed the corresponding bugs
          > (BTW should I add the compilable examples to the distribution and to boost
          > CVS?).

          I think yes -- `examples' in `test?' and `bench?' directories are
          not very useful as introductory examples.

          BTW, IMHO examples for hermitian matrix/adaptor are wrong
          -- if A=A^H, then diagonal elements must have imaginary part
          equal to zero [a + bi == a - bi]. Initialization loop can be:

          for (int i = 0; i < ha.size1 (); ++ i)
          for (int j = 0; j <= i; ++ j)
          hal (i, j) = std::complex<double> (3 * i + j, 3 * i - 3 * j);

          Regards,

          fres
        • jhr.walter@t-online.de
          Hi Kresimir, you wrote: [snip] ... documentation ... documentation ... bugs ... boost ... Done, new download available. ... Corrected. Thanks, Joerg
          Message 4 of 5 , Feb 2, 2003
            Hi Kresimir,

            you wrote:

            [snip]

            > > P.S.:I've just uploaded my current version containing mostly
            documentation
            > > updates. Due to some bug reports regarding the examples in the
            documentation
            > > I've now compiled every example and hopefully fixed the corresponding
            bugs
            > > (BTW should I add the compilable examples to the distribution and to
            boost
            > > CVS?).
            >
            > I think yes -- `examples' in `test?' and `bench?' directories are
            > not very useful as introductory examples.

            Done, new download available.

            > BTW, IMHO examples for hermitian matrix/adaptor are wrong
            > -- if A=A^H, then diagonal elements must have imaginary part
            > equal to zero [a + bi == a - bi]. Initialization loop can be:
            >
            > for (int i = 0; i < ha.size1 (); ++ i)
            > for (int j = 0; j <= i; ++ j)
            > hal (i, j) = std::complex<double> (3 * i + j, 3 * i - 3 * j);

            Corrected.

            Thanks,

            Joerg
          Your message has been successfully submitted and would be delivered to recipients shortly.