🌊

JavaScriptで数式微分を作る(2)

2022/12/24に公開約12,500字

前回

前回はボトムアップ型自動微分を作りました。動くところまで作れました。
https://zenn.dev/ihashiguchi/articles/35f2a3ab8fbf7b

今回やること

微分したものをテキストで書きだすことで微分した式を得られるか実験します。

前回との変更点

BUAD.js diff
@@ -15,15 +15,23 @@ class BUAD {
      * @param {numbers} v
      * @param {numbers} dv
      * @param {Decimal.Config|undefined} decimalConfig
+     * @param {string|undefined} valString
+     * @param {string|undefined} dvalString
      */
-    constructor(v = 0, dv = 0, decimalConfig) {
+    constructor(v = 0, dv = 0, decimalConfig, valString, dvalString) {
         this.setBUADDecimal(decimalConfig);
         this.val = v;
         this.dval = dv;
+        /** @type {string} */
+        this.valString = valString ? valString : v.toString();
+        /** @type {string} */
+        this.dvalString = dvalString ? dvalString : '0';
     }

     selectX() {
         this.dval = new this.BUADDecimal("1");
+        this.valString = 'x';
+        this.dvalString = '1';
     }

     /** @returns {Decimal} */
@@ -78,14 +86,16 @@ class BUAD {
      *
      * @param {numbers} v
      * @param {numbers} dv
+     * @param {string} valString
+     * @param {string} dvalString
      * @returns {BUAD}
      */
-    create(v, dv) {
-        return new BUAD(v, dv, this.decimalConfig);
+    create(v, dv, valString, dvalString) {
+        return new BUAD(v, dv, this.decimalConfig, valString, dvalString);
     }

     copy() {
-        return this.create(new this.BUADDecimal(this.val), new this.BUADDecimal(this.dval));
+        return this.create(new this.BUADDecimal(this.val), new this.BUADDecimal(this.dval), this.decimalConfig, this.valString, this.dvalString);
     }

     toString() {
@@ -100,7 +110,7 @@ class BUAD {
      */
     add(x = undefined) {
         if (x === undefined) return this.copy();
-        return this.create(this.val.add(x.val), this.dval.add(x.dval));
+        return this.create(this.val.add(x.val), this.dval.add(x.dval), `(${this.valString})+(${x.valString})`, `(${this.dvalString})+(${x.dvalString})`);
     }

     /**
@@ -109,8 +119,8 @@ class BUAD {
      * @returns {BUAD}
      */
     sub(x = undefined) {
-        if (x === undefined) return this.create(this.val.negated(), this.dval.negated());
-        return this.create(this.val.sub(x.val), this.dval.sub(x.dval));
+        if (x === undefined) return this.create(this.val.negated(), this.dval.negated(), `-(${this.valString})`, `-(${this.dvalString})`);
+        return this.create(this.val.sub(x.val), this.dval.sub(x.dval), `(${this.valString})-(${x.valString})`, `(${this.dvalString})-(${x.dvalString})`);
     }

     /**
@@ -127,7 +137,7 @@ class BUAD {
      * @returns {BUAD}
      */
     mul(x) {
-        return this.create(this.val.mul(x.val), this.dval.mul(x.val).add(this.val.mul(x.dval)));
+        return this.create(this.val.mul(x.val), this.dval.mul(x.val).add(this.val.mul(x.dval)), `(${this.valString})*(${x.valString})`, `(${this.dvalString})*(${x.valString})+(${this.valString})*(${x.dvalString})`);
     }

     /**
@@ -136,7 +146,7 @@ class BUAD {
      * @returns {BUAD}
      */
     div(x) {
-        return this.create(this.val.div(x.val), (this.dval.mul(x.val).sub(this.val.mul(x.dval))).div(x.val.pow(2)));
+        return this.create(this.val.div(x.val), (this.dval.mul(x.val).sub(this.val.mul(x.dval))).div(x.val.pow(2)), `(${this.valString})/(${x.valString})`, `((${this.dvalString})*(${x.valString})-(${this.valString})*(${x.dvalString}))/(${x.valString})^2`);
     }

     // Comparison functions
@@ -183,7 +193,7 @@ class BUAD {
      */
     sqrt() {
         const sqrtX = this.val.sqrt();
-        return this.create(sqrtX, this.dval.div(new this.BUADDecimal("2").mul(sqrtX)));
+        return this.create(sqrtX, this.dval.div(new this.BUADDecimal("2").mul(sqrtX)), `sqrt(${this.valString})`, `(${this.dvalString})/(2*sqrt(${this.valString}))`);
     }

     /**
@@ -194,9 +204,13 @@ class BUAD {
      */
     pow(y) {
         const xPym1 = this.val.pow(y.val.sub(new this.BUADDecimal(1)));
+        const xPym1Str = `(${this.valString})^{(${y.valString})-1}`;
         const yxP = y.val.mul(this.dval);
-        const xyPlnx = (this.val.comparedTo(new Decimal(0)) === 0 || y.dval.comparedTo(new Decimal(0)) === 0) ? new Decimal(0) : this.val.mul(y.dval).mul(this.val.ln());
-        return this.create(this.val.pow(y.val), xPym1.mul(yxP.add(xyPlnx)));
+        const yxPStr = `(${y.valString})*(${this.dvalString})`;
+        const xyPlnx = (this.val.comparedTo(new this.BUADDecimal(0)) === 0 || y.dval.comparedTo(new this.BUADDecimal(0)) === 0) ? new this.BUADDecimal(0) : this.val.mul(y.dval).mul(this.val.ln());
+        const xyPlnxStr = y.dvalString === '0' ? '0' : `(${this.valString})*(${y.dvalString})*(ln(${this.valString}))`;
+        const string = `(${xPym1Str})*((${yxPStr})+(${xyPlnxStr}))`;
+        return this.create(this.val.pow(y.val), xPym1.mul(yxP.add(xyPlnx)), `(${this.valString})^{${y.valString}}`, string);
     }

     /**
@@ -205,7 +219,7 @@ class BUAD {
      */
     exp() {
         const expX = this.val.exp();
-        return this.create(expX, this.dval.mul(expX));
+        return this.create(expX, this.dval.mul(expX), `e^{${this.valString}}`, `(${this.dvalString})*e^{${this.valString}}`);
     }

     /**
@@ -213,7 +227,7 @@ class BUAD {
      * @returns {BUAD}
      */
     log() {
-        return this.create(this.val.ln(), this.dval.div(this.val));
+        return this.create(this.val.ln(), this.dval.div(this.val), `ln(${this.valString})`, `(${this.dvalString})/(${this.valString})`);
     }

     /**
@@ -221,7 +235,7 @@ class BUAD {
      * @returns {BUAD}
      */
     sin() {
-        return this.create(this.val.sin(), this.dval.mul(this.val.cos()));
+        return this.create(this.val.sin(), this.dval.mul(this.val.cos()), `sin(${this.valString})`, `(${this.dvalString})*cos(${this.valString})`);
     }

     /**
@@ -229,7 +243,7 @@ class BUAD {
      * @returns {BUAD}
      */
     cos() {
-        return this.create(this.val.cos(), this.dval.mul(this.val.sin().negated()));
+        return this.create(this.val.cos(), this.dval.mul(this.val.sin().negated()), `cos(${this.valString})`, `(${this.dvalString})*-sin(${this.valString})`);
     }

     /**
@@ -237,7 +251,7 @@ class BUAD {
      * @returns {BUAD}
      */
     asin() {
-        return this.create(this.val.asin(), this.dval.div(new this.BUADDecimal(1).sub(this.val.pow(2))));
+        return this.create(this.val.asin(), this.dval.div(new this.BUADDecimal(1).sub(this.val.pow(2))), `asin(${this.valString})`, `(${this.dvalString})/sqrt(1-(${this.valString})^2)`);
     }

     /**
@@ -245,7 +259,7 @@ class BUAD {
      * @returns {BUAD}
      */
     acos() {
-        return this.create(this.val.acos(), this.dval.negated().div(new this.BUADDecimal(1).sub(this.val.pow(2))));
+        return this.create(this.val.acos(), this.dval.negated().div(new this.BUADDecimal(1).sub(this.val.pow(2))), `acos(${this.valString})`, `-(${this.dvalString})/sqrt(1-(${this.valString})^2)`);
     }

     /**
@@ -253,6 +267,6 @@ class BUAD {
      * @returns {BUAD}
      */
     atan() {
-        return this.create(this.val.atan(), this.dval.div(new this.BUADDecimal(1).add(this.val.pow(2))));
+        return this.create(this.val.atan(), this.dval.div(new this.BUADDecimal(1).add(this.val.pow(2))), `atan(${this.valString})`, `(${this.dvalString})/(1+(${this.valString})^2)`);
     }
 }
\ No newline at end of file

解説

constructor
      * @param {numbers} v
      * @param {numbers} dv
      * @param {Decimal.Config|undefined} decimalConfig
+     * @param {string|undefined} valString
+     * @param {string|undefined} dvalString
      */
-    constructor(v = 0, dv = 0, decimalConfig) {
+    constructor(v = 0, dv = 0, decimalConfig, valString, dvalString) {
         this.setBUADDecimal(decimalConfig);
         this.val = v;
         this.dval = dv;
+        /** @type {string} */
+        this.valString = valString ? valString : v.toString();
+        /** @type {string} */
+        this.dvalString = dvalString ? dvalString : '0';
     }

元の式を保存するvalStringと微分した式を保存するdvalStringを追加しました。
例えば定数を宣言した時は(new BUAD(2))valString='2',dvalString='0'となります。

selectX
     selectX() {
         this.dval = new this.BUADDecimal("1");
+        this.valString = 'x';
+        this.dvalString = '1';
     }

変数として設定するとvalStringが数値から'x',dvalString'0'から'1'に変更されます。

create,copy
      *
      * @param {numbers} v
      * @param {numbers} dv
+     * @param {string} valString
+     * @param {string} dvalString
      * @returns {BUAD}
      */
-    create(v, dv) {
-        return new BUAD(v, dv, this.decimalConfig);
+    create(v, dv, valString, dvalString) {
+        return new BUAD(v, dv, this.decimalConfig, valString, dvalString);
     }

     copy() {
-        return this.create(new this.BUADDecimal(this.val), new this.BUADDecimal(this.dval));
+        return this.create(new this.BUADDecimal(this.val), new this.BUADDecimal(this.dval), this.decimalConfig, this.valString, this.dvalString);
     }

constructor の変更に伴って、今の式の情報と微分した式の情報が渡されるようになりました。

四則演算と基本的な関数

それぞれについては省略しますが、それぞれの式と微分した式をテキストにしているだけです。
例として掛け算だけ紹介しておきます。

f\left(x\right)=g\left(x\right)h\left(x\right)\\ f^\prime\left(x\right)=g^\prime\left(x\right)h\left(x\right)+g\left(x\right)h^\prime\left(x\right)

これに従って

multiply text
valString = `(${this.valString})*(${x.valString})`;
dvalString = `(${this.dvalString})*(${x.valString})+(${this.valString})*(${x.dvalString})`;

となります。
式で微分を表しているものはdvalString、そのままのものはvalStringを持ってきて当てはめています。

べき乗

前回書いたとおり、条件によっては正しい答えが出ないので一部条件分岐してあります。
テキストで出すときはそれに合わせてあります。

const xyPlnx =
  this.val.comparedTo(new this.BUADDecimal(0)) === 0 ||
  y.dval.comparedTo(new this.BUADDecimal(0)) === 0
    ? new this.BUADDecimal(0)
    : this.val.mul(y.dval).mul(this.val.ln());
const xyPlnxStr =
  y.dvalString === "0"
    ? "0"
    : `(${this.valString})*(${y.dvalString})*(ln(${this.valString}))`;
f\left(x\right)=g\left(x\right)^{h\left(x\right)}\\ f^\prime\left(x\right)=g\left(x\right)^{h\left(x\right)-1}\left(g^\prime\left(x\right)h\left(x\right)+g\left(x\right)h^\prime\left(x\right)\ln{g\left(x\right)}\right)

xyPlnxg\left(x\right)h^\prime\left(x\right)\ln{g\left(x\right)}の部分を表しています。

実際に動かす

微分した式を表示するための場所を追加します。

sketch.js diff
@@ -2,9 +2,9 @@ const canvasSize = 600;

 function setup() {
     createCanvas(canvasSize, canvasSize);
+    drawFunction();
     textAlign(RIGHT, TOP);
     textSize(20);
-    drawFunction();
 }
 function drawFunction() {
     clear();
@@ -14,6 +14,8 @@ function drawFunction() {
     x.val = 5;
     x.selectX();
     strokeWeight(2);
+    const derivativeFormula = document.createElement('p');
+    const formula = document.createElement('p');
     for (let x = -500; x <= 500; x += 1) {
         const decimalX = new BUAD(Decimal(x.toString()).div(100));
         decimalX.selectX();
@@ -34,6 +36,13 @@ function drawFunction() {
         //atan(asin(x))
         //const y = decimalX.asin().atan();

+        const yStr = y.dvalString.replace(/(?<!\*)[+-]?\(0\)[+-]?(?!\*)/g, '').replace(/\*?\(1\)\*?/g, '').replace(/\*?\(x\)'\*?/g, '').replace(/\(x\)/g, 'x').replace(/\(-x\)'\*/g, '-');
+        console.log(yStr);
+        formula.textContent = y.valString;
+        document.querySelector('body').appendChild(formula);
+        derivativeFormula.textContent = yStr;
+        document.querySelector('body').appendChild(derivativeFormula);
+
         stroke('black');
         point((5 + x / 100) * canvasSize / 10, -(-5 + y.val.toNumber()) * canvasSize / 10);
         stroke('green');

微分前の式と微分後の式を表示する場所を作ってそれぞれに表示します。

const yStr = y.dvalString
  .replace(/(?<!\*)[+-]?\(0\)[+-]?(?!\*)/g, "")
  .replace(/\*?\(1\)\*?/g, "")
  .replace(/\*?\(x\)'\*?/g, "")
  .replace(/\(x\)/g, "x")
  .replace(/\(-x\)'\*/g, "-");

これは上から順に

f\left(x\right)=f\left(x\right)\pm0\\ f\left(x\right)=f\left(x\right)\cdot1\\ f\left(x\right)=f\left(x\right)\cdot x^\prime=f\left(x\right)\cdot1\\ f\left(x\right)\cdot x=f\left(x\right)\cdot\left(x\right)\\ -f\left(x\right)=\left(-x\right)^\prime\cdot f\left(x\right)

を表しています。

cosh,sinh

\cosh x=\frac{e^x+e^{-x}}{2}\\ (((e^{x})+((-)*e^{-x}))*(2)-((e^{x})+(e^{-x}))*(0))/(2)^2=\frac{e^x-e^{-x}}{2}=\sinh x\\ (\cosh x)^\prime=\sinh x

となって正しいことが確認できます。

まとめ

微分した式をテキストで書きだせることを確認できました。次回はいよいよ数式微分を作っていきます。

Discussion

ログインするとコメントできます