Commit b26540ea authored by Laurent Belcour's avatar Laurent Belcour

Correction of a parallel issue in the rational fitting.

parent c85342d1
......@@ -214,7 +214,8 @@ function* plugins_manager::get_function(const std::string& filename)
// Parameters of the function object
int nX, nY;
params::input param_in; params::output param_out;
params::input param_in = params::UNKNOWN_INPUT;
params::output param_out = params::UNKNOWN_OUTPUT;
arguments args;
// Test for the first line of the file. Should be a ALTA FUNC HEADER
......@@ -331,7 +332,7 @@ function* plugins_manager::get_function(const arguments& args)
}
else
{
std::string filename = args["func"];
const std::string filename = args["func"];
FunctionPrototype myFunction = open_library<FunctionPrototype>(filename, "provide_function");
if(myFunction != NULL)
{
......
......@@ -122,14 +122,16 @@ class quadratic_program
QuadProgPP::Vector<double> ci(0.0, n) ;
QuadProgPP::Matrix<double> CE(0.0, m, 0) ;
QuadProgPP::Vector<double> ce(0.0, 0) ;
/*
// Update the ci column with the delta parameter
// (See Celis et al. 2007 p.12)
Eigen::JacobiSVD<Eigen::MatrixXd, Eigen::HouseholderQRPreconditioner> svd(Eigen::MatrixXd::Map(&CI[0][0], m, n));
const double sigma_m = svd.singularValues()(std::min(m, n)-1) ;
const double sigma_M = svd.singularValues()(0) ;
delta = sigma_M / sigma_m ;
/*/
delta = 1.0;
//*/
// Select the size of the result vector to
// be equal to the dimension of p + q
for(int i=0; i<m; ++i)
......
......@@ -77,7 +77,7 @@ bool rational_fitter_parallel::fit_data(const data* dat, function* fit, const ar
int nb_sol_found = 0;
int nb_sol_tested = 0;
#pragma omp parallel for shared(nb_sol_found, nb_sol_tested, min_delta, mean_delta), schedule(dynamic,1)
#pragma omp parallel for shared(args, nb_sol_found, nb_sol_tested, min_delta, mean_delta), schedule(dynamic,1)
for(int j=1; j<i; ++j)
{
// Compute the number of coefficients in the numerator and in the denominator
......@@ -90,18 +90,22 @@ bool rational_fitter_parallel::fit_data(const data* dat, function* fit, const ar
// Allocate a rational function and set it to the correct size, dimensions
// and parametrizations.
rational_function* rk = dynamic_cast<rational_function*>(plugins_manager::get_function(args));
rational_function* rk = NULL;
#pragma omp critical (args)
{
rk = dynamic_cast<rational_function*>(plugins_manager::get_function(args));
}
if(rk == NULL)
{
std::cerr << "<<ERROR>> unable to obtain a rational function from the plugins manager" << std::endl;
throw;
}
rk->setParametrization(r->input_parametrization());
rk->setParametrization(r->output_parametrization());
rk->setDimX(r->dimX()) ;
rk->setDimY(r->dimY()) ;
rk->setMin(r->min()) ;
rk->setMax(r->max()) ;
if(rk == NULL)
{
std::cerr << "<<ERROR>> unable to obtain a rational function from the plugins manager" << std::endl;
throw;
}
// Set the rational function size
rk->setSize(temp_np, temp_nq);
......@@ -135,7 +139,8 @@ bool rational_fitter_parallel::fit_data(const data* dat, function* fit, const ar
}
}
delete rk; // memory clean
if(rk != NULL)
delete rk; // memory clean
#pragma omp critical (nb_sol_tested)
{
......
......@@ -29,11 +29,11 @@ rational_function_chebychev_1d::rational_function_chebychev_1d()
setDimY(0);
}
rational_function_chebychev_1d::rational_function_chebychev_1d(int np, int nq) :
rational_function_chebychev_1d::rational_function_chebychev_1d(int nX, int nY, int np, int nq) :
rational_function_1d(np, nq)
{
setDimX(0);
setDimY(0);
setDimX(nX);
setDimY(nY);
this->resize(np, nq);
}
......@@ -115,9 +115,7 @@ rational_function_1d* rational_function_chebychev::get(int i)
{
if(rs[i] == NULL)
{
rs[i] = new rational_function_chebychev_1d(np, nq);
rs[i]->setDimX(dimX());
rs[i]->setDimY(dimY());
rs[i] = new rational_function_chebychev_1d(dimX(), dimY(), np, nq);
// Test if the input domain is not empty. If one dimension of
// the input domain is a point, I manually inflate this dimension
......
......@@ -17,7 +17,7 @@ class rational_function_chebychev_1d : public rational_function_1d
public: // methods
rational_function_chebychev_1d() ;
rational_function_chebychev_1d(int np, int nq) ;
rational_function_chebychev_1d(int nX, int nY, int np, int nq) ;
rational_function_chebychev_1d(const vec& a, const vec& b) ;
virtual ~rational_function_chebychev_1d() {}
......@@ -52,17 +52,19 @@ public: // methods
// Update the function
virtual void update(const vec& in_a,
const vec& in_b)
const vec& in_b)
{
rational_function_1d::update(in_a, in_b);
// Get the size of the input vector
const int np = in_a.size();
const int nq = in_b.size();
for(int k=0; k<np; ++k)
{
_p_coeffs[k].a = in_a[k];
}
const int nq = in_b.size();
for(int k=0; k<nq; ++k)
{
_q_coeffs[k].a = in_b[k];
......
Markdown is supported
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