Playgroundで数値アルゴリズムを学ぶ

中学校では、数学の描画ほど恐ろしいものはありませんでした。

多くの問題にはすぐに利用できる解析的解法がなく、解決できるとしても膨大な計算作業を必要とする問題もあります。この場合、アルゴリズムの助けが必要になります。

アルゴリズムにうんざりしなかったことを願います。もし吐いてしまったら、明るい面を見ましょう。もう1回ランチを食べられますよ :]

このチュートリアルでは、アルゴリズムと、それを使用して分析的に難しい問題を解決する方法について学習します。プレイグラウンドを使用すると、ソリューションを簡単に視覚化してレビューしやすくなります。

数学が好きでなかったり、物理学やコンピューターサイエンスにあまり興味がなかったりする場合でも、このチュートリアルには興味深い内容が含まれているはずです。必要なのは微積分と基礎物理学の知識だけです。

数値アルゴリズムを使用して 2 つの難しい問題を解決する方法を学びます。しかし、より明確に言えば、これら 2 つの問題は解析的に解決することもできます。分析方法が機能しない場合はアルゴリズムが最適ですが、それでも 2 つの方法を比較することで、その本質を理解しやすくなります。

数値アルゴリズムとは何ですか?

簡単に言えば、数値アルゴリズムは、閉じた形式の解析解に依存しない数学的問題を解決するための方法のクラスです。

問題は、閉じた形式の解析解とは何かということです。

既知の数値と有限数のステップを使用して正確な結果を得ることができる式がある場合、それは閉形式の解析解と呼ばれます。

もっと簡単に言えば、代数を使って未知の量の問題を解く式を見つけ、特定の既知の数値を代入して結果を得ることができる場合、それは閉じた形式の解析手法を取得したことを意味します。

数値アルゴリズムはいつ使用すればよいですか?

多くの問題にはすぐに利用できる解析的解法がなく、解決できるとしても膨大な計算作業を必要とする問題もあります。この場合、アルゴリズムの助けが必要になります。

たとえば、限られた時間内に多数のオブジェクトの動作を計算する物理エンジンを作成する必要があるとします。このとき、数値アルゴリズムを使用してより速く結果を得ることができます。

これには、計算が高速になると結果がそれほど正確ではなくなるという副作用があります。ただし、ほとんどの場合、この結果で十分です。

天気予報は数値計算の恩恵を受けます。天候が変化するスピードと影響を与える要因の数はどちらも驚異的です。これは非常に複雑なシステムであり、将来を予測するタスクを実行できるのは数値シミュレーションのみです。

外は太陽が輝いているのに iPhone がいつも雨が降りそうだと知らせてくるのは、こうしたアルゴリズムが欠如しているからかもしれません。

[[142780]]

始める

ウォーミングアップとして、与えられた数値の平方根を計算するゲームをしてみましょう。どちらの方法でも二分法が使用されます。驚くべきことに、この方法は既にご存知かもしれませんが、その名前は知らないかもしれません。

子供の頃を思い出してください。100 以内の数字を選んで、他の人にそれを当ててもらうというゲームをしたことがあるかもしれません。彼が推測した数字が高すぎるか低すぎるかをほのめかすだけです。

ゲームが始まります。シャオミンはシャオチアンに推測を始めるように言いました。シャオチアンは「1だと思います」と言いました。シャオミンはそれは小さすぎると言いました。シャオチアンは再び2と推測しましたが、シャオミンはそれは小さすぎると言いました。その後、シャオチアンは3、4…と推測し、最終的に5を選択し、それが正解でした。

5ステップで正解しました。悪くないですね。しかし、シャオミンが78を選択した場合、正しく推測するにはしばらく時間がかかります。

このゲームを二分法でプレイすると、推測速度がはるかに速くなります。

二分法

数字は 1 から 100 の間であることがわかっているので、1 つずつ推測したり、ランダムに推測したりする以外に方法はありません。数字を2つの区間に分けます: a = [1,50]、b = [51, 100]。

次に、その数字がどの区間にあるかを判断する必要があります。これは非常に簡単です。数字が 50 かどうかを尋ねるだけです。数値が 50 未満の場合、間隔 b は考慮されません。次に、区間 a を 2 つに分割し、この手順を繰り返します。

例えば：

数は60です。間隔は、a = [1,50]、b = [51, 100]です。

最初のステップでは、Xiaoqiang は 50 (間隔 a の上限) と言いましたが、Xiaoming は小さすぎると言いました。この時点で、Xiaoqiang はその数字が区間 b 内にあるはずだとわかっていました。

2番目のステップは、区間bをc=[51,75]、d=[76,100]に分解することです。彼は依然として範囲 C の上限である 75 を選択しました。Xiaoqiang は彼に、それは大きすぎると言いました。これは、数値が区間 c 内にある必要があることを意味します。

分解を続けます。 。 。

この方法を使用すると、7 回の試行で正解を得ることができますが、1 つずつ試行すると 60 回の試行が必要になります。

1. 50 -> 小さい

2. 75 -> 大きい

3. 62 -> 大きい

4. 56 -> 小さい

5. 59 -> 小さい

6. 61 -> 大きい

7. 60 -> そうです！

x の平方根を計算する場合も、手順は同様です。数 x の平方根は 0 から x までの範囲にあります。これは、(0, x] と書くこともできます。数値 x が 1 以上の場合は、[1, x] と書くこともできます。

区間を分解すると、a=(0, x/2]、b=(x/2, x]が得られます。

数値xが9の場合、区間は[1, 9]、分解された区間はa=[1, 5]、b=(5, 9]、中央の値mは(1+9)/2 = 5です。

次に、m * m - x が目的の精度より大きいかどうかを確認します。ここでは、m*m が x より大きいか小さいか等しいかを確認します。大きい場合は、区間 (0, x/2] を使用してチェックを続行し、そうでない場合は、区間 (x/2, x] を使用します。

実際の手順を見てみましょう。m=5 に初期化し、精度値が 0.1 になることを期待します。

1. m * m - xを計算します: 5 * 5 - 9 = 11

2. 結果が期待値以下かどうかを確認します: 25 - 11 <= 0.1?

3. 明らかに満たされていない場合は、次の間隔に進みます。m * m は 9 より大きいですか?はい。次に、区間[1, 5]を使用すると、テスト値は(1+5)/2=3になります。

4. m * m - xを計算します: ３ * 3 - 9 = 0

5. 期待値が満たされているかどうかを確認します: 9 - 9 <= 0.1?

6. 完了。

注: 括弧の意味は何だろうと疑問に思うかもしれません。間隔には固定の形式（下限と上限）があります。異なる括弧は異なる間隔の境界を表します。丸括弧は境界が区間外にあること (開区間) を示し、角括弧は境界が区間内にあること (閉区間) を示します。たとえば、(0, a] には a が含まれますが、0 は含まれません。上記の例では、区間 (0, a/2] には a/2 が含まれますが、0 は含まれません。区間 (a/2, a] は、a/2 より大きく a より小さいすべての数値を表します。

遊び場で二分法を使う

さて、学んだ理論を適用する時が来ました。バイナリ検索アルゴリズムを自分で実装してみましょう。新しいプレイグラウンド ファイルを作成し、次のコードを追加します。

関数 bisection(x: Double) -> Double {
 //1  
    var 下限 = 1.0  
 //2  
    var 上 = x
 //3  
    var m = (下限 + 上限) / 2  
    var イプシロン = 1e- 10  
 //4  
 （fabs（m * m - x）>イプシロン） {
 //5  
        m = (下限 + 上限) / 2  
 m * m > x {の場合
            上 = m
        }それ以外{
            下 = m
        }
    } 
   
戻るm
 }

1. 間隔の下限を1に設定する

2. 区間の上限をxに設定する

3. 中間値を見つけ、期待される精度を定義する

4. 操作が精度要件を満たしているかどうかを確認する

5. 満足できない場合は、新しい中間値を見つけ、新しい間隔を定義して、検索を続けます。

関数をテストするには、次のコードを追加します。

 bis = bisection( 2.5 )とします。

m = (lower + upper) / 2 の行の右側を見ると、コードが 35 回実行されていることがわかります。つまり、答えを見つけるのに 35 ステップかかるということです。

#p#

すぐに、遊び場の素敵な機能の利点がわかります。値の完全な履歴が利用可能です。

バイナリ検索法では、実際の結果の近似値を正確に計算できます。数値履歴グラフを使用すると、数値アルゴリズムがどのように正しい解に近づくかを確認できます。

option+cmd+enter を押してアシスタント エディターを開き、m = (lower + upper) / 2 行の右側にある丸いボタンをクリックして、アシスタント エディターに履歴を追加します。

メソッドが少しずつ正しい結果にジャンプしていくのがわかります。

古典数学もできる

次のアルゴリズムは古代にまで遡り、紀元前 1750 年のバビロンで生まれ、アレクサンドリアのヘロンの著書『Metrica 100 AD』に記載されています。それが「ヒーローメソッド」と呼ばれる理由です。このリンクからHeroについて詳しく知ることができます。

このメソッドは関数を使用します。ここでは、a の平方根を計算します。関数曲線の「ゼロ点」、つまり関数値が 0 である点を見つけることができれば、x 値を解くことで a の平方根が得られます。

これを行うには、まず開始値として任意の x 値を選択し、その値での接線値 (関数の接線) を計算し、次に接線上の 0 点 (つまり、線が x 軸と交差する点) を確認します。次に、この 0 ポイントを再び開始値として使用し、精度が満たされるまでこのプロセスを複数回繰り返します。

なぜなら、正接値が計算されるたびに、真の値である 0 に近づくからです。

このプロセスにより、実際の結果に徐々に近づいていきます。次の図は、a = 9 で開始値が 1 の場合の a の平方根の解を示しています。

開始点 x0 = 1 では赤い接線が生成され、次に x1 点が生成され、次に紫色の線が生成されます。次に x2 点が生成され、青い線と接続され、最後に正解が見つかります。

古典数学と遊び場が出会うとき

私たちには、遊び場など、古代バビロニア人にはなかった良いものが沢山あります。次のコードをプレイグラウンドに追加して確認します。

関数ヘロン(x: Double) -> Double {
 //1  
  var xOld = 0.0  
  var xNew = (x + 1.0 ) / 2.0  
  var イプシロン = 1e- 10   
   
 //2  
 (fabs(xNew - xOld) > イプシロン)の場合
//3  
    x古い = x新しい
    x新 = (x旧 + x / x旧) / 2  
  } 
   
戻るxNew
 }

このコードは何をするのでしょうか?

1. 現在の結果を保存する変数を作成します。 xOld は最後の計算の結果であり、xNew は実際の結果です。

• xNewを初期化する方法を見つけることは良い出発点です

• イプシロンは望ましい精度である

• この例では、平方根を小数点以下 10 桁まで計算します。

2. while ループを使用して、目的の精度が達成されているかどうかを確認します。

3. 精度が達成できない場合は、xNew を xOld に設定し、次の反復に進みます。

次のコードを使用して関数を検証します。

      
彼女を = ヘロン( 2.5 )

Hero メソッドでは、正しい結果を見つけるのに 5 回の反復のみが必要です。

xNew = (xOld + x/xOld)/2 の行の右側にある丸いボタンをクリックして値の履歴を追加すると、最初の反復で適切な近似値が見つかることがわかります。

調和振動子の運動をシミュレートする

このセクションでは、数値積分アルゴリズムを使用して、基本的なタイプの動的システムである単純な調和システムの動きをシミュレートする方法について説明します。

このシステムは、振り子の揺れやバネの振動など、多くの現象を記述できます。特に、オブジェクトが一定時間移動するシナリオを記述できます。

この例では、バネに質量が取り付けられた物体があると仮定します。問題を単純化するために、抗力と重力を無視します。したがって、作用する唯一の力は、物体を元の位置に戻すバネの力です。

これらの仮定の下では、問題に対処するには 2 つの力学法則を使用するだけで済みます。

• ニュートンの運動の第二法則。物体に作用する力とその加速度の関係を説明します。

• フックの法則は、物体の弾性力とたわみの比例関係を説明します。ここで、k はバネ係数、x はオブジェクトのオフセットです。

単純な振動子方程式

物体に作用する力は弾性力のみなので、上記の 2 つの式を次のように変形できます。

簡略化した書き方:

これは共振周波数の二乗として記録されることもあります。

より正確な式は次のようになります。

ここで、A は振幅であり、これはここではオブジェクトのオフセット、つまり位相差として知られています。どちらの値も初期化時に固定値に設定されます。

時間を t = 0 に設定すると、振幅と位相の差を簡単に計算できます。

例を見てみましょう。質量 2kg の物体がバネに接続されており、バネの弾性係数は k=196N/m です。初期時間 t=0 では、バネは 0.1 m 移動しています。振幅、位相差、共振周波数は数式で計算できます。

プレイグラウンドでこれを試してみましょう:

この数式を使用して任意の時点での値を計算する前に、コードを追加する必要があります。

プレイグラウンド ファイルに戻り、最後に次のコードを追加します。

 //1  
 typealias ソルバー = (Double, Double, Double, Double, Double) -> Void 
   
 //2  
構造体ハーモニックオシレータ{
  var kSpring = 0.0  
  var 質量 = 0.0  
  var位相 = 0.0  
  var振幅 = 0.0  
  var デルタT = 0.0   
   
  init(kSpring: Double、質量: Double、位相: Double、振幅: Double、デルタT: Double) {
    self.kSpring = kSpring
    自己質量 = 質量
    self.phase = フェーズ
    自己振幅 = 振幅
    自己.デルタT = デルタT
  } 
   
 //3  
  関数solveUsingSolver(solver: Solver) {
    ソルバー(kSpring、質量、位相、振幅、デルタT)
  }
 }

これらのコードブロックでは何が起こっているのでしょうか?

1. 関数型の typealias を定義します。関数には Double 型の 5 つのパラメータがあり、void を返します。

2. 共振器を定義する構造体を作成します。

3. この関数は、振動問題を解決するために使用される Closure を作成します。 (実際の計算コードはありません)

正確な解決策を見てみましょう。

コードは次のとおりです。

 funcsolveExact(振幅: Double、位相: Double、kSpring: Double、質量: Double、t: Double) {
  変数x = 0.0  
 //1  
  オメガ = sqrt(kSpring / mass) とします。 
   
  変数 i = 0.0   
   
 i < 100.0の場合{
 //2  
    x = 振幅 * sin(ω * i + 位相)
    私 += t
  }
 }

この方法には、運動方程式を解くために必要なすべてのパラメータが含まれています。

1.共振周波数を計算する

2. while ループ内でオブジェクトの現在の位置を計算し、次の計算のために i を自己増分変数として使用します。

次のテスト コードを追加します。

 osci = HarmonicOscillator(kSpring: 0.5 、質量: 10 、位相: 10 、振幅: 50 、デルタT: 0.1 )とします。
 osci.solveSolver を使って(solveExact)

このソリューションのメソッドは少し奇妙です。パラメータを受け取りますが、データを返さず、何も表示しません。

#p#

これを実行するとどのような利点がありますか?

この関数を使用する目的は、while ループ内の振動プロセス中の特定の動的値を計算することです。プレイグラウンドでは、これらの動的な値の履歴を観察できます。

x = 振幅 * sin(ω * i + 位相) の値履歴を追加すると、動きの軌跡を確認できます。

[[142803]]

最初の正確なアルゴリズムが検証されたので、数値計算スキームを見てみましょう。

オイラー法

オイラー法は数値積分を計算するために使用される方法です。この方法は、1768 年にレオンハルト オイラーの著書『積分法の原理』で初めて提案されました。この方法の詳細については、http://en.wikipedia.org/wiki/Euler_method を参照してください。

オイラー法の中心的な考え方は、破線を使用して曲線を近似することです。

具体的な方法は、与えられた点の傾きを計算し、同じ傾きの折れ線を描くことです。この破線の終わりで、傾斜を計算し続けて別の線を描きます。ご覧のとおり、アルゴリズムの精度はポリラインの長さに依存します。

deltaT が何をするのか知りたいですか?

この数値アルゴリズムではステップ サイズの使用が必要であり、これはアルゴリズムの精度にとって重要です。ステップ サイズが大きいと精度は低下しますが、実行は高速になります。逆に、精度は上がり、速度は低下します。

deltaT 変数はステップ サイズを表します。この値を 0.1 に初期化します。つまり、オブジェクトの位置を 0.1 秒ごとに計算します。オイラーのアルゴリズムでは、これはポリラインの x 軸への投影の長さが 0.1 であることを意味します。

コードを書く前に、この式をもう一度確認する必要があります。

2 階微分方程式は 2 つの 1 階微分方程式に変換できます。次のように書くことができる。

差の商を使用して変換を完了できます。

以下のものも入手できます:

これらの方程式がわかれば、オイラー法を直接実装できます。

solveExact メソッドの後に次のコードを追加します。

 funcsolveEuler(振幅: Double、位相: Double、kSpring: Double、質量: Double、t: Double) {
 //1  
  var x = 振幅 * sin(位相)
  オメガ = sqrt(kSpring / mass) とします。
  変数 i = 0.0  
 //2  
  var v = 振幅 * オメガ * cos(位相);
  var vold = v
  var xoldオイラー = x 
   
 i < 100 の場合{
 //3  
    v -= オメガ * オメガ * x * t
 //4  
    x += ボルト * t
    xoldオイラー = x 
   
    vold = v
    私+=t
    }
 }

コードの内容:

1. 現在の位置とオメガ値を設定します。

2. 現在の速度を計算します。

3. while ループで、一次関数を使用して新しい速度を計算します。

4. 速度を使用して新しい位置を計算し、最後に i の値をステップ サイズ t だけ増加します。

次に、次のテスト コードを追加します。

 osci.solveSolver(オイラーを解析する)

xoldEuler = x の位置に値履歴を追加して確認します。この方法では、振幅が増加する正弦関数が表示されることがわかります。オイラー法は正確なアルゴリズムではないため、ここでのステ​​ップ サイズが 0.1 と大きすぎると、計算結果も不正確になります。

[[142811]]

これは、ステップ サイズが 0.01 の別の関数グラフです。明らかにこちらの方が優れています。したがって、オイラー法を使用する場合は、ステップ サイズが小さいほど結果が良くなることに留意する必要があります。ただし、妥協的なステップ サイズを使用すると、実装が簡単になります。

[[142812]]

ベロシティ・ヴェルレ

最後に説明するアルゴリズムは、Velocity Verlet と呼ばれます。これはオイラー法と同じ考え方に基づいていますが、新しい位置の計算方法が若干異なります。

新しい位置を計算するとき、オイラー法では実際の加速度を無視します。式は次のとおりです。一方、速度のヴェルレ アルゴリズムでは、計算時に加速度を考慮します。これにより、ステップの長さが等しい場合に、より良い結果が得られます。

solveEuler メソッドの後に新しいコードを追加します。

 funcsolveVerlet(振幅: Double、位相: Double、kSpring: Double、質量: Double、t: Double) {
 //1  
  var x = 振幅 * sin(位相)
  var xoldVerlet = x
  オメガ = sqrt(kSpring / mass) とします。 
   
  var v = 振幅 * オメガ * cos(位相)
  var vold = v
  変数 i = 0.0   
   
 i < 100 の場合{
 //2  
    x = xoldVerlet + v * t + 0.5 * オメガ * オメガ * t * t
    v -= オメガ * オメガ * x * t
    xoldVerlet = x
    私+=t
  }
 }

コードの内容:

1. ループ前のすべてのコードはオイラー法と同じです。

2. 速度に基づいて新しい位置を計算し、i の値を増やします。

プレイグラウンドでコードをテストします。

 osci.solveSolver(solveVerlet) を使用します

[[142815]]

次に何をすればいいでしょうか?

完成したプロジェクトはこのリンクからダウンロードできます。

デジタルの世界を巡るこの旅を楽しんでいただければ幸いです。いくつかのアルゴリズム、あるいは古典数学の興味深い歴史を知っているだけでも役に立つことがあります。