0

The code line: gsl_blas_daxpy(-a,&gsl_matrix_column(D, q).vector,y); cause the error

error C2102: '&' requires l-value

, now the problem is that I have no control of the GSL functions so I don't know how to figure this out (removing the "&" didn't work)

afterwards i get

error C2198: 'gsl_blas_daxpy' : too few arguments for call

I'm using Visual studio 2010.

GSL_EXPORT int  gsl_blas_daxpy (double alpha,
                                const gsl_vector * X,
                                gsl_vector * Y);


#include <stdio.h>
#include <math.h>
#include <time.h>
#include <gsl/gsl_vector.h>
#include <gsl/gsl_matrix.h>
#include <gsl/gsl_blas.h>

#define M (10) // Number of columns in dictionary */
#define N ((int)(M/2)) // Number of rows in dictionary */
int K = 0.07*M; //Number of non-zero elements in signal - the sparsity 
int P=1; //number of signals
double epsilon = 1.0e-7; // Residual error 
int numOfIterations = N; /* Max num of iterations - same as num of elements in signal */

double sign(double x){return (x>=0) - (x<0);} // Sign function 

int main(int argc, char** argv)
{
int n, m, k, iter, q;
double normi, normf, tmp , norm=sqrt(N), htime;
gsl_matrix *D;  // A random dictionary used for encoding the sparse signal  NxM
gsl_vector *x;  // Sparse info signal (encoder input) MxP
gsl_vector *z;  // Evaluated Sparse info signal (decoder output)  MxP
gsl_vector *r;  // Residual error vector MxP
gsl_vector *y;  // Sparse representation of signal (encoder output) NxP
gsl_vector_view v;
clock_t start; //for measuring performance

printf("\nDictionary is:NxM=%dx%d,and the signal sparsity is K=%d", N, M, K);


srand(time(NULL)); //Initialize srand
start =clock(); //Initialize  clock

/* Initiallize D as a Bernoulli random dictionary */
D = gsl_matrix_alloc (N, M);
for(m=0; m<M; m++)
    {
    for(n=0; n<N; n++)
     {
        tmp=sign(2.0*rand()/(double)RAND_MAX-1.0)/norm;
        gsl_matrix_set (D, n, m, tmp);    //D[n,m]=tmp
     }
   }

/* Create a random K-sparse info signal */
x = gsl_vector_alloc(M);
for(k=0; k<K; k++)
 {
    gsl_vector_set(x, rand()%M, 2.0*rand()/(float)RAND_MAX - 1.0); //put random values at k random positions
 }

/* Allocate memory for solution (evaluated signal) */
z = gsl_vector_calloc(M); 

/* Allocate memory for residual vector */
r = gsl_vector_calloc(M);

/* Allocate memory for the encoded signal vector (its representation) */
y = gsl_vector_alloc(N);

htime=((double)clock()-start)/CLOCKS_PER_SEC;
printf("\nTime data allocation: %f", htime);

/* Encoding the signal (x to y) */
start = clock();
gsl_blas_dgemv(CblasNoTrans, 1, D, x, 0, y); // y = Dx
htime=((double)clock()-start)/CLOCKS_PER_SEC;
printf("\nTime for encoding: %f", htime);

/* Decoding the signal */
start = clock();
normi = gsl_blas_dnrm2(y); // ||y|| (L2 norm)
epsilon = sqrt(epsilon * normi); 
normf = normi;
iter = 0;

/*iterate till the computational error is small enough*/
while(normf > epsilon && iter < numOfIterations)
    {
        gsl_blas_dgemv(CblasTrans, 1, D, y, 0, r); // r=D'*y
        q = gsl_blas_idamax(r); //index of max element in residual vector
        tmp = gsl_vector_get(r, q); //the  max element in r
        gsl_vector_set(z, q, gsl_vector_get(z, q)+tmp); // z[q]=z[q]+ tmp
        v=gsl_matrix_column(D, q); // choose the dictrionary's atom (coloum) with the index of largest element in r
        gsl_blas_daxpy(-tmp,&v.vector,y); // y = y-tmp*v
        normf = gsl_blas_dnrm2(y); // ||y|| (L2 norm)
        iter++;
    }

htime = ((double)clock()-start)/CLOCKS_PER_SEC;
printf("\nTime for decoding: %f", htime);
tmp = 100.0*(normf*normf)/(normi*normi); // the error at end of algorithm
printf("\nComputation residual error: %f",tmp);

/* Check the solution (evaluated signal) against the original signal */
printf("\nSolution (first column),Reference (second column):");
getchar(); // wait for pressing a key
for(m=0; m<M; m++)
    {   
        printf("\n%.3f\t%.3f", gsl_vector_get(x, m),gsl_vector_get(z, m));
    }

normi = gsl_blas_dnrm2(x);
gsl_blas_daxpy(-1.0, x, z); // z = z-x
normf = gsl_blas_dnrm2(z); // ||z|| (L2 norm)
tmp = 100.0*(normf*normf)/(normi*normi); //final error
printf("\nSolution residual error: %f\n",tmp);

/* Memory clean up and shutdown*/
gsl_vector_free(y); gsl_vector_free(r);
gsl_vector_free(z); gsl_vector_free(x);
gsl_matrix_free(D);
getchar();
return EXIT_SUCCESS;
}
Idan Banani
  • 131
  • 2
  • 13
  • 7
    `#define N ((int)(M/2))` No. Please, **NO**. Also, are you writing in C and C++ at the same time? – Bartek Banachewicz Mar 04 '13 at 14:44
  • How is `gsl_blas_daxpy` declared? What are the argument types? What is `gsl_matrix_column(D, q).vector`? – Some programmer dude Mar 04 '13 at 14:45
  • 2
    Your compiler tells you exactly what is the problem. `gsl_matrix_column(D, q).vector` isn't an addressable slot in memory. – SomeWittyUsername Mar 04 '13 at 14:46
  • This is pure c. The tag c++ is an error – qPCR4vir Mar 04 '13 at 15:08
  • 1
    If your original problem is solved, then accept an answer (click the check mark next to the answer which you feel answered the question best). If you solved it yourself, make an answer, and accept it. If you have another question, make another post. Don't tack on additional questions to a solved post. – Benjamin Lindley Mar 05 '13 at 00:08

