Provided by: librsb-dev_1.2.0-rc7-6_amd64 bug

NAME

       librsb - rsb-examples - Examples of usage of librsb.

DESCRIPTION

SYNOPSIS

Detailed Description

       Examples of usage of librsb.

       The following fully working example programs illustrate correct ways of using the library.
       The script displayed here should be sufficient to build them.

       #!/bin/bash
       # Script to build the librsb example programs.

       LIBRSB_CONFIG=${LIBRSB_CONFIG:-librsb-config}

       for s in *.c
       do
               p=${s/.c/}
               rm -f $p
               CFLAGS=`${LIBRSB_CONFIG} --I_opts`
               LDFLAGS=`${LIBRSB_CONFIG} --ldflags --extra_libs`
               CC=`${LIBRSB_CONFIG} --cc`
               cmd="$CC $CFLAGS $s $LDFLAGS -o $p"
               echo $cmd
               $cmd
       done

       if test x"yes" = x"yes" ; then
       # activated if you have built the Fortran modules and installed them in the right path.
       for s in *.F90
       do
               p=${s/.F90/}
               rm -f $p
               CFLAGS=`${LIBRSB_CONFIG} --I_opts`
               LDFLAGS=`${LIBRSB_CONFIG} --ldflags --extra_libs`
               FC=`${LIBRSB_CONFIG} --fc`
               cmd="$FC $CFLAGS $s $LDFLAGS -o $p"
               echo $cmd
               $cmd
       done
       fi

       /*

       Copyright (C) 2008-2015 Michele Martone

       This file is part of librsb.

       librsb is free software; you can redistribute it and/or modify it
       under the terms of the GNU Lesser General Public License as published
       by the Free Software Foundation; either version 3 of the License, or
       (at your option) any later version.

       librsb is distributed in the hope that it will be useful, but WITHOUT
       ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
       FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
       License for more details.

       You should have received a copy of the GNU Lesser General Public
       License along with librsb; see the file COPYING.
       If not, see <http://www.gnu.org/licenses/>.

       */
       /*!
        ingroup rsb-examples
        @file
        @author Michele Martone
        @brief This is a first "hello RSB" example program.

        include hello.c
       */
       #include <rsb.h> /* librsb header to include */
       #include <stdio.h>       /* printf() */

       int main(const int argc, char * const argv[])
       {
               /*!
                 A Hello-RSB program.

                 This program shows how to use the rsb.h interface correctly to:

                 - initialize the library using #rsb_lib_init()
                 - set library options using #rsb_lib_set_opt()
                 - revert such changes
                 - allocate (build) a single sparse matrix in the RSB format
                   using #rsb_mtx_alloc_from_coo_const()
                 - prints information obtained via #rsb_mtx_get_info_str()
                 - multiply the matrix times a vector using #rsb_spmv()
                 - deallocate the matrix using #rsb_mtx_free()
                 - finalize the library using #rsb_lib_exit(RSB_NULL_EXIT_OPTIONS)

                 In this example, we use #RSB_DEFAULT_TYPE as matrix type.
                 This type depends on what was configured at library build time.
                * */
               struct rsb_mtx_t *mtxAp = NULL; /* matrix structure pointer */
               const int bs = RSB_DEFAULT_BLOCKING;
               const int brA = bs, bcA = bs;
               const RSB_DEFAULT_TYPE one = 1;
               rsb_type_t typecode = RSB_NUMERICAL_TYPE_DEFAULT;
               rsb_err_t errval = RSB_ERR_NO_ERROR;
               const rsb_nnz_idx_t nnzA = 4;           /* matrix nonzeroes count */
               const rsb_coo_idx_t nrA = 3;            /* matrix rows count */
               const rsb_coo_idx_t ncA = 3;            /* matrix columns count */
               /* nonzero row indices coordinates: */
               rsb_coo_idx_t IA[] = {0,1,2,2};
               /* nonzero column indices coordinates: */
               rsb_coo_idx_t JA[] = {0,1,2,2};
               RSB_DEFAULT_TYPE VA[] = {11,22,32,1};/* values of nonzeroes */
               RSB_DEFAULT_TYPE X[] = { 0, 0, 0 };     /* X vector's array */
               const RSB_DEFAULT_TYPE B[] = { -1, -2, -5 }; /* B vector's array */
               char ib[200];

               printf("Hello, RSB!0);
               printf("Initializing the library...0);
               if((errval = rsb_lib_init(RSB_NULL_INIT_OPTIONS)) !=
                               RSB_ERR_NO_ERROR)
               {
                       printf("Error initializing the library!0);
                       goto err;
               }
               printf("Correctly initialized the library.0);

               printf("Attempting to set the"
                      " RSB_IO_WANT_EXTRA_VERBOSE_INTERFACE library option.0);
               {
                       rsb_int_t evi=1;
                       /* Setting a single optional library parameter. */
                       errval = rsb_lib_set_opt(
                               RSB_IO_WANT_EXTRA_VERBOSE_INTERFACE, &evi);
                       if(errval != RSB_ERR_NO_ERROR)
                       {
                               char errbuf[256];
                               rsb_strerror_r(errval,&errbuf[0],sizeof(errbuf));
                               printf("Failed setting the"
                               " RSB_IO_WANT_EXTRA_VERBOSE_INTERFACE"
                               " library option (reason string:3s).0,errbuf);
                               if(errval&RSB_ERRS_UNSUPPORTED_FEATURES)
                               {
                                 printf("This error may be safely ignored.0);
                               }
                               else
                               {
                                 printf("Some unexpected error occurred!0);
                                 goto err;
                               }
                       }
                       else
                       {
                               printf("Setting back the "
                                       "RSB_IO_WANT_EXTRA_VERBOSE_INTERFACE"
                                       " library option.0);
                               evi = 0;
                               errval = rsb_lib_set_opt(RSB_IO_WANT_EXTRA_VERBOSE_INTERFACE,
                                               &evi);
                               errval = RSB_ERR_NO_ERROR;
                       }
               }

               mtxAp = rsb_mtx_alloc_from_coo_const(
                       VA,IA,JA,nnzA,typecode,nrA,ncA,brA,bcA,
                       RSB_FLAG_NOFLAGS    /* default format will be chosen */
                       |RSB_FLAG_DUPLICATES_SUM/* duplicates will be summed */
                               ,&errval);
               if((!mtxAp) || (errval != RSB_ERR_NO_ERROR))
               {
                       printf("Error while allocating the matrix!0);
                       goto err;
               }
               printf("Correctly allocated a matrix.0);
               printf("Summary information of the matrix:0);
               /* print out the matrix summary information  */
               rsb_mtx_get_info_str(mtxAp,"RSB_MIF_MATRIX_INFO__TO__CHAR_P",
                               ib,sizeof(ib));
               printf("%s",ib);
               printf("0);

               if((errval =
                       rsb_spmv(RSB_TRANSPOSITION_N,&one,mtxAp,B,1,&one,X,1))
                               != RSB_ERR_NO_ERROR )
               {
                       printf("Error performing a multiplication!0);
                       goto err;
               }
               printf("Correctly performed a SPMV.0);
               rsb_mtx_free(mtxAp);
               printf("Correctly freed the matrix.0);
               if((errval = rsb_lib_exit(RSB_NULL_EXIT_OPTIONS))
                               != RSB_ERR_NO_ERROR)
               {
                       printf("Error finalizing the library!0);
                       goto err;
               }
               printf("Correctly finalized the library.0);
               printf("Program terminating with no error.0);
               return 0;
       err:
               rsb_perror(NULL,errval);
               printf("Program terminating with error.0);
               return -1;
       }

       /*

       Copyright (C) 2008-2015 Michele Martone

       This file is part of librsb.

       librsb is free software; you can redistribute it and/or modify it
       under the terms of the GNU Lesser General Public License as published
       by the Free Software Foundation; either version 3 of the License, or
       (at your option) any later version.

       librsb is distributed in the hope that it will be useful, but WITHOUT
       ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
       FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
       License for more details.

       You should have received a copy of the GNU Lesser General Public
       License along with librsb; see the file COPYING.
       If not, see <http://www.gnu.org/licenses/>.

       */
       /*!
        ingroup rsb-examples
        @file
        @author Michele Martone
        @brief This is a first "hello RSB" example program using
               a Sparse BLAS interface.

        include hello-spblas.c
       */
       #include <rsb.h> /* for rsb_lib_init */
       #include <blas_sparse.h> /* Sparse BLAS on the top of librsb */
       #include <stdio.h>       /* printf */

       int main(const int argc, char * const argv[])
       {
               /*!
                * A Hello/Sparse BLAS program.
                *
                * This program shows how to use the blas_sparse.h
                * interface correctly to:
                *
                * - initialize the library using #rsb_lib_init()
                * - allocate (build) a single sparse matrix in the RSB
                *   format using #BLAS_duscr_begin()/#BLAS_duscr_insert_entries()
                *   /#BLAS_duscr_end()
                * - extract one matrix element with #BLAS_dusget_element()
                * - multiply the matrix times a vector using #BLAS_dusmv()
                * - deallocate the matrix using #BLAS_usds()
                * - finalize the library using
                *   #rsb_lib_exit(#RSB_NULL_EXIT_OPTIONS)
               */
       #ifndef RSB_NUMERICAL_TYPE_DOUBLE
               printf("'double' type configured out."
               " Please reconfigure the library with it and recompile.0);
               return 0;
       #else /* RSB_NUMERICAL_TYPE_DOUBLE */
               blas_sparse_matrix A = blas_invalid_handle; /* handle for A */
               const int nnz = 4;      /* number of nonzeroes of matrix A */
               const int  nr = 3;      /* number of A's rows */
               const int  nc = 3;      /* number of A's columns */
               /* A's nonzero elements row indices (coordinates): */
               int   IA[] = { 0, 1, 2, 2 };
               /* A's nonzero elements column indices (coordinates): */
               int   JA[] = { 0, 1, 0, 2 };
               /* A's nonzero values (matrix coefficients): */
               double VA[] = { 11.0, 22.0, 13.0, 33.0  };
               /* the X vector's array: */
               double X[] = { 0.0, 0.0, 0.0 };
               /* the B vector's array: */
               double B[] = { -1.0, -2.0, -2.0 };
               /* the (known) result array: */
               double AB[] = { 11.0+26.0, 44.0, 66.0+13.0 };
               /* rsb error variable: */
               rsb_err_t errval = RSB_ERR_NO_ERROR;
               int i;

               printf("Hello, RSB!0);
               /* initialize the library */
               if((errval = rsb_lib_init(RSB_NULL_INIT_OPTIONS))
                               != RSB_ERR_NO_ERROR)
               {
                       goto err;
               }
               printf("Correctly initialized the library.0);

               /* initialize a matrix descriptor */
               A = BLAS_duscr_begin(nr,nc);
               if( A == blas_invalid_handle )
               {
                       goto err;
               }

               /* specify properties (e.g.: symmetry)*/
               if( BLAS_ussp(A,blas_lower_symmetric) != 0 )
               {
                       goto err;
               }

               /* get properties (e.g.: symmetry) */
               if( BLAS_usgp(A,blas_lower_symmetric) != 1 )
               {
                       printf("Symmetry property non set ?!0);
                       goto err;
               }

               /* insert the nonzeroes (here, all at once) */
               if( BLAS_duscr_insert_entries(A, nnz, VA, IA, JA)
                               == blas_invalid_handle)
               {
                       goto err;
               }

               /* finalize (allocate) the matrix build  */
               if( BLAS_duscr_end(A) == blas_invalid_handle )
               {
                       goto err;
               }
               printf("Correctly allocated a matrix.0);

               VA[0] = 0.0;
               if( BLAS_dusget_element(A, IA[0], JA[0], &VA[0]) )
               {
                       goto err;
               }

               /* a check */
               if( VA[0] != 11.0 )
               {
                       goto err;
               }

               /* compute X = X + (-1) * A * B   */
               if(BLAS_dusmv(blas_no_trans,-1,A,B,1,X,1))
               {
                       goto err;
               }

               for( i = 0 ; i < nc; ++i )
                       if( X[i] != AB[i] )
                       {
                               printf("Computed SPMV result seems wrong. Terminating.0);
                               goto err;
                       }
               printf("Correctly performed a SPMV.0);

               /* deallocate matrix A */
               if( BLAS_usds(A) )
               {
                       goto err;
               }
               printf("Correctly freed the matrix.0);

               /* finalize the library */
               if((errval = rsb_lib_exit(RSB_NULL_EXIT_OPTIONS))
                               != RSB_ERR_NO_ERROR)
               {
                       goto err;
               }
               printf("Correctly finalized the library.0);
               printf("Program terminating with no error.0);

               return 0;
       err:
               rsb_perror(NULL,errval);
               printf("Program terminating with error.0);
               return -1;
       #endif /* RSB_NUMERICAL_TYPE_DOUBLE */
       }

       /*

       Copyright (C) 2008-2015 Michele Martone

       This file is part of librsb.

       librsb is free software; you can redistribute it and/or modify it
       under the terms of the GNU Lesser General Public License as published
       by the Free Software Foundation; either version 3 of the License, or
       (at your option) any later version.

       librsb is distributed in the hope that it will be useful, but WITHOUT
       ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
       FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
       License for more details.

       You should have received a copy of the GNU Lesser General Public
       License along with librsb; see the file COPYING.
       If not, see <http://www.gnu.org/licenses/>.

       */
       /*!
        ingroup rsb-examples
        @file
        @author Michele Martone
        @brief This is a first "hello RSB" example program using
               a Sparse BLAS interface.

        include hello-spblas.c
       */
       #include <rsb.h> /* for rsb_lib_init */
       #include <blas_sparse.h> /* Sparse BLAS on the top of librsb */
       #include <stdio.h>       /* printf */

       int main(const int argc, char * const argv[])
       {
               /*!
                * A Hello/Sparse BLAS program.
                *
                * This program shows how to use the blas_sparse.h
                * interface correctly to:
                *
                * - initialize the library using #rsb_lib_init()
                * - allocate (build) a single sparse matrix in the RSB
                *   format using #BLAS_duscr_begin()/#BLAS_duscr_insert_entries()
                *   /#BLAS_duscr_end()
                * - extract one matrix element with #BLAS_dusget_element()
                * - multiply the matrix times a vector using #BLAS_dusmv()
                * - deallocate the matrix using #BLAS_usds()
                * - finalize the library using
                *   #rsb_lib_exit(#RSB_NULL_EXIT_OPTIONS)
               */
       #ifndef RSB_NUMERICAL_TYPE_DOUBLE
               printf("'double' type configured out."
               " Please reconfigure the library with it and recompile.0);
               return 0;
       #else /* RSB_NUMERICAL_TYPE_DOUBLE */
               blas_sparse_matrix A = blas_invalid_handle; /* handle for A */
               const int nnz = 4;      /* number of nonzeroes of matrix A */
               const int  nr = 3;      /* number of A's rows */
               const int  nc = 3;      /* number of A's columns */
               /* A's nonzero elements row indices (coordinates): */
               int   IA[] = { 0, 1, 2, 2 };
               /* A's nonzero elements column indices (coordinates): */
               int   JA[] = { 0, 1, 0, 2 };
               /* A's nonzero values (matrix coefficients): */
               double VA[] = { 11.0, 22.0, 13.0, 33.0  };
               /* the X vector's array: */
               double X[] = { 0.0, 0.0, 0.0 };
               /* the B vector's array: */
               double B[] = { -1.0, -2.0, -2.0 };
               /* the (known) result array: */
               double AB[] = { 11.0+26.0, 44.0, 66.0+13.0 };
               /* rsb error variable: */
               rsb_err_t errval = RSB_ERR_NO_ERROR;
               int i;

               printf("Hello, RSB!0);
               /* initialize the library */
               if((errval = rsb_lib_init(RSB_NULL_INIT_OPTIONS))
                               != RSB_ERR_NO_ERROR)
               {
                       goto err;
               }
               printf("Correctly initialized the library.0);

               /* initialize a matrix descriptor */
               A = BLAS_duscr_begin(nr,nc);
               if( A == blas_invalid_handle )
               {
                       goto err;
               }

               /* specify properties (e.g.: symmetry)*/
               if( BLAS_ussp(A,blas_lower_symmetric) != 0 )
               {
                       goto err;
               }

               /* get properties (e.g.: symmetry) */
               if( BLAS_usgp(A,blas_lower_symmetric) != 1 )
               {
                       printf("Symmetry property non set ?!0);
                       goto err;
               }

               /* insert the nonzeroes (here, all at once) */
               if( BLAS_duscr_insert_entries(A, nnz, VA, IA, JA)
                               == blas_invalid_handle)
               {
                       goto err;
               }

               /* finalize (allocate) the matrix build  */
               if( BLAS_duscr_end(A) == blas_invalid_handle )
               {
                       goto err;
               }
               printf("Correctly allocated a matrix.0);

               VA[0] = 0.0;
               if( BLAS_dusget_element(A, IA[0], JA[0], &VA[0]) )
               {
                       goto err;
               }

               /* a check */
               if( VA[0] != 11.0 )
               {
                       goto err;
               }

               /* compute X = X + (-1) * A * B   */
               if(BLAS_dusmv(blas_no_trans,-1,A,B,1,X,1))
               {
                       goto err;
               }

               for( i = 0 ; i < nc; ++i )
                       if( X[i] != AB[i] )
                       {
                               printf("Computed SPMV result seems wrong. Terminating.0);
                               goto err;
                       }
               printf("Correctly performed a SPMV.0);

               /* deallocate matrix A */
               if( BLAS_usds(A) )
               {
                       goto err;
               }
               printf("Correctly freed the matrix.0);

               /* finalize the library */
               if((errval = rsb_lib_exit(RSB_NULL_EXIT_OPTIONS))
                               != RSB_ERR_NO_ERROR)
               {
                       goto err;
               }
               printf("Correctly finalized the library.0);
               printf("Program terminating with no error.0);

               return 0;
       err:
               rsb_perror(NULL,errval);
               printf("Program terminating with error.0);
               return -1;
       #endif /* RSB_NUMERICAL_TYPE_DOUBLE */
       }

       /*

       Copyright (C) 2008-2015 Michele Martone

       This file is part of librsb.

       librsb is free software; you can redistribute it and/or modify it
       under the terms of the GNU Lesser General Public License as published
       by the Free Software Foundation; either version 3 of the License, or
       (at your option) any later version.

       librsb is distributed in the hope that it will be useful, but WITHOUT
       ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
       FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
       License for more details.

       You should have received a copy of the GNU Lesser General Public
       License along with librsb; see the file COPYING.
       If not, see <http://www.gnu.org/licenses/>.

       */
       /*!
        ingroup rsb-examples
        @file
        @author Michele Martone
        @brief This is a first "RSB autotuning" example program.

        include autotuning.c
       */
       #include <rsb.h> /* librsb header to include */
       #include <stdio.h>       /* printf() */
       #include <ctype.h>       /* isdigit() */
       #include <stdlib.h>      /* atoi() */
       /* #include "rsb_internals.h" */

       int tune_from_file(char * const filename, rsb_int_t wvat)
       {
               struct rsb_mtx_t *mtxMp = NULL;
               /* spmv specific variables */
               const RSB_DEFAULT_TYPE alpha = 1;
               const RSB_DEFAULT_TYPE beta = 1;
               rsb_flags_t order = RSB_FLAG_WANT_COLUMN_MAJOR_ORDER;
               const rsb_coo_idx_t nrhs = 2;  /* number of right hand sides */
               rsb_trans_t transA = RSB_TRANSPOSITION_N; /* transposition */
               rsb_nnz_idx_t ldB = 0;
               rsb_nnz_idx_t ldC = 0;
               /* misc variables */
               rsb_err_t errval = RSB_ERR_NO_ERROR;
               rsb_time_t dt;
               char ib[200];
               const char*is = "RSB_MIF_MATRIX_INFO__TO__CHAR_P";
               /* misc variables */
               /* input autotuning variables */
               rsb_int_t oitmax = 1 /*15*/;    /* auto-tune iterations */
               rsb_time_t tmax = 0.1;   /* time per autotune operation */
               /* output autotuning variables */
               rsb_flags_t flagsA = RSB_FLAG_NOFLAGS;
               int ione = 1;
               rsb_type_t typecodea [] = RSB_MATRIX_SPBLAS_TYPE_CODES_ARRAY;
               int typecodei;

               errval = rsb_lib_init(RSB_NULL_INIT_OPTIONS);

               if( (errval) != RSB_ERR_NO_ERROR )
                       goto err;

               errval = rsb_lib_set_opt(RSB_IO_WANT_VERBOSE_TUNING, &wvat );

               /*
               errval = rsb_lib_set_opt(RSB_IO_WANT_EXTRA_VERBOSE_INTERFACE, &ione);
               */

               if( (errval) != RSB_ERR_NO_ERROR )
                       goto err;

               printf("Loading matrix from file

               mtxMp = rsb_file_mtx_load(filename, flagsA, typecodea[0], &errval);

               if( (errval) != RSB_ERR_NO_ERROR )
                       goto err;

               for( typecodei = 0 ; typecodei < RSB_IMPLEMENTED_TYPES; ++typecodei )
               {
                       rsb_type_t typecode = typecodea[typecodei];
                       struct rsb_mtx_t *mtxAp = NULL;
                       struct rsb_mtx_t *mtxOp = NULL;
                       rsb_real_t sf = 0.0;
                       rsb_int_t tn = 0;

                       sf = 0.0;
                       tn = 0;

                       printf("Considering %c clone.0,typecode);

                       errval = rsb_mtx_clone(&mtxAp, typecode, transA, NULL, mtxMp,
                                       flagsA);

                       if( (errval) != RSB_ERR_NO_ERROR )
                               goto err;

                       printf("Base matrix:0);
                       rsb_mtx_get_info_str(mtxAp,is,ib,sizeof(ib));
                       printf("%s0,ib);

                       dt = -rsb_time();
                       errval = rsb_tune_spmm(NULL, &sf, &tn, oitmax, tmax, transA,
                            &alpha, mtxAp, nrhs, order, NULL, ldB, &beta, NULL, ldC);

                       dt += rsb_time();
                       if(tn == 0)
                       printf("After %lfs, autotuning routine did not find a better"
                               " threads count configuration.0,dt);
                       else
                       printf("After %lfs, thread autotuning declared speedup of %lg x,"
                               " when using threads count of %d.0,dt,sf,tn);
                       printf("0);

                       dt = -rsb_time();

                       mtxOp = mtxAp;
                       errval = rsb_tune_spmm(&mtxAp, &sf, &tn, oitmax, tmax, transA,
                               &alpha, NULL, nrhs, order, NULL, ldB, &beta, NULL, ldC);
                       if( (errval) != RSB_ERR_NO_ERROR )
                               goto err;

                       dt += rsb_time();
                       if( mtxOp == mtxAp )
                       {
                               printf("After %lfs, global autotuning found old matrix optimal,"
                               " with declared speedup %lg x when using %d threads0,dt,sf,tn);
                       }
                       else
                       {
                               printf("After %lfs, global autotuning declared speedup of %lg x,"
                               " when using threads count of %d and a new matrix:0,dt,sf,tn);
                               rsb_mtx_get_info_str(mtxAp,is,ib,sizeof(ib));
                               printf("%s0,ib);
                       }
                       printf("0);

                       /* user is expected to:
                       errval = rsb_lib_set_opt(RSB_IO_WANT_EXECUTING_THREADS,&tn);
                       and use mtxAp in SpMV.
                       */
                       rsb_mtx_free(mtxAp);
                       mtxAp = NULL;
               }
               rsb_mtx_free(mtxMp);
               mtxMp = NULL;

               goto ret;
       ret:
               return 0;
       err:
               rsb_perror(NULL,errval);
               printf("Program terminating with error.0);
               return -1;
       }

       int main(const int argc, char * const argv[])
       {
               /*!
                Autotuning example.
                */
               /* matrix variables */
               struct rsb_mtx_t *mtxAp = NULL; /* matrix structure pointer */
               const int bs = RSB_DEFAULT_BLOCKING;
               rsb_coo_idx_t nrA = 500; /* number of rows */
               rsb_coo_idx_t ncA = 500; /* number of cols */
               rsb_type_t typecode = RSB_NUMERICAL_TYPE_DEFAULT;
               rsb_coo_idx_t rd = 1; /* every rd rows one is non empty */
               rsb_coo_idx_t cd = 4; /* every cd cols one is non empty */
               rsb_nnz_idx_t nnzA = (nrA/rd)*(ncA/cd); /* nonzeroes */
               rsb_coo_idx_t*IA = NULL;
               rsb_coo_idx_t*JA = NULL;
               RSB_DEFAULT_TYPE*VA = NULL;
               /* spmv specific variables */
               const RSB_DEFAULT_TYPE alpha = 1;
               const RSB_DEFAULT_TYPE beta = 1;
               RSB_DEFAULT_TYPE*Cp = NULL;
               RSB_DEFAULT_TYPE*Bp = NULL;
               rsb_flags_t order = RSB_FLAG_WANT_COLUMN_MAJOR_ORDER;
               const rsb_coo_idx_t nrhs = 2;  /* number of right hand sides */
               rsb_trans_t transA = RSB_TRANSPOSITION_N; /* transposition */
               rsb_nnz_idx_t ldB = nrA;
               rsb_nnz_idx_t ldC = ncA;
               /* misc variables */
               rsb_err_t errval = RSB_ERR_NO_ERROR;
               size_t so = sizeof(RSB_DEFAULT_TYPE);
               size_t si = sizeof(rsb_coo_idx_t);
               rsb_time_t dt,odt;
               rsb_int_t t,tt = 100;   /* will repeat spmv tt times */
               char ib[200];
               const char*is = "RSB_MIF_MATRIX_INFO__TO__CHAR_P";
               /* misc counters */
               rsb_coo_idx_t ci;
               rsb_coo_idx_t ri;
               rsb_coo_idx_t ni;
               rsb_int_t nrhsi;
               /* misc variables */
               rsb_time_t etime = 0.0;
               /* input autotuning variables */
               rsb_int_t oitmax = 15;  /* auto-tune iterations */
               rsb_time_t tmax = 0.1;   /* time per autotune operation */
               /* input/output autotuning variables */
               rsb_int_t tn = 0;       /* threads number */
               /* output autotuning variables */
               rsb_real_t sf = 0.0;     /* speedup factor obtained from auto tuning */
               rsb_int_t wvat = 1;     /* want verbose autotuning; see documentation
                                          of RSB_IO_WANT_VERBOSE_TUNING */

               if(argc > 1 && !isdigit(argv[1][0]) )
                       return tune_from_file(argv[1],wvat);

               if(argc > 1)
               {
                       nrA = ncA = atoi(argv[1]);
                       if ( nrA < RSB_MIN_MATRIX_DIM || (nrA > (RSB_MAX_MATRIX_DIM) ))
                               goto err;

                       nnzA = (nrA/rd)*(ncA/cd);
                       ldB = nrA;
                       ldC = ncA;
               }

               printf("Creating %d x %d matrix with %d nonzeroes.0,nrA,ncA,nnzA);

               IA = calloc(nnzA, si);
               JA = calloc(nnzA, si);
               VA = calloc(nnzA, so);
               Bp = calloc(nrhs*ncA ,so);
               Cp = calloc(nrhs*nrA ,so);

               if( ! ( VA && IA && JA && Bp && Cp ) )
                       goto err;

               for(nrhsi=0;nrhsi<nrhs;++nrhsi)
                       for(ci=0;ci<ncA/cd;++ci)
                               Bp[nrhsi*ldC+ci] = 1.0;

               for(nrhsi=0;nrhsi<nrhs;++nrhsi)
                       for(ri=0;ri<nrA/rd;++ri)
                               Cp[nrhsi*ldC+ri] = 1.0;

               ni = 0;

               for(ci=0;ci<ncA/cd;++ci)
                       for(ri=0;ri<nrA/rd;++ri)
                       {
                               VA[ni] = nrA * ri + ci,
                               IA[ni] = ri;
                               JA[ni] = ci;
                               ni++;
                       }

               if((errval = rsb_lib_init(RSB_NULL_INIT_OPTIONS))
                               != RSB_ERR_NO_ERROR) goto err;

               errval = rsb_lib_set_opt(RSB_IO_WANT_VERBOSE_TUNING, &wvat );

               mtxAp = rsb_mtx_alloc_from_coo_const(
                       VA,IA,JA,nnzA,typecode,nrA,ncA,bs,bs,
                       RSB_FLAG_NOFLAGS,&errval);

               /* VA, IA, JA are not necessary anymore */
               free(VA);
               free(IA);
               free(JA);
               VA = NULL;
               IA = NULL;
               JA = NULL;

               if((!mtxAp) || (errval != RSB_ERR_NO_ERROR))
                       goto err;

               printf("Allocated matrix of %zd nonzeroes:0,(size_t)nnzA);
               rsb_mtx_get_info_str(mtxAp,is,ib,sizeof(ib));
               printf("%s0,ib);

               dt = - rsb_time();
               for(t=0;t<tt;++t)
                       /*
                          If nrhs == 1, the following is equivalent to
                          rsb_spmv(transA,&alpha,mtxAp,Bp,1,&beta,Cp,1);
                       */
                       rsb_spmm(transA,&alpha,mtxAp,nrhs,order,Bp,ldB,&beta,Cp,ldC);
               dt += rsb_time();
               odt = dt;
               printf("Before auto-tuning, %d multiplications took %lfs.0,tt,dt);

               printf("Threads autotuning (may take more than %lfs)...0,
                               oitmax*tmax);
               dt = -rsb_time();
               errval = rsb_tune_spmm(NULL, &sf, &tn, oitmax, tmax, transA,
                               &alpha, mtxAp, nrhs, order, Bp, ldB, &beta, Cp, ldC);
               dt += rsb_time();
               if(errval != RSB_ERR_NO_ERROR)
                       goto err;

               if(tn == 0)
               printf("After %lfs, autotuning routine did not find a better"
                               " threads count configuration.0,dt);
               else
               printf("After %lfs, autotuning routine declared speedup of %lg x,"
                               " when using threads count of %d.0,dt,sf,tn);

               errval = rsb_lib_set_opt(RSB_IO_WANT_EXECUTING_THREADS,&tn);
               if(errval != RSB_ERR_NO_ERROR)
                       goto err;

               rsb_mtx_get_info_str(mtxAp,is,ib,sizeof(ib));
               printf("%s0,ib);

               dt = -rsb_time();
               for(t=0;t<tt;++t)
                       /*rsb_spmv(transA,&alpha,mtxAp,Bp,1,&beta,Cp,1);*/
                       rsb_spmm(transA,&alpha,mtxAp,nrhs,order,Bp,ldB,&beta,Cp,ldC);
               dt += rsb_time();
               printf("After threads auto-tuning, %d multiplications took %lfs"
                               "  --  effective speedup of %lg x0,tt,dt,odt/dt);
               odt = dt;

               tn = 0; /* this will restore default threads count */
               errval = rsb_lib_set_opt(RSB_IO_WANT_EXECUTING_THREADS,&tn);
               if(errval != RSB_ERR_NO_ERROR)
                       goto err;
               errval = rsb_lib_get_opt(RSB_IO_WANT_EXECUTING_THREADS,&tn);
               if(errval != RSB_ERR_NO_ERROR)
                       goto err;

               printf("Matrix autotuning (may take more than %lfs; using %d"
                               " threads )...0, oitmax*tmax, tn);

               /* A negative tn will request also threads autotuning: */
               /* tn = -tn; */

               dt = -rsb_time();
               errval = rsb_tune_spmm(&mtxAp, &sf, &tn, oitmax, tmax, transA,
                               &alpha,  NULL, nrhs, order, Bp, ldB, &beta, Cp, ldC);
               dt += rsb_time();

               if(errval != RSB_ERR_NO_ERROR)
                       goto err;

               if(tn == 0)
               printf("After %lfs, autotuning routine did not find a better"
                               " threads count configuration.0,dt);
               else
               printf("After %lfs, autotuning routine declared speedup of %lg x,"
                               " when using threads count of %d.0,dt,sf,tn);

               rsb_mtx_get_info_str(mtxAp,is,ib,sizeof(ib));
               printf("%s0,ib);

               dt = -rsb_time();
               for(t=0;t<tt;++t)
                       /*rsb_spmv(transA,&alpha,mtxAp,Bp,1,&beta,Cp,1);*/
                       rsb_spmm(transA,&alpha,mtxAp,nrhs,order,Bp,ldB,&beta,Cp,ldC);
               dt += rsb_time();
               printf("After threads auto-tuning, %d multiplications took %lfs"
                               "  --  further speedup of %lg x0,tt,dt,odt/dt);

               rsb_mtx_free(mtxAp);
               free(Cp);
               free(Bp);

               errval = rsb_lib_get_opt(RSB_IO_WANT_LIBRSB_ETIME,&etime);
               if(errval == RSB_ERR_UNSUPPORTED_FEATURE)
               {
                       printf("librsb timer-based profiling is not supported in "
                       "this build. If you wish to have it, re-configure librsb "
                       "with its support. So you can safely ignore the error you"
                       " might just have seen printed out on screen.0);
                       errval = RSB_ERR_NO_ERROR;
               }
               else
               if(etime) /* This will only work if enabled at configure time. */
                       printf("Elapsed program time is %5.2lfs0,etime);

               if((errval = rsb_lib_exit(RSB_NULL_EXIT_OPTIONS))
                               !=RSB_ERR_NO_ERROR)
                       goto err;
               return 0;
       err:
               rsb_perror(NULL,errval);
               printf("Program terminating with error.0);
               return -1;
       }

       /*

       Copyright (C) 2008-2015 Michele Martone

       This file is part of librsb.

       librsb is free software; you can redistribute it and/or modify it
       under the terms of the GNU Lesser General Public License as published
       by the Free Software Foundation; either version 3 of the License, or
       (at your option) any later version.

       librsb is distributed in the hope that it will be useful, but WITHOUT
       ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
       FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
       License for more details.

       You should have received a copy of the GNU Lesser General Public
       License along with librsb; see the file COPYING.
       If not, see <http://www.gnu.org/licenses/>.

       */
       /*!
        ingroup rsb-examples
        @file
        @author Michele Martone
        @brief This is an example program using a Sparse BLAS interface
               and reading from file using the RSB library.

        include io-spblas.c
       */
       #include <rsb.h> /* for rsb_lib_init */
       #include <blas_sparse.h>
       #include <stdio.h>

       int main(const int argc, char * const argv[])
       {
       #ifndef RSB_NUMERICAL_TYPE_DOUBLE
               printf("Skipping a test because of 'double' type opted out.0);
               return 0;
       #else /* RSB_NUMERICAL_TYPE_DOUBLE */
               blas_sparse_matrix A = blas_invalid_handle;
               rsb_type_t typecode = RSB_NUMERICAL_TYPE_DOUBLE;
               rsb_char_t * filename = argc > 1 ? argv[1] : "pd.mtx";

               printf("Hello, RSB!0);
               if((rsb_perror(NULL,
                       rsb_lib_init(RSB_NULL_INIT_OPTIONS)))!=RSB_ERR_NO_ERROR)
               {
                       printf("Error while initializing the library.0);
                       goto err;
               }

               printf("Correctly initialized the library.0);

               A = rsb_load_spblas_matrix_file_as_matrix_market(filename,
                               typecode );
               if( A == blas_invalid_handle )
               {
                       printf("Error while loading matrix %s from file.0,
                                       filename);
                       goto err;
               }

               printf("Correctly loaded and allocated a matrix"
                               " from file %s.0,filename);

               if( BLAS_usgp(A,blas_symmetric) == 1 )
                       printf("Matrix is symmetric0);

               if( BLAS_usgp(A,blas_hermitian) == 1 )
                       printf("Matrix is hermitian0);

               printf("Now SPMV with NULL vectors will be attempted,"
                               " resulting in an error (so don't worry).0);

               if(BLAS_dusmv(blas_no_trans,-1,A,NULL,1,NULL,1))
               {
                       printf("Correctly detected an error condition.0);
                       goto okerr;
               }

               printf("No error detected ?0f you see this line printed out,"
                       " please report as a bug, because the above NULL pointers"
                       " should have been detected0);
               return -1;

       okerr:
               printf("Program correctly recovered from intentional"
                               " error condition.0);
               if(BLAS_usds(A))
               {
                       printf("Error while freeing the matrix!0);
                       goto err;
               }

               printf("Correctly freed the matrix.0);
       err:
               if(rsb_perror(NULL,
                       rsb_lib_exit(RSB_NULL_EXIT_OPTIONS))!=RSB_ERR_NO_ERROR)
               {
                       printf("Failed finalizing the library.0);
                       goto ferr;
               }

               printf("Correctly finalized the library.0);
               return 0;
       ferr:
               return -1;
       #endif /* RSB_NUMERICAL_TYPE_DOUBLE */
       }

       /*

       Copyright (C) 2008-2015 Michele Martone

       This file is part of librsb.

       librsb is free software; you can redistribute it and/or modify it
       under the terms of the GNU Lesser General Public License as published
       by the Free Software Foundation; either version 3 of the License, or
       (at your option) any later version.

       librsb is distributed in the hope that it will be useful, but WITHOUT
       ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
       FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
       License for more details.

       You should have received a copy of the GNU Lesser General Public
       License along with librsb; see the file COPYING.
       If not, see <http://www.gnu.org/licenses/>.

       */
       /*!
        @file
        @author Michele Martone
        @brief A toy program showing instantiation, transposition and other
        operations on a single matrix.
        ingroup rsb-examples

        include transpose.c
       */
       #include <rsb.h>
       #include <stdio.h>       /* printf */

       int main(const int argc, char * const argv[])
       {
               struct rsb_mtx_t *mtxAp = NULL;
               rsb_blk_idx_t brA = RSB_DEFAULT_BLOCKING, bcA=RSB_DEFAULT_BLOCKING;
               rsb_err_t errval = RSB_ERR_NO_ERROR;
               rsb_nnz_idx_t nnzA = 4;
               rsb_coo_idx_t  nrA = 3;
               rsb_coo_idx_t  ncA = 3;
               rsb_coo_idx_t    IA[] = { 0, 1, 2, 0 };
               rsb_coo_idx_t    JA[] = { 0, 1, 2, 2 };
               RSB_DEFAULT_TYPE VA[] = { 11, 22, 33, 13 };
               RSB_DEFAULT_TYPE XV[] = { 0,0,0,0,0,0 };
               rsb_coo_idx_t  vl = 0;
               rsb_type_t typecode = RSB_NUMERICAL_TYPE_DEFAULT;

               /* library initialization */
               if(rsb_lib_init(RSB_NULL_INIT_OPTIONS)!=RSB_ERR_NO_ERROR)
               {
                       return -1;
               }

               /* allocation */
               mtxAp = rsb_mtx_alloc_from_coo_const(
                               VA,IA,JA,nnzA,typecode,nrA,ncA,
                               brA,bcA,RSB_FLAG_NOFLAGS,NULL);
               if(!mtxAp)
               {
                       return -1;
               }

               /* printout */
               if(RSB_ERR_NO_ERROR!=(errval = rsb_file_mtx_save(mtxAp,NULL)))
               {
                       if(errval != RSB_ERR_UNSUPPORTED_FEATURE)
                               goto err;
               }

               /* matrix transposition */
               if( RSB_ERR_NO_ERROR != (errval =
                       rsb_mtx_clone(&mtxAp,RSB_NUMERICAL_TYPE_SAME_TYPE,
                       RSB_TRANSPOSITION_T,NULL,mtxAp,RSB_FLAG_IDENTICAL_FLAGS)))
               {
                       goto err;
               }

               /* printout */
               if(RSB_ERR_NO_ERROR!=(errval = rsb_file_mtx_save(mtxAp,NULL)))
               {
                       if(errval != RSB_ERR_UNSUPPORTED_FEATURE)
                               goto err;
               }

               rsb_mtx_free(mtxAp);

               /* doing the same after load from file */
               mtxAp = rsb_file_mtx_load("pd.mtx",
                       RSB_FLAG_NOFLAGS,typecode,NULL);
               if(!mtxAp)
               {
                       return -1;
               }

               /* printout */
               if(RSB_ERR_NO_ERROR!=(errval = rsb_file_mtx_save(mtxAp,NULL)))
               {
                       if(errval != RSB_ERR_UNSUPPORTED_FEATURE)
                               goto err;
               }

               /* one can see dimensions in advance, also */
               if(RSB_ERR_NO_ERROR!=(errval =
                       rsb_file_mtx_get_dims("pd.mtx",&nrA,&ncA,&nnzA,NULL)))
               {
                       if(errval != RSB_ERR_UNSUPPORTED_FEATURE)
                               goto err;
               }

               /* A matrix can be rendered to Postscript. */
               {
                       if(RSB_ERR_NO_ERROR!=(errval =
                       rsb_mtx_rndr("pd.eps",mtxAp,512,512,RSB_MARF_EPS_B)))
                               goto err;
               }

               rsb_mtx_free(mtxAp);

               /* also vectors can be loaded */
               if(RSB_ERR_NO_ERROR!=(errval =
                       rsb_file_vec_load("vf.mtx",typecode,NULL,&vl )))
                       goto err;
               /* we expecy vf.mtx to be 6 rows long */
               if( vl != 6 )
               {
                       goto err;
               }

               if(RSB_ERR_NO_ERROR!=(errval =
                       rsb_file_vec_load("vf.mtx",typecode,XV, NULL )))
                       goto err;

               /* matrices can be rendered from file to a pixelmap as well */
               {
                       unsigned char pixmap[3*2*2];

                       if(RSB_ERR_NO_ERROR!=(errval =
                       rsb_file_mtx_rndr(pixmap,"pd.mtx",2,2,2,RSB_MARF_RGB)))
                               goto err;
               }

               if(RSB_ERR_NO_ERROR != rsb_lib_exit(RSB_NULL_EXIT_OPTIONS))
               {
                       goto err;
               }
               return 0;
       err:
               rsb_perror(NULL,errval);
               return -1;
       }

       /*

       Copyright (C) 2008-2015 Michele Martone

       This file is part of librsb.

       librsb is free software; you can redistribute it and/or modify it
       under the terms of the GNU Lesser General Public License as published
       by the Free Software Foundation; either version 3 of the License, or
       (at your option) any later version.

       librsb is distributed in the hope that it will be useful, but WITHOUT
       ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
       FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
       License for more details.

       You should have received a copy of the GNU Lesser General Public
       License along with librsb; see the file COPYING.
       If not, see <http://www.gnu.org/licenses/>.

       */
       /*!
        @file
        @author Michele Martone
        @brief A toy program implementing the power method
               for computing matrix eigenvalues.
        ingroup rsb-examples

        include power.c
       */

       #include <stdio.h>       // printf
       #include <math.h>        // sqrt
       #include <stdlib.h>      // calloc
       #include <rsb.h>

       int main(const int argc, char * const argv[])
       {
               int WANT_VERBOSE = 0;
               struct rsb_mtx_t *mtxAp = NULL;
               const int bs = RSB_DEFAULT_BLOCKING;
               int i;
               const int br = bs, bc = bs; /* bs x bs blocked */
               rsb_err_t errval = 0;
               rsb_nnz_idx_t nnzA = 4;
               rsb_coo_idx_t  nrA = 3;
               rsb_coo_idx_t  ncA = 3;
               rsb_int_t it = 0, maxit = 100;
               const rsb_coo_idx_t    IA[] = { 0, 1, 2, 0 };
               const rsb_coo_idx_t    JA[] = { 0, 1, 2, 2 };
               const RSB_DEFAULT_POSSIBLY_FIRST_BLAS_TYPE VA[] = { 11, 22, 33, 13 };
               const RSB_DEFAULT_POSSIBLY_FIRST_BLAS_TYPE ZERO = 0;

               RSB_DEFAULT_POSSIBLY_FIRST_BLAS_TYPE norm = 0.0, /* nu */
               oldnorm = 1.0, /* oldnorm */
               *b1 = NULL, *b2 = NULL,
               *bnow = NULL, *bnext = NULL;/* b1 and b2 aliases */
               rsb_type_t typecode = RSB_NUMERICAL_TYPE_FIRST_BLAS;
               size_t ds = 0;
               /* tolerance */
               const RSB_DEFAULT_POSSIBLY_FIRST_BLAS_TYPE tol = 1e-14;

               /* library initialization */
               if(rsb_lib_init(RSB_NULL_INIT_OPTIONS)!=RSB_ERR_NO_ERROR)
                       return -1;

               /* allocation */
               mtxAp = rsb_mtx_alloc_from_coo_const(VA,IA,JA,nnzA,
                               typecode,nrA,ncA,br,bc,RSB_FLAG_NOFLAGS,NULL);
               if(!mtxAp)
                       return -1;

               ds = (nrA)*sizeof(RSB_DEFAULT_POSSIBLY_FIRST_BLAS_TYPE);
               b1 = calloc(1,ds);
               b2 = calloc(1,ds);

               if(! (b1 && b2))
               {
                       errval = RSB_ERR_ENOMEM;
                       goto err;
               }

               for( i = 0; i < nrA; ++i )
                       b1[i] = 1;

               bnow = b1, bnext = b2;/* b,b' */

               while( fabs(norm-oldnorm) > tol && it<maxit )
               {
                       ++ it;
                       oldnorm = norm;
                       /* b'<-Ab */
                       if(( rsb_spmv(RSB_TRANSPOSITION_N,NULL,mtxAp,bnow,
                               1,&ZERO,bnext,1)) != RSB_ERR_NO_ERROR )
                               goto err;
                       /* nu<-||Ab||^2 */
                       norm = 0;
                       for(i=0;i<nrA;++i)
                               norm += bnext[i]*bnext[i];
                       /* nu<-||Ab|| */
                       norm = sqrt(norm);
                       norm = 1.0/norm;
                       /* b'<- Ab / ||Ab|| */
                       for(i=0;i<nrA;++i)
                               bnext[i] *= norm;
                       norm = 1.0/norm;
                       printf("it:%d norm:%lg norm diff:%lg0,it,norm,norm-oldnorm);

                       {void *tmp=bnow;bnow=bnext;bnext=tmp;/* pointers swap */}
                       if(WANT_VERBOSE)
                       {
                               printf("norm:%lg0,norm);
                               if(isinf(norm))
                               /* isinf is a C99 feature (need correct
                                * compilation flags) */
                                       goto err;

                               for(i=0;i<2;++i)
                                       printf("x[%d]=%lg0,i,((double*)bnext)[i]);
                       }
               }
               /* the biggest eigenvalue should be in bnow */

               rsb_mtx_free(mtxAp);
               free(b1);
               free(b2);
               if(rsb_lib_exit(RSB_NULL_EXIT_OPTIONS)!=RSB_ERR_NO_ERROR)
                       goto err;
               if( it == maxit )
               {
                       printf("ERROR: hit iterations limit without convergence!");
                       errval=RSB_ERR_GENERIC_ERROR;
               }
               return 0;
       err:
               rsb_perror(NULL,errval);
               return -1;
       }

       !
       ! Copyright (C) 2008-2016 Michele Martone
       !
       ! This file is part of librsb.
       !
       ! librsb is free software; you can redistribute it and/or modify it
       ! under the terms of the GNU Lesser General Public License as published
       ! by the Free Software Foundation; either version 3 of the License, or
       ! (at your option) any later version.
       !
       ! librsb is distributed in the hope that it will be useful, but WITHOUT
       ! ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
       ! FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
       ! License for more details.
       !
       ! You should have received a copy of the GNU Lesser General Public
       ! License along with librsb; see the file COPYING.
       ! If not, see <http://www.gnu.org/licenses/>.
       !

             SUBROUTINE blas_sparse_mod_example(res)
             USE blas_sparse
             USE rsb ! For the second part of the example
             IMPLICIT NONE
             INTEGER :: res, istat = 0, i
             TYPE(c_ptr),TARGET :: mtxap = c_null_ptr ! matrix pointer
             INTEGER :: a
             INTEGER,PARAMETER :: transn = blas_no_trans
             INTEGER,PARAMETER :: incx = 1
             INTEGER,PARAMETER :: incy = 1
             REAL(KIND=8),PARAMETER :: alpha = 3
       ! Symmetric (declared via lower triangle) matrix based example, e.g.:
       ! 1 0
       ! 1 1
             ! declaration of VA,IA,JA
             !INTEGER,PARAMETER :: nr = 100
             INTEGER,PARAMETER :: nr = 20
             INTEGER,PARAMETER :: nc = nr
             INTEGER,PARAMETER :: nnz = (nr*(nr+1))/2 ! half the square
             INTEGER :: nt = 0
             INTEGER :: ic, ir
             INTEGER,PARAMETER :: ia(nnz) = (/ (((ir), ic=1,ir), ir=1,nr ) /) ! (/1, 2, 2/)
             INTEGER,PARAMETER :: ja(nnz) = (/ (((ic), ic=1,ir), ir=1,nr ) /) ! (/1, 1, 2/)
             REAL(KIND=8),PARAMETER :: va(nnz) = (/ ((1, ic=1,ir), ir=1,nr ) /) ! (/1, 1, 1/)
             REAL(KIND=8) :: x(nc) = (/((1), ir=1,nc)/) ! reference x ! (/1, 1/)
             REAL(KIND=8),PARAMETER :: cy(nr) = (/((alpha+alpha*nr), ir=1,nr)/) ! reference cy after ! (/9, 9/)
             REAL(KIND=8) :: y(nr) = (/((alpha), ir=1,nr)/) ! y will be overwritten ! (/3, 3/)
             ! First example part: pure blas_sparse code.
             res = 0
             CALL duscr_begin(nr,nc,a,res)
             IF (res.NE.0) GOTO 9999
             CALL ussp(a,blas_lower_symmetric,istat)
             IF (istat.NE.0) GOTO 9997
             CALL ussp(a,blas_rsb_spmv_autotuning_on,istat) ! (experimental) turns auto-tuning + thread setting on
             IF (istat.NE.0) print *,"autotuning returned nonzero:", istat &
              &," ...did you enable autotuning ?"
             !
             ! First style example
             CALL uscr_insert_entries(a,nnz,va,ia,ja,istat)
             IF (istat.NE.0) GOTO 9997
             CALL uscr_end(a,istat)
             IF (istat.NE.0) GOTO 9997
             ! CALL ussp(A,blas_rsb_duplicates_sum,istat)
             ! CALL uscr_insert_entries(A,nnz,VA,IA,JA,istat) ! uncomment this to activate add of coefficients to pattern
             CALL usgp(a,blas_rsb_spmv_autotuning_on,nt)  ! (experimental)
             IF (nt.NE.0) print*,"autotuner chose ",nt," threads"
             CALL ussp(a,blas_rsb_spmv_autotuning_off,istat) ! (experimental) turns auto-tuning + thread setting off
             IF (istat.NE.0) GOTO 9997

             CALL usmv(transn,alpha,a,x,incx,y,incy,istat)
             IF (istat.NE.0) GOTO 9997
             !
             DO i = 1, nr
                   IF (y(i).NE.cy(i)) print *, "first check results are not ok"
                   IF (y(i).NE.cy(i)) GOTO 9997
             END DO
             !
             y(:) = alpha ! reset
             !
             ! Second style example
             CALL ussp(a,blas_rsb_autotune_next_operation,istat) ! (experimental) turns auto-tuning + thread setting on
             IF (istat.NE.0) GOTO 9997
             CALL usmv(transn,alpha,a,x,incx,y,incy,istat)
             CALL usmm(blas_colmajor,transn,1, alpha,a,x,nr,y,nc,istat) ! Equivalent to the above (as long as incx=incy=1).
             CALL usmm(blas_colmajor,transn,1,-alpha,a,x,nr,y,nc,istat) ! Subtract the last usmm call contribution.
             IF (istat.NE.0) GOTO 9997
             !
             DO i = 1, nr
                   IF (y(i).NE.cy(i)) print *,"second check results are not ok"
                   IF (y(i).NE.cy(i)) GOTO 9997
             END DO
             !
             print *, "check results are ok"

             ! Second part of the example: access to the rsb.h interface via
             ! the ISO C Binding interface.
             mtxap = rsb_blas_get_mtx(a) ! get pointer to rsb structure (as in the rsb.h API)
             IF(nr.LT.5) istat = rsb_file_mtx_save(mtxap,c_null_ptr) ! write to stdout (only if matrix small enough)

             GOTO 9998
       9997      res = -1
       9998      CONTINUE
             CALL usds(a,istat)
             IF (istat.NE.0) res = -1
       9999      CONTINUE
             end SUBROUTINE blas_sparse_mod_example

             PROGRAM main
             USE rsb, ONLY: rsb_lib_init, rsb_lib_exit, c_ptr, c_null_ptr,&
              & rsb_io_want_extra_verbose_interface,rsb_io_want_verbose_tuning,&
              & rsb_lib_set_opt
             USE iso_c_binding
             IMPLICIT NONE
             INTEGER :: res = 0, passed = 0, failed = 0
             !TYPE(C_PTR),PARAMETER :: EO = RSB_NULL_EXIT_OPTIONS
             !TYPE(C_PTR),PARAMETER :: IO = RSB_NULL_INIT_OPTIONS
             ! Note: using C_NULL_PTR instead of the previous lines because of http://gcc.gnu.org/bugzilla/show_bug.cgi?id=59411
             TYPE(c_ptr),PARAMETER :: eo = c_null_ptr
             TYPE(c_ptr),PARAMETER :: io = c_null_ptr
             INTEGER,TARGET::ione=1
             res = rsb_lib_init(io)
             res = rsb_lib_set_opt(rsb_io_want_verbose_tuning,c_loc(ione))

             CALL blas_sparse_mod_example(res)
             IF (res.LT.0) failed = failed + 1
             IF (res.EQ.0) passed = passed + 1

             res = rsb_lib_exit(eo)

             print *, "FAILED:", failed
             print *, "PASSED:", passed
             IF (failed .GT. 0) THEN
              stop 1
             END IF
             END PROGRAM

       !
       ! Copyright (C) 2008-2016 Michele Martone
       !
       ! This file is part of librsb.
       !
       ! librsb is free software; you can redistribute it and/or modify it
       ! under the terms of the GNU Lesser General Public License as published
       ! by the Free Software Foundation; either version 3 of the License, or
       ! (at your option) any later version.
       !
       ! librsb is distributed in the hope that it will be useful, but WITHOUT
       ! ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
       ! FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
       ! License for more details.
       !
       ! You should have received a copy of the GNU Lesser General Public
       ! License along with librsb; see the file COPYING.
       ! If not, see <http://www.gnu.org/licenses/>.
       !
             SUBROUTINE rsb_mod_example1(res)
             USE rsb
             USE iso_c_binding
             IMPLICIT NONE
             INTEGER ::res
             INTEGER,TARGET :: istat = 0, i
             INTEGER :: transt = rsb_transposition_n ! Please note that this interface is unfinished
             INTEGER :: incx = 1, incy = 1
             REAL(KIND=8),TARGET :: alpha = 3, beta = 1
       ! 1 1
       ! 1 1
             ! declaration of VA,IA,JA
             INTEGER :: nnz = 4
             INTEGER :: nr = 2
             INTEGER :: nc = 2
             INTEGER :: nrhs = 1
             INTEGER :: order = rsb_flag_want_column_major_order ! rhs layout
             INTEGER :: flags = rsb_flag_noflags
             INTEGER,TARGET :: ia(4) = (/0, 1, 1,0/)
             INTEGER,TARGET :: ja(4) = (/0, 0, 1,1/)
             REAL(KIND=8),TARGET :: va(4) = (/1,1,1,1/)
             REAL(KIND=8),TARGET :: x(2) = (/1, 1/)! reference x
             REAL(KIND=8),TARGET :: cy(2) = (/9, 9/)! reference cy after
             REAL(KIND=8),TARGET :: y(2) = (/3, 3/)! y will be overwritten
             TYPE(c_ptr),TARGET :: mtxap = c_null_ptr ! matrix pointer
             REAL(KIND=8) :: tmax = 2.0 ! tuning max time
             INTEGER :: titmax = 2 ! tuning max iterations
             INTEGER,TARGET :: ont = 0     ! optimal number of threads

             res = 0
             mtxap = rsb_mtx_alloc_from_coo_const(c_loc(va),c_loc(ia),c_loc(ja)&
              &,nnz,&
              & rsb_numerical_type_double,nr,nc,1,1,flags,c_loc(istat))

             IF (istat.NE.rsb_err_no_error) GOTO 9997

             istat = rsb_file_mtx_save(mtxap,c_null_ptr)

             ! Structure autotuning:
             istat = rsb_tune_spmm(c_loc(mtxap),c_null_ptr,c_null_ptr,titmax,&
              & tmax,&
              & transt,c_loc(alpha),c_null_ptr,nrhs,order,c_loc(x),nr,&
              & c_loc(beta),c_loc(y),nc)

             IF (istat.NE.rsb_err_no_error) GOTO 9997

             ! Thread count autotuning:
             istat = rsb_tune_spmm(c_null_ptr,c_null_ptr,c_loc(ont),titmax,&
              & tmax,&
              & transt,c_loc(alpha),mtxap,nrhs,order,c_loc(x),nr,c_loc(beta),&
              & c_loc(y),nc)
             print *, "Optimal number of threads:", ont

             y(:) = (/3, 3/)! reference y
             IF (istat.NE.rsb_err_no_error) GOTO 9997

             istat = rsb_file_mtx_save(mtxap,c_null_ptr)
             IF (istat.NE.rsb_err_no_error) GOTO 9997

             istat = rsb_spmv(transt,c_loc(alpha),mtxap,c_loc(x),incx,&
              & c_loc(beta),c_loc(y),incy)
             IF (istat.NE.rsb_err_no_error) GOTO 9997
             DO i = 1, 2
                   IF (y(i).NE.cy(i)) print *, "type=d dims=2x2 sym=g diag=g &
             &blocks=1x1 usmv alpha= 3 beta= 1 incx=1 incy=1 trans=n is not ok"
                   IF (y(i).NE.cy(i)) GOTO 9997
             END DO
             print*,"type=d dims=2x2 sym=g diag=g blocks=1x1 usmv alpha= 3&
              & beta= 1 incx=1 incy=1 trans=n is ok"
             GOTO 9998
       9997      res = -1
       9998      CONTINUE
             mtxap = rsb_mtx_free(mtxap)
             IF (istat.NE.rsb_err_no_error) res = -1
       ! 9999      CONTINUE
             istat = rsb_perror(c_null_ptr,istat)
             end SUBROUTINE rsb_mod_example1

             SUBROUTINE rsb_mod_example2(res)
             USE rsb
             USE iso_c_binding
             IMPLICIT NONE
             INTEGER,TARGET :: errval
             INTEGER :: res
             INTEGER :: transt = rsb_transposition_n  ! no transposition
             INTEGER :: incx = 1, incb = 1        ! X, B vectors increment
             REAL(KIND=8),TARGET :: alpha = 3,beta = 1
             INTEGER :: nnza = 4, nra = 3, nca = 3     ! nonzeroes, rows, columns of matrix A
             INTEGER,TARGET :: ia(4) = (/1, 2, 3, 3/)  ! row    indices
             INTEGER,TARGET :: ja(4) = (/1, 2, 1, 3/)  ! column indices
             INTEGER(C_SIGNED_CHAR) :: typecode = rsb_numerical_type_double
             INTEGER :: flags =rsb_flag_default_matrix_flags+rsb_flag_symmetric
             REAL(KIND=8),TARGET :: va(4) = (/11.0, 22.0, 13.0, 33.0/) ! coefficients
             REAL(KIND=8),TARGET :: x(3) = (/   0,    0,    0/)
             REAL(KIND=8),TARGET :: b(3) = (/-1.0, -2.0, -2.0/)
             TYPE(c_ptr),TARGET  :: mtxap = c_null_ptr
             TYPE(c_ptr)  :: mtxapp = c_null_ptr
             REAL(KIND=8),TARGET :: etime = 0.0
             !TYPE(C_PTR),PARAMETER :: EO = RSB_NULL_EXIT_OPTIONS
             !TYPE(C_PTR),PARAMETER :: IO = RSB_NULL_INIT_OPTIONS
             ! Note: using C_NULL_PTR instead of the previous lines because of http://gcc.gnu.org/bugzilla/show_bug.cgi?id=59411
             TYPE(c_ptr),PARAMETER :: eo = c_null_ptr
             TYPE(c_ptr),PARAMETER :: io = c_null_ptr

             errval = rsb_lib_init(io)                ! librsb initialization
             IF (errval.NE.rsb_err_no_error) &
              & stop "error calling rsb_lib_init"
       #if defined(__GNUC__) && (__GNUC__ == 4) && (__GNUC_MINOR__ < 5)
       #define RSB_SKIP_BECAUSE_OLD_COMPILER 1
       #endif
       #ifndef RSB_SKIP_BECAUSE_OLD_COMPILER
             mtxap = rsb_mtx_alloc_from_coo_begin(nnza,typecode,nra,nca,flags,&
              & c_loc(errval)) ! begin matrix creation
             errval = rsb_mtx_set_vals(mtxap,&
              & c_loc(va),c_loc(ia),c_loc(ja),nnza,flags) ! insert some nonzeroes
             mtxapp = c_loc(mtxap) ! Old compilers like e.g.: Gfortran 4.4.7 will NOT compile this.
             IF (errval.NE.rsb_err_no_error) &
              & stop "error calling rsb_mtx_set_vals"
             errval = rsb_mtx_alloc_from_coo_end(mtxapp)                   ! end matrix creation
             IF (errval.NE.rsb_err_no_error) &
              & stop "error calling rsb_mtx_alloc_from_coo_end"
             errval = rsb_spmv(transt,c_loc(alpha),mtxap,c_loc(x),&
              & incx,c_loc(beta),c_loc(b),incb) ! X := X + (3) * A * B
             IF (errval.NE.rsb_err_no_error)&
              & stop "error calling rsb_spmv"
             mtxap = rsb_mtx_free(mtxap)                                 ! destroy matrix

             ! The following is optional and depends on configure options, so it is allowed to fail
             errval = rsb_lib_get_opt(rsb_io_want_librsb_etime,c_loc(etime))
             IF (errval.EQ.rsb_err_no_error)&
              & print*,"Time spent in librsb is:",etime
             ! IF (errval.NE.0)STOP "error calling rsb_lib_get_opt"
             errval = rsb_err_no_error

             IF (errval.NE.rsb_err_no_error) &
              & stop "error calling rsb_mtx_free"
       #else
             print*,"You have an old Fortran compiler not supporting C_LOC."
             print*,"Skipping a part of the test"
       #endif
             errval=rsb_lib_exit(eo)                 ! librsb finalization
             IF (errval.NE.rsb_err_no_error)&
              & stop "error calling rsb_lib_exit"
             print *, "rsb module fortran test is ok"
             res = errval
             end SUBROUTINE rsb_mod_example2

             PROGRAM main
             USE rsb
             IMPLICIT NONE
             INTEGER :: res = rsb_err_no_error, passed = 0, failed = 0
             !TYPE(C_PTR),PARAMETER :: EO = RSB_NULL_EXIT_OPTIONS
             !TYPE(C_PTR),PARAMETER :: IO = RSB_NULL_INIT_OPTIONS
             ! Note: using C_NULL_PTR instead of the previous lines because of http://gcc.gnu.org/bugzilla/show_bug.cgi?id=59411
             TYPE(c_ptr),PARAMETER :: eo = c_null_ptr
             TYPE(c_ptr),PARAMETER :: io = c_null_ptr

             res = rsb_lib_init(io)

             CALL rsb_mod_example1(res)
             IF (res.LT.0) failed = failed + 1
             IF (res.EQ.0) passed = passed + 1

             res = rsb_lib_exit(eo)

             CALL rsb_mod_example2(res)
             IF (res.LT.0) failed = failed + 1
             IF (res.EQ.0) passed = passed + 1

             print *, "FAILED:", failed
             print *, "PASSED:", passed
             IF (failed.GT.0) THEN
              stop 1
             END IF
             END PROGRAM

Author

       librsb was written by Michele Martone; this documentation has been generated by Doxygen.

SEE ALSO

       rsb-examples rsb-spblas.h rsb.h