https://www.beiz.jp/download_P/sea-ocean/00047.html（一部改変）
その５ 波：波頭は尖る！
前章でsin波を発生させて自然の水面に起こる波の元を作れるようになりました。ではここで、おもむろにですが実際の海面を見てみましょう：
https://www.beiz.jp/download_P/sea-ocean/00047.html（一部改変）
波の高い所（なみがしら：波頭）に注目です。そう、波って波頭が尖るんですよね。これは水の特徴的な性質の一つで、波の頂上付近と底とで水粒子の動き方が違うからなるのだそうです。さらに言うと、波の谷の所は結構平べったくもなります。つまり、波の断面を模式的に描くとこんな感じになっているわけです：
より正確には波の形は「トロコイド」という曲線でより良く表されます。トロコイドというのは、円が何らかの別の図形にピッタリ沿って滑らないように転がった時の、円の外周及び内部点のある点が描く軌道の事を言います。例えば平らな地面を走る自転車のタイヤにある空気を入れる所は「サイクロイド」というトロコイド曲線になります：
サイクロイドは転がる円の半径と、その円の中心から内部点までの距離が半径に等しい場合の特殊ケースです。この内部点までの距離が円の半径より短くなると、波の曲線が現れます：
縦軸の向きを入れ替えています
前章のsin波はこれに比べると大分と丸みを帯びた波であり、波頭と底が相似な形をしているため、そのままだと上のトロコイド曲線とは少しかけ離れていて、水面の波っぽくない見た目になってしまいます。そこで、sin波の波頭を尖らせて波の形に近付ける方法を見ていく事にしましょう。
① べき乗の性質を利用
例えば次のような0～1の間に定義された青色の直線について、値をべき乗してみます：
aの値が大きくなる程直線はどんどん凸型になっていくのが分かります。でも0と1は必ず通ります。0のべき乗は0、1のべき乗は1なので当然です。またYの値が小さい所はよりゼロに近付く傾向が近いのも見て取れます。
このべき乗の性質をsin波にそのまま適用してみます。ただその為にはsin波を0～1の間にする変換が必要です。sin波は-1～1の間で振動するので、まず+1すると0～2の間での振動に変わります。さらに0.5を掛ける（半分にする）と上の直線同様に0～１の間の振動に変換できます：
で、この右辺をべき乗すると、
a = 1, 2, 3の時のグラフ
こんな感じでべき乗の値が大きくなると頂上部は尖り、底面部は平たくなります。少なくとも元のsin波よりは水の波っぽい雰囲気が出ました。もちろん理論的に正しい訳ではないので、正解のaがある訳ではありません。あれこれ触りながらいい感じに見える値を模索する事になります。この波便宜上「べき乗sin波」と名前を付けておきます。
前章の波の式を用いるとべき乗環状波、べき乗直進波それぞれは次の式のようになります：
ハイトマップを作成するとこんな感じになります：
左がa=1、つまり元のsin波で発生させた環状波、右はa=4で同じ位置と波長で発生させたべき乗sin波です。右のほうが黒い部分が長く、白い部分がきゅっと狭まっていますね。少なくとも元のsin波よりは「波頭が尖る」に寄っているハイトマップになっています。
② トロコイド曲線をそのまま使用できないのか？
波がトロコイド曲線である事実があるのなら、それをそのまま使えば良いのでは？と思うかもしれません。トロコイド曲線は中心点の回転角度θを媒介変数とした時に以下の式で表現されます：
rcは転がる円の半径、rdはその円の中心から軌道を描く点までの距離で、波の形にするにはrc>rdにします。θは中心円の転がる角度で、2πまで行くと一回転となります。1回転するとsinθはゼロになるのでxの値はrc2πとなります。これが波の1周期、すなわち波長となります：
このθに値を入れると軌道を描く点の位置が算出されます。ただ、上の式だと上下が逆さまの波になってしまったり色々と波として扱いにくいので、式を変形して波として操作しやすい形に整え直しましょう。
まず波の高さをHと表し直し、上下逆さまを直すため上yの式の左辺をマイナスします：
さらに波が一番低くなるcosθ=-1の時に波の高さHを0にしたいのでHの値をオフセットします：
この状態で波の高さは0～2rdまで変化します。今rdは振幅として扱いたいので上のHをマイナスrdすると-rd～rdの範囲に収まります：
という事で、最終的には[振幅rd]×cosθという、言われてみれば当たり前の式になりました(^-^;。
さて、このrdはrcより小さいという制約があります。rdがrcに近い程波は尖り、同じ値になると完全に尖ります。rdを決める際にいちいちrcの値を気にするのも面倒なので、rdをrcに対する割合に変更してしまいます：
Rは0～1までの尖り係数で、0だと波は立たず、1に近付くにつれてどんどん波が高くなり、また尖っていきます。さらにrcは波長Lで置き換えらえるので、
と波長と尖り係数でパラメータを算出できるようになります。ただ、波長に対して波高は最大でも約1/7程度（0.143）、振幅だと1/14程度であるという事を前章でお伝えしましたが、上の式でR/2πはR=1の時でおよそ1/6.28なのでちょっと高過ぎです。そこでR=1の時に比率が1/14になるようにrdの右辺を半分にします：
これでrcとrdも整いました。
波の位相は座標xの値をずらす事で実現できます。前章同様に位相速度vは重力gとの関連式、
で表せるので、これに時間を掛け算すればある時間経過した時の波の移動距離が算出出来ます。
以上の式変形の結果を全部盛り込んだシェーダ用のトロコロイド波の式は次のように記述出来ます：
上の式で実際に描いてみたグラフがこちら：
横軸が距離x、縦軸が波の高さHで、青いのがトロコイド波、オレンジは波長と振幅を同じにしたsin波です。波長L=1で尖り係数R=1にしています。しっかり尖っていて、実際の波の形により近づいていますね(^-^)。
この実際の水の特性を反映させた何かいい感じのトロコイド波なのですが、実は実用上の大きな問題が一つあります。それはこの式が媒介変数表現になっている事です。前章までのsin波は右辺に座標（距離）があって、結果として波の高さHが出てきました。所が上の波の高さHの右辺にはθしかありません。このθを座標（距離）であるxで表現できれば話は終わりなのですが、残念ながら上のxの式をθの式にする、つまり逆関数を解析的に求める事は出来ません。
求める事が出来ないなら多項式などで近似すれば…と思うかもしれませんが、それも収束が良くない事が知られています（トロコイドのxの式はKeplerの方程式で、その逆関数はBessel関数を係数とするフーリエ級数の式で表され、それが収束が遅い事は分かっている…そうです（参照：特殊関数グラフィックスライブラリー、Keplerの逆関数））。その為、このままだとある座標x（＝距離d）の波の高さHを求める事が出来ないのです。
さぁ困ったぞとなるのですが、リアルタイム性を犠牲にするなら、解決策はあります。θを探索すれば良いのです。
③ 距離dになるθをニュートン法で求める
②の最後に上げた位置xの式。ちょっと色々ごちゃごちゃしているので固定的な係数をA,B,Cで表現し直します：
すっきりさせた所で、xを右辺に移項します：
この式を満たすθを求めるというのは「解の探索」に他なりません。解を探索する一番分かりやすい方法の一つが「ニュートン法」です。ニュートン法の細かな導出などはWikipedia等が参考になるかと思いますが、検索の動きをイメージすると以下のようになります：
青いラインが関数f(x)で、赤い所がX軸と交わっている所、すなわち解の一つです。この解を探すために、まず適当なXの値x0に対応した値f(0)を求めます（①）。そのf(0)の位置の接線を次に求め、接線とX軸とが交わる所x1を計算します（②）。x1から同じように対応するf(x1)を求め（③）、接線を頼ってx2を算出(④）。この作業を繰り返していくとどんどん求めたい解の方へ近づいて行きます。とあるxnの時のf(xn)が十分にゼロに近い値になったらそのxnを解として採用します。
接線の方へ降りて行った先のxを求める式は以下の通りです：
式の中に微分式があるので、関数は少なくとも1階微分可能である必要があります。
③冒頭の距離の式の左辺をf(θ)とすると、その微分式は、
です。よって、解を探す漸化式は、
となります。
後は簡単で、適当な初期値θ0を右辺に放り込むと左辺に次のθ1が算出されるので、それをまた右辺に入れて…を繰り返すだけです。
スタートとなる初期値θ0ですが、これは求めたい解に近い程収束が早くなります。ここは目星が欲しい所です。そこでxの式のグラフをご覧ください：
横軸はθ、縦軸はf(θ)です。青のがそのグラフになります。これを見ると青いグラフはふにゃふにゃしてはいますが、一定の角度で増加しているのが分かります。これは式を見ると良く分かります：
f(θ)は左側の大かっこにある傾きがAで(C-x)を通る直線の式に、右の大かっこにある振幅Bのsin関数がまとわりついたような式になっています。実際オレンジの直線が直線部分だけを描画した物です。実際に求めたい解θはθ軸と青いグラフの交点です。グラフのだとおよそ3.7くらいでしょうか。それに対してオレンジの直線のθの解は4.0位で非常に近いですよね。ですから、このオレンジの直線の解θを初期値θ0として使うと、求めたいf(θ)の解に非常に近い所からスタートできる事になります。直線のθの解は極めて簡単で、
このθ0を初期値として使うと、かなりな速度で収束します。実際表計算でテストしてみた所、何とθ2位でもうf(θ)が小数点以下5桁以下の精度でゼロに漸近していました！これは極めて速い収束と言えます(^-^)。
もしシェーダ内にこれを書くとするなら、ループ回数を3回と決め打ちしてしまっても良いかもしれません。そうすれば条件文を排除出来るので、その点でのパフォーマンス劣化は防げます。とは言っても計算量は多くなるので全体的なパフォーマンスとの兼ね合いを見て採用を判断するのが賢明かなと思います。
この波トロコイドで作成したハイトマップは次のような感じになりました：
左がべき乗sin波で3乗にセットした場合、そして右側が波トロコイドのハイトマップです。sin波も白い高い所がググっと尖ってはいるのですが、波トロコイドの方は一筋ピカーンと光っていますよね。それだけしっかりと尖っているという事です。
前章と本章で結構な自由度を持った波を発生させる事が出来るようになりました。ただ一つの波だけでは水面を再現する事は出来ません。水面の波はあらゆる要因で発生した無限の波の重なりによって出来ています。そこで次の章では「波の合成」について見ていきます。
