Mentions légales du service

Skip to content
Snippets Groups Projects
Commit 99beb4d7 authored by Alban Bellot's avatar Alban Bellot
Browse files

add symetric norm

parent fb5753a9
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
spm.c
+ 39
− 0
View file @ 99beb4d7
......@@ -68,6 +68,45 @@ static int (*conversionTable[3][3][6])(pastix_spm_t*) = {
};
/**
*******************************************************************************
*
* @ingroup pastix_spm
*
* spmInit - Init the spm structure given as parameter
*
*******************************************************************************
*
* @param[in,out] spm
* The sparse matrix to init.
*
*******************************************************************************/
void
spmInit( pastix_spm_t *spm )
{
spm->mtxtype = PastixGeneral;
spm->flttype = PastixComplex64;
spm->fmttype = PastixCSC;
spm->gN = 0;
spm->n = 0;
spm->gnnz = 0;
spm->nnz = 0;
spm->gNexp = 0;
spm->nexp = 0;
spm->gnnzexp = 0;
spm->nnzexp = 0;
spm->dof = 1;
spm->dofs = NULL;
spm->colptr = NULL;
spm->rowptr = NULL;
spm->loc2glob = NULL;
spm->values = NULL;
}
/**
*******************************************************************************
......
spm.h
+ 1
− 0
View file @ 99beb4d7
......@@ -41,6 +41,7 @@ struct pastix_spm_s {
if > 0, constant degree of freedom
otherwise, irregular degree of freedom (refer to dofs) */
pastix_int_t *dofs; /*< Number of degrees of freedom per unknown (NULL, if dof > 0) */
pastix_int_t *colptr; /*< List of indirections to rows for each vertex */
pastix_int_t *rowptr; /*< List of edges for each vertex */
pastix_int_t *loc2glob; /*< Corresponding numbering from local to global */
......
z_spm_norm.c
+ 108
− 17
View file @ 99beb4d7
......@@ -139,11 +139,10 @@ z_spmMaxNorm( const pastix_spm_t *spm )
pastix_complex64_t *valptr = (pastix_complex64_t*)spm->values;
double tmp, norm = 0.;
for(i=0; i <spm->nnz; i++, valptr++) {
for(i=0; i < spm->nnzexp; i++, valptr++) {
tmp = cabs( *valptr );
norm = norm > tmp ? norm : tmp;
}
return norm;
}
......@@ -171,17 +170,22 @@ z_spmMaxNorm( const pastix_spm_t *spm )
double
z_spmInfNorm( const pastix_spm_t *spm )
{
pastix_int_t col, row, i, baseval;
pastix_int_t col, row, i, j, k, ii, jj, baseval;
pastix_int_t dofi, dofj;
pastix_int_t *dofs=spm->dofs;
pastix_complex64_t *valptr = (pastix_complex64_t*)spm->values;
double norm = 0.;
double *sumcol;
MALLOC_INTERN( sumcol, spm->gN, double );
memset( sumcol, 0, spm->gN * sizeof(double) );
MALLOC_INTERN( sumcol, spm->gNexp, double );
memset( sumcol, 0, spm->gNexp * sizeof(double) );
baseval = spmFindBase( spm );
switch( spm->fmttype ){
case PastixCSC:
/* Original */
/*
for( col=0; col < spm->gN; col++ )
{
for( i=spm->colptr[col]-baseval; i<spm->colptr[col+1]-baseval; i++ )
......@@ -189,8 +193,10 @@ z_spmInfNorm( const pastix_spm_t *spm )
sumcol[spm->rowptr[i]-baseval] += cabs( valptr[i] );
}
}
*/
/* Add the symmetric/hermitian part */
/*
if ( (spm->mtxtype == PastixHermitian) ||
(spm->mtxtype == PastixSymmetric) )
{
......@@ -201,6 +207,55 @@ z_spmInfNorm( const pastix_spm_t *spm )
sumcol[col] += cabs( valptr[i] );
}
}
}
*/
/* Dofs */
for( i=0; i < spm->n; i++ )
{
dofi = ( spm->dof > 0 ) ? spm->dof : dofs[i+1] - dofs[i];
for(k=spm->colptr[i]; k < spm->colptr[i+1]; k++)
{
j = spm->rowptr[k - baseval] - baseval;
dofj = ( spm->dof > 0 ) ? spm->dof : dofs[j+1] - dofs[j];
row = dofs[j];
for(ii=0; ii < dofi; ii++)
{
for(jj=0; jj < dofj; jj++, valptr++)
{
{
sumcol[row + jj] += cabs( *valptr );
}
}
}
}
}
/* Add the symmetric/hermitian part */
if ( (spm->mtxtype == PastixHermitian) ||
(spm->mtxtype == PastixSymmetric) )
{
valptr = (pastix_complex64_t*)spm->values;
for(i=0; i < spm->n; i++)
{
col = dofs[i];
dofi = ( spm->dof > 0 ) ? spm->dof : dofs[i+1] - dofs[i];
for(k=spm->colptr[i]; k < spm->colptr[i+1]; k++)
{
j = spm->rowptr[k - baseval] - baseval;
dofj = ( spm->dof > 0 ) ? spm->dof : dofs[j+1] - dofs[j];
for(ii=0; ii < dofi; ii++)
for(jj=0; jj < dofj; jj++,valptr++)
{
{
if(i != j)
{
sumcol[col + ii] += cabs( *valptr );
}
}
}
}
}
}
break;
......@@ -250,7 +305,7 @@ z_spmInfNorm( const pastix_spm_t *spm )
return PASTIX_ERR_BADPARAMETER;
}
for( i=0; i<spm->gN; i++)
for( i=0; i<spm->gNexp; i++) //gNexp
{
if(norm < sumcol[i])
{
......@@ -286,34 +341,70 @@ z_spmInfNorm( const pastix_spm_t *spm )
double
z_spmOneNorm( const pastix_spm_t *spm )
{
pastix_int_t col, row, i, baseval;
pastix_int_t col, row, i, j, k, ii, jj, baseval;
pastix_int_t dofi, dofj;
pastix_complex64_t *valptr = (pastix_complex64_t*)spm->values;
double norm = 0.;
double *sumrow;
pastix_int_t* dofs=spm->dofs;
MALLOC_INTERN( sumrow, spm->gN, double );
memset( sumrow, 0, spm->gN * sizeof(double) );
MALLOC_INTERN( sumrow, spm->gNexp, double );
memset( sumrow, 0, spm->gNexp * sizeof(double) );
baseval = spmFindBase( spm );
switch( spm->fmttype ){
case PastixCSC:
for( col=0; col < spm->gN; col++ )
/*
for( col=0; col < spm->gN; col++ )
{
for( i=spm->colptr[col]-baseval; i<spm->colptr[col+1]-baseval; i++ )
{
sumrow[col] += cabs( valptr[i] );
}
}
*/
for(i=0; i < spm->n; i++)
{
for( i=spm->colptr[col]-baseval; i<spm->colptr[col+1]-baseval; i++ )
col = dofs[i];
dofi = ( spm->dof > 0 ) ? spm->dof : dofs[i+1] - spm->dofs[i];
for(k=spm->colptr[i]; k < spm->colptr[i+1]; k++)
{
sumrow[col] += cabs( valptr[i] );
j = spm->rowptr[k - baseval] - baseval;
dofj = ( spm->dof > 0 ) ? spm->dof : dofs[j+1] - spm->dofs[j];
for(ii=0; ii < dofi; ii++)
for(jj=0; jj < dofj; jj++, valptr++)
{
{
sumrow[col + ii] += cabs( *valptr );
}
}
}
}
/* Add the symmetric/hermitian part */
if ( (spm->mtxtype == PastixHermitian) ||
(spm->mtxtype == PastixSymmetric) )
{
for( col=0; col < spm->gN; col++ )
valptr = (pastix_complex64_t*)spm->values;
for(i=0; i < spm->n; i++)
{
for( i=spm->colptr[col]-baseval+1; i<spm->colptr[col+1]-baseval; i++ )
dofi = ( spm->dof > 0 ) ? spm->dof : dofs[i+1] - spm->dofs[i];
for(k=spm->colptr[i]; k < spm->colptr[i+1]; k++)
{
sumrow[spm->rowptr[i]-baseval] += cabs( valptr[i] );
j = spm->rowptr[k - baseval] - baseval;
row = dofs[j];
dofj = ( spm->dof > 0 ) ? spm->dof : dofs[j+1] - spm->dofs[j];
for(ii=0; ii < dofi; ii++)
for(jj=0; jj < dofj; jj++, valptr++)
{
{
if(i != j)
{
sumrow[row + jj] += cabs( *valptr );
}
}
}
}
}
}
......@@ -365,7 +456,7 @@ z_spmOneNorm( const pastix_spm_t *spm )
return PASTIX_ERR_BADPARAMETER;
}
for( i=0; i<spm->gN; i++)
for( i=0; i<spm->gNexp; i++)
{
if(norm < sumrow[i])
{
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment