人気ブログランキング |

ブログトップ

電子工作やってみたよ

最近、毎日がジタバタしてシヌーという感じ。
仕事では客先で新型商品が出るというのでそれに対応した試験機のハード、ソフトを請け負ってアップアップ
個人的には村の役を引き受けたので、人間関係を含めてガチャガチャ
わたしの脳みそは、ボケかかり。  どうしたらいいんだ。

そうそう ブログでは 「人の悪口と愚痴は書くな」と言われたけどつい出ちゃう。  
やばいやばい。

c0335218_10495804.jpg
パイを計算するのに級数展開というのがあるそうです。
今回 比較的わかり易そうな下図二つ目のライプニッツ級数の式でやってみました。
左辺のパイの下にある4を右辺に移して、左辺は パイ= の形にしました。
しかし こうするとパイがでる、なんてどうして考え付いたのでしょうね。 すごい頭ですね。


c0335218_19244143.jpg



全プログラムリスト 
DM42(HP42Sと同じ) 動作はUSBで電源供給なので電池の3倍速で動いています。

パイを求めるライプニッツ級数
使用変数
00  現時点でのパイ
01 現在の項の分母
02 演算回数 2000 を減算して=> 0
03 1回前のパイ
04 1回前の分母

LBL "AA" ソフトスタート
LCRG   変数オールクリア
1
STO 01  最初の分母 1を設定
4     最初のパイの値を設定     
STO 00
LBL "BB"
2000   途中経過表示のカウンタ(演算2000回ごとに表示)
STO 02
4      項の分子の読み込み
ENTER
2      項の分母の階差の読み込み
RCL + 01 今までの分母に階差を加算する
STO 01  現在の分母をしまう
÷     分子(4)÷分母
STO - 00  過去のパイの演算結果から今の値を引く
RCL 00     現時点のパイの値を一時退避
STO 03
RCL 01  現時点の項の分母を一時退避
STO 04
LBL "CC"
4      項の分子の読み込み
ENT
2      項の分母の階差の読み込み
RCL + 01 今までの分母に階差を加算する
STO 01  現在の分母をしまう
÷     分子(4)÷分母
STO + 00 過去のパイの演算結果から今の値を引く
RCL 03      1回前のパイ
RCL 04  1回前の項の分母
PSE    2000回ごとに引き算の時の値を表示する
RCL 00  現時点のパイ
RCL 01  現時点の項の分母
PSE    2000回ごとに足し算の時の値を表示する
RCL 02  演算回数を読み込む
1
-     演算回数を
X キ 0   演算回数がゼロかチェック
GTO "CC"  ゼロでないなら"CC"へ行く
GTO "BB"  ゼロならば"BB"へ行き2000をセットし直し
END    ソフト 終わり

結果
R/S で演算を停止させたときのレジスタ X,Y,Z,Tを書きます。

T 3.14159254855
Z 19,039,999 元のミスの値 1903999
Y 3.14159275863
X 19,040,001 元のミスの値 1904001

TとYの値は加算する項と減算する項の違いで正解を挟んで上下の値。
計算ごとに段々近づくのですね。
Xはこの時の項の分母となります。
すなわち 演算回数はXの値の約半分の950万回(元ミス95万回)となります。
これだけ計算しても有効桁数は 3.141592 の 7桁しかありません。
とても手計算で出来るやり方ではないですね。

=====================================
追記 2019-04-30
計算を続けてやってみました。

T 3.14159258965
Z 31,279,999
Y 3.14159271753
X 31,280,001

精度は上がらず上下の位置がずれました。
ゆえに 演算桁数を上げずに計算項数を増やしても精度は上がらないと思います。


=====================================
=====================================
追記 2019-04-30
すがわらさんからコメントいただいたのでここに返事を書いてみます。

すがわら氏
横入りすいません。VS入れたのでテストとしてC#で計算させました。FORTRANと違って桁数が少ないので大して精度が出ないのですね。ライプニッツ式で1000万項計算させて、 3.1415926035897882384621752928が得られました。100万項だと、 3.1415921535892932379938932832でした。 無限級数の和は積分に置換できるので、積分からの逆変換じゃないですかね。π/4だからarctanの計算が怪しいかな。

高井
すがわらさんの ライプニッツ式で1000万項計算させて、
3.1415926035897882384621752928
これの9桁目の0は5の移し間違えでないですか、
このままだと3.1415926の計8桁しか有効でありませんね。
ここに5を置くと 3.1415926535897 までの計14桁まで正しくなりますね。

私思うのですがこのライプニッツ式ではいくら演算の項を増やしても精度は上がらないきがします。
右辺2項目の1/3(あるいは4/3)の所で0.333333333333で最後切り捨てられるので演算ソフトの有効桁数以上の精度は出ないでしょう。
HP42Sの場合 12桁(内部演算15桁)でやってますから良くてそこどまりでしょう。

外国のサイトでHP15CやHP32SiiのRPN法でパイを級数を解いているところがありました。
そこ調べてどんなことやってるのか調べて載せてみますね。
ひとの作ったソフトの解析するのを楽しいと感じているうちはまだボケていないのかな。


参考サイト 2個
https://akatuki-724.blogspot.com/2014/05/hp-prime.html

www.hpmuseum.org/cgi-sys/cgiwrap/hpmuseum/articles.cgi?read=899

========


もう 疲れたので とりあえずここまで。 あとで追加します。
****************************************************************


電子工作ランキング

にほんブログ村 その他趣味ブログ 電子工作へ
にほんブログ村

by telmic-gunma | 2019-04-29 19:28 | HP電卓 | Comments(6)

DM42で乱数を遊ぶ (3)

前回モンテカルロ法を使って遊んでみましたが、2017年11月16日 HP42S でやったのと同じことをDM42 でサンプリング数を増やしてやってみました。
モンテカルロ法でランダムにやるのでなく等間隔で均一にやってみようというわけです。
DM42ならHP42Sの750倍高速ですし、何より電池を消耗しないという大きなメリットがあります。
図のNo1はモンテカルロ法で、ひとつの真四角の中に乱数でランダムに点を置いて半円(1/4円)の内側にあるかを調べて数の比率で円周率を求めます。

c0335218_10495804.jpg


c0335218_10462042.jpg



全プログラムリスト 
DM42(HP42Sと同じ) 動作はUSBで電源供給なので電池の3倍速で動いています。

01 LBL "AA"
02 0
03 STO X * X座標のポイント
04 STO Y   * Y座標のポイント
05 STO T   * トータルポイント数
06 STO P   * 円内部の点の数 (x^2 + Y^2)が1以下
07   LBL  "BB"     * ループ先頭ラベル
08 XEQ "HH" * 円内部の点の数を数えるサブルーチンへ行く(x^2 + Y^2)が1以下
09 0.001 * 1ステップのY移動値
10 STO +Y   * Yに1ステップ分を加算する。
11 RCL Y   * 1ステップ分移動したY値を読み込む
12 1      * 円の内側かいなかチェックする値(1)をセット
13 X>Y ?    * 円の内側ならば
14 GTO "BB"  * ラベルBBへ行く
15 0      * 0を読み込む
16 STO Y   * Yの座標位置をゼロとする
17 0.001    * 座標 1ステップ分読み込み
18 STO +X   * Y座標 1ステップ分加算
19 RCL X   * 1ステップ分移動したX値を読み込む
20 1      * 円の内側かいなかチェックする値(1)をセット
21 X>Y ?    * 円の内側ならば
22 GTO "BB"  * ラベルBBへ行く
23 RCL T   * トータルポイント数
24 RCL P   * 円内部の点の数 (x^2 + Y^2)が1以下
25 STOP    *


26   LBL  "HH"  *  円内部の点の数を数えるサブルーチンへ行く(x^2 + Y^2)が1以下
27 1      *
28 STO +T   * 測定トータル数を+1
29 RCL X   *
30 X^2 * 現在のX座標値を2乗する
31 RCL Y   *
32 X^2     * 現在のY座標値を2乗する
33 +      * 現在のX座標値を2乗とY座標値の2乗を加算する
34 1       *
35 X<Y ?    * 結果が1より大きければ円の外なので何もせずにリターン
36 RTN     * 結果が1より大きければ円の外なので何もせずにリターン
37 1      *
38 STO +P   * 結果が1以下ならば円の内側なので円内部の点の数に1加算してもどる
39 RTN     *

使用した変数
X * X座標のポイント
Y   * Y座標のポイント
T   * トータルポイント数
P   * 円内部の点の数 (x^2 + Y^2)が1以下



結果A     100dig * 100dig = 10,000点の時
0.01ごとに0~1まで
円の内側の点  7,953
真四角の全点 10,000
円周率  ( 7,953 * 4 )  /  10,000   =  3.1812
*有効桁は2桁


結果B 1,000dig * 1,000dig = 1,000,000点の時
DM42 USB電源にて  8分55秒
0.001ごとに0~1まで
円の内側の点  786,338
真四角の全点 1,000,000
円周率  ( 786,338 * 4 ) / 1,000,000 = 3.14554
*有効桁は3桁


結果C 10,000dig * 10,000dig = 100,000,000点の時
DM42 USB電源にて  約15時間
0.0001ごとに0~1まで
円の内側の点  78,549,762
真四角の全点 100,000,000
円周率  ( 78,549,762 * 4 ) / 100,000,000 = 3.14199048
*有効桁は4桁

========


結果A サンプリング数 10,000点の時   3.1812
結果B サンプリング数 1,000,000点の時 3.14554
結果C サンプリング数 100,000,000点の時 3.14199048

乱数ではなく縦横規則正しいポイントならばサンプリング数に比例して精度が上がりますね。
しかしこれはまったくの遊びですね。高速のDM42を使って15時間もかけたのに有効数字4桁なんて。
DM42でUSB電源はいいですね、今までのように長時間でも電池代金を気にしなくていいのですから。 
これHP42Sでやったならば途中で電池切れをおこして答えは出ないでしょう。
15時間の750倍とは11250時間 これは469日 なんと1.28年  
まあ 電池がもったいし、一年もやっている人絶対いないでしょう。


ふと思ったのですが、乱数を使ったこのモンテカルロ法は、
最近はやりの選挙で行われている「出口調査」やNHKの世論調査と同じですね。
精度はせいぜい2桁ぐらいだけれど、少ないサンプリング数で大まかな動向を出せると言う事でしょう。


このHP電卓のRPN言語は、本当にアセンブラと同じですね。
コメントの書いてないリストは、自分で作った物でさえ後から見ると全く動作が解りませんでした。
これからは、ブログラムの一行ごとに真面目にコメントを書いていきましょう。

****************************************************************
つい先日雪が降っていたのに、もう庭に沢山の花が咲きだしました。
花は綺麗ですけど、地面の映りがうまくないですね。 どうしたらいいかな。
c0335218_07412113.jpg


白いチューリップ
c0335218_07414935.jpg


これなんだろう。

c0335218_07421112.jpg


これ スイセンですよね。
色や形 たくさんの種類があるみたいですね。
c0335218_07430123.jpg



c0335218_07433259.jpg


c0335218_07434999.jpg


タマネギ 300本くらい植えました。 みんな元気に育っています。
向こう側はニンニクです。 去年食べられなかったものを植えておいたら育ってきました。

c0335218_07441168.jpg


長ネギ

c0335218_07482767.jpg




電子工作ランキング

にほんブログ村 その他趣味ブログ 電子工作へ
にほんブログ村

by telmic-gunma | 2019-04-17 20:37 | HP電卓 | Comments(6)

DM42で乱数を遊ぶ (2)

DM42の高速性を生かして乱数の数を増やしたら精度が上がるのか遊んでみました。
前回乱数の発生のばらつきは,サンプリング数を増やすと小さくなるようなので、
以前やったように モンテカルロ法で円周率を求めてみます。
これは x,y 二つの値を乱数でだして座標としそれが円の内部にあるか外部なのかで円周率を求めるものです。 


c0335218_10080924.jpg


全プログラムリスト 
DM42(HP42Sと同じ) 動作はUSBで電源供給なので電池の3倍速で動いています。

01 LBL "CC" **************************
02 123.4567  乱数発生器の初期値 何でもよい
03 SEED    乱数初期値セット
04 CLRG   すべての変数をゼロにする

05 LBL "DD" **************************
06 1
07 STO + 3  全体の回数 カウントアップ
08 RAN  乱数作成 X
09 X^2    2乗して01にしまう
10 STO 01
11 RAN  乱数作成 Y
12 X^2    2乗する
13 RCL 01
14 +     乱数のXとYを足す
15 1     乱数を2乗して足したものが1を超えたか "D
16 X < Y ?
17 GTO "DD"   超えていたらデータ無視して次の乱数の計算に入る。
18 1
19 STO + 04  1以内ならばDをカウントアップする
20 GTO "DD"   次の乱数の計算に戻る


測定は10分間やりました。
変数名  データー
 01  Xのランダム数
 02  未使用
 03  全体の回数  開始から1分後 54880 10分後 555760
 04  1以下の回数 開始から1分後 43129 10分後 436878

1分後円周率   (  43129 X 4 ) / 54880 = 3.14351
10分後円周率  ( 436878 X 4 ) / 555760 = 3.14436


========


乱数のサンプリング数を増やせば乱数のバラつきの比率は減る方向でしたが、乱数を使って何かやるのは、サンプリング数を増やしても精度は上がらないようですね。
乱数の使い方自体 まだ勉強する必要があるな。

HP電卓(含むDM42)でプログラムを組むのは もうゲームですね。 役立つかどうかは二の次。

****************************************************************

4月も半ば近いというのに雪景色。
自動車のタイヤ、まだ冬用のままで良かったです。

今年は、村の役を引き受けたので、なんだかんだでてんてこ舞い。
この前書いた「ジタバタ」しまくりです。
先日県議会議員の選挙で立会人を初めてやりました。
朝6時から夜7時まで投票場に缶詰になり(投票場から外に出るのは禁止されてます)
一日中 投票箱をにらんでいました。
そして次は、中学校の入学式の来賓、
まあ 端っこに座っているだけでしたがこれも初めて。
今年は、この後何があるのかな。

c0335218_10584652.jpg



c0335218_10585790.jpg




電子工作ランキング

にほんブログ村 その他趣味ブログ 電子工作へ
にほんブログ村

by telmic-gunma | 2019-04-11 10:36 | HP電卓 | Comments(4)

消費電流を測るためにDM42の内部を覗いたのでその時の写真を載せます。

c0335218_22254080.jpg


c0335218_22255602.jpg


c0335218_22260412.jpg

DM42の消費電流を測ったときの様子です。
測定機は三和計器のPC5000です。
c0335218_22261211.jpg



うわー  また雪だ
c0335218_07235924.jpg



by telmic-gunma | 2019-04-04 07:26 | HP電卓 | Comments(0)

DM42で乱数を遊ぶ

DM42になって計算が非常に高速になり、かつUSB外部電源を使用すれば、
電池の消耗を気にしないで済むというメリットがうまれました。
乱数を扱うと、どうしても計算量が大きくなりますが、今までのHP電卓の遅さに悩んでいたことが、
このDM42を使うとかなり改善されます。

c0335218_10083064.jpg

( データ 1 )
HP32Sで乱数を発生させました。
30分で約40000個 バラつきは、4%でした。

( データ 2 )
DM42(電池)で乱数を発生させました。
1分で約30000個 バラつきは、8%でした。

( データ 3 )
DM42 (USB外部電源)で乱数を発生させました。
10分で約1000000個 バラつきは、0.8%でした。


大まかに言ってサンプル数を30倍するとバラつきは 1/10になりました。
いつも思うのですが、CPUで乱数を発生させているなら、現在の状況を監視させて、乱数のバラつきを
自分で補正をかけるようにすればバラつきのない乱数が出来るはずです。
次はこの自分で補正をかけながら乱数を発生するソフトを作って見たいですね。
このDM42の高速性があればかなり面白いこと出来るとおもいます。




乱数精度のデータ ( 表示文字サイズが統一できずバラついてます)


( データ 1 )
HP32S 乱数の精度 (2017-10-31 HP電卓で遊ぶ=>パイの計算 より)

出力されたデータ (30分間測定)
変数名  データー範囲   出力回数
 A     0~1     3920
 B     1~2     3913
 C     2~3     4075
 D     3~4     3937
 E     4~5     3937
 F     5~6     3965
 G     6~7     3925
 H     7~8     4007
 I     8~9      4017
 J     9~10     3983

最大値  4075
最小値  3913
差    162

生データの合計  39677 
データの平均値  39677 / 10 = 3967.7
差 / 平均値   162 / 3867.7 = 0.040829
乱数の発生頻度のバラツキ  4.08%


**********************************************

DM42での乱数精度テストプログラム

*01 LBL "AA" ■ [PGM.FNC] LBL "AA" ENT 先頭ラベルをAAとする

*02 123.456                 乱数の種データ

*03 SEED ■ [PROB] SEED      乱数の種データを設定する
*04 CLRG ■ [CLEAR] CLRG        数値変数 00~99をクリアする

*05 LBL "BB" ■ [PGM.FNC] LBL "BB" ENT  ループ先頭ラベルをBBとする

*06 RAN ■ [PROB] RAN 乱数を出力する  0~1

*07 10                 10をセット

*08 ×                  乱数を10倍  0~ 10

*09 STO II                IIを配列ポインターとして乱数をセット

*10 1                  乱数カウンターに加算する数値

*11 STO . IND + II           ポインターで振り分けられた所に1を加算する

*12 GTO "BB"             ループの先頭 BBへ行く


注1、 STO II まだ未定義なら  STO ENT II ENT の手順で定義する。
注2、 ソフトを終了させるには 時計で時間を計り [ R/S ] を押して停止させる。


( データ 2 )

DM42(電池) 1分間

0~1 3120
1~2 3108
2~3 3182
3~4 3119
4~5 3163
5~6 3053
6~7 3094
7~8 3078
8~9 3309
9~10 3050

最大値  3309
最小値  3050
差    259

生データの合計  31276 

データの平均値計算  31276 / 10 = 3127.6

差 / 平均値   259 / 3127.6= 0.0828

乱数の発生頻度のバラツキ  8.28%



**********************************************


( データ 3 )
DM42 ( USB 外部電源 10分間 )

0~1 101804
1~2 101900
2~3 101896
3~4 102420
4~5 101843
5~6 102249
6~7 102119
7~8 102431
8~9 102341
9~10 101548

最大値  102431
最小値  101541
差    883

生データの合計  1020551 

データの平均値計算  1020551 / 10 = 102055.1

差 / 平均値   883 / 102055.1= 0.0086521

乱数の発生頻度のバラツキ  0.865%


これで終わりです。




by telmic-gunma | 2019-04-03 22:06 | HP電卓 | Comments(0)