bcsc_cnorm.c

Go to the documentation of this file.

/**
 *
 * @file bcsc_cnorm.c
 *
 * @copyright 2004-2023 Bordeaux INP, CNRS (LaBRI UMR 5800), Inria,
 *                      Univ. Bordeaux. All rights reserved.
 *
 * @version 6.3.2
 * @author Mathieu Faverge
 * @author Theophile Terraz
 * @date 2023-07-21
 * @generated from /builds/solverstack/pastix/bcsc/bcsc_znorm.c, normal z -> c, Wed Dec 13 12:09:45 2023
 *
 **/
#include "common.h"
#include "bcsc/bcsc.h"
#include <math.h>
#include "frobeniusupdate.h"
 
/**
 *******************************************************************************
 *
 * @ingroup bcsc_internal
 *
 * @brief Compute the max norm of a bcsc matrix.
 *
 *******************************************************************************
 *
 * @param[in] bcsc
 *          TODO
 *
 *******************************************************************************
 *
 * @retval The norm of the matrix.
 *
 *******************************************************************************/
float
bcsc_cnorm_max( const pastix_bcsc_t *bcsc )
{
    float temp;
    float norm = 0.;
    pastix_complex32_t *valptr = (pastix_complex32_t*)bcsc->Lvalues;
    pastix_int_t i, j, bloc;
 
    for( bloc=0; bloc < bcsc->cscfnbr; bloc++ )
    {
        for( j=0; j < bcsc->cscftab[bloc].colnbr; j++ )
        {
            for( i = bcsc->cscftab[bloc].coltab[j]; i < bcsc->cscftab[bloc].coltab[j+1]; i++ )
            {
                temp = cabsf(valptr[i]);
                if(norm < temp)
                {
                    norm = temp;
                }
            }
        }
    }
 
    return norm;
}
 
/**
 *******************************************************************************
 *
 * @ingroup bcsc_internal
 *
 * @brief Compute the infinity norm of a bcsc matrix.
 * The infinity norm is equal to the maximum value of the sum of the
 * absolute values of the elements of each rows.
 *
 *******************************************************************************
 *
 * @param[in] bcsc
 *          The Pastix bcsc.
 *
 *******************************************************************************
 *
 * @retval The norm of the matrix.
 *
 *******************************************************************************/
float
bcsc_cnorm_inf( const pastix_bcsc_t *bcsc )
{
    float norm = 0.;
    pastix_complex32_t *valptr;
    int i, j, bloc;
 
    if( bcsc->Uvalues != NULL )
    {
        float sum;
 
        valptr = (pastix_complex32_t*)(bcsc->Uvalues);
        for( bloc=0; bloc<bcsc->cscfnbr; bloc++ )
        {
            for( j=0; j<bcsc->cscftab[bloc].colnbr; j++ )
            {
                sum = 0.;
                for( i = bcsc->cscftab[bloc].coltab[j]; i < bcsc->cscftab[bloc].coltab[j+1]; i++ )
                {
                    sum += cabsf(valptr[i]);
                }
                if( sum > norm ) {
                    norm = sum;
                }
            }
        }
    }
    else {
        float *sumrow;
        valptr = (pastix_complex32_t*)bcsc->Lvalues;
 
        MALLOC_INTERN( sumrow, bcsc->gN, float);
        memset( sumrow, 0, bcsc->gN * sizeof(float) );
 
        for( bloc=0; bloc < bcsc->cscfnbr; bloc++ )
        {
            for( j=0; j < bcsc->cscftab[bloc].colnbr; j++ )
            {
                for( i = bcsc->cscftab[bloc].coltab[j]; i < bcsc->cscftab[bloc].coltab[j+1]; i++ )
                {
                    sumrow[ bcsc->rowtab[i] ] += cabsf(valptr[i]);
                }
            }
        }
 
        for( i=0; i<bcsc->gN; i++)
        {
            if(norm < sumrow[i])
            {
                norm = sumrow[i];
            }
        }
        memFree_null( sumrow );
    }
 
    return norm;
}
 
/**
 *******************************************************************************
 *
 * @ingroup bcsc_internal
 *
 * @brief Compute the norm 1 of a bcsc matrix.
 * Norm 1 is equal to the maximum value of the sum of the
 * absolute values of the elements of each columns.
 *
 *******************************************************************************
 *
 * @param[in] bcsc
 *          The Pastix bcsc.
 *
 *******************************************************************************
 *
 * @retval The norm of the matrix.
 *
 *******************************************************************************/
float
bcsc_cnorm_one( const pastix_bcsc_t *bcsc )
{
    pastix_complex32_t *valptr = (pastix_complex32_t*)bcsc->Lvalues;
    float sum, norm = 0.;
    int i, j, bloc;
 
    for( bloc=0; bloc<bcsc->cscfnbr; bloc++ )
    {
        for( j=0; j<bcsc->cscftab[bloc].colnbr; j++ )
        {
            sum = 0.;
            for( i = bcsc->cscftab[bloc].coltab[j]; i < bcsc->cscftab[bloc].coltab[j+1]; i++ )
            {
                sum += cabsf(valptr[i]);
            }
            if( sum > norm ) {
                norm = sum;
            }
        }
    }
 
    return norm;
}
 
/**
 *******************************************************************************
 *
 * @ingroup bcsc_internal
 *
 * @brief Compute the frobenius norm of a bcsc matrix.
 *
 *******************************************************************************
 *
 * @param[in] bcsc
 *          The Pastix bcsc.
 *
 *******************************************************************************
 *
 * @retval The norm of the matrix
 *
 *******************************************************************************/
float
bcsc_cnorm_frobenius( const pastix_bcsc_t *bcsc )
{
    float scale = 0.;
    float sum = 1.;
    float norm;
    float *valptr = (float*)bcsc->Lvalues;
    pastix_int_t i, j, bloc;
 
    for( bloc=0; bloc < bcsc->cscfnbr; bloc++ )
    {
        for( j=0; j < bcsc->cscftab[bloc].colnbr; j++ )
        {
            for( i = bcsc->cscftab[bloc].coltab[j]; i < bcsc->cscftab[bloc].coltab[j+1]; i++, valptr++ )
            {
                frobenius_update( 1, &scale, &sum, valptr);
#if defined(PRECISION_z) || defined(PRECISION_c)
                valptr++;
                frobenius_update( 1, &scale, &sum, valptr);
#endif
            }
        }
    }
 
    norm = scale*sqrtf(sum);
 
    return norm;
}
 
/**
 *******************************************************************************
 *
 * @brief Compute the norm of an bcsc matrix
 *
 *******************************************************************************
 *
 * @param[in] ntype
 *          = PastixMaxNorm: Max norm
 *          = PastixOneNorm: One norm
 *          = PastixInfNorm: Infinity norm
 *          = PastixFrobeniusNorm: Frobenius norm
 *
 * @param[in] bcsc
 *          The bcsc structure describing the matrix.
 *
 *******************************************************************************
 *
 * @retval The norm of the matrix.
 *
 *******************************************************************************/
float
bcsc_cnorm( pastix_normtype_t    ntype,
            const pastix_bcsc_t *bcsc )
{
    float norm = 0.;
 
    if(bcsc == NULL)
    {
        return -1.;
    }
 
    switch( ntype ) {
    case PastixMaxNorm:
        norm = bcsc_cnorm_max( bcsc );
        break;
 
    case PastixInfNorm:
        norm = bcsc_cnorm_inf( bcsc );
        break;
 
    case PastixOneNorm:
        norm = bcsc_cnorm_one( bcsc );
        break;
 
    case PastixFrobeniusNorm:
        norm = bcsc_cnorm_frobenius( bcsc );
        break;
 
    default:
        fprintf(stderr, "bcsc_cnorm: invalid norm type\n");
        return -1.;
    }
 
    return norm;
}