Mentions légales du service

Skip to content
Snippets Groups Projects
Commit b47e770c authored by hhakim's avatar hhakim
Browse files

Fix a segfault in EigTJParallelSparse C++ test. 

The bug was due to a buffer over-read in `choices' MATIO variable because of a wrong hard-coded value of J (used next as size of choices).
Now the dimensions are dynamically obtained from the .mat file.
Other change: use the proper function Mat_VarFree instead of system free() for freeing the variable.
parent 3d23639a
No related branches found
No related tags found
No related merge requests found
Checking pipeline status
Hide whitespace changes
Inline Side-by-side
misc/test/src/C++/EigTJParallelSparse.cpp.in
+ 117
− 114
View file @ b47e770c
#include "faust_EigTJParallel.h"
#include "faust_MatDense.h"
#include "faust_init_from_matio_params.h"
#include "faust_init_from_matio_core.h"
#include <complex>
using namespace Faust;
typedef @TEST_FPP@ FPP;
typedef @TEST_FPP2@ FPP2;
/***
* This test results have to be compared with what misc/test/src/Matlab/test_GivensDiagPar.m outputs.
* That's a way to validate the C++ impl. of Parallel EigTJ.
*/
int main()
{
string configFilename = "@FAUST_DATA_MAT_DIR@/test_GivensDiagParallel_Lap_U_J_choices.mat";
Faust::MatDense<FPP,Cpu> Lap;
Faust::MatSparse<FPP,Cpu> sLap;
init_faust_mat_from_matio(Lap,configFilename.c_str(),"Lap");
sLap = Lap;
int J = 5472; //TODO: should be retrieved from .mat file (not be a hard-coded value)
EigTJParallel<FPP,Cpu, FPP> algo(sLap,J,/*t=*/ Lap.getNbRow()/2,
/*verbosity */ 0, /* stoppingError */ 0.0, /* errIsRel */ true, /* enable_large_Faust */ true);
algo.compute_facts();
vector<pair<int,int>> coord_choices = algo.get_coord_choices();
mat_t *matfp;
matfp = Mat_Open(configFilename.c_str(),MAT_ACC_RDONLY);
if ( NULL == matfp ) {
fprintf(stderr,"Error opening MAT file %s", configFilename.c_str());
return EXIT_FAILURE;
}
matvar_t* cell_choices = Mat_VarRead(matfp, "choices");
int* p_choices = new int[J];
int* q_choices = new int[J];
cout << "=======================================" << endl;
if ( NULL != cell_choices) {
for (int i = 0; i < J*2; i+=2 )
{
p_choices[i/2] = (int)((double*) (cell_choices->data))[i];
q_choices[i/2] = (int)((double*) (cell_choices->data))[i+1];
}
free(cell_choices);
}
Mat_Close(matfp);
for(int j=0;j<J;j++)
printf("ite=%d (1-base index) ref. p=%d q=%d, algo. p=%d q=%d eq=%d\n", j, p_choices[j], q_choices[j], coord_choices[j].first+1, coord_choices[j].second+1, p_choices[j] == coord_choices[j].first+1 && q_choices[j] == coord_choices[j].second+1);
Faust::MatDense<FPP,Cpu> fourier_diag(Lap.getNbRow(), Lap.getNbCol());
fourier_diag.setEyes();
const vector<Faust::MatSparse<FPP,Cpu>>& givens_facts = algo.get_facts();
for(int i=givens_facts.size()-1;i>=0;i--)
givens_facts[i].multiply(fourier_diag, 'N');
cout << "fourier_diag fro norm: " << fourier_diag.norm() << endl;
const Faust::MatSparse<FPP,Cpu> D = algo.get_Dspm();
Faust::MatDense<FPP,Cpu> full_D = Faust::MatDense<FPP,Cpu>(D);
cout << "D fro norm:" << D.norm() << endl;
Faust::MatDense<FPP,Cpu> fourier_diag2 = algo.compute_fourier();
Faust::MatDense<FPP,Cpu> ordered_fourier_diag2 = algo.compute_fourier(true);
Faust::MatDense<FPP,Cpu> * ordered_fourier_diag = fourier_diag.get_cols(algo.get_ord_indices());
fourier_diag2 -= fourier_diag;
ordered_fourier_diag2 -= *ordered_fourier_diag;
cout << "norm(fourier_diag2-fourier_diag): " << fourier_diag2.norm() << endl;
cout << "norm(ordered_fourier_diag2-ordered_fourier_diag): " << ordered_fourier_diag2.norm() << endl;
//verify approx. fourier is properly computed (when ordered or not)
assert(fourier_diag2.norm()==0);
assert(ordered_fourier_diag2.norm()==0);
Faust::MatDense<FPP,Cpu> ordered_D(algo.get_Dspm(true));
cout << "orderded D fro norm: " << ordered_D.norm() << endl;
cout << "ordered_D:" << endl;
ordered_D.Display();
cout << "ordered D eles" << endl;
for(int i=0;i<ordered_D.getNbRow();i++)
{
cout << ordered_D.getData()[i*ordered_D.getNbRow()+i] << " " << full_D(i,i) << endl;
}
cout << "ordered fourier_diag fro norm: " << ordered_fourier_diag->norm() << endl;
cout << "Lap fro norm: " << Lap.norm() << endl;
Faust::MatDense<FPP,Cpu> tmp = *ordered_fourier_diag;
tmp.multiplyRight(ordered_D);
ordered_fourier_diag->transpose();
tmp.multiplyRight(*ordered_fourier_diag);
ordered_fourier_diag->transpose();
tmp -= Lap;
cout << tmp.norm()/Lap.norm() << endl;
delete []p_choices;
delete []q_choices;
return 0;
}
#include "faust_EigTJParallel.h"
#include "faust_MatDense.h"
#include "faust_init_from_matio_params.h"
#include "faust_init_from_matio_core.h"
#include <complex>
using namespace Faust;
typedef @TEST_FPP@ FPP;
typedef @TEST_FPP2@ FPP2;
/***
* This test results have to be compared with what misc/test/src/Matlab/test_GivensDiagPar.m outputs.
* That's a way to validate the C++ impl. of Parallel EigTJ.
*/
int main()
{
string configFilename = "@FAUST_DATA_MAT_DIR@/test_GivensDiagParallel_Lap_U_J_choices.mat";
Faust::MatDense<FPP,Cpu> Lap;
Faust::MatSparse<FPP,Cpu> sLap;
init_faust_mat_from_matio(Lap,configFilename.c_str(),"Lap");
sLap = Lap;
mat_t *matfp;
matfp = Mat_Open(configFilename.c_str(),MAT_ACC_RDONLY);
if ( NULL == matfp ) {
fprintf(stderr,"Error opening MAT file %s", configFilename.c_str());
return EXIT_FAILURE;
}
auto mv_info = Mat_VarReadInfo(matfp, "choices");
matvar_t* cell_choices = Mat_VarRead(matfp, "choices");
// Mat_VarPrint(cell_choices, 0);
int J = cell_choices->dims[1];
EigTJParallel<FPP,Cpu, FPP> algo(sLap,J,/*t=*/ Lap.getNbRow()/2,
/*verbosity */ 0, /* stoppingError */ 0.0, /* errIsRel */ true, /* enable_large_Faust */ true);
algo.compute_facts();
vector<pair<int,int>> coord_choices = algo.get_coord_choices();
int* p_choices = new int[J];
int* q_choices = new int[J];
cout << "=======================================" << endl;
if ( NULL != cell_choices) {
for (int i = 0; i < J*2; i+=2 )
{
p_choices[i/2] = (int)((double*) (cell_choices->data))[i];
q_choices[i/2] = (int)((double*) (cell_choices->data))[i+1];
}
Mat_VarFree(cell_choices);
}
Mat_Close(matfp);
for(int j=0;j<J;j++)
printf("ite=%d (1-base index) ref. p=%d q=%d, algo. p=%d q=%d eq=%d\n", j, p_choices[j], q_choices[j], coord_choices[j].first+1, coord_choices[j].second+1, p_choices[j] == coord_choices[j].first+1 && q_choices[j] == coord_choices[j].second+1);
Faust::MatDense<FPP,Cpu> fourier_diag(Lap.getNbRow(), Lap.getNbCol());
fourier_diag.setEyes();
const vector<Faust::MatSparse<FPP,Cpu>>& givens_facts = algo.get_facts();
for(int i=givens_facts.size()-1;i>=0;i--)
givens_facts[i].multiply(fourier_diag, 'N');
cout << "fourier_diag fro norm: " << fourier_diag.norm() << endl;
const Faust::MatSparse<FPP,Cpu> D = algo.get_Dspm();
Faust::MatDense<FPP,Cpu> full_D = Faust::MatDense<FPP,Cpu>(D);
cout << "D fro norm:" << D.norm() << endl;
Faust::MatDense<FPP,Cpu> fourier_diag2 = algo.compute_fourier();
Faust::MatDense<FPP,Cpu> ordered_fourier_diag2 = algo.compute_fourier(true);
Faust::MatDense<FPP,Cpu> * ordered_fourier_diag = fourier_diag.get_cols(algo.get_ord_indices());
fourier_diag2 -= fourier_diag;
ordered_fourier_diag2 -= *ordered_fourier_diag;
cout << "norm(fourier_diag2-fourier_diag): " << fourier_diag2.norm() << endl;
cout << "norm(ordered_fourier_diag2-ordered_fourier_diag): " << ordered_fourier_diag2.norm() << endl;
//verify approx. fourier is properly computed (when ordered or not)
assert(fourier_diag2.norm()==0);
assert(ordered_fourier_diag2.norm()==0);
Faust::MatDense<FPP,Cpu> ordered_D(algo.get_Dspm(true));
cout << "orderded D fro norm: " << ordered_D.norm() << endl;
cout << "ordered_D:" << endl;
ordered_D.Display();
cout << "ordered D eles" << endl;
for(int i=0;i<ordered_D.getNbRow();i++)
{
cout << ordered_D.getData()[i*ordered_D.getNbRow()+i] << " " << full_D(i,i) << endl;
}
cout << "ordered fourier_diag fro norm: " << ordered_fourier_diag->norm() << endl;
cout << "Lap fro norm: " << Lap.norm() << endl;
Faust::MatDense<FPP,Cpu> tmp = *ordered_fourier_diag;
tmp.multiplyRight(ordered_D);
ordered_fourier_diag->transpose();
tmp.multiplyRight(*ordered_fourier_diag);
ordered_fourier_diag->transpose();
tmp -= Lap;
cout << tmp.norm()/Lap.norm() << endl;
delete []p_choices;
delete []q_choices;
return 0;
}
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment