Sparse multiplication Matrice comme (maxmin) en C ++ en utilisant les bibliothèques Octave

Question

J'implémenter une fonction de maxmin, il fonctionne comme la multiplication de la matrice, mais au lieu de produits sommateurs il obtient max min entre deux nombres pointwise. Un exemple de mise en œuvre naïve est

double mx = 0;
double mn = 0;
for (i = 0; i < rowsC;i++)
{
    for(j = 0; j < colsC;j++)
    {
        mx = 0;
        for(k = 0; k < colsA; k++)
        { 
            if (a(i, k) < b(k, j))
                mn = a(i,k);
            else
                mn = b(k,j);

            if (mn > mx)
                mx = mn;
        } 
        c(i, j) = mx;
    }
}

Je le coder en tant que fichier-oct Octave donc je dois utiliser oct.h structure de données. Le problème est que je veux mettre en œuvre une version rare, mais généralement vous besoin d'une référence à l'élément suivant non nul dans une ligne ou dans une colonne comme dans cet exemple (voir 4.3 algorithme): http://www.eecs.harvard.edu/ ~ Ellard / Q-97 / HTML / root / node20.html

Il fait row_p-> donne ensuite l'élément suivant de la ligne non nulle (même pour la colonne). Est-il possible de faire la même chose avec la classe d'octave Matrice Creuse? Ou est-il une autre façon de mettre en œuvre la multiplication de matrices creuses je peux adopter pour ma fonction maxmin?

Answer 1

Je ne sais pas si anyoe serait jamais intéressé, mais je réussi à trouver une solution. Le code de la solution est une partie de fl-core1.0 un ensemble de base de logique floue pour Octave et il est distribué sous LGPL. (Le code se fonde sur certaines fonctions d'octave)

// Calculate the S-Norm/T-Norm composition of sparse matrices (single thread)
void sparse_compose(octave_value_list args)
{
    // Create constant versions of the input matrices to prevent them to be filled by zeros on reading.
    // a is the const reference to the transpose of a because octave sparse matrices are column compressed
    // (to cycle on the rows, we cycle on the columns of the transpose).
    SparseMatrix atmp = args(0).sparse_matrix_value();
    const SparseMatrix a = atmp.transpose();
    const SparseMatrix b = args(1).sparse_matrix_value();

    // Declare variables for the T-Norm and S-Norm values 
    float snorm_val;    
    float tnorm_val;    

    // Initialize the result sparse matrix
    sparseC = SparseMatrix((int)colsB, (int)rowsA, (int)(colsB*rowsA));

    // Initialize the number of nonzero elements in the sparse matrix c
    int nel = 0;
    sparseC.xcidx(0) = 0;

    // Calculate the composition for each element
    for (int i = 0; i < rowsC; i++)
    {
        for(int j = 0; j < colsC; j++)
        {

            // Get the index of the first element of the i-th column of a transpose (i-th row of a)
            // and the index of the first element of the j-th column of b
            int ka = a.cidx(i);
            int kb = b.cidx(j);
            snorm_val = 0;

            // Check if the values of the matrix are really not 0 (it happens if the column of a or b hasn't any value)
            // because otherwise the cidx(i) or cidx(j) returns the first nonzero element of the previous column
            if(a(a.ridx(ka),i)!=0 && b(b.ridx(kb),j)!=0)
            {
                // Cicle on the i-th column of a transpose (i-th row of a) and j-th column of b
                // From a.cidx(i) to a.cidx(i+1)-1 there are all the nonzero elements of the column i of a transpose (i-th row of a)
                // From b.cidx(j) to b.cidx(j+1)-1 there are all the nonzero elements of the column j of b
                while ((ka <= (a.cidx(i+1)-1)) && (kb <= (b.cidx(j+1)-1)))
                {

                    // If a.ridx(ka) == b.ridx(kb) is true, then there's a nonzero value on the same row
                    // so there's a k for that a'(k, i) (equals to a(i, k)) and b(k, j) are both nonzero
                    if (a.ridx(ka) == b.ridx(kb))
                    {
                        tnorm_val = calc_tnorm(a.data(ka), b.data(kb)); 
                        snorm_val = calc_snorm(snorm_val, tnorm_val);
                        ka++;
                        kb++;
                    }

                    // If a.ridx(ka) == b.ridx(kb) ka should become the index of the next nonzero element on the i column of a 
                    // transpose (i row of a)
                    else if (a.ridx(ka) < b.ridx(kb))           
                        ka++;
                    // If a.ridx(ka) > b.ridx(kb) kb should become the index of the next nonzero element on the j column of b
                    else
                        kb++;
                }
            }

            if (snorm_val != 0)
            {
                // Equivalent to sparseC(i, j) = snorm_val;
                sparseC.xridx(nel) = j;
                sparseC.xdata(nel++) = snorm_val;
            }
        }
        sparseC.xcidx(i+1) = nel;
    }

    // Compress the result sparse matrix because it is initialized with a number of nonzero element probably greater than the real one
    sparseC.maybe_compress();

    // Transpose the result
    sparseC = sparseC.transpose();
}