Anyone can show me an example of the function dspr in openblas

Question

I am trying to use the function dspr in Openblas with Rcpp. the purpose of dspr is:

A := alpha*x*x**T + A

In my case, firstly I define A as a matrix with all the elements are 0, alpha=1, x=(1,3), so, the final matrix A should be {(1,3),(3,9)}, but I never get the right result, I set the parameters as follow:

cblas_dspr(CblasColMajor,CblasUpper,2, 1, &(x[0]),1, &(A[0]));

Can anyone tell me how to set the right parameters of dspr? Thanks.

score 0 · Answer 1 · answered Dec 28 '13 at 17:40

0

The file /usr/include/cblas,h on my machine shows the following signature for the C interface of the BLAS:

void cblas_dspr(const enum CBLAS_ORDER order, const enum CBLAS_UPLO Uplo,
                const int N, const double alpha, const double *X,
                const int incX, double *Ap);

Try that. You get the beginning of Rcpp vectors via x.begin() or via &(x[0]) as you did.

There is nothing pertinent to Rcpp here though.

answered Dec 28 '13 at 17:40

Dirk Eddelbuettel

360,940
56
644
725

Thanks. I checked netlib for dspr and put the parameters into dspr as mentioned, however, I always get a matrix A finally {(1,9),(3,0)}, so I wonder whether I put the wrong parameters into the function dspr. – Jason_C Dec 28 '13 at 18:33

Markus-Hermann · Answer 2 · 2014-01-31T08:53:02.353

Repeated from your own post: The BLAS dyadic product performs A := alpha*x*x' + A

So A would need to be initialized with zero values. In addition do not forget that A is an upper triangular matrix.

For further reading I recommend these links:

https://developer.apple.com/library/ios/documentation/Accelerate/Reference/BLAS_Ref/Reference/reference.html

https://github.com/foundintranslation/Kaldi/blob/master/tools/ATLAS/interfaces/blas/C/src/cblas_dspr.c

However, you wanted an example. Here goes:

/** dspr_demo.cpp */

#include <cblas.h>
#include <iostream>
#include <cstdlib>

using namespace std;

int main(int argc, char** argv)
{
  int n=2;
  double* x = (double*)malloc(n*sizeof(double));
  double* upperTriangleResult = (double*)malloc(n*(n+1)*sizeof(double)/2);

  for (int j=0;j<n*(n+1)/2;j++) upperTriangleResult[j] = 0;
  x[0] = 1; x[1] = 3;

  cblas_dspr(CblasRowMajor,CblasUpper,n,1,x,1,upperTriangleResult);
  double*& A = upperTriangleResult;
  cout << A[0] << "\t" << A[1] << endl << "*\t" << A[2] << endl;

  free(upperTriangleResult); free(x);
  return EXIT_SUCCESS;
}

Anyone can show me an example of the function dspr in openblas

2 Answers2