潜降2015年09月13日 16:37

潜降
潜降


宿題になっている、自作減圧計算ソフトの作成だが、もちろん、全く進んでいない(エンリッチやディープのお勉強に勤しんでいるので)。

しかし、まあ、何もやらないというのも気になるので、昔懐かしいベーシックでやってみようかと、少し手を付けてみた。

(Chipmunk Basic)
http://www.nicholson.com/rhn/basic/

マックでもできるという、軽いベーシックで、ちょっと癖がある感じはするが、インストールとかしなくても動くので、環境を変えずに使うことが出来る。

(Chipmunk BASICプログラミング入門)
http://pharm.ph.sojo-u.ac.jp/~kumayaku/KH/CMBasic/CMBasic.html

後ろの方に、リファレンスが付いていて、一部日本語訳されていたりするので、若干抵抗なく始められるかもしれない。

次はジャバスクリプトとかで書くつもりで、最終目標はアプリとしてアンドロイドに突っ込めるようにしようと思っている。

何度も断るつもりだが、ビュールマンZHーL16Aを、何の加工もせずに実際の潜水に使ったりしてはいけない。

間違いなく減圧症を引き起こす。

もちろん、M値をどう設定するかということにもよるが、潜水士のテキストにも「これを守らなければ罰則が科せられる最低限のものである」と、明記されている。

減圧停止を必要とする潜水は、レクリエーショナルレベルでは、全世界の全ての指導団体で禁じられており、テクニカルレベルの講習を受け、必要な器材を整えない限り行ってはならない。

浮沈子が減圧ソフトを書いてみようと思ったのは、あくまでプログラミングの勉強と、厚生労働省が示している潜水士の基準を確認するためだ。

一人一人の潜水士が、自ら国家の定める基準を確認できないでどうする?。

潜水士テキストの「改訂にあたって」には、こう書かれている。

「改正高圧則の趣旨があくまで減圧表は潜水者が自分で作成することとされています・・・」

どこにそんな趣旨が書いてあるかは知らないが、事業者のいいなりで潜ってはならないということなんだろう。

じゃあ、自分でどうやって確認すればいいのかということになれば、電卓叩くか表計算ソフトに式を入れて、仕掛けを作ってやるしかない。

表計算ソフトは、フリーのものも充実しているので、浮沈子はそれで作ってみた。

個人的には(潜水士の受験勉強としても)これで十分だったが、だれもが簡単に確認できる仕掛けはない。

減圧を伴う潜水を奨励するわけではないし、法に基づく最低限の基準を、方程式ではなく、実際の減圧表として確認するためにはどこかにその基準が明記されていればいい。

今のところ、高圧則を改正するに当たって検討会が設けられ、その報告の資料として示されているのが、公式(?)の減圧表だ。

(高気圧作業安全衛生規則改正検討会報告書)
http://www.mhlw.go.jp/stf/shingi/0000037880.html

「・報告書
・報告書資料」

報告資料の方である(ファイル、重いです)。

まあ、これを見ればいいんだろうが、ヘリウム混合とかになると、素の減圧表ではなく、安全率が掛かっていたりするので少し不明確になる。

違法かどうかの判定には使えない。

もちろん、巷の減圧表(NSK減圧表)に比べたら、驚くほどおおらかなので、とてもじゃないがこれに従って潜る気にはなれないな。

浮沈子は、指導団体の講習を受けて、もちっと安心できる減圧スケジュールを使うつもりだ。

だから、あくまでも、お勉強の範囲である。

さて、ベーシックといっても、弄っていたのは40年近くも昔の話である・・・。

ちゃんと勉強したわけではなく、雑誌の中に書いてあるプログラムを真似して、若干変更を加えては、プログラミングの楽しさを味わっていたに過ぎない。

データの保存に、カセットテープを使っていた時代である・・・。

当然、記憶は薄れ、いや、殆どなく、覚えているコマンドは、cls(クリアスクリーン)くらいで、初学者と変わらない。

手ほどきを受けるのにちょうどいいページがあった。

(BASIC入門 (全11回))
http://dotinstall.com/lessons/basic_basic

「・#01 BASICとはなにか? (02:29)
・#02 「hello world」を表示させてみよう (02:57)
・#03 saveとloadを使ってみよう (01:34)
・#04 gotoを使ってみよう (02:41)
・#05 変数を使ってみよう (02:25)
・#06 ifを使ってみよう (02:29)
・#07 ifとgotoを組み合わせてみよう (02:50)
・#08 for、whileを使ってみよう (02:56)
・#09 gosubでサブルーチンに飛ばそう (02:30)
・#10 inputでユーザー入力を受け付けよう (02:39)
・#11 数あてゲームを作ってみよう (02:54)」

このページを見つけたので、チップマンクベーシックを選択したといってもいい。

ベーシックには、いろいろ種類があるが、どれを選んでも大して変わりはないだろう。

とりあえず、ロジカルな部分をどうやってコーディングするかというところを、とっつきのいいプログラミング言語で実装してみて、インターフェースや他言語への移植は、それから試みようと考えている。

元より、プログラムで飯を食う気はないので、いろはの「い」位で書ければいい。

小難しい仕組みを使ったり、玄人の使う記法を真似するつもりもない。

そうはいっても、コンパートメントは16もあり、それに対して演算式を深度と時間に応じて使っていくわけだから、配列位は使わないと話にならない。

再帰的な関数の記法や、グローバル変数(ワケワカ)の活用などは、今回はパスだ。

とはいっても、ベーシックの変数は、自動的にグローバル変数になってしまうわけで、誤って弄ってしまわないように気を付けなければならない。

本当は、フローチャートとかを書いて、その中で、変数がどう変わっていくかというのも記述して、正確に扱わなければいけないんだろうが、プロの開発ではないので、いきなりぶっつけで書きだしたら、案の定、躓いた!。

減圧の際に、計算の根幹となる「log2」の記法を間違えてしまった(log(2)と書くのが正しい)。

値が出なくて、何時間もさまよった挙句、リファレンスから読み直して分かった。

とりあえず、潜降するところは書けたので、以下に示す。

「10 rem (Pa+Pb)*NN2+R*NN2*(t-1/k)-((Pa+Pb)*NN2-QN2-R*NN2/k)*e^(-k*t)
20 pa = 100
30 pb = 0
40 nn2 = 0.79
50 r = 100
60 t = 4
70 qn2 = 74.5207
80 e = exp(1)
90 dim harf_time(16)
100 dim k(16)
110 dim bunatu(16)
120 for i = 1 to 16
130 read harf_time(i)
140 k(i) = log(2)/harf_time(i)
150 bunatu(i) = (pa+pb)*nn2+r*nn2*(t-1/k(i))-((pa+pb)*nn2-qn2-r*nn2/k(i))*e^(-k(i)*t)
160 print i,harf_time(i),k(i),bunatu(i)
170 next i
180 data 5,8,12.5,18.5,27,38.3,54.3,77,109,146,187,239,305,390,498,635
190 end」

初めの10行目は、窒素分圧を求める式を、参考までに書いている。

毎分10mで、40mまで潜降するというパターンで計算している。

70行目までは、値の入力をそのまま書いているが、この所は対話式に入力できるようにしようと思っている。

80行目の自然対数の底は、関数の値で持っていたのでそれを突っ込んでいる。

90行目から110行目は、配列を定義している。

kの値は、ハーフタイムごとに決まっているので、ダイレクトに入力してもいいんだが、ハーフタイムを明確にしておくために、一手間かけている。

120行目から170行目までは、配列の中の値を順番に入れたり出したりするための繰り返し処理をしている。

130行目で、ハーフタイムを呼び出し、140行目でkの値を出している。

150行目では、10行目の式そのものに、kだけ配列に直して代入しているのがわかる。

160行目で、念のためにプリント文で画面に表示させて、正しく計算されているかを確認する。

表計算ソフトと見比べてみると、正しい値になっているようだ(そっちは、大丈夫なのかあ?)。

ハーフタイムの一覧とかは、最終的にはファイルから読み込ませようと考えているが、とりあえずプログラム中に突っ込んでおく(180行目)。

190行目で、省略可能なエンド文を書いておく。

配列の計算方法を勉強するために、事例的に書いただけで、さっきのグローバル変数とローカル変数の扱いをどうしようかとか、全く考えていない。

実行するとこうなる。

「>run
1 5 0.138629 149.864007
2 8 0.086643 124.777521
3 12.5 0.055452 108.00371
4 18.5 0.037467 97.683619
5 27 0.025672 90.641226
6 38.3 0.018098 86.000261
7 54.3 0.012765 82.675667
8 77 9.001911E-03 0.300653
9 109 6.359148E-03 78.618322
10 146 4.747583E-03 77.58653
11 187 3.706669E-03 76.917703
12 239 2.900197E-03 76.398221
13 305 2.272614E-03 75.993184
14 390 1.777300E-03 75.673028
15 498 1.391862E-03 75.423595
16 635 1.091570E-03 75.229085」

滞底(正確には、潜降時間も含む概念だが、ここでは海底にいるという意味)については、ここで出したbunatu(16)という変数セット(配列)を、QN2に代入していく。

各プロセス内では、QN2を上書きしてしまってもいいのだが、どこのステージで計算した値であるかが分からなくなってしまっては困るので、何らかの変数に代入しておいて、区別できるようにしておいた方がいいかもしれない(全然検討していません!)。

今後の作業としては、以下の通りになる。

・潜降(済)
・滞底(正確には、滞底時間ー潜降時間)
・第1回停止までの浮上(ここで、M値と比較して判定)
・第1回減圧停止(次のステージに浮上できるできるまで減圧停止)
・第2回浮上(正確には、この浮上時間の間の窒素排出も含めて、ここで判定する)
・第2回減圧停止(浮上まで、以下同じ)

プログラミング的には、第1回減圧停止深度が分かった段階で、区切りをつけてもいいような気がする。

これが分かれば、後はひたすら同じ処理を繰り返して浮上まで持ってくればいいだけだ。

窒素だけ考えていればいい空気潜水とか、ナイトロックス潜水は、深度が分かればM値が決まるのでシンプルだが、ヘリウムが入ってくると、途端にややっこしくなる。

そのM値に安全率を与えるということになれば、さらに複雑だ。

実際のことを考えれば、酸素減圧やエアブレイクにも配慮しなければならない。

そこまでやるかどうかは未定だ(表計算でもやってません)。

どこまでやるかは、気分次第。

どっかで、誰かがやってくれれば、即、中断しよう。

だから、当然、焦ってやることはない。

さて、飯でも食いに行こうかな・・・。

滞底2015年09月13日 22:52

滞底
滞底


引き続き、チップマンクベーシックによる減圧ソフトの研究(?)である。

(潜降)
http://kfujito2.asablo.jp/blog/2015/09/13/7790890

潜降に引き続き、海底での作業時間を入れて、各コンパートメント(半飽和組織)の窒素分圧がどうなっているかを計算していく。

「10 cls
20 rem (Pa+Pb)*NN2+R*NN2*(t-1/k)-((Pa+Pb)*NN2-QN2-R*NN2/k)*e^(-k*t)
30 pa = 100
40 pb = 0
50 nn2 = 0.79
60 input "潜行速度(m/分) ? ",r
70 r = 10*r
80 input "作業深度(m) ? ",m
90 t = 10*m/r
100 qn2 = 74.5207
110 e = exp(1)
120 dim harf_time(16)
130 dim k(16)
140 dim bunatu(16)
150 for i = 1 to 16
160 read harf_time(i)
170 k(i) = log(2)/harf_time(i)
180 bunatu(i) = (pa+pb)*nn2+r*nn2*(t-1/k(i))-((pa+pb)*nn2-qn2-r*nn2/k(i))*e^(-k(i)*t)
190 next i
200 data 5,8,12.5,18.5,27,38.3,54.3,77,109,146,187,239,305,390,498,635
210 pb = 10*m
220 dim qn2(16)
230 r = 0
240 input "作業時間(分) ? ",t
250 for i = 1 to 16
260 qn2(i) = bunatu(i)
270 bunatu(i) = (pa+pb)*nn2+r*nn2*(t-1/k(i))-((pa+pb)*nn2-qn2(i)-r*nn2/k(i))*e^(-k(i)*t)
280 print i,harf_time(i),bunatu(i)
290 next i
300 end」

