ARM-software / CMSIS-DSP

CMSIS-DSP embedded compute library for Cortex-M and Cortex-A
https://arm-software.github.io/CMSIS-DSP
Apache License 2.0
454 stars 122 forks source link

QR Decomposition seems to coverge very quickly. How to adjust parameters. #168

Closed DominikRzecki closed 2 months ago

DominikRzecki commented 3 months ago

Hello,

I have a question regarding the implementation of the qr algorithm.

I am trying to compute the eigenvalues and eigenvectors of a matrix.

for (size_t i = 0; i < N; ++i) {
            arm_mat_qr(tmpD, Q, R);
            tmpD = R * Q;
            tmpV = tmpV * Q;
}

where arm_mat_qr is:

template<class T1, size_t Rows1, size_t Columns1>
void arm_mat_qr(matrix& A, matrix& Q, matrix& R) {
      T1 threshold = DEFAULT_HOUSEHOLDER_THRESHOLD_F32;
      T1 tau[Columns1] {0};
      T1 tmpA[Rows1] {0};
      T1 tmpB[Columns1] {0};

      arm_mat_qr_f32(&A.instance_, threshold, &R.instance_, &Q.instance_, tau, tmpA, tmpB);

     // Zeroing lower triangle
      for(size_t i = 0; i < 8; ++i) {
          for(size_t j = 0; j < i; ++j) {
              R.at(i, j) = 0;
          }
      }
}

Matlab:

    D = A;
    V = eye(size(A));
    for i = 0:30
       [Q, R] = qr(D);
       D = R * Q;
       V = V * Q;
    end

Matlab implementation provides me with the correct eigenvalues and eigenvectors. In CMSIS during convergence the eigenvalues on the diagonal of R seem to miss the correct eigenvalues and converge towards the Identity Matrix rapidly during the first 2 iterations. Eigenvectors also disappear quickly, but when they are still there only the first column matches with Matlab.

I think that the algorithm with my current parameters is to aggressive. Do I maybe have to adjust the threshold, TmpA and TmpB? If so, then how?

With regards Dominik Rzecki

christophe0606 commented 3 months ago

@DominikRzecki Hi, Can you share the matrix so that I can try on my side ?

christophe0606 commented 2 months ago

Without additional information nothing can be done. So, I close it.

DominikRzecki commented 2 weeks ago

Since I haven't responded on time:

The arm_mat_qr_f32 function works flawlessly indeed.

It seems that gcc with optimizations enabled has reordered my code in a way that breaks it's logic, because some class didn't have a properly implemented constructor. It even made the qr algorithm behave differently.

christophe0606 commented 2 weeks ago

@DominikRzecki Glad to know. Thank you.