## triangularize a matrix in one go with drotg & dlasr

Open discussion regarding features, bugs, issues, vendors, etc.

### triangularize a matrix in one go with drotg & dlasr

Hello,

I have been trying to use Givens rotations via a combination of drotg and dlasr but have not been able to get it right. Can anyone please help me out. Please find below the self contained code of what I am trying to acomplish.

Best regards,
Giovanni

Code: Select all
`/* * test.cc * *  Created on: May 31, 2012 *      Author: Giovanni Azua */#include <assert.h>#include <math.h>#include <stdio.h>#include <stdlib.h>#include <mkl_lapack.h>#include <mkl_lapacke.h>#include <mkl_cblas.h>static const int ROWS = 3;static const int COLS = 3;static double matrix[ROWS][COLS];static void print_matrix() {   int rows = ROWS;   int cols = COLS;   printf("matrix of dimensions %dx%d: \n", rows, cols);   for (int i=0; i < rows; ++i) {      for (int j=0; j < cols; ++j) {         printf("%f ", matrix[j][i]);      }      printf("\n");   }}int main(int argc, char** argv) {   // create a simple matrix in "column major" with 1's in the   // upper triangule and 2s in the band below the diagonal   // now I want to get rid of the 2s using Givens rotations   matrix = 1; matrix = 1; matrix = 1;   matrix = 2; matrix = 1; matrix = 1;   matrix = 0; matrix = 2; matrix = 1;   print_matrix();   double c[COLS];   double s[COLS];   // create the rotations   for (int j = 0; j < (COLS - 1); j++) {      double a = matrix[j][j];      double b = matrix[j + 1][j];      cblas_drotg(&a, &b, &c[j], &s[j]);   }   // apply the rotations along the diagonal   char side      = 'L';   char pivot     = 'V';   char direct    = 'F';   lapack_int m   = ROWS;   lapack_int n   = COLS;   lapack_int lda = ROWS;   dlasr(&side, &pivot, &direct, &m, &n, c, s, (double*) matrix, &lda);   print_matrix();   return 0;}`
bravegag

Posts: 12
Joined: Sat Mar 24, 2012 3:23 pm

### Re: triangularize a matrix in one go with drotg & dlasr

Hello,

I glanced at the code quickly. I assume that you take a Hessenberg matrix and then you want to triangularize it through Givens rotations.
If this is the case, the rotation generation and the rotation application needs to be interlaced.
You need to create the first rotation, then apply it, then create the second, then apply it, etc.
Because the first rotation affects the second and the third row and so you need to apply it before creating the second rotation which needs elements in the third row.
Not sure this helps though. As I said as just glanced at the code.

Julien.
Julien Langou

Posts: 835
Joined: Thu Dec 09, 2004 12:32 pm
Location: Denver, CO, USA

### Re: triangularize a matrix in one go with drotg & dlasr

Hi Julien,

Thank you. Indeed, it is a upper Hessenberg matrix. I was actually aiming to do the triangulation in a "BLAS level 2" fashion that's why I was trying to use dlasr. Do you have example code for using this function dlasr?

From the top of your head can you point me to one example of this use-case of trying to triangularize in "one go" and using Given rotations a matrix in upper Hessenberg form i.e. when one column in the middle has been deleted?

Many TIA,
Best regards,
Giovanni
bravegag

Posts: 12
Joined: Sat Mar 24, 2012 3:23 pm

### Re: triangularize a matrix in one go with drotg & dlasr

Never mind, figured it out already. Please see below the working example self-contained C++ code:

Code: Select all
`/* * test.cc * *  Created on: May 31, 2012 *      Author: Giovanni Azua */#include <assert.h>#include <math.h>#include <stdio.h>#include <stdlib.h>#include <mkl_lapack.h>#include <mkl_lapacke.h>#include <mkl_cblas.h>static const int ROWS = 3;static const int COLS = 3;static double matrix[ROWS][COLS];static void print_matrix() {   int rows = ROWS;   int cols = COLS;   printf("matrix of dimensions %dx%d: \n", rows, cols);   for (int i=0; i < rows; ++i) {      for (int j=0; j < cols; ++j) {         printf("%f ", matrix[j][i]);      }      printf("\n");   }}int main(int argc, char** argv) {   // create a simple matrix in "column major" with 1's in the   // upper triangule and 2s in the band below the diagonal   // now I want to get rid of the 2s using Givens rotations   matrix = 1; matrix = 1; matrix = 1;   matrix = 2; matrix = 1; matrix = 1;   matrix = 0; matrix = 2; matrix = 1;   print_matrix();   // Givens's sins and cos   double cc[COLS];   double ss[COLS];   // create the rotations   for (int j = 0; j < COLS; j++) {      double a = matrix[j][j];      double b = matrix[j][j + 1];      double c = 0.0;      double s = 0.0;      double r = 0.0;      dlartg(&a, &b, &c, &s, &r);      // could also use cblas_drotg       //cblas_drotg(&a, &b, &c, &s);      for (int jj = j; jj < COLS; ++jj) {         cc[jj] = c;         ss[jj] = s;      }      // apply the rotations to trailing columns on the right      char side      = 'L';      char pivot     = 'V';      char direct    = 'F';      lapack_int m   = ROWS - j;      lapack_int n   = COLS - j;      lapack_int lda = ROWS;      dlasr(&side, &pivot, &direct, &m, &n, &cc[j], &ss[j], ((double*) matrix) + j*ROWS + j, &lda);      print_matrix();   }   return 0;}`
bravegag

Posts: 12
Joined: Sat Mar 24, 2012 3:23 pm