第３７回セーリングヨット研究会向け資料 2009. 7. 4
Excel ソルバーを用いない定常帆走性能の推定計算（ＶＰＰ）
金沢工業大学 増山 豊１ はじめに
第33 回セーリングヨット研究会において“Excelを用いた定常帆走性能の推定計算(VPP)”1)を 報告した。この報告におけるVPP計算の手段として、Microsoft Excelに組み込まれているソル バー機能を用いた。また、真風向を自動的に変化させて自動計算するために、やはりExcelに組 み込まれているVBA機能を用いた。これらの機能を用いることによって、面倒なプログラミン グをしなくても手軽にVPP計算を行うことが可能であった。 しかしながら、Excel2007 のバージョンになって、ソルバーの格納フォルダー名が変化するな どしたため、円滑な使用ができなくなった。今後、さらなるバージョンアップ（？）などによっ てまた使い勝手が変わるのも困ると思い、この際原点に帰ってソルバーを用いずに全てプログラ ミングによってVPP 計算するものに改めることにした。なお、Excel の VBA 機能だけは変わら ないものと考えられるので、プログラミング言語として使わせてもらうことにする。２ 変更した点
VPP は基本的に、非線形の連立方程式を解くプログラムである。このような連立方程式を解 く方法として、ニュートン・ラプソン法と呼ばれる手法が有名であり、これを用いることにした。 （なお、Excel ソルバーでは準ニュートン法という手法を用いている。）これらのプログラミン グは、かってVPP 解析を始めた頃に Fortran や Basic で用いたものであり、今回、Excel VBA に用いるために、古びたBasic のテキストをあらためて引っ張り出してきたという次第である。 また、以前の計算では、船体抵抗やセール流体力係数のグラフを、多項式で近似した式で求め ていた。これについては、必ずしも厳密に元の曲線と一致しない場合もあったため、多項式近似 をやめ、スプライン曲線で補間して求めることにした。さらに、以前のやり方ではExcel シート のセルの中に計算式を埋め込んでいたが、これも全てプログラムの中で計算することにした。こ れによって、うっかりセル内の式を消してしまうという恐れも無くなったものと考えている。 これらをまとめると、前報告からの変更点は以下のとおりである。 （１） Excel ソルバーを用いないで、ニュートン・ラプソン法で連立方程式を解く。 （２） 船体抵抗やセール流体力係数のグラフの値を多項式近似で求めず、数表データを読み込 んでこれからスプライン補間して求める。 （３） Excel シートのセル内には式は埋め込まず、全てプログラム内で計算する。３ 解析対象艇や定式化
解析対象とするのは前報告と同じで、筆者が基本設計したKIT-34 級艇“FAIR V”で、全長10.6m、排水量 3.4ton の外洋セーリングヨットである。この諸元や帆装図などは、前報告を参 照頂きたい。（前報告は、セーリングヨット研究会ホームページに掲載しています。）また、VPP の釣合式や、船体に作用する流体力の計算式などについても、前報告を参照頂きたい。
４
Excel ファイルの使い方
VPP 計算に用いる Excel ファイルの“釣合計算”シートを図１に示す。この表示内容は、前 報告のシート内容とまったく同じである。違う点は、上記のように前報告では、各々のセルの中 に計算式が埋め込まれていたのに対し、今回のものは単に数値が表示されているだけである。し たがって、前報告のシートではVBA を用いないで単独でソルバー計算を行うことが可能であっ たが、今回のものではソルバーは使えず、プログラム計算しかできない。以下、図中の番号と対 応させながら説明する。 ―――――――――― 入力 ―――――――――――― ① 真風速UTの入力。（なお、下段のγTの欄は真風向を表すが、入力する欄ではなくプログラ ム計算の終了時の値が出力される。） ② （これは入力欄ではなく、４つの解：u（x方向速度）、β（横流れ角）、δ（舵角）、φ（ヒ ール角）が出力される欄。プログラム計算終了時は①の下段に表示されたγT時の解が出 力されている。） ③ 水密度、空気密度、船体諸元、セール面積などの入力。 ④ 船体の流体力微係数の入力。 「計算条件」A において、真風向を指定する条件（計算開始、終了、変更間隔の角度）と、4 つの解のための初期値を入力する。 B の「VPP 実行」ボタンを押すことによって計算がスタートする。 ―――――――――― 出力 ―――――――――――― 以下、⑤～⑭の欄には計算結果が出力されるが、プログラム計算終了時は①の下段に表示さ れたγT時の結果が出力されている。あるγT時のこれらの値を知りたい場合は、Aの計算終了 の真風向の角度をそのγTにすればよい。 ⑤ 船体直立抵抗計算結果。 ⑥ 船体流体力係数。 ⑦ 船体に作用する流体力。 ⑧ ラダーに作用する流体力。 ⑨ 相対風速、相対風向の計算結果。 ⑩ セール流体力係数。 ⑪ セールに作用する流体力 ⑫ 船体復原モーメント ⑬ 船体とセールに働くX、Y 力と、K、N モーメントの合計値。計算が収束（釣合う）して いればこれらの値はほぼ０となる。⑭ 上記の合計残差。計算が収束していればほぼ０となる。なお、許容誤差はプログラム中で 指示している。 以上の計算結果のまとめは真風向ごとにC 部分に出力される。ここで、D の「残差」は⑭の 値が示されているが、これが十分小さくなっていない場合は、解が収束していないことを意味し ているので、初期値を変えるなどして再度求める必要がある。 A の初期値の目安はほぼ次のとおりである。 u0 =4m/s、 β0 =5°、 δ0 =－10°、 φ0 =－30° 計算では右舷側から風を受けた状態（スターボードタック）を考えているので、通常βは正の値、 δとφは負の値となる。（ヨットは通常、若干ウエザーヘルムとなる（風上へ切り上がろうとす る、すなわちスターボードタックの場合は右回頭する）ように設定されている。これを抑えて釣 合うように舵を切るため、δが負となる。図１内の図の定義参照） なお、プログラム計算において A の初期値を用いるのは、計算開始時の真風向の時だけであ る。２回目以降の繰り返し計算では、その前の釣合状態の値を初期値として用いている。このた め、まれに計算途中の真風向で収束しない場合があるが、このような場合は、その収束しなかっ た真風向から、初期値を変えて計算を開始してみるとよい。 図２のE～G部分は計算結果のグラフを表示している。（なお、これは図１のシートの右側部 分である。）Ｆはポーラーダイアグラムであり、Ｇは、真風向に対する帆走状態量（リーウェイ 角、ヒール角、舵角、艇速）の変化を示している。なお、ポーラーダイアグラムＦはExcelのレ ーダーチャートグラフを流用しているので、ここでは 10°毎（360°表示のためにデータが 36 個必要）しか表示できない。プログラム計算において真風向の刻みを 10°以下にすることは可 能であるが、このレーダーチャートでは 10°毎の点のみプロットするようにしている。また同 じ理由で、ポーラーダイアグラムＦにはリーウエイ角は含まれていないので注意を要する。な お、Ｇのグラフの横軸（=γT＋β）には、リーウエイ角を含めている。 図３は“係数”シートの内容を示しており、Ｈに直立直進時の船体抵抗、Ｉにセール流体力 係数の表とグラフを示している。前報告では、これらの値は多項式近似した式で求めていたが、 ここではこれらの表の値を読み込んで、スプライン補間して任意の相対風向における値を求めて いる。したがって船体抵抗やセール流体力を変更して計算を行いたい場合は、これらの表の値を 変更すればよい。なおこの時データ数を変更する場合は、“データ数”欄の数値を変更する必要 がある。またこの表ではセール流体力係数中のＣＹとＣＫは正の値で表示しているが、スターボ ードタックの計算を行うため、プログラム中で負にしている。（入力ではこれらを正の値で入力 のこと）
５ プログラム内容
このプログラム計算は、このシートに含まれているVBA プログラム（図４）によって行っている。プログラムは次の手順で表示させることができる。（なお、以下はExcel2003 以前のバー ジョンの場合である。） （１） 図１を表示させる。 （２） 「ツール」→「マクロ」→「マクロ」→「マクロ」ウインドウの表示。 （３） 「マクロ」ウインドウの中に「FAIRVPPbyNR」が表示され、すでに選択されている ので、「編集」ボタンを押すと、「コードウインドウ」と呼ばれる画面が現れ、プログ ラムが書き込まれている。 （４） なお、プログラムの実行は「マクロ」ウインドウの中の、「実行」ボタンを押すことに よってもできる。 以下、図4 のプログラム中の番号と対応させながら説明する。 （１） 計算条件および初期値の入力 １．“係数”シートの H から船体直立抵抗データと、I からセール流体力係数データを読み 込む。なお、ここでは確認のためプリントアウトするようにしている。各々のデータ数 もD3、N3 から読み込むことに注意。 ２．“釣合計算”シートへ戻って、初期値をA から読み込む。（なお、ここで読み込んだ値を 1 回目の初期値とするが、これらのセルの値は変更しないので、A の値は変化しない。 また、実行開始後の2 回目以降の初期値は、前回の収束結果を初期値として用いている。） ３．③、④から船体データ、船体流体力微係数などを読み込む。 （２）ポーラーダイアグラム入力用 E 表のクリア （３）繰り返し計算の設定と開始 （４）各回の計算結果の出力 （５）ポーラーダイアグラム入力用 E 表への出力 （６）メインプログラムの終了 （７）ニュートンラプソン法 （８）ガウス・ジョーダン法（連立方程式の解法） （９）釣合状態の計算 （１０）船体に作用する X 力の計算 （１１）船体に作用する Y 力の計算 （１２）船体に作用する K モーメントの計算 （１３）船体に作用する N モーメントの計算 （１４）セールに作用する力とモーメントの計算 （１５）スプライン補間
６ 結言
以前の手法に比べてVPP 計算を全プログラム化しただけで、船体データやセール流体力係数 は同じものを使っているので、当然計算結果にほとんど変化はない。しかしながら、以前はこれらのデータを多項式近似して求めていたのに対し、今回はスプライン補間して求めているので、 多少の違いが見られる。このようにすることによって、船体抵抗やセール流体力係数を変化させ たい場合は、Excel シート上の表の値を変えるだけでよいので、より汎用性が高まったものと考 えている。 なお、船体データの表現には未だ水槽実験から求められた流体力微係数を用いており、水槽実 験を行わないで性能推定を行うレベルには至っていない。今後、Delft シリーズのデータなどを 参考に、船体形状データさえあれば実験を行わなくても性能推定計算ができるような形へ発展さ せたいと考えている。
参考文献
（１） 増山 豊：Excel を用いた定常帆走性能の推定計算(VPP)、第 33 回セーリングヨット研 究会配布資料、(2007) セーリングヨットの釣合計算 by Masuyama 2009.5.3艇名：FAIR V 改良バラスト計算条件 水密度 rhow＝ 1025 kg/m^3 空気密度 rho＝ 1.2 kg/m^3
新抵抗グラフによる 計算開 40 deg 諸元
UT= 5 m/s 計算終了 180 deg 排水量 disp= 3775 kg セール面積 S= 56.4 m2
γT= 180 deg 計算間隔 10 deg 喫水線長さ L= 8.55 m セール風圧中 zce= -6.4 m （重心より上が負） 上記γT時の解 自動計算の初期値 喫水深さ D= 1.94 m メタセンタ高さ GM= 1.31 m
ｘ方向速度 u= 2.24 m/s u0= 4 m/s 船体流体力微係数 微係数変更 2007/1/18
リーウェイ角 β= -0.10 deg β0= 5 deg Xvv= 3.37E-01 Yv= -6.04E-01 Kv= 3.02E-01 Nv= -2.27E-02 舵角 δ= 0.09 deg δ0= -15 deg Xvφ= 0.00E+00 Yφ= 1.61E-02 Kφ= -8.08E-03 Nφ= -2.36E-02 ヒール角 φ= 0.44 deg Φ0= -35 deg Xφφ= -9.83E-04 Yvvv= 3.25E+00 Kvvv= -1.33E+00 Nvvv= -3.41E-03 φ= 0.0076 rad Xvvvv= -1.88E+00 Yvvφ= -4.68E-01 Kvvφ= 1.33E-01 Nvvφ= -9.20E-04 Yvφφ= 6.61E-01 Kvφφ= -3.33E-01 Nvφφ= -2.78E-02 横流れ速度 v= 0.003724377 m/s Yφφφ= -4.03E-02 Kφφφ= 2.08E-02 Nφφφ= 3.07E-02 基準横流れ v0= 0.001659315 舵流体力係数
合成速度 VB= 2.24 m/s Cxδ= -6.210E-02 Cyδ= -1.78E-01 Ckδ= 9.32E-02 CNδ= 9.10E-02
VB= 4.36 kt 図の方向を正とする。スターボードタック
直立直進抵抗 相対風 で考えるので、ヒール角（Φ）は負になる。
合計 RT= 28.85 kgf UA= 2.76 m/s また、舵角（δ）も一般に負になる。
Ｘ力 0.000000 282.99 N γA= 179.73 deg 0.27
Y力 0.000000 2009.5.3 waterco= 42826.37 airco= 256.94 deg Kモーメント 0.000000 ソルバーを使わないで 船体に働く流体力係数 セール流体力係数 Ｎモーメン 0.000000 Newton-Raphson法による XH0= 0.00000 Xs0= 1.101 2009.6.13 YH0= -0.00088 Ys0= -0.199 残差 1.626E-11 セール力を補間で求める KH0= 0.00044 Ks0= -0.165 NH0= -0.00022 Ns0= 0.018 船体に働く力 セールに働く力 XH= -282.95 N XS= 282.96 N YH= -37.69 N YS= 51.05 N KH= 36.58 N-m KS= 318.22 N-m NH= -79.46 N-m NS= 21.19 N-m 舵に働く力 舵流入角減少率 Xδ= -0.01 N γｒ= 0.066 Ｙδ= -13.36 N 船体復原モーメント 有効舵角 Kδ= 13.55 N-m KHheel= -368.35 N-m αR= 0.100 Nδ= 58.28 N-m 真風向 X軸方向速度 リーウェイ角 ヒール角 舵角 艇速 VB 艇速 γT+β 相対風速 相対風向 VB/UT セール推進力セール横押 セールヨーモセールヒール船体直進抵 残差 (deg) (m/s) (deg) (deg) (deg) (m/s) (knot) (deg) (m/s) (deg) (N) (N) (Nm) (Nm) (kgf) (N&Nm)
40 2.79 2.13 -19.58 -4.65 2.79 5.42 42.13 7.31 25.18 0.56 479.97 -2063.28 57.08 -14147.80 48.12 8.73296E-10 50 3.07 1.62 -21.96 -5.51 3.08 5.98 51.62 7.32 30.76 0.62 657.19 -2216.68 430.12 -15826.41 62.94 7.64545E-07 60 3.27 1.18 -21.76 -5.82 3.28 6.37 61.18 7.18 36.43 0.66 830.90 -2199.85 848.63 -15657.07 78.49 1.65596E-08 70 3.45 0.81 -20.55 -5.78 3.45 6.71 70.81 6.95 42.00 0.69 978.14 -2095.38 1169.53 -14802.55 92.40 2.72341E-08 80 3.55 0.57 -18.06 -5.35 3.55 6.90 80.57 6.59 47.90 0.71 1086.85 -1887.99 1319.94 -13042.64 104.09 7.04743E-08 90 3.59 0.41 -14.75 -4.59 3.59 6.97 90.41 6.13 54.20 0.72 1131.81 -1587.89 1290.02 -10694.28 110.34 5.48224E-09 100 3.58 0.30 -11.38 -3.69 3.58 6.97 100.30 5.61 61.03 0.72 1107.97 -1252.03 1126.10 -8279.82 109.73 2.10424E-09 110 3.53 0.22 -8.11 -2.74 3.53 6.86 110.22 5.03 68.77 0.71 1006.47 -910.53 894.34 -5916.83 100.89 5.55558E-10 120 3.36 0.16 -5.11 -1.87 3.36 6.53 120.16 4.41 78.73 0.67 842.14 -588.67 678.52 -3724.85 85.14 1.58877E-07 130 3.11 0.08 -2.75 -1.22 3.11 6.04 130.08 3.83 91.73 0.62 645.34 -324.09 551.13 -1996.14 65.53 4.83838E-09 140 2.78 0.07 -1.29 -0.76 2.78 5.40 140.07 3.38 108.29 0.56 468.16 -177.08 362.61 -916.60 47.65 1.2311E-07 150 2.59 0.02 -0.38 -0.34 2.59 5.04 150.02 3.04 124.87 0.52 390.93 -66.14 182.30 -258.91 39.84 1.21404E-10 160 2.44 -0.03 0.13 -0.10 2.44 4.74 159.97 2.83 142.79 0.49 336.39 1.26 93.34 106.92 34.29 3.48235E-09 170 2.32 -0.08 0.38 0.05 2.32 4.51 169.92 2.74 161.33 0.46 304.14 39.52 40.92 280.23 31.01 3.4213E-11 180 2.24 -0.10 0.44 0.09 2.24 4.36 179.90 2.76 179.73 0.45 282.96 51.05 21.19 318.22 28.85 1.626E-11 ＶＰＰ実行 ⑨ ⑩ ⑪ ① ② ③ ④ ⑤ ⑥ ⑦ ⑫ ⑬ A B C ⑧ ⑭ D γ γ U V u v x y A B T A T U -δ WIND HB HB β 図1 “釣合計算”シートの左側（入力データと計算結果）
VBAによる帆走性能表 艇名：FAIR V UT= 5 m/s
GammaT VB [m/s VB/UT Polar Diagram には、リーウェイ角を含んでいないことに注意 下記の図には、リーウェイ角を含んでいる
0 0.000 0.000 10 0.000 0.000 20 0.000 0.000 30 0.000 0.000 40 2.789 0.558 50 3.076 0.615 60 3.276 0.655 70 3.454 0.691 80 3.549 0.710 90 3.587 0.717 100 3.584 0.717 110 3.530 0.706 120 3.361 0.672 130 3.109 0.622 140 2.777 0.555 150 2.591 0.518 160 2.439 0.488 170 2.322 0.464 180 2.245 0.449 190 2.322 0.464 200 2.439 0.488 210 2.591 0.518 220 2.777 0.555 230 3.109 0.622 240 3.361 0.672 250 3.530 0.706 260 3.584 0.717 270 3.587 0.717 280 3.549 0.710 290 3.454 0.691 300 3.276 0.655 310 3.076 0.615 320 2.789 0.558 330 0.000 0.000 340 0.000 0.000 350 0.000 0.000 360 0.000 0.000 Polar Diagram (VB [m/s]) 0.0 1.0 2.0 3.0 4.0 90 180 -90 ヒール角、リーウェイ角、舵角 -25 -20 -15 -10 -5 0 5 10 0 30 60 90 120 150 180 真風向＋リーウェイ角 [deg] リーウェイ角 (deg) ヒール角 (deg) 舵角 (deg) 艇速 (knot)
Polar Diagram (VB/UT)
0.0 0.2 0.4 0.6 0.8 90 -90 180 E F G 図２ “釣合計算”シートの右側（ポーラーダイアグラムと帆走パラメータのグラフ） 艇名：FAIR V 新抵抗グラフ データ数 直立時 データ数 船体直立抵抗 Dragdatan= 12 セール流体力係数 Saildatan= 22 ＶB(m/s) D(kgf) gammaa CX CY CK CN 0.5 2.0 1 10 0.034 0.260 0.200 0.000 1 1.0 6.0 2 15 0.107 0.700 0.600 -0.015 2 1.5 12.0 3 20 0.195 1.060 0.930 -0.040 3 2.0 22.0 4 25 0.295 1.280 1.100 -0.080 4 2.5 36.0 5 30 0.404 1.415 1.248 -0.097 5 3.0 57.0 6 35 0.519 1.460 1.285 -0.100 6 3.5 96.0 7 40 0.636 1.470 1.294 -0.098 7 4.0 178.0 8 45 0.753 1.452 1.280 -0.092 8 4.5 300.0 9 50 0.866 1.400 1.217 -0.078 9 5.0 430.0 10 60 1.068 1.243 1.074 -0.040 10 5.5 570.0 11 70 1.220 1.065 0.911 0.003 11 6.0 700.0 12 80 1.303 0.880 0.736 0.046 12 90 1.312 0.688 0.565 0.092 13 100 1.272 0.500 0.403 0.110 14 110 1.200 0.450 0.298 0.100 15 120 1.243 0.280 0.162 0.080 16 130 1.250 0.138 0.055 0.060 17 140 1.240 0.023 -0.032 0.050 18 150 1.230 -0.076 -0.105 0.043 19 160 1.200 -0.150 -0.144 0.029 20 170 1.150 -0.188 -0.162 0.022 21 180 1.100 -0.199 -0.165 0.018 22 0.5 2 10 0.034 0.26 0.2 0 1 6 15 0.107 0.7 0.6 -0.015 1.5 12 20 0.195 1.06 0.93 -0.04 2 22 25 0.295 1.28 1.1 -0.08 2.5 36 30 0.404 1.415 1.248 -0.0972 3 57 35 0.519 1.46 1.285 -0.1 3.5 96 40 0.636 1.47 1.294 -0.098 4 178 45 0.753 1.452 1.28 -0.092 4.5 300 50 0.866 1.4 1.217 -0.078 5 430 60 1.068 1.243 1.074 -0.04 5.5 570 70 1.22 1.065 0.911 0.0029 6 700 80 1.303 0.88 0.736 0.046 0 100 200 300 400 500 600 700 800 0.0 2.0 4.0 6.0 8.0 VB(m/s) D(kg f) 0 0.2 0.4 0.6 0.8 1 1.2 1.4 0 50 GammaA100 150 CX -0.4 -0.2 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 0 50 100 150 GammaA CY -0.4 -0.2 0 0.2 0.4 0.6 0.8 1 1.2 1.4 0 50 100 150 GammaA CK -1 -0.8 -0.6 -0.4 -0.2 0 0.2 0.4 0.6 0.8 1 0 50 100 150 200 GammaA CN H I 図３ “係数”シート（船体抵抗とセール流体力係数の表とグラフ）