4 Answers4

4

gsl_matrix_column(D, q).vector is an R-value. You can't take its address. You need an L-value, so assign it to a named variable first, then pass the address of that variable to the function.

Benjamin Lindley
  • 101,917
  • 9
  • 204
  • 274
1
    gsl_vector_view c=gsl_matrix_column(D, q);
    gsl_blas_daxpy(-a,&c.vector,y);

I think, introducing a temporal variable led you pass a pointer to it to the function.

EDIT: Well, trying to understand the problem, I wanted to know what the function expect:

int gsl_blas_daxpy (double alpha, const gsl_vector * x, gsl_vector * y)

and

gsl_vector_view gsl_matrix_column (gsl_matrix * m, size_t j)

witj some explanation:

A vector view can be passed to any subroutine which takes a vector argument just as a directly allocated vector would be, using &view.vector.

and an example:

   for (j = 0; j < 10; j++)
     {
       gsl_vector_view column = gsl_matrix_column (m, j);
       double d;

       d = gsl_blas_dnrm2 (&column.vector);

       printf ("matrix column %d, norm = %g\n", j, d);
     }
qPCR4vir
  • 3,521
  • 1
  • 22
  • 32
1

If you make a more permanent home for the return value of gsl_matrix_column, (this particular) problem will go away.

Here is some simplified code that illustrates how one might capture a return value in an addressable slot:

struct _foo {
  int i;
};

struct _foo bar () {
  struct _foo result = { 5 };
  return result;
}

/* won't compile; 'lvalue required as unary & operand */
void qux () {
  int *j = &bar().i;
}

/* compiles OK */
void qal () {
  struct _foo result = bar();
  int* j = &result.i;
}
gcbenison
  • 11,723
  • 4
  • 44
  • 82
1

Now we have another problem:

Here another answer:

Are you aware that int K= 0.7 is K=0 ??

#define M (10) // Number of columns in dictionary */
int K = 0.07*M; //Number of non-zero elements in signal - the sparsity 

alloc do not initialice the vector x. x will contain garbage values, not 0. Did you meant x = gsl_vector_calloc(M); with c? It will set x to 0.

/* Create a random K-sparse info signal */
x = gsl_vector_alloc(M);
for(k=0; k<K; k++)      // K=0, for get skiped and x not modified.
 {
    gsl_vector_set(x, rand()%M, 2.0*rand()/(float)RAND_MAX - 1.0); //put random values at k random positions
 }

(And here you will have at most K random values, but possible lest)

qPCR4vir
  • 3,521
  • 1
  • 22
  • 32
  • What is your advice for figuring out the casting/rounding problem? I can't pass the build error... #define M (10) #define N ((int)(M/2)) double sparsity= 0.07*M; int K = (int)floor(sparsity); – Idan Banani Mar 05 '13 at 08:22
  • @IdanBanani You have NO PROBLEM with casting/rounding. It is simple that 7% of 10 is very small. Use a bigger M (100?) or more non-zero elements 30% ?? M=10, K=0.3*M=3 for example. – qPCR4vir Mar 05 '13 at 08:58