実行すると、こんな感じ。

「潜行速度(m/分) ? 10
作業深度(m) ? 40
作業時間(分) ? 60
1 5 394.940152
2 8 393.507218
3 12.5 384.697745
4 18.5 363.601758
5 27 329.772568
6 38.3 290.679262
7 54.3 249.796866
8 77 211.630537
9 109 178.973644
10 146 156.265934
11 187 140.344926
12 239 127.283084
13 305 116.656588
14 390 107.972641
15 498 101.027956
16 635 95.501135
>」

今回は、スクリーンをクリアにしてから実行するとか、インプット文を使って対話的な入力をするなど、見かけにも気を使っている。

えっ、入力間違えたらどうするのかってえ?。

全部チャラにして、もう一度やり直すに決まってんじゃん!。

そういうところは、いいかげんでもいい。

アマチュアだからな。

コードも汚いが、気にしない気にしない・・・。

200行目までは、概ね前回のコードをそのまま使っている。

補足すると、60行目から90行目が、対話的に値を取り込んで、換算しているところだ。

潜降速度はkPa単位に、潜降時間は潜降速度と作業深度から導いている。

コーディングの度に、この数字を入れることになるので、もう少し直接書き込んだままにしておいてもいい。

210行目は、作業深度のゲージ圧(kPa)を入力値から導いている。

220行目は、前回の潜降ステージの分圧を受け入れるために、配列を定義している。

本当は、ここでfor文を回して、配列を作っておくのがきれいなんだろうが、どうせ計算するのにfor文を回すので、260行目に押し込んでいる。

汚いやり口だな。

この場所に入れないと、正しく計算してくれないので注意だ。

また、さっきのpbとか、rとかtとかも、迷うことなく上書きしてしまっている。

本来なら、異なる値を取るわけだから、ステージ毎に異なる変数名を与えるのが筋だろう。

グローバル変数のまま、使いまわしている。

肝心のbunatu(16)も同じだ。

270行目で、上書きしてしまっている(あーあ・・・)。

280行目で確認のため、kの値を除いて出力している(だんだん、手抜きになって来てるのが分かるな)。

290行目でfor文を閉じて300行目で、ちゃんと省略せずにend書いている。

ざっとこんな感じで、浮上開始時点の窒素分圧を得ることが出来た。

これで、浮上開始となるが、浮上時には浮上速度を入力させることにする。

潜降速度と異なる場合も想定されるし、何より向きが逆だから、符号を反転させる必要がある。

速度というのは、ベクトルなので、犬が西向きゃ尾は東なのだ(ワケワカ・・・)。

M値の計算もしなければならないし、ここから一気に難しくなるんだが、まあ、なんとかなるだろう。

for文を回して配列を使えるようになったので、嬉しくて仕方ない(たんじゅん・・・)。

初回減圧停止深度の判定は、表計算のマクロでも、一番苦労したところだ。

水面から3m置きに逆に計算してきて、M値を下回るところでループを抜けるようにしたいのだが、16個の値のそれぞれの深度における判定方法、16個全部がその深度でクリアしたことの判定方法、ループからのうまい抜け方、初回減圧停止深度が確定して以降の処理に、うまく繋げるための工夫など、課題が多いセクションだ。

下から計算していく方法もないではないが、浮上した時のM値と比較しないといけない。

結局、減圧を終えて浮上していく間の分圧の変化も織り込むことになるので、初回停止時間1分毎に、浮上を含めて計算して、3m上に上がった時の分圧を出しておかなければならないのだ。

めんどくせー!。

まあ、パソコンが面倒くさいだけなんだがな。

そういうことを考えると、水面から逆算して初回停止深度を求め、そこから3m置きに1分刻みで停止と浮上を組み合わせて計算を繰り返すという戦略が妥当だろう。

この手順は表計算ソフトで実証済みなので、途中の値の検証も可能だ(精神的には、極めて楽ですな)。

配列の全てがクリアする時の深度を、きちんと捉えないといけないので、そこのところをどう書くかで悩んでいる。

分圧の配列も上書きしてしまって、必要なら、その時間の値を使って再計算する方法でもいい。

手作業を、そのまま置き換える形だが、うまいやり方を知らないのでそうするしかないだろう。

素人なので、ゴリゴリ書くしかない。

分圧の式も、サブルーチンにして使いまわせばいいんだろうが、大した式ではないので、そういうところには拘らない。

次回は、M値を扱う。

M値再登場!2015年09月13日 23:51

M値再登場!


(潜降)
http://kfujito2.asablo.jp/blog/2015/09/13/7790890

(滞底)
http://kfujito2.asablo.jp/blog/2015/09/13/7791360

今回は、M値を取り込むところだけ。

「10 rem cls
20 rem
30 rem <<<<潜降>>>>
40 rem
50 pa = 100
60 pb = 0
70 nn2 = 0.79
80 r = 10 : rem input "潜行速度(m/分) ? ",r
90 r = 10*r
100 m = 40 : rem input "作業深度(m) ? ",m
110 t = 10*m/r
120 qn2 = 74.5207
130 e = exp(1)
140 dim harf_time(16)
150 dim k(16)
160 dim bunatu(16)
170 for i = 1 to 16
180 read harf_time(i)
190 k(i) = log(2)/harf_time(i)
200 bunatu(i) = (pa+pb)*nn2+r*nn2*(t-1/k(i))-((pa+pb)*nn2-qn2-r*nn2/k(i))*e^(-k(i)*t)
210 next i
220 data 5,8,12.5,18.5,27,38.3,54.3,77,109,146,187,239,305,390,498,635
230 rem
240 rem <<<<作業>>>>
250 rem
260 pb = 10*m
270 dim qn2(16)
280 r = 0
290 t = 60 : rem input "作業時間(分) ? ",t
300 for i = 1 to 16
310 qn2(i) = bunatu(i)
320 bunatu(i) = (pa+pb)*nn2+r*nn2*(t-1/k(i))-((pa+pb)*nn2-qn2(i)-r*nn2/k(i))*e^(-k(i)*t)
330 next i
340 rem
350 rem <<<<M値>>>>
360 rem
370 dim mvn2a(16)
380 for i = 1 to 16
390 read mvn2a(i)
400 next i
410 data 126.885,109.185,94.381,82.446,73.918,63.153,56.483,51.133,48.246,43.709,40.774,38.68,34.463,33.161,30.765,29.284
420 dim mvn2b(16)
430 for i = 1 to 16
440 read mvn2b(i)
450 next i
460 data 0.5578,0.6514,0.7222,0.7825,0.8126,0.8434,0.8693,0.891,0.9092,0.9222,0.9319,0.9403,0.9477,0.9544,0.9602,0.9653
470 for i = 1 to 16
480 print i,mvn2a(i),mvn2b(i)
490 next i
500 end」

クリアスクリーン(10行目)、インプット文によるデータ入力(80、100、290行目)をコメントアウトして規定値を直書きしていたり、内容ごとに見出しを付けたり(30、240、350行目辺り)、少し手を入れているが、基本は同じだ。

M値(の係数)については、例のページからコピペする。

(改正高気圧作業安全衛生規則が施行されます: 高気圧作業安全衛生規則第八条第二項等の規定に基づく厚生労働大臣が定める方法等(平成26年12月1日厚生労働省告示第457号)参照)
http://www.mhlw.go.jp/stf/seisakunitsuite/bunya/0000071152.html

神の告示である。

(神の告示)
http://kfujito2.asablo.jp/blog/2015/08/26/7750625

まあいい。

こっちは、それどころじゃない。

それを、配列(mvn2a(16)、mvn2b(16))に取り込んでいる。

実行するとこんな感じ。

「>run
1 126.885 0.5578
2 109.185 0.6514
3 94.381 0.7222
4 82.446 0.7825
5 73.918 0.8126
6 63.153 0.8434
7 56.483 0.8693
8 51.133 0.891
9 48.246 0.9092
10 43.709 0.9222
11 40.774 0.9319
12 38.68 0.9403
13 34.463 0.9477
14 33.161 0.9544
15 30.765 0.9602
16 29.284 0.9653
>」

第一と第二半飽和組織のところが、画面ではうまくプリントできないんだが、細かいことは気にしない。

もう、配列のハンドリングは慣れたものだ。

どうやら、read文とdata文は、出てきた順番に対応するようで、心配していたバッティングとかはない。

さすがに、その辺りは言語環境を作る時に配慮しているようだ。

M値の係数とかは、どこで取り込んでもいいので、本来ならデータを読み込むディビジョン作って、纏めて書いておくのが筋だろう(コボルかよ?)。

これに、深度のデータを入れれば、M値が求まる。

ちなみに、窒素だけだとこうなる。

M=((Pa+Pc)/mvn2b)-mvn2a

海中で水深15mの半飽和時間5分のコンパートメントのM値は、上式によりこうなる。

((100+150)/0.5578)-126.885=321.3043151667264

この数字を、やれ水深だとか、圧力だとか(圧力ですが)といって、実際のダイビングの際の目安としようとしても、おそらく何の意味もない。

減圧モデルにおける、計算途中のただの数字だ。

どのコンパートメントを重視するかは、ダイビングのプロファイルが関係してくるし、そもそも、M値の係数がさじ加減によるものなので、当てにはならない。

しかし、この係数は臨床的に裏付けがあるので、まるっきりのデタラメではない。

というよりも、この係数こそが減圧表の命なのだ。

神の数字だな。

どっかの出版社が掲出していたホームページ上の数字が間違ってたことがあったが、浮沈子の指摘で差し替えられたりもした。

ったく、もう・・・。

まあ、どうでもいいんですが。

さすがに、この数字については、パブリックコメントでも疑義が寄せられていたが、厚労省は実績に基づくとして譲らなかった。

(「高気圧作業安全衛生規則の一部を改正する省令及び高気圧作業安全衛生規則第八条第二項等の規定に基づく厚生労働大臣が定める方法等に係る異見募集について」に対して寄せられた御意見等について:「異見募集」って、何なんだあ?)
http://search.e-gov.go.jp/servlet/Public?CLASSNAME=PCMMSTDETAIL&id=495130288&Mode=2

「Q:減圧表の計算の根拠について、a 値や安全率をどのようにして導いているのか。」

「A:今般新たに規定した減圧管理方法は検討会報告書を踏まえたものです。具体的な導出方法については検討会報告書第3(2)や参考文献である厚生労働科学研究「新しい標準減圧表作成に伴う実地調査および検証調査研究」の報告書第4章から第5章、その他専門図書で御確認ください。」

間違いないから、安心しろと言われてもなあ・・・。

とにかく、係数の取り込みは成功した。

当該深度が分かれば、16のコンパートメント全てのM値のセットを出せる。

その深度の窒素分圧は、コンパートメント毎に算出可能である。

あとは、海面からの深度を3m刻みにして、分圧とM値を突合して、判定ロジックを書いて実装するだけだ。

そこの所がキモなのは分かっている。

だけど、今日は限界だな。

明日は、インスピの点検とか、石垣ダイビングの準備で忙しいしな。

帰って来てからやることになるだろう。

外堀は埋まってきた。

石垣を登れば、天守閣はすぐそこである。

うん?、石垣繋がりかあ?。

まあ、そういうことにしておこうか・・・。