Mentions légales du service

Skip to content
Snippets Groups Projects
Commit 0abc7722 authored by FELŠÖCI Marek's avatar FELŠÖCI Marek
Browse files

Use native C data types

parent b5c5d947
No related branches found
No related tags found
No related merge requests found
...@@ -32,42 +32,6 @@ ...@@ -32,42 +32,6 @@
typedef enum { MPF_FALSE, MPF_TRUE } Logical; typedef enum { MPF_FALSE, MPF_TRUE } Logical;
/**************************************************************************************/
/*** Scalar Types Definitions *****************************************************/
/**************************************************************************************/
/*! \brief Simple real type.
This scalar type corresponds to the SIMPLE_PRECISION value in the ScalarType enumeration.
*/
typedef float S_type ;
/**************************************************************************************/
/*! \brief Double real type.
This scalar type corresponds to the DOUBLE_PRECISION value in the ScalarType enumeration.
*/
typedef double D_type ;
/**************************************************************************************/
/*! \brief Simple complex type.
This scalar type corresponds to the SIMPLE_COMPLEX value in the ScalarType enumeration.
*/
typedef struct {
/*! \brief real part */
float r ;
/*! \brief imaginary part */
float i ;
} C_type;
/**************************************************************************************/
/*! \brief Double complex type.
This scalar type corresponds to the DOUBLE_COMPLEX value in the ScalarType enumeration.
*/
typedef struct {
/*! \brief real part */
double r ;
/*! \brief imaginary part */
double i ;
} Z_type;
/**************************************************************************************/ /**************************************************************************************/
/*! \brief All the scalar types /*! \brief All the scalar types
...@@ -75,224 +39,24 @@ typedef struct { ...@@ -75,224 +39,24 @@ typedef struct {
They are transmitted to the vectors and matrices creation routine such as MatCreate(). They are transmitted to the vectors and matrices creation routine such as MatCreate().
*/ */
typedef enum { typedef enum {
/*! \brief Simple real type (float in C, REAL in fortran, S_type in MPF) */ /*! \brief Simple real type (float in C, REAL in fortran, float in MPF) */
SIMPLE_PRECISION=0, SIMPLE_PRECISION=0,
/*! \brief Double real type (double in C, REAL*8 in fortran, D_type in MPF) */ /*! \brief Double real type (double in C, REAL*8 in fortran, double in MPF) */
DOUBLE_PRECISION=1, DOUBLE_PRECISION=1,
/*! \brief Simple complex type (doesn't exist in C, COMPLEX in fortran, C_type in MPF) */ /*! \brief Simple complex type (doesn't exist in C, COMPLEX in fortran, float _Complex in MPF) */
SIMPLE_COMPLEX=2, SIMPLE_COMPLEX=2,
/*! \brief Double complex type (doesn't exist in C, DOUBLE COMPLEX in fortran, Z_type in MPF) */ /*! \brief Double complex type (doesn't exist in C, DOUBLE COMPLEX in fortran, double _Complex in MPF) */
DOUBLE_COMPLEX=3, DOUBLE_COMPLEX=3,
/*! \brief Number of scalar types available. */ /*! \brief Number of scalar types available. */
nbScalarType=4 nbScalarType=4
} ScalarType; } ScalarType;
/**************************************************************************************/
/*** Sparse Scalar Types Definitions **********************************************/
/**************************************************************************************/
/*! \brief Simple real sparse type.
This sparse scalar type corresponds to the position (row, column) of a non-zero element
followed by its value (S_type).
*/
typedef struct {
/*! \brief Row position */
int i ;
/*! \brief Column position */
int j ;
/*! \brief Value */
S_type v;
} SnonZero;
/**************************************************************************************/
/*! \brief Double real sparse type.
This sparse scalar type corresponds to the position (row, column) of a non-zero element
followed by its value (D_type).
*/
typedef struct {
/*! \brief Row position */
int i ;
/*! \brief Column position */
int j ;
/*! \brief Value */
D_type v;
} DnonZero;
/**************************************************************************************/
/*! \brief Simple complex sparse type.
This sparse scalar type corresponds to the position (row, column) of a non-zero element
followed by its value (C_type).
*/
typedef struct {
/*! \brief Row position */
int i ;
/*! \brief Column position */
int j ;
union {
float _Complex cv;
/*! \brief Value */
C_type v;
};
} CnonZero;
/**************************************************************************************/
/*! \brief Double complex sparse type.
This sparse scalar type corresponds to the position (row, column) of a non-zero element
followed by its value (Z_type).
*/
typedef struct {
/*! \brief Row position */
int i ;
/*! \brief Column position */
int j ;
union {
double _Complex cv;
/*! \brief Value */
Z_type v;
};
} ZnonZero;
/**************************************************************************************/
/*! \brief Name of the scalar types */ /*! \brief Name of the scalar types */
_EXTERN_TEST_FEMBEM_ char *MPF_scalName[4] _EXTERN_TEST_FEMBEM_ char *MPF_scalName[4]
#ifdef ___TEST_FEMBEM___ #ifdef ___TEST_FEMBEM___
={"SIMPLE REAL", "DOUBLE REAL", "SIMPLE COMPLEX", "DOUBLE COMPLEX"} ={"SIMPLE REAL", "DOUBLE REAL", "SIMPLE COMPLEX", "DOUBLE COMPLEX"}
#endif #endif
; ;
/* ================================================================================== */
/*! \brief Multi-arithmetics constant.
This structure is used to store any floating point constants in all of the four possible
arithmetics. This is an array of 4 void* pointers, pointing on the constant in S, D, C and Z precision
(in that order). For instance, if X is an MpfConst, then X[SIMPLE_COMPLEX] is a pointer on a C_type.
*/
typedef void* MpfConst [nbScalarType] ;
/* ================================================================================== */
/* +1 in Floating point */
/*! \brief +1 in simple real precision. */
_EXTERN_TEST_FEMBEM_ S_type Mpf_s_pone
#ifdef ___TEST_FEMBEM___
=1.
#endif
;
/*! \brief +1 in double real precision. */
_EXTERN_TEST_FEMBEM_ D_type Mpf_d_pone
#ifdef ___TEST_FEMBEM___
=1.
#endif
;
/*! \brief +1 in simple complex precision. */
_EXTERN_TEST_FEMBEM_ C_type Mpf_c_pone
#ifdef ___TEST_FEMBEM___
={1., 0.}
#endif
;
/*! \brief +1 in double complex precision. */
_EXTERN_TEST_FEMBEM_ Z_type Mpf_z_pone
#ifdef ___TEST_FEMBEM___
={1., 0.}
#endif
;
/*! \brief +1 in all precisions. */
_EXTERN_TEST_FEMBEM_ MpfConst Mpf_pone
#ifdef ___TEST_FEMBEM___
={&Mpf_s_pone, &Mpf_d_pone, &Mpf_c_pone, &Mpf_z_pone}
#endif
;
/* ================================================================================== */
/* 0 in Floating point */
/*! \brief Zero in simple real precision. */
_EXTERN_TEST_FEMBEM_ S_type Mpf_s_zero
#ifdef ___TEST_FEMBEM___
=0.
#endif
;
/*! \brief Zero in double real precision. */
_EXTERN_TEST_FEMBEM_ D_type Mpf_d_zero
#ifdef ___TEST_FEMBEM___
=0.
#endif
;
/*! \brief Zero in simple complex precision. */
_EXTERN_TEST_FEMBEM_ C_type Mpf_c_zero
#ifdef ___TEST_FEMBEM___
={0., 0.}
#endif
;
/*! \brief Zero in double complex precision. */
_EXTERN_TEST_FEMBEM_ Z_type Mpf_z_zero
#ifdef ___TEST_FEMBEM___
={0., 0.}
#endif
;
/*! \brief Zero in all precisions. */
_EXTERN_TEST_FEMBEM_ MpfConst Mpf_zero
#ifdef ___TEST_FEMBEM___
={&Mpf_s_zero, &Mpf_d_zero, &Mpf_c_zero, &Mpf_z_zero}
#endif
;
/* ================================================================================== */
/* -1 in Floating point */
/*! \brief -1 in simple real precision. */
_EXTERN_TEST_FEMBEM_ S_type Mpf_s_mone
#ifdef ___TEST_FEMBEM___
=-1.
#endif
;
/*! \brief -1 in double real precision. */
_EXTERN_TEST_FEMBEM_ D_type Mpf_d_mone
#ifdef ___TEST_FEMBEM___
=-1.
#endif
;
/*! \brief -1 in simple complex precision. */
_EXTERN_TEST_FEMBEM_ C_type Mpf_c_mone
#ifdef ___TEST_FEMBEM___
={-1., 0.}
#endif
;
/*! \brief -1 in double complex precision. */
_EXTERN_TEST_FEMBEM_ Z_type Mpf_z_mone
#ifdef ___TEST_FEMBEM___
={-1., 0.}
#endif
;
/*! \brief -1 in all precisions. */
_EXTERN_TEST_FEMBEM_ MpfConst Mpf_mone
#ifdef ___TEST_FEMBEM___
={&Mpf_s_mone, &Mpf_d_mone, &Mpf_c_mone, &Mpf_z_mone}
#endif
;
/* ================================================================================== */
/*! \brief 0 in integer type. */
_EXTERN_TEST_FEMBEM_ int Mpf_i_zero
#ifdef ___TEST_FEMBEM___
=0
#endif
;
/*! \brief 1 in integer type. */
_EXTERN_TEST_FEMBEM_ int Mpf_i_one
#ifdef ___TEST_FEMBEM___
=1
#endif
;
/* ================================================================================== */
/*! \brief Size (in bytes) of the different scalar types. */
_EXTERN_TEST_FEMBEM_ size_t Mpf_stype_size[nbScalarType]
#ifdef ___TEST_FEMBEM___
={sizeof(S_type), sizeof(D_type), sizeof(C_type), sizeof(Z_type)}
#endif
;
/*! \brief Size (in bytes) of the different non-zero scalar types. */
_EXTERN_TEST_FEMBEM_ size_t Mpf_nonzerostype_size[nbScalarType]
#ifdef ___TEST_FEMBEM___
={sizeof(SnonZero), sizeof(DnonZero), sizeof(CnonZero), sizeof(ZnonZero)}
#endif
;
#ifdef HAVE_HMAT #ifdef HAVE_HMAT
enum mpf_hmat_engine { mpf_hmat_seq, mpf_hmat_starpu, mpf_hmat_toyrt }; enum mpf_hmat_engine { mpf_hmat_seq, mpf_hmat_starpu, mpf_hmat_toyrt };
struct mpf_hmat_settings_t { struct mpf_hmat_settings_t {
...@@ -411,15 +175,6 @@ int Mpf_progressBar(int currentValue, int maximumValue) ; ...@@ -411,15 +175,6 @@ int Mpf_progressBar(int currentValue, int maximumValue) ;
int MpfArgGetDouble( int *Argc, char **argv, int rflag, const char *name, double *val ) ; int MpfArgGetDouble( int *Argc, char **argv, int rflag, const char *name, double *val ) ;
int MpfArgHasName( int *Argc, char **argv, int rflag, const char *name ); int MpfArgHasName( int *Argc, char **argv, int rflag, const char *name );
int MpfArgGetInt( int *Argc, char **argv, int rflag, const char *name, int *val ); int MpfArgGetInt( int *Argc, char **argv, int rflag, const char *name, int *val );
int S_sortAndCompact(SnonZero *mat, int *size);
int sortAndCompact(DnonZero *mat, int *size);
int C_sortAndCompact(CnonZero *mat, int *size);
int Z_sortAndCompact(ZnonZero *mat, int *size);
int Z_sortAndCompactL(ZnonZero *mat, size_t *size);
int compare_SnonZero_Line(const void *pt1, const void *pt2);
int compare_DnonZero_Line(const void *pt1, const void *pt2);
int compare_CnonZero_Line(const void *pt1, const void *pt2);
int compare_ZnonZero_Line(const void *pt1, const void *pt2);
Time get_time(void); Time get_time(void);
Time time_interval(Time t1, Time t2); Time time_interval(Time t1, Time t2);
Time add_times(Time t1, Time t2); Time add_times(Time t1, Time t2);
......
...@@ -23,9 +23,9 @@ int produitClassiqueBEM(int ffar, void *sol) { ...@@ -23,9 +23,9 @@ int produitClassiqueBEM(int ffar, void *sol) {
for (int j=0 ; j<nbPts ; j+=sparseRHS) { for (int j=0 ; j<nbPts ; j+=sparseRHS) {
int i, k, ierr; int i, k, ierr;
double coord_i[4], coord_j[3]; double coord_i[4], coord_j[3];
Z_type Aij ; double _Complex Aij ;
Z_type *z_sol=sol, *z_rhs=rhs ; double _Complex *z_sol=sol, *z_rhs=rhs ;
D_type *d_sol=sol, *d_rhs=rhs ; double *d_sol=sol, *d_rhs=rhs ;
ierr=computeCoordCylinder(j, &(coord_j[0])) ; CHKERRA(ierr); ierr=computeCoordCylinder(j, &(coord_j[0])) ; CHKERRA(ierr);
/* Boucle sur les lignes */ /* Boucle sur les lignes */
...@@ -35,12 +35,11 @@ int produitClassiqueBEM(int ffar, void *sol) { ...@@ -35,12 +35,11 @@ int produitClassiqueBEM(int ffar, void *sol) {
ierr=computeKernelBEM(&(coord_i[0]), &(coord_j[0]), i==j, (double _Complex *)&Aij) ; ierr=computeKernelBEM(&(coord_i[0]), &(coord_j[0]), i==j, (double _Complex *)&Aij) ;
if (complexALGO) { if (complexALGO) {
for (k=0 ; k<nbRHS ; k++) { for (k=0 ; k<nbRHS ; k++) {
z_sol[(size_t)k*nbPts+i].r += Aij.r * z_rhs[(size_t)k*nbPts+j].r - Aij.i * z_rhs[(size_t)k*nbPts+j].i ; z_sol[(size_t)k*nbPts+i] += Aij * z_rhs[(size_t)k*nbPts+j];
z_sol[(size_t)k*nbPts+i].i += Aij.r * z_rhs[(size_t)k*nbPts+j].i + Aij.i * z_rhs[(size_t)k*nbPts+j].r ;
} }
} else { } else {
for (k=0 ; k<nbRHS ; k++) { for (k=0 ; k<nbRHS ; k++) {
d_sol[(size_t)k*nbPts+i] += Aij.r * d_rhs[(size_t)k*nbPts+j] ; d_sol[(size_t)k*nbPts+i] += creal(Aij) * d_rhs[(size_t)k*nbPts+j] ;
} }
} }
} /* if testDofNeighbour(.. */ } /* if testDofNeighbour(.. */
...@@ -62,7 +61,7 @@ int produitClassique(void **sol) { ...@@ -62,7 +61,7 @@ int produitClassique(void **sol) {
int ierr; int ierr;
double temps_initial, temps_final, temps_cpu; double temps_initial, temps_final, temps_cpu;
void *solCLA = MpfCalloc((size_t)nbPts*nbRHS, complexALGO ? sizeof(Z_type) : sizeof(D_type)) ; CHKPTRQ(solCLA) ; void *solCLA = MpfCalloc((size_t)nbPts*nbRHS, complexALGO ? sizeof(double _Complex) : sizeof(double)) ; CHKPTRQ(solCLA) ;
temps_initial = getTime (); temps_initial = getTime ();
ierr = produitClassiqueBEM(2, solCLA) ; CHKERRQ(ierr) ; ierr = produitClassiqueBEM(2, solCLA) ; CHKERRQ(ierr) ;
......
...@@ -19,8 +19,8 @@ frobenius_update( double *scale, double *sumsq, const double value ) ...@@ -19,8 +19,8 @@ frobenius_update( double *scale, double *sumsq, const double value )
int computeRelativeError(void *sol, void *ref, double *eps) { int computeRelativeError(void *sol, void *ref, double *eps) {
double normRef[2], normDelta[2]; double normRef[2], normDelta[2];
Z_type *z_sol=sol, *z_ref=ref ; double _Complex *z_sol=sol, *z_ref=ref ;
D_type *d_sol=sol, *d_ref=ref ; double *d_sol=sol, *d_ref=ref ;
*eps=0.; *eps=0.;
normDelta[0] = 1.; normDelta[0] = 1.;
...@@ -29,11 +29,11 @@ int computeRelativeError(void *sol, void *ref, double *eps) { ...@@ -29,11 +29,11 @@ int computeRelativeError(void *sol, void *ref, double *eps) {
normRef[1] = 0.; normRef[1] = 0.;
for (size_t i=0 ; i<(size_t)nbPts*nbRHS ; i++) { for (size_t i=0 ; i<(size_t)nbPts*nbRHS ; i++) {
if (complexALGO) { if (complexALGO) {
frobenius_update( normDelta, normDelta+1, (z_sol[i].r - z_ref[i].r) ); frobenius_update( normDelta, normDelta+1, (creal(z_sol[i]) - creal(z_ref[i])) );
frobenius_update( normDelta, normDelta+1, (z_sol[i].i - z_ref[i].i) ); frobenius_update( normDelta, normDelta+1, (cimag(z_sol[i]) - cimag(z_ref[i])) );
frobenius_update( normRef, normRef+1, z_ref[i].r ); frobenius_update( normRef, normRef+1, creal(z_ref[i]) );
frobenius_update( normRef, normRef+1, z_ref[i].i ); frobenius_update( normRef, normRef+1, cimag(z_ref[i]) );
} else { } else {
frobenius_update( normDelta, normDelta+1, (d_sol[i] - d_ref[i]) ); frobenius_update( normDelta, normDelta+1, (d_sol[i] - d_ref[i]) );
frobenius_update( normRef, normRef+1, d_ref[i] ); frobenius_update( normRef, normRef+1, d_ref[i] );
...@@ -52,14 +52,14 @@ int computeRelativeError(void *sol, void *ref, double *eps) { ...@@ -52,14 +52,14 @@ int computeRelativeError(void *sol, void *ref, double *eps) {
int computeVecNorm(void *ref, double *norm) { int computeVecNorm(void *ref, double *norm) {
double normRef ; double normRef ;
Z_type *z_ref=ref ; double _Complex *z_ref=ref ;
D_type *d_ref=ref ; double *d_ref=ref ;
*norm=0.; *norm=0.;
normRef=0. ; normRef=0. ;
for (size_t i=0 ; i<(size_t)nbPts*nbRHS ; i++) { for (size_t i=0 ; i<(size_t)nbPts*nbRHS ; i++) {
if (complexALGO) { if (complexALGO) {
normRef += z_ref[i].r*z_ref[i].r+z_ref[i].i*z_ref[i].i ; normRef += z_ref[i]*z_ref[i];
} else { } else {
normRef += d_ref[i]*d_ref[i] ; normRef += d_ref[i]*d_ref[i] ;
} }
......
...@@ -300,8 +300,7 @@ void computeDenseBlockFEMBEM_Cprec(int *LigInf, int *LigSup, int *ColInf, int *C ...@@ -300,8 +300,7 @@ void computeDenseBlockFEMBEM_Cprec(int *LigInf, int *LigSup, int *ColInf, int *C
inxBloc = (*tailleL)*(ic-ColInf[0] +1) + il - LigInf[0] + 1 + (*iloc - 1) + (*jloc - 1) * (*tailleL) + inxBloc = (*tailleL)*(ic-ColInf[0] +1) + il - LigInf[0] + 1 + (*iloc - 1) + (*jloc - 1) * (*tailleL) +
(*tailleS)*id; (*tailleS)*id;
if (stype == SIMPLE_COMPLEX) { if (stype == SIMPLE_COMPLEX) {
((C_type*)bloc)[inxBloc].r = (float)(((Z_type*)zbloc)[inxBloc].r) ; ((float _Complex*)bloc)[inxBloc] = (float)(((double _Complex*)zbloc)[inxBloc]) ;
((C_type*)bloc)[inxBloc].i = (float)(((Z_type*)zbloc)[inxBloc].i) ;
} else { } else {
((float*)bloc)[inxBloc] = (float)(((double*)zbloc)[inxBloc]) ; ((float*)bloc)[inxBloc] = (float)(((double*)zbloc)[inxBloc]) ;
} }
......
...@@ -17,7 +17,7 @@ int prepareTEST(void) { ...@@ -17,7 +17,7 @@ int prepareTEST(void) {
/* Cree le second membre du produit mat-vec */ /* Cree le second membre du produit mat-vec */
printf("Computing RHS...\n") ; printf("Computing RHS...\n") ;
rhs = MpfCalloc((size_t)nbPts*nbRHS, complexALGO ? sizeof(Z_type) : sizeof(D_type)) ; CHKPTRQ(rhs) ; rhs = MpfCalloc((size_t)nbPts*nbRHS, complexALGO ? sizeof(double _Complex) : sizeof(double)) ; CHKPTRQ(rhs) ;
ierr = computeRhs() ; CHKERRQ(ierr) ; ierr = computeRhs() ; CHKERRQ(ierr) ;
leave_context(); leave_context();
...@@ -30,14 +30,14 @@ double getTime() { ...@@ -30,14 +30,14 @@ double getTime() {
int displayArray(char *s, void *f, int m, int n) { int displayArray(char *s, void *f, int m, int n) {
int i, j ; int i, j ;
Z_type *z_f=f ; //double _Complex *z_f=f ;
D_type *d_f=f ; double *d_f=f ;
for (j=0 ; j<n ; j++) for (j=0 ; j<n ; j++)
for (i=0 ; i<m ; i++) for (i=0 ; i<m ; i++)
if (complexALGO) //if (complexALGO)
printf("%s[%d, %d]=(%e, %e)\n", s, i, j, z_f[(size_t)j*m+i].r, z_f[(size_t)j*m+i].i) ; //printf("%s[%d, %d]=(%e, %e)\n", s, i, j, z_f[(size_t)j*m+i].r, z_f[(size_t)j*m+i].i) ;
else //else
printf("%s[%d, %d]=%e\n", s, i, j, d_f[(size_t)j*m+i]) ; printf("%s[%d, %d]=%e\n", s, i, j, d_f[(size_t)j*m+i]) ;
return 0 ; return 0 ;
} }
......
...@@ -6,7 +6,7 @@ int computeRhs(void) { ...@@ -6,7 +6,7 @@ int computeRhs(void) {
double coord[3], OP[3] ; double coord[3], OP[3] ;
int i, j, ierr ; int i, j, ierr ;
double complex *z_rhs=rhs ; double complex *z_rhs=rhs ;
D_type *d_rhs=rhs ; double *d_rhs=rhs ;
double small = simplePrec ? 1.e-20 : 1.e-100; double small = simplePrec ? 1.e-20 : 1.e-100;
double k=2.*M_PI/lambda; double k=2.*M_PI/lambda;
...@@ -46,14 +46,14 @@ void computeDenseBlockData (int *LigInf, int *LigSup, int *ColInf, int *ColSup, ...@@ -46,14 +46,14 @@ void computeDenseBlockData (int *LigInf, int *LigSup, int *ColInf, int *ColSup,
int jloc_____ ; int jloc_____ ;
int il,ic,inxBloc; int il,ic,inxBloc;
size_t inxGlob ; size_t inxGlob ;
S_type *s_bloc=bloc ; float *s_bloc=bloc ;
D_type *d_bloc=bloc ; double *d_bloc=bloc ;
C_type *c_bloc=bloc ; float _Complex *c_bloc=bloc ;
Z_type *z_bloc=bloc ; double _Complex *z_bloc=bloc ;
ASSERTA(context); ASSERTA(context);
contextTestFEMBEM *myCtx = (contextTestFEMBEM *)context; contextTestFEMBEM *myCtx = (contextTestFEMBEM *)context;
Z_type *z_data=(Z_type *)myCtx->dataVec ; double _Complex *z_data=(double _Complex *)myCtx->dataVec ;
D_type *d_data=(D_type *)myCtx->dataVec ; double *d_data=(double *)myCtx->dataVec ;
tailleL__ = *tailleL ; tailleL__ = *tailleL ;
iloc_____ = *iloc - 1; iloc_____ = *iloc - 1;
...@@ -66,14 +66,13 @@ void computeDenseBlockData (int *LigInf, int *LigSup, int *ColInf, int *ColSup, ...@@ -66,14 +66,13 @@ void computeDenseBlockData (int *LigInf, int *LigSup, int *ColInf, int *ColSup,
switch (stype) { switch (stype) {
case SIMPLE_PRECISION: case SIMPLE_PRECISION:
s_bloc[inxBloc] = (S_type)d_data[inxGlob]; s_bloc[inxBloc] = (float)d_data[inxGlob];
break; break;
case DOUBLE_PRECISION: case DOUBLE_PRECISION:
d_bloc[inxBloc] = d_data[inxGlob]; d_bloc[inxBloc] = d_data[inxGlob];
break; break;
case SIMPLE_COMPLEX: case SIMPLE_COMPLEX:
c_bloc[inxBloc].r = (S_type)z_data[inxGlob].r; c_bloc[inxBloc] = (float _Complex)z_data[inxGlob];
c_bloc[inxBloc].i = (S_type)z_data[inxGlob].i;
break; break;
case DOUBLE_COMPLEX: case DOUBLE_COMPLEX:
z_bloc[inxBloc] = z_data[inxGlob]; z_bloc[inxBloc] = z_data[inxGlob];
......
...@@ -48,7 +48,7 @@ int testHMAT(double * relative_error) { ...@@ -48,7 +48,7 @@ int testHMAT(double * relative_error) {
hmat_info_t info; hmat_info_t info;
hi->get_info(hmatrix, &info); hi->get_info(hmatrix, &info);
printf("Compressed size: %ld\n", info.compressed_size); printf("Compressed size: %ld\n", info.compressed_size);
printf("<PERFTESTS> AssembledSizeMb = %f \n", (double)info.compressed_size*Mpf_stype_size[stype]/1024./1024.) ; printf("<PERFTESTS> AssembledSizeMb = %f \n", (double)info.compressed_size*sizeof(double)/1024./1024.) ;
printf("Uncompressed size: %ld\n", info.uncompressed_size); printf("Uncompressed size: %ld\n", info.uncompressed_size);
printf("Ratio: %g\n", (double)info.compressed_size/info.uncompressed_size); printf("Ratio: %g\n", (double)info.compressed_size/info.uncompressed_size);
} }
...@@ -77,7 +77,7 @@ int testHMAT(double * relative_error) { ...@@ -77,7 +77,7 @@ int testHMAT(double * relative_error) {
hmat_info_t info; hmat_info_t info;
hi->get_info(hmatrix, &info); hi->get_info(hmatrix, &info);
printf("Compressed size: %ld\n", info.compressed_size); printf("Compressed size: %ld\n", info.compressed_size);
printf("<PERFTESTS> FactorizedSizeMb = %f \n", (double)info.compressed_size*Mpf_stype_size[stype]/1024./1024.) ; printf("<PERFTESTS> FactorizedSizeMb = %f \n", (double)info.compressed_size*sizeof(double)/1024./1024.) ;
printf("Uncompressed size: %ld\n", info.uncompressed_size); printf("Uncompressed size: %ld\n", info.uncompressed_size);
printf("Ratio: %g\n", (double)info.compressed_size/info.uncompressed_size); printf("Ratio: %g\n", (double)info.compressed_size/info.uncompressed_size);
} }
......
...@@ -391,241 +391,6 @@ int MpfArgGetString( int *Argc, char **argv, int rflag, const char *name, char * ...@@ -391,241 +391,6 @@ int MpfArgGetString( int *Argc, char **argv, int rflag, const char *name, char *
return 1; return 1;
} }
/* ================================================================================== */ /* ================================================================================== */
/*! \brief Sort the array mat, remove the zeros and sum the elements placed at the same position
For simple precision only (S)
\param mat the sparse matrix to sort
\param size the size (modified on output)
\return 0 for success
*/
int S_sortAndCompact(SnonZero *mat, int *size) {
int ie, il ;
/* Trie le tableau */
qsort(mat, *size, sizeof(SnonZero), compare_SnonZero_Line) ;
/* somme les valeurs situees a la meme position et enleve les zeros */
ie=-1 ; /* ie = position d'ecriture */
for (il=0 ; il<*size ; il++) /* il = position de lecture */
if (mat[il].v) {
if (ie>=0 && compare_SnonZero_Line( &(mat[ie]), &(mat[il]) ) == 0) /* s'ils sont a la meme position, on somme */
mat[ie].v += mat[il].v ;
else { /* sinon, on garde l'element 'il' en le recopiant dans une nouvelle case */
ie ++ ;
if (il>ie) mat[ie] = mat[il] ;
}
}
*size=ie+1 ;
return 0 ;
}
/* ================================================================================== */
/*! \brief Sort the array mat, remove the zeros and sum the elements placed at the same position
For double precision only (D)
\param mat the sparse matrix to sort
\param size the size (modified on output)
\return 0 for success
*/
int sortAndCompact(DnonZero *mat, int *size) {
int ie, il ;
/* Trie le tableau */
qsort(mat, *size, sizeof(DnonZero), compare_DnonZero_Line) ;
/* somme les valeurs situees a la meme position et enleve les zeros */
ie=-1 ; /* ie = position d'ecriture */
for (il=0 ; il<*size ; il++) /* il = position de lecture */
if (mat[il].v) {
if (ie>=0 && compare_DnonZero_Line( &(mat[ie]), &(mat[il]) ) == 0) /* s'ils sont a la meme position, on somme */
mat[ie].v += mat[il].v ;
else { /* sinon, on garde l'element 'il' en le recopiant dans une nouvelle case */
ie ++ ;
if (il>ie) mat[ie] = mat[il] ;
}
}
*size=ie+1 ;
return 0 ;
}
/* ================================================================================== */
/*! \brief Sort the array mat, remove the zeros and sum the elements placed at the same position
For simple complex only (C)
\param mat the sparse matrix to sort
\param size the size (modified on output)
\return 0 for success
*/
int C_sortAndCompact(CnonZero *mat, int *size) {
int ie, il ;
/* Trie le tableau */
qsort(mat, *size, sizeof(CnonZero), compare_CnonZero_Line) ;
/* somme les valeurs situees a la meme position et enleve les zeros */
ie=-1 ; /* ie = position d'ecriture */
for (il=0 ; il<*size ; il++) /* il = position de lecture */
if (mat[il].v.r || mat[il].v.i) {
if (ie>=0 && compare_CnonZero_Line( &(mat[ie]), &(mat[il]) ) == 0) { /* s'ils sont a la meme position, on somme */
mat[ie].v.r += mat[il].v.r ;
mat[ie].v.i += mat[il].v.i ;
} else { /* sinon, on garde l'element 'il' en le recopiant dans une nouvelle case */
ie ++ ;
if (il>ie) mat[ie] = mat[il] ;
}
}
*size=ie+1 ;
return 0 ;
}
/* ================================================================================== */
/*! \brief Sort the array mat, remove the zeros and sum the elements placed at the same position
For double complex only (Z)
\param mat the sparse matrix to sort
\param size the size (modified on output)
\return 0 for success
*/
int Z_sortAndCompact(ZnonZero *mat, int *size) {
int ie, il ;
/* Trie le tableau */
qsort(mat, *size, sizeof(ZnonZero), compare_ZnonZero_Line) ;
/* somme les valeurs situees a la meme position et enleve les zeros */
ie=-1 ; /* ie = position d'ecriture */
for (il=0 ; il<*size ; il++) /* il = position de lecture */
if (mat[il].v.r || mat[il].v.i) {
if (ie>=0 && compare_ZnonZero_Line( &(mat[ie]), &(mat[il]) ) == 0) { /* s'ils sont a la meme position, on somme */
mat[ie].v.r += mat[il].v.r ;
mat[ie].v.i += mat[il].v.i ;
} else { /* sinon, on garde l'element 'il' en le recopiant dans une nouvelle case */
ie ++ ;
if (il>ie) mat[ie] = mat[il] ;
}
}
*size=ie+1 ;
return 0 ;
}
/* ================================================================================== */
/*! \brief Sort the array mat, remove the zeros and sum the elements placed at the same position
For double complex only (Z) with size_t indices
\param mat the sparse matrix to sort
\param size the size (modified on output)
\return 0 for success
*/
int Z_sortAndCompactL(ZnonZero *mat, size_t *size) {
ssize_t ie, il ;
/* Trie le tableau */
qsort(mat, *size, sizeof(ZnonZero), compare_ZnonZero_Line) ;
/* somme les valeurs situees a la meme position et enleve les zeros */
ie=-1 ; /* ie = position d'ecriture */
for (il=0 ; il<*size ; il++) /* il = position de lecture */
if (mat[il].v.r || mat[il].v.i) {
if (ie>=0 && compare_ZnonZero_Line( &(mat[ie]), &(mat[il]) ) == 0) { /* s'ils sont a la meme position, on somme */
mat[ie].v.r += mat[il].v.r ;
mat[ie].v.i += mat[il].v.i ;
} else { /* sinon, on garde l'element 'il' en le recopiant dans une nouvelle case */
ie ++ ;
if (il>ie) mat[ie] = mat[il] ;
}
}
*size=ie+1 ;
return 0 ;
}
/* ================================================================================== */
/*! \brief compares two SnonZeros by lines
\param pt1 pointer on the first SnonZero
\param pt2 pointer on the second SnonZero
\return -1, 0 or +1 if \a pt1 is before, equal or after \a pt2.
*/
int compare_SnonZero_Line(const void *pt1, const void *pt2) {
SnonZero *s1=(SnonZero*)pt1 ;
SnonZero *s2=(SnonZero*)pt2 ;
if (s1->i > s2->i)
return(+1) ;
if (s1->i < s2->i)
return(-1) ;
if (s1->j > s2->j)
return(+1) ;
if (s1->j < s2->j)
return(-1) ;
return(0) ;
}
/* ================================================================================== */
/*! \brief compares two DnonZeros by lines
\param pt1 pointer on the first DnonZero
\param pt2 pointer on the second DnonZero
\return -1, 0 or +1 if \a pt1 is before, equal or after \a pt2.
*/
int compare_DnonZero_Line(const void *pt1, const void *pt2) {
DnonZero *d1=(DnonZero*)pt1 ;
DnonZero *d2=(DnonZero*)pt2 ;
if (d1->i > d2->i)
return(+1) ;
if (d1->i < d2->i)
return(-1) ;
if (d1->j > d2->j)
return(+1) ;
if (d1->j < d2->j)
return(-1) ;
return(0) ;
}
/* ================================================================================== */
/*! \brief compares two CnonZeros by lines
\param pt1 pointer on the first CnonZero
\param pt2 pointer on the second CnonZero
\return -1, 0 or +1 if \a pt1 is before, equal or after \a pt2.
*/
int compare_CnonZero_Line(const void *pt1, const void *pt2) {
CnonZero *c1=(CnonZero*)pt1 ;
CnonZero *c2=(CnonZero*)pt2 ;
if (c1->i > c2->i)
return(+1) ;
if (c1->i < c2->i)
return(-1) ;
if (c1->j > c2->j)
return(+1) ;
if (c1->j < c2->j)
return(-1) ;
return(0) ;
}
/* ================================================================================== */
/*! \brief compares two ZnonZeros by lines
\param pt1 pointer on the first ZnonZero
\param pt2 pointer on the second ZnonZero
\return -1, 0 or +1 if \a pt1 is before, equal or after \a pt2.
*/
int compare_ZnonZero_Line(const void *pt1, const void *pt2) {
ZnonZero *z1=(ZnonZero*)pt1 ;
ZnonZero *z2=(ZnonZero*)pt2 ;
if (z1->i > z2->i)
return(+1) ;
if (z1->i < z2->i)
return(-1) ;
if (z1->j > z2->j)
return(+1) ;
if (z1->j < z2->j)
return(-1) ;
return(0) ;
}
/* ================================================================================== */
Time get_time() { Time get_time() {
Time result ; Time result ;
#ifdef _WIN32 #ifdef _WIN32
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment