more missing implicits

[jimka2001]

I'm trying to use the eig function to diagonalize an nxn matrix, and verify that it worked. I.e., I want to compute V and inv(V) where V is the matrix of eigenvectors. And multiply V, D, inv(V) to verify that I get the original matrix (or something very close thereto).

The code I'm trying to compile is given below.

There are two problems. 1. I don't know how to make the build-in norm function work, so I've written my own maxabsnorm, which seems to work, although it is ugly. and 2, there seem to be some missing implicits which I don't see how to create or include.

Is there already a function which will give me the classical eigendecomposition? Or do I really need to write this myself?

/Users/jnewton/Repos/scalain_e/scalain-e-course-code/src/test/scala/TranscendentalsTestSuite.scala:526:19
could not find implicit value for parameter op: breeze.linalg.operators.OpMulMatrix.Impl2[breeze.linalg.DenseMatrix[Double],breeze.linalg.DenseMatrix[breeze.math.Complex],That]
      val dm2 = v * d * inv(v)

Here is the code.

  test("breeze verify eig"){
    // verify the behavior of eigen decomposition
    import breeze.linalg._
    import breeze.math._
    import scala.util.Random
    def maxabsnorm[T](dim:Int, m:DenseMatrix[T],abs:T=>Double):Double = {
      (0 until dim).foldLeft(abs(m(0,0))){(acc:Double,row:Int) =>
        val x = (0 until dim).foldLeft(abs(m(0,0))){
          (acc:Double,col:Int) => math.max(acc,abs(m(row,col)))}
        math.max(acc,x)
      }
    }
    for{dim <- 1 to 4} {
      val dm1: DenseMatrix[Double] = DenseMatrix.tabulate(dim, dim)((_, _)=>Random.between(-10.0, 10.0))
      val eig.Eig(er,ei,v) = eig(dm1)
      val zero:Complex = Complex(0.0,0.0)
      def dia(k:Int,j:Int):Complex = if (k==j) Complex(er(j),ei(j)) else zero
      def cabs(c:Complex):Double = c.abs
      val d:DenseMatrix[Complex] = DenseMatrix.tabulate(dim,dim)(dia)
      val dm2 = v * d * inv(v)

      assert(maxabsnorm(dim,dm1 - dm2, cabs) < 0.001 )
    }

[darrenjw]

I suspect that you need to import breeze.numerics._ But more generally, you might like the Breeze chapter in my notes on Scala for statistical computing and data science, as they potentially answer some of your questions: https://github.com/darrenjw/scala-course/raw/master/scscala.pdf

[dlwh]

The issue is that d is a DenseMatrix[Complex] and in general Breeze refuses to operate on vectors/matrices of different types. (The error message is telling you that, though it's obviously pretty hard to tell in the noise.) If you convert v to complex it should work.

On Wed, Apr 7, 2021 at 3:20 AM Jim Newton @.***> wrote:

I'm trying to use the eig function to diagonalize an nxn matrix, and verify that it worked. I.e., I want to compute V and inv(V) where V is the matrix of eigenvectors. And multiply V, D, inv(V) to verify that I get the original matrix (or something very close thereto).

The code I'm trying to compile is given below.

There are two problems. 1. I don't know how to make the build-in norm function work, so I've written my own maxabsnorm, which seems to work, although it is ugly. and 2, there seem to be some missing implicits which I don't see how to create or include.

Is there already a function which will give me the classical eigendecomposition? Or do I really need to write this myself?

/Users/jnewton/Repos/scalain_e/scalain-e-course-code/src/test/scala/TranscendentalsTestSuite.scala:526:19 could not find implicit value for parameter op: breeze.linalg.operators.OpMulMatrix.Impl2[breeze.linalg.DenseMatrix[Double],breeze.linalg.DenseMatrix[breeze.math.Complex],That] val dm2 = v d inv(v)

Here is the code.

test("breeze verify eig"){ // verify the behavior of eigen decomposition import breeze.linalg. import breeze.math. import scala.util.Random def maxabsnorm[T](dim:Int, m:DenseMatrix[T],abs:T=>Double):Double = { (0 until dim).foldLeft(abs(m(0,0))){(acc:Double,row:Int) => val x = (0 until dim).foldLeft(abs(m(0,0))){ (acc:Double,col:Int) => math.max(acc,abs(m(row,col)))} math.max(acc,x) } } for{dim <- 1 to 4} { val dm1: DenseMatrix[Double] = DenseMatrix.tabulate(dim, dim)((, )=>Random.between(-10.0, 10.0)) val eig.Eig(er,ei,v) = eig(dm1) val zero:Complex = Complex(0.0,0.0) def dia(k:Int,j:Int):Complex = if (k==j) Complex(er(j),ei(j)) else zero def cabs(c:Complex):Double = c.abs val d:DenseMatrix[Complex] = DenseMatrix.tabulate(dim,dim)(dia) val dm2 = v d inv(v)

  assert(maxabsnorm(dim,dm1 - dm2, cabs) < 0.001 )
}

— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub https://github.com/scalanlp/breeze/issues/804, or unsubscribe https://github.com/notifications/unsubscribe-auth/AAACLIKE7HCMB2ZSIQBQTGTTHQWW7ANCNFSM42QP4BHA .

[jimka2001]

I think that gets me further. Thanks for the advice. I've made sure the matrices in question are declared as Complex and contain Complex entries. But the inv function doesn't seem to work. How can I find the inverse of a Complex matrix?

    for{dim <- 1 to 4} {
      val dm1: DenseMatrix[Double] = DenseMatrix.tabulate(dim, dim)((_, _)=>Random.between(-10.0, 10.0))
      val eig.Eig(er,ei,v) = eig(dm1)
      val vc:DenseMatrix[Complex] = DenseMatrix.tabulate(dim, dim)((j,k)=>Complex(v(j,k),0.0))
      val zero:Complex = Complex(0.0,0.0)
      def dia(k:Int,j:Int):Complex = if (k==j) Complex(er(j),ei(j)) else zero
      def cabs(c:Complex):Double = c.abs
      val d:DenseMatrix[Complex] = DenseMatrix.tabulate(dim,dim)(dia)
      val dm2:DenseMatrix[Complex] = vc * d * inv(vc)

      assert(maxabsnorm(dim,dm1 - dm2, cabs) < 0.001 )
    }

I see the following error

/Users/jnewton/Repos/scalain_e/scalain-e-course-code/src/test/scala/TranscendentalsTestSuite.scala:527:50
could not find implicit value for parameter impl: breeze.linalg.inv.Impl[breeze.linalg.DenseMatrix[breeze.math.Complex],VR]
      val dm2:DenseMatrix[Complex] = vc * d * inv(vc)

[jimka2001]

BTW, I tried to declare dm1 as DenseMatrix[Complex] containing complex numbers each of whose imaginary part is 0.0. In this case the eig function call refuses to compile.

Also if dm1 has type DenseMatrix[Double] and dm2 has type DenseMatrix[Complex] then dm1-dm2 fails to compile.

Scala implicits are a big mystery to me, but it seems like some implicit conversions from Double to Complex may be missing? I don't know if you just need to convert Double to Complex (which is easy), or whether for every type H and for every arity you also need to convert H[...,Double,...] to H[...,Complex,...] which sounds like a real nightmare.

Otherwise working with eigenvalues is going to be really tedious. Right? I noticed some comments in the eig function inline header, that allude to this problem. The comment is:

  // TODO: probably we should just return an eigenValues: DV[Complex] ?

[jimka2001]

I don't know whether breeze is supposed to be able to compute the inverse of a complex matrix. But I wrote a function to do it (perhaps just as a workaround). If M has dimension nxn, we can create a matrix of (2n)x(2n), putting the real parts in quadrant 2 and 4, the imaginary parts in quadrant 1, and the negatives of the imaginary parts in quadrant 3. Then we have a matrix of reals which is invertible iff the original matrix is invertible. So I can call inv on that, and then construct the complex matrix (which is the inverse) by reading the reals from quadrant 2 and the imaginaries from quadrant 1.

    def invc(dim:Int, vc:DenseMatrix[Complex]):DenseMatrix[Complex] = {
      val w:DenseMatrix[Double] = DenseMatrix.tabulate(dim*2, dim*2) {
        case (j,k) if j < dim && k < dim => vc(j,k).real
        case (j,k) if j < dim => vc(j,k-dim).imag
        case (j,k) if k < dim => -1 * vc(j-dim,k).imag
        case (j,k) => vc(j-dim,k-dim).real
      }
      val winv:DenseMatrix[Double]  = inv(w)
      val vinvc:DenseMatrix[Complex] = DenseMatrix.tabulate(dim, dim)((j,k) => Complex(winv(j,k),winv(j,k+dim)))

      vinvc
    }

[jimka2001]

I suspect that you need to import breeze.numerics._ But more generally, you might like the Breeze chapter in my notes on Scala for statistical computing and data science, as they potentially answer some of your questions: https://github.com/darrenjw/scala-course/raw/master/scscala.pdf

That's stupendous! thanks.