Eigen::unaryExpr and const for random number engine

I have the following class:

class TestRandom
{

    public:
        struct scalar_op
          {
            std::mt19937 rng;                  
            std::normal_distribution<double> norm;  
            scalar_op(){}
            scalar_op(unsigned long seed){rng = std::mt19937(seed);}
            void reset(unsigned long seed=std::mt19937::default_seed){ rng142.seed(seed);norm.reset(); }
            double operator()(void) const
            { 
                    double n_noise =norm(rng);
                    return n_noise;
            }
          };

        TestRandom(){}
        void reset(unsigned long seed=std::mt19937::default_seed){rand.reset(seed);}

        /// Draw nn samples from the gaussian and return them
        Eigen::MatrixXd samples(int n)
        {    
             Eigen::MatrixXd noise(3,n);
             noise.unaryExpr(rand);
             return corr_noise;
        }
    private:
        scalar_op rand; 
}; //TestRandom  

What i am trying to achieve is first reset the random number generator as

TestRandom r1;
 r1.reset(1)

Then fill up a matrix with random samples as

  Eigen::MatrixXd normal_random = r1.samples(1); // just a sample of random numbers

But I get a compilation error of

 error C3848: expression having type 'const std::mt19937' would lose some const-volatile qualifiers in order to call 'unsigned int std::mersenne_twister<_Ty,_Wx,_Nx,_Mx,_Rx,_Px,_Ux,_Sx,_Bx,_Tx,_Cx,_Lx>::operator ()(void)

If I remove the const of the () operator, code compiles but the computation inside operator() is not executed (I can see this in the debugger).

It would be very kind of someone to point out my mistake and possibly provide way to resolve this. Thanks in advance