From 8b05e840e80de29f5fcf917dd6b248440f30c1c3 Mon Sep 17 00:00:00 2001
From: Amaury Johnan <amjohnen@gmail.com>
Date: Thu, 10 Apr 2014 10:50:13 +0000
Subject: [PATCH] clean up
---
Numeric/MetricBasis.cpp | 255 +---------------------------------------
Numeric/MetricBasis.h | 3 +-
2 files changed, 3 insertions(+), 255 deletions(-)
diff --git a/Numeric/MetricBasis.cpp b/Numeric/MetricBasis.cpp
index 1f3b277870..1e2e22e192 100644
--- a/Numeric/MetricBasis.cpp
+++ b/Numeric/MetricBasis.cpp
@@ -10,14 +10,6 @@
#include <cmath>
#include <queue>
#include "OS.h"
-#if defined(HAVE_SCIP)
-#include "scip/scip.h"
-#include "scip/scipdefplugins.h"
-#include "soplex/soplex.h"
-#include <sstream>
-#include "scip/cons.h"
-#include "scip/struct_scip.h"
-#endif
namespace {
double cubicCardanoRoot(double p, double q)
@@ -122,7 +114,6 @@ void MetricCoefficient::_fillInequalities(int metricOrder)
Msg::Info("ncoeff %d", ncoeff);
int countP3 = 0, countJ2 = 0;
double tol = .99;
- double tol2 = .001;
for (int i = 0; i < ncoeff; i++) {
for (int j = i; j < ncoeff; j++) {
double A = 1, B = 1;
@@ -251,19 +242,15 @@ void MetricCoefficient::_fillInequalities(int metricOrder)
}
}
- int ret = _lightenInequalities(tol, tol2, countJ2, countP3);
- if (ret) {
- Msg::Warning("Heavy Inequalities (error %d)", ret);
- }
+ _lightenInequalities(tol, countJ2, countP3);
Msg::Info("A : %d / %d", v.size(), ncoeff*(ncoeff+1)/2);
Msg::Info("J2 : %d / %d", countJ2, njac*(njac+1)/2);
Msg::Info("P3 : %d / %d", countP3, ncoeff*(ncoeff+1)*(ncoeff+2)/6);
}
-int MetricCoefficient::_lightenInequalities(double tol, double tol2, int &countj, int &countp)
+void MetricCoefficient::_lightenInequalities(double tol, int &countj, int &countp)
{
-#ifndef BEFORE
std::map<int, std::vector<IneqData> >::iterator it, itbeg[3], itend[3];
int cnt[3] = {0,0,0};
@@ -293,206 +280,6 @@ int MetricCoefficient::_lightenInequalities(double tol, double tol2, int &countj
}
countj -= cnt[0];
countp -= cnt[1];
- return 0;
-
-#else
- std::map<int, std::vector<IneqData> >::iterator itJ, itP;
- itJ = _ineqJ2.begin();
- itP = _ineqP3.begin();
-
- while (itJ != _ineqJ2.end() && itP != _ineqP3.end()) {
- double verif = 0;
- for (unsigned int i = 0; i < itJ->second.size(); ++i) {
- verif += itJ->second[i].val;
- }
- if (verif < 1. - 1e-7 || verif > 1. + 1e-7) Msg::Warning("verif J %g", verif);
- verif = 0;
- for (unsigned int i = 0; i < itP->second.size(); ++i) {
- verif += itP->second[i].val;
- }
- if (verif < 1.-1e-7 || verif > 1.+1e-7) Msg::Warning("verif P %g", verif);
- // TODO : remove
-
- std::sort(itJ->second.begin(), itJ->second.end(), gterIneq());
- std::sort(itP->second.begin(), itP->second.end(), gterIneq());
-
- bool print = false;
- if (true || itJ->second.back().val < tol && itP->second.back().val < tol) {
- print = true;
- Msg::Info(" ");
- //Msg::Info("sizes %d %d", itJ->second.size(), itP->second.size());
- Msg::Info("JJ-before-JJJJJJJJJJJJJJJJJJJ");
- for (unsigned int i = 0; i < itJ->second.size(); ++i) {
- Msg::Info("%g", itJ->second[i].val / tol);
- }
- Msg::Info("PP-before-PPPPPPPPPPPPPPPPPPP");
- for (unsigned int i = 0; i < itP->second.size(); ++i) {
- Msg::Info("%g", itP->second[i].val / tol);
- }
- }
-
-#if defined(HAVE_SCIP)
- SCIP* scip;
- SCIP_CALL(SCIPcreate(&scip));
- SCIP_CALL(SCIPincludeDefaultPlugins(scip));
- //SCIP_CALL(SCIPsetMessagehdlr(NULL, NULL));
- SCIP_CALL(SCIPcreateProb(scip, "inequalities", NULL, NULL,
- NULL, NULL, NULL, NULL, NULL));
- SCIP_CALL(SCIPsetObjsense(scip, SCIP_OBJSENSE_MAXIMIZE));
-
- int nj = itJ->second.size();
- int np = itP->second.size();
-
- SCIP_VAR** varj = new SCIP_VAR*[nj];
- SCIP_VAR** varp = new SCIP_VAR*[np];
- double* cj = new double[nj];
- double* cp = new double[np];
-
- for (int i = 0; i < nj; ++i) {
- std::ostringstream namebuf;
- namebuf.str("");
- namebuf << "j" << i;
-
- SCIP_CALL(SCIPcreateVar(scip, &varj[i], namebuf.str().c_str(), .0, 1.,
- 1., SCIP_VARTYPE_BINARY, TRUE, FALSE,
- NULL, NULL, NULL, NULL, NULL));
- SCIP_CALL(SCIPaddVar(scip, varj[i]));
- cj[i] = itJ->second[i].val;
- }
- for (int i = 0; i < np; ++i) {
- std::ostringstream namebuf;
- namebuf.str("");
- namebuf << "p" << i;
-
- SCIP_CALL(SCIPcreateVar(scip, &varp[i], namebuf.str().c_str(), .0, 1.,
- 1., SCIP_VARTYPE_BINARY, TRUE, FALSE,
- NULL, NULL, NULL, NULL, NULL));
- SCIP_CALL(SCIPaddVar(scip, varp[i]));
- cp[i] = itP->second[i].val;
- }
-
- SCIP_CONS* cons[3];
- for (int i = 0; i < nj; ++i) {
- Msg::Info("+cj%d (%g)", i, cj[i]);
- }
-
- Msg::Info("<= tol %g", tol);
- SCIP_CALL(SCIPcreateConsLinear(scip, &cons[0], "bound J",
- nj, varj, cj, .0, tol, TRUE,
- TRUE, TRUE, TRUE, TRUE, FALSE,
- FALSE, FALSE, FALSE, FALSE));
- SCIPconsPrint(cons[0], scip->set, NULL, NULL);
- SCIP_CALL(SCIPcreateConsLinear(scip, &cons[1], "bound P",
- np, varp, cp, .0, tol, TRUE,
- TRUE, TRUE, TRUE, TRUE, FALSE,
- FALSE, FALSE, FALSE, FALSE));
- SCIP_CALL(SCIPcreateConsLinear(scip, &cons[2], "difference",
- nj, varj, cj, -tol2, tol2, TRUE,
- TRUE, TRUE, TRUE, TRUE, FALSE,
- FALSE, FALSE, FALSE, FALSE));
- for (int i = 0; i < np; ++i) {
- //SCIP_CALL(SCIPaddCoefLinear(scip, cons[2], varp[i], -cp[i]));
- }
-
- SCIP_CALL(SCIPsolve(scip));
- SCIPconsPrint(cons[0], scip->set, NULL, NULL);
- SCIP_SOL* sol = SCIPgetBestSol(scip);
- if (!sol) {
- Msg::Error("Scip did not find a solution");
- }
- else {
- int ind = 0;
- for (int i = 0; i < nj; ++i) {
- Msg::Info("j%d = %g", i, SCIPgetSolVal(scip, sol, varj[i]));
- if (SCIPgetSolVal(scip, sol, varj[i]) > .5) {
- itJ->second[ind] = itJ->second.back();
- itJ->second.pop_back();
- }
- else ++ind;
- }
- ind = 0;
- for (int i = 0; i < np; ++i) {
- if (SCIPgetSolVal(scip, sol, varp[i]) > .5) {
- itP->second[ind] = itP->second.back();
- itP->second.pop_back();
- }
- else ++ind;
- }
- }
-
- for (int i = 0; i < nj; ++i)
- SCIP_CALL(SCIPreleaseVar(scip, &varj[i]));
- for (int i = 0; i < np; ++i)
- SCIP_CALL(SCIPreleaseVar(scip, &varp[i]));
- //SCIP_CALL(SCIPreleaseCons(scip, &cons[0]));
- //SCIP_CALL(SCIPreleaseCons(scip, &cons[1]));
- //SCIP_CALL(SCIPreleaseCons(scip, &cons[2]));
- SCIP_CALL(SCIPfree(&scip));
- delete[] varj;
- delete[] varp;
- delete[] cj;
- delete[] cp;
-
-#else
-
- double rmvedJ = 0, rmvedP = 0;
- bool changed = true;
-
- while (changed) {
- changed = false;
- double potentialJ = rmvedJ, potentialP = rmvedP;
- int j = itJ->second.size();
- int p = itP->second.size();
- bool found = false;
-
- while (j > 1 && !found && potentialJ <= tol) {
- double tmp = itP->second[p-1].val;
- //if (potentialP - potentialJ > -tol2) found = true;
- while (potentialP + tmp <= tol && potentialP + tmp - potentialJ <= tol2) {
- potentialP += tmp;
- --p;
- if (p < 1) Msg::Fatal("me");
- tmp = itP->second[p-1].val;
- if (potentialP - potentialJ >= -tol2) found = true;
- }
- if (!found) {
- potentialJ += itJ->second[j-1].val;
- --j;
- if (j < 1) Msg::Fatal("ma");
- potentialP = rmvedP;
- p = itP->second.size();
- }
- }
- if (found) {
- rmvedJ = potentialJ, rmvedP = potentialP;
- //countJ2 -= itJ->second.size() - j;
- //countP3 -= itP->second.size() - p; // TODO
- itJ->second.erase(itJ->second.begin() + j, itJ->second.end());
- //Msg::Info("%d %d", p, itP->second.size());
- itP->second.erase(itP->second.begin() + p, itP->second.end());
- changed = true;
- }
- }
-#endif
-
-
- if (print) {
- //Msg::Info("sizes %d %d", itJ->second.size(), itP->second.size());
- Msg::Info("JJ-after-JJJJJJJJJJJJJJJJJJJ");
- for (unsigned int i = 0; i < itJ->second.size(); ++i) {
- Msg::Info("%g", itJ->second[i].val / tol);
- }
- Msg::Info("PP-after-PPPPPPPPPPPPPPPPPPP");
- for (unsigned int i = 0; i < itP->second.size(); ++i) {
- Msg::Info("%g", itP->second[i].val / tol);
- }
- }
-
- ++itJ;
- ++itP;
- }
- return 0;
-#endif
}
void MetricCoefficient::getCoefficients(fullMatrix<double> &coeff, bool bezier)
@@ -1170,44 +957,6 @@ double MetricCoefficient::_maxq(const fullMatrix<double> &coeff) const
return ans;
}
-void MetricCoefficient::_minMaxA(
- const fullMatrix<double> &coeff, double &min, double &max) const
-{
- int i = _inequations(0, 0);
- int j = _inequations(0, 1);
-
- // TODO : si coeff(0, m) < 0
-
- min = max = 0;
- for (int l = 1; l < 7; ++l) {
- min += coeff(i, l) * coeff(j, l);
- }
- min /= coeff(i, 0) * coeff(j, 0);
- // TODO : si coeff(0, i|j) == 0
- max = min;
-
- for (int k = 1; k < _inequations.size1(); ++k) {
- double tmp = 0;
- i = _inequations(k, 0);
- j = _inequations(k, 1);
- for (int l = 1; l < 7; ++l) {
- tmp += coeff(i, l) * coeff(j, l);
- }
- tmp /= coeff(i, 0) * coeff(j, 0);
- min = std::min(tmp, min);
- max = std::max(tmp, max);
- }
-
- if (min < 0)
- Msg::Fatal("minA ne devrait pas etre negatif");
-
- min = std::sqrt(min);
- max = std::sqrt(max);
- double tmp = min;
- min = 1./max;
- max = 1./tmp;
-}
-
void MetricCoefficient::_minMaxA2(
const fullMatrix<double> &coeff, double &min, double &max) const
{
diff --git a/Numeric/MetricBasis.h b/Numeric/MetricBasis.h
index 2ea5729fdb..56f9229dc9 100644
--- a/Numeric/MetricBasis.h
+++ b/Numeric/MetricBasis.h
@@ -72,7 +72,7 @@ private:
private:
void _computeBezCoeff();
void _fillInequalities(int order);
- int _lightenInequalities(double, double, int&, int&); //TODO change
+ void _lightenInequalities(double, int&, int&); //TODO change
void _interpolateBezierPyramid(const double *uvw, double *minmaxQ);
void _computeRmin1(fullMatrix<double>&, fullVector<double>&,
double &RminLag, double &RminBez, int) const;
@@ -83,7 +83,6 @@ private:
double _minq(const fullMatrix<double>&) const;
double _maxp(const fullMatrix<double>&) const;
double _maxq(const fullMatrix<double>&) const;
- void _minMaxA(const fullMatrix<double>&, double &min, double &max) const;
void _minMaxA2(const fullMatrix<double>&, double &min, double &max) const;
void _minMaxJacobianSqr(const fullVector<double>&, double &min, double &max) const;
void _minJ2P3(const fullMatrix<double>&, const fullVector<double>&, double &min) const;
--
GitLab