From b8896bd2a65e676e81a2782409cba02dba538f55 Mon Sep 17 00:00:00 2001
From: "francois.grand" <francois.grand@irstea.fr>
Date: Thu, 19 Oct 2017 12:03:16 +0200
Subject: [PATCH] =?UTF-8?q?dichotomie=20:=20-=20ajout=20du=20nom=20du=20pa?=
 =?UTF-8?q?ram=C3=A8tre=20analytique=20=C3=A0=20utiliser=20(cas=20o=C3=B9?=
 =?UTF-8?q?=20Equation()=20peut=20calculer=20plusieurs=20param=C3=A8tres)?=
 =?UTF-8?q?=20-=20modif=20de=20la=20m=C3=A9thode=20de=20calcul=20du=20sens?=
 =?UTF-8?q?=20de=20variation=20de=20la=20fonction=20-=20d=C3=A9termination?=
 =?UTF-8?q?=20de=20l'intervalle=20de=20recherche=20initial=20:=20si=20on?=
 =?UTF-8?q?=20ne=20trouve=20pas=20de=20solution=20avec=20le=20sens=20de=20?=
 =?UTF-8?q?variation=20de=20la=20fonction,=20on=20cherche=20dans=20l'autre?=
 =?UTF-8?q?=20direction=20(cas=20de=20fonctions=20non=20monotones)?=
MIME-Version: 1.0
Content-Type: text/plain; charset=UTF-8
Content-Transfer-Encoding: 8bit

---
 src/dichotomie.ts | 226 ++++++++++++++++++++--------------------------
 1 file changed, 100 insertions(+), 126 deletions(-)

diff --git a/src/dichotomie.ts b/src/dichotomie.ts
index 7ab30e06..f25f9f67 100644
--- a/src/dichotomie.ts
+++ b/src/dichotomie.ts
@@ -175,15 +175,16 @@ class SearchInterval extends Interval {
         return this.targets.intervalHasValue(targ);
     }
 
-    // isIncreasingFunction(): boolean {
-    // non fiable : certaines fonctions (coef section puiss par ex) sont globalement monotones mais peuvent inverser localement leur sens de variation
-    //     return this.targetUpper > this.targetLower;
-    // }
+    toString(): string {
+        return super.toString() + " step=" + this._step;
+    }
 }
 
 
 /**
  * calcul par dichotomie
+ * principe : y=f(x1,x2,...), on connait une valeur de y, on cherche par ex le x2 correspondant, les autres xi étant fixés.
+ * la méthode Equation() calcule analytiquement y=f(xi)
  */
 export class Dichotomie extends Debug {
     /**
@@ -196,29 +197,35 @@ export class Dichotomie extends Debug {
      */
     private _startIntervalMaxSteps = 100;
 
+    /**
+     * nom du paramètre calculé analytiquement par Equation()
+     */
+    private analyticalSymbol: string;
+
     /**
     * Construction de la classe.
     * @param nub Noeud de calcul contenant la méthode de calcul Equation
     * @param sVarCalc Nom de la variable à calculer
+    * @param targetSymbol nom du paramètre calculé analytiquement par Equation()
     */
-    constructor(private nub: Nub, private sVarCalc: string, dbg: boolean = false) {
+    constructor(private nub: Nub, private sVarCalc: string, dbg: boolean = false, targetSymbol: string = undefined) {
         super(dbg);
         this._paramX = this.nub.getParameter(this.sVarCalc);
+        if (targetSymbol == undefined)
+            this.analyticalSymbol = this.nub.getFirstAnalyticalParameter().symbol;
+        else
+            this.analyticalSymbol = targetSymbol;
     }
 
     /**
      * Valeur inconnue à rechercher
      */
-    get vX(): number {
-        return this.paramX.v;
+    private get vX(): number {
+        return this._paramX.v;
     }
 
-    set vX(vCalc: number) {
-        this.paramX.v = vCalc;
-    }
-
-    get paramX() {
-        return this._paramX;
+    private set vX(vCalc: number) {
+        this._paramX.v = vCalc;
     }
 
     /**
@@ -239,9 +246,8 @@ export class Dichotomie extends Debug {
      * Il faudra s'assurer que cette première variable correspond à la méthode de calcul la plus rapide
      */
     private Calcul(): Result {
-        let targetParam: ParamDefinition = this.nub.getFirstAnalyticalParameter();
-        let r: Result = this.nub.Equation(targetParam.symbol);
-        this.debug("dicho : Calcul(vX=" + this.sVarCalc + "=" + this.vX + ") -> " + targetParam.symbol + "=" + r.vCalc);
+        let r: Result = this.nub.Equation(this.analyticalSymbol);
+        this.debug("dicho : Calcul(vX=" + this.sVarCalc + "=" + this.vX + ") -> " + this.analyticalSymbol + "=" + r.vCalc);
         return r;
     }
 
@@ -252,68 +258,44 @@ export class Dichotomie extends Debug {
 
     /**
      * Détermine si une fonction est croissante ou décroissante.
-     * On fait ça en sondant plusieurs points sur un intervalle : on ne peut pas
-     * prendre 2 valeurs de la fonction pour 2 valeurs successives de la variable car certaines
-     * fonctions (Q=f(coef de la section puissance) par ex) ne sont pas strictement monotones (globalement
-     * monotones mais peuvent inverser localement leur sens de variation)
-     * @param d domaine de définition de la variable
+     * @param x point auquel on calcul la dérivée
+     * @param dom domaine de définition de la variable
      */
-    private isIncreasingFunction(d: ParamDomain): boolean {
-        let n = 20;
-        let variation = 0;
-        let min, max, step;
+    private isIncreasingFunction(x: number, dom: Interval): boolean {
+        let epsilon = 1e-8;
+        let bounds = new Interval(x - epsilon, x + epsilon);
+        bounds.setBounds(bounds.intersect(dom)); // au cas où l'on sorte du domaine de la variable de la fonction
 
-        switch (d.domain) {
-            case ParamDomainValue.INTERVAL:
-                min = d.minValue;
-                max = d.maxValue;
-                break;
+        let y1 = this.CalculX(bounds.min).vCalc;
+        let y2 = this.CalculX(bounds.max).vCalc;
+        return y2 > y1;
+    }
 
-            case ParamDomainValue.ANY:
-            case ParamDomainValue.NOT_NULL:
-                min = -100;
-                max = 100;
-                break;
+    /**
+     * recherche l'intervalle contenant la valeur cible
+     * @param rTarget valeur cible
+     * @param intSearch intervalle de départ
+     * @param intMax intervalle maximum (domaine de définition de la variable)
+     */
+    private searchTarget(rTarget: number, intSearch: SearchInterval, intMax: Interval) {
+        this.debug("searchTarget : Target " + rTarget + " dans " + intSearch.toString());
 
-            case ParamDomainValue.POS:
-                min = 1e-8;
-                max = 200;
+        let n = 0;
+        let ok: boolean = false;
+        do {
+            var inters: Interval = intSearch.intersect(intMax);
+            if (inters.length == 0)
                 break;
-
-            case ParamDomainValue.POS_NULL:
-                min = 0;
-                max = 200;
+            intSearch.setBounds(inters);
+            if (intSearch.hasTargetValue(rTarget)) {
+                ok = true;
                 break;
-
-            default:
-                throw "unsupported parameter domain value"
-        }
-
-        let x = min;
-        step = (max - min) / n;
-        let lastY;
-        for (let i = 0; i < n; i++) {
-            if (x == 0 && d.domain == ParamDomainValue.NOT_NULL)
-                x = 1e-8;
-            let y = this.CalculX(x).vCalc;
-
-            if (i > 0) {
-                if (y > lastY)
-                    variation++;
-                else
-                    variation--;
             }
-            lastY = y;
-
-            x += step;
-        }
-
-        if (variation == 0) {
-            let e = new ErrorMessage(ErrorCode.ERROR_DICHO_FUNCTION_VARIATION);
-            throw e;
-        }
+            intSearch.growStep(2);
+            intSearch.next();
+        } while (n++ < this._startIntervalMaxSteps);
 
-        return variation > 0;
+        return { ok, intSearch };
     }
 
     /**
@@ -322,83 +304,65 @@ export class Dichotomie extends Debug {
      * @param rInit  valeur initiale de la variable
      */
     private getStartInterval(rTarget: number, rInit: number): any {
-        let prmDom: ParamDomain = this.paramX.getDomain();
-        let min, max, step;
+        this.debug("intervalle de départ");
+
+        let prmDom: ParamDomain = this._paramX.getDomain();
+        let step = 0.001;
 
+        let min, max;
         switch (prmDom.domain) {
             case ParamDomainValue.INTERVAL:
                 min = prmDom.minValue;
                 max = prmDom.maxValue;
-                step = (max - min) / 100.0;
-                break;
-
-            case ParamDomainValue.ANY:
-                min = -Infinity;
-                max = Infinity;
-                step = 1;
                 break;
 
             case ParamDomainValue.NOT_NULL:
-                min = -Infinity;
-                max = Infinity;
-                step = 1;
                 if (rInit == 0)
                     rInit = 1e-8;
-                break;
-
-            case ParamDomainValue.POS:
-                min = 1e-8;
-                max = Infinity;
-                step = 1;
-                break;
-
-            case ParamDomainValue.POS_NULL:
-                min = 0;
-                max = Infinity;
-                step = 1;
-                break;
+            // pas de break
 
             default:
-                throw "unsupported parameter domain value";
+                let a = ParamDomain.getDefaultBounds(prmDom.domain);
+                min = a.min;
+                max = a.max;
         }
 
         let intMax: Interval = new Interval(min, max);
         intMax.checkValue(rInit);
 
-        let intSearch: SearchInterval = new SearchInterval(this, rInit, rInit + step, step);
+        // initialisation de l'intervalle de recherche
+        let intSearch: SearchInterval = new SearchInterval(this, rInit - step / 2, rInit + step / 2, step);
+        intSearch.setBounds(intSearch.intersect(intMax)); // au cas où l'on sorte du domaine de la variable de la fonction
 
-        let inc = this.isIncreasingFunction(prmDom);
+        // sens de variation de la fonction
+        let inc = this.isIncreasingFunction(rInit, intMax);
 
         let initTarget = this.CalculX(rInit).vCalc;
         // if (initTarget > rTarget && inc || initTarget < rTarget && !inc)
         if (BoolIdentity(initTarget > rTarget, inc))
-            intSearch.reverse();
+            intSearch.reverse(); // par ex, la fonction est croissante et la valeur initiale de la variable a une image par la fonction > valeur cible
 
-        let n = 0;
-        let ok: boolean = false;
-        do {
-            var inters: Interval = intSearch.intersect(intMax);
-            if (inters.length == 0)
-                break;
-            intSearch.setBounds(inters);
-            if (intSearch.hasTargetValue(rTarget)) {
-                ok = true;
-                break;
-            }
-            intSearch.growStep(2);
-            intSearch.next();
-        } while (n++ < this._startIntervalMaxSteps);
+        // on cherche dans une première direction
 
-        // intervalle trouvé
+        let a = this.searchTarget(rTarget, intSearch, intMax);
+        if (a.ok)
+            return a;
 
-        if (ok)
-            return { ok, intSearch };
+        // il se peut que la fonction ne soit pas monotone et qu'il faille chercher dans l'autre direction
+
+        let oldStepSign = intSearch.step > 0 ? 1 : -1;
+        intSearch = new SearchInterval(this, rInit + step / 2, rInit + step, step * -oldStepSign);
+        intSearch.setBounds(intSearch.intersect(intMax));  // au cas où l'on sorte du domaine de la variable de la fonction
+
+        a = this.searchTarget(rTarget, intSearch, intMax);
+        if (a.ok)
+            return a;
+
+        // gestion de l'erreur
 
         // la valeur cible de la fonction est elle trouvable ?    
 
-        let aps: string = this.nub.getFirstAnalyticalParameter().symbol;
         let m;
-
         let res: ErrorMessage;
         let errDomain = false;
         switch (prmDom.domain) {
@@ -422,28 +386,30 @@ export class Dichotomie extends Debug {
 
         if (errDomain) {
             res = new ErrorMessage(ErrorCode.ERROR_DICHO_INIT_DOMAIN);
-            res.extraVar["targetSymbol"] = aps; // symbole de la variable calculée par la fonction
+            res.extraVar["targetSymbol"] = this.analyticalSymbol; // symbole de la variable calculée par la fonction
             res.extraVar["targetValue"] = rTarget; // valeur cible pour la fonction
             res.extraVar["variableInterval"] = intMax.toString(); // intervalle de valeurs pour la variable d'entrée de la fonction
-            res.extraVar["variableSymbol"] = this.paramX.symbol; // symbole de la variable d'entrée de la fonction
+            res.extraVar["variableSymbol"] = this._paramX.symbol; // symbole de la variable d'entrée de la fonction
         }
         else {
-            if (intSearch.step < 0)
+            // if (intSearch.step < 0)
+            if (BoolIdentity(initTarget > rTarget, inc))
                 res = new ErrorMessage(ErrorCode.ERROR_DICHO_INITVALUE_HIGH);
             else
                 res = new ErrorMessage(ErrorCode.ERROR_DICHO_INITVALUE_LOW);
 
-            res.extraVar["variableSymbol"] = this.paramX.symbol; // symbole de la variable de la fonction
+            res.extraVar["variableSymbol"] = this._paramX.symbol; // symbole de la variable de la fonction
             res.extraVar["variableInitValue"] = rInit; // valeur initiale de la variable
-            res.extraVar["targetSymbol"] = aps; // symbole de la variable calculée par la fonction
+            res.extraVar["targetSymbol"] = this.analyticalSymbol; // symbole de la variable calculée par la fonction
             res.extraVar["targetValue"] = rTarget; // valeur cible pour la fonction
             res.extraVar["initTarget"] = initTarget; // valeur de la fonction pour la valeur initiale de la variable
         }
 
-        return { ok, res };
+        return { ok: false, res };
     }
 
     /**
+     * effectue récursivement la recherche dichotomique
      * @param Xmin valeur min de l'intervalle
      * @param Xmax valeur max de l'intervalle
      * @param Y valeur cible de la fonction
@@ -451,6 +417,8 @@ export class Dichotomie extends Debug {
      * @param Ymax valeur de la fonction pour Xmax
      */
     private dichoSearch(Ytarg: number, rTol: number, Xmin: number, Xmax: number, Ymin: number, Ymax: number): Result {
+        this.debug("dichoSearch : yTarget " + Ytarg + " X " + Xmin + "->" + Xmax + " tol " + rTol);
+
         let Xmid = (Xmin + Xmax) / 2;
 
         if (Math.abs(Xmin - Xmax) <= rTol)
@@ -473,6 +441,12 @@ export class Dichotomie extends Debug {
         }
     }
 
+    /**
+     * trouve la valeur x pour y=f(x) avec un y donné
+     * @param rTarget valeur de y
+     * @param rTol tolérance pour l'arrêt de la recherche
+     * @param rInit valeur initiale approximative de x
+     */
     Dichotomie(rTarget: number, rTol: number, rInit: number): Result {
         // console.log("-----");
         // for (let x = 0; x <= 1; x += 0.1)
@@ -481,7 +455,8 @@ export class Dichotomie extends Debug {
 
         // recherche de l'intervalle de départ
 
-        this.debug("intervalle de départ");
+        this.debug("Dichotomie : valeur cible de " + this.analyticalSymbol + "=" + rTarget + ", valeur initiale de " + this.sVarCalc + "=" + rInit);
+
         let r = this.getStartInterval(rTarget, rInit);
         if (r.ok)
             var interv: SearchInterval = r.intSearch;
@@ -492,7 +467,6 @@ export class Dichotomie extends Debug {
 
         // Dichotomie
 
-        this.debug("start dicho");
         return this.dichoSearch(rTarget, rTol, interv.min, interv.max, interv.targetLower, interv.targetUpper);
     }
 }
-- 
GitLab