Sorry about my ignorance for the following question! It is probably one you get all the time, but I couldn't figure it out.

When writing a template function so that it might work either on a matrix or matrix expression, I understand that the type for the argument to the function should not be a Matrix, but rather MatrixBase.

I would like if the function returned an expression, rather than a matrix, so that the evaluation is lazy (and might possibly be rewritten later). However, the compiler yelled at me when I tried to make the return type be MatrixBase.

My code compiles fine if the return type is Matrix, but I am concerned that this forces an evaluation and creates a temporary product. Is this the case? Or is it accepting the matrix expression as type matrix?

Here is a snippet of my code:

Code: Select all

typedef std::complex<double> Complex;
typedef Eigen::Matrix<Complex, Eigen::Dynamic, Eigen::Dynamic> fieldType;
...
template<typename Derived>
  inline const fieldType Btimes(const Eigen::MatrixBase<Derived>& field) const {
    return field*r2plusa22 - sinSquaredTheta*field*a2Delta ;
  };
...
  inline const fieldType dissipationTerms() const {
    return
      BPhidot() * twosMr2minusa2minusDeltar
      + Btimes( dByPhi(Phi) * twosarminusM );
  };


In the second function, I am passing a GeneralProduct to Btimes(). This way compiles just fine. I thought that the return types of dissipationTerms() and Btimes() would have been something like Eigen::MatrixBase< fieldType >, but the compiler complains about not being able to convert between types when I do that.

How do I do this right, to avoid evaluation and creation of a temporary?

Also: when would one return a reference to a matrix or matrix expression, as compared to just returning an object? I'm not using the results of these functions as lvalues, but I am still not clear on when you whether I should be returning a reference (yes, I know, I'm a newbie, sorry!).

Cheers
Leo

Those are not newbie questions at all but concern rather advanced uses of the library. To answer you I need more information: which Eigen version ? average size of the matrices ? the type of each variable (r2plusa22? sinSquaredTheta?, etc.)

Thanks for taking your time to answer my questions!

This is part of a spectral PDE solver I'm working on. Things of type fieldType are values on a one dimensional grid, and each grid point has some number of angular harmonics. The number of harmonics of interest is probably not very large, something like 10 or 30. The number of radial grid points is larger, probably something like 3000. So field is something like (10 x 3000).
Acting on the left are angular operators like sinSquaredTheta, which is sparse in the angular decomposition; this is a SparseMatrix of size (10 x 10), and is pentadiagonal or tridiagonal.
Acting on the right are radial operators like multiplying by r or any other purely radial function. These are stored as vectors, so typical length 3000. To act as an operator it should look like a diagonal matrix -- I dropped this in my code snippet. It should really look like r2plusa22.asDiagonal() (this is only a compile time penalty, right?)

I'm using the devel version of eigen.