triangularize a matrix in one go i.e. drotg & dlasr

triangularize a matrix in one go i.e. 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.

Many thanks in advance,
Best regards,
Giovanni

/*

 * test.cc

 *

 *  Created on: May 31, 2012

 *      Author: Giovanni Azua

 */

#include

#include

#include

#include 
#include

#include

#include 
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[0][0] = 1; matrix[1][0] = 1; matrix[2][0] = 1;

    matrix[0][1] = 2; matrix[1][1] = 1; matrix[2][1] = 1;

    matrix[0][2] = 0; matrix[1][2] = 2; matrix[2][2] = 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;

}

2 posts / 0 new
Last post
For more complete information about compiler optimizations, see our Optimization Notice.

Oh, never mind, figured it out. Kind of hate that I need to copy over the c and s values at line 53.

/*

 * rotations_sandbox.cc

 *

 *  Created on: May 31, 2012

 *      Author: Giovanni Azua

 */

#include

#include

#include

#include 
#include

#include

#include 
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[0][0] = 1; matrix[1][0] = 1; matrix[2][0] = 1;

	matrix[0][1] = 2; matrix[1][1] = 1; matrix[2][1] = 1;

	matrix[0][2] = 0; matrix[1][2] = 2; matrix[2][2] = 1;

	print_matrix();
	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);
		for (int jj = j; jj < COLS; ++jj) {

			cc[jj] = c;

			ss[jj] = s;

		}
		// apply the rotations along the diagonal

		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;

}

Leave a Comment

Please sign in to add a comment. Not a member? Join today