Casio Python のポテンシャル

Python Casio Python
 Casioグラフ関数電卓の Python を使ってみる
     - Casio Python のポテンシャル
目次


初版:2020/06/07
更新:2020/06/29
更新:2020/11/03
更新:2021/03/25

 次の記事 - 1. はじめに:電卓で作る初めてのスクリプト



カシオ・グラフ関数電卓にある Pythonモード - Casio Python

[2020/11/03 更新]
カシオグラフ関数電卓に Pythonモード <パイソン・モード> が搭載されるようになり、casioplot モジュールが追加されたので、それで遊んでみようと思います。Casio Basic やPCでのプログラミング経験があれば、Casio Python は意外に敷居が低いかも知れません。

管理人はネットでチョコチョコっと調べながら、結構使えています。そこで、 プログラミング経験はあるが Python 初心者の人向けに Casio Python のサンプルプログラム作成の連載をしばらくやってみようと思います。Casio Basic プログラムを Casio Python へ移植することからスタートし、Python らしいスクリプトの書き方の紹介へ進んでゆく予定です。併せて、Casio Python のリファレンス情報、特に Casio Python 特有の仕様があれば、それを公開してゆこうと思います。

本連載の全体像は 目 次 をご覧ください。

さて、Casio Python の連載は、管理人がそのポテンシャルの高さに驚いたことが全ての始まりなので、今回はそのあたりを紹介します。
==== 2020/11/03 更新ここまで ====


0. Casio Python のポテンシャル

Python の準備

Casio Python 搭載モデルの準備
Casio Python
搭載電卓を入手して、必要に応じてOSのアップデートをします。
現時点では casioplotモジュールが追加された fx-CG50 OS3.40 以降を強くお勧めします。

ちなみに、Casio Python が搭載されているのは、以下のモデルです (2020/06/07 現在);
fx-CG50, OS3.20 以降
fx-9750GIII, OS3.21 以降 (北米限定発売モデル)
fx-9860GIII, OS3.21 以降 (ヨーロッパ限定発売モデル)
GRAPH 90+E, OS3.20 以降/ GRAPH 35+EII, OS3.21 以降 (フランス専用モデル)
 それぞれが fx-CG50fx-9860GIII と同じ

[2020/10/22 追記]
Pythonモードで使い物になるスクリプトを書くにはグラフィックモジュール casioplot が不可欠です。casioplot は OSアップデートで追加するしかありませんが、それが可能なのは下記のモデルとOSバージョンです (2020/10/15 現在)
 fx-CG50, OS3.40以降     このモデルの概要 
fx-9750GIII, OS3.40以降    このモデルの概要 
fx-9860GIII, OS3.40以降    このモデルの概要 
 GRAPH 90+E, OS3.40以降
 GRAPH 35+EII, OS3.30以降

fx-CG50 のOSは、2020/10/15 時点で OS3.50 にアップデート可能です。OSアップデートは、それぞれのモデルが発売されている地域のカシオサイトからアップデートファイル と アップデートされたマニュアルが入手できます。具体的なアップデートサイトについては上記の "このモデルの概要" ページで紹介しています。

Pythonモードでは、casioplotモジュールを使わないと、まともに動作するスクリプトになりません。本連載は主に fx-CG50 OS3.40以降を前提にして進めてゆくつもりです。但し、モノクロ液晶搭載モデル fx-9750GIII や fx-9860GIII の Pythonモードに casioplot を追加するOSアップデートファイルが 2020/10/15 に公開されたので、これらをアップデートした Casio Python についても紹介する予定です (目次参照)。

 fx-CG50, OS3.40以降: 日本語マニュアルがあり国内どこでも入手できます。
 fx-9750GIII, OS3.40以降: 最安値のモデルです。国内でもAmazon USAで購入可で安価(英文マニュアル)。
fx-9860GIII, OS3.40以降: ヨーロッパ専用モデル。ほぼ同じ機能の fx-9750GIII をお勧めする。

モノクロ液晶搭載モデルのうち、現状最も安価に入週できる fx-9750GIII も無視できないとことがあります。これらモデルの Pythonモードはメモリの制限が強く、使い勝手が fx-CG50 よりも悪いのは残念です。

さて、Pythonモードのアップデートにより新機能や新モジュールが追加された場合は、随時 対応することにします。


Python スクリプト保存フォルダの作成
工場出荷状態では、ストレージメモリのルートにスクリプトファイルが保存されるようになっています。

これから多くのスクリプトファイルを作ってゆくので、ストレージメモリ直下に専用フォルダ、例えば @Python フォルダを作ると良いと思います。電卓内に新たなフォルダを作るには、PCとリンクしてからエクスプローラで電卓の内部を表示して、そこでフォルダを作ります。もし C.Basic を既にお使いの場合は、C.Basicプログラム専用フォルダを作り、そこへ移動しておくことも合わせてお勧めします。これで Python スクリプトと C.Basicプログラムが混在することを防げます。


Python スクリプトファイルの正体
スクリプトファイルは、アスキー文字で書かれたテキストファイルです。Python の構造制御に不可欠なインデントは半角スペースです。電卓上でスクリプトを書く場合、自動インデント機能によりスペース2個が付加されます。また、改行は Windows標準の CR LF (0x0d 0x0a)になっています。

PCに転送すれば、テキストエディタを用いて、スクリプトの確認や編集ができます。PCで編集する場合は、アスキーコード以外を入力しないように留意しましょう。



Casio Python とは? [2020/06/09 の夜に大幅改訂]

Casio Python は、カシオグラフ関数電卓にある Pythonモードを指します。ある程度の関数があらかじめビルトインされていますが、追加ビルトインするモジュールはカシオから提供されていません。自分で作った関数をモジュールとして使うことはできます。Casio Basic では自分で関数を作って使うことはできないので、その分 Casio Pythonは使い勝手が良いと言えます。

Pythonモードが公開された当初は、mathモジュールとrandimモジュールが追加された状態でした。その後 casioplotモジュールが追加されて、ようやくグラフィック描画が可能になりました (fx-CG50, fx-9750GIII, fx-9860GIII)。但し、ピクセル単位でドットを描画すること、ピクセルの色を得ること、テキスト文字を表示することのみができる程度で、Casio Basic のようにグラフを自在に描く機能はまた公開されていません。今後カシオからさらに新しいモジュールが提供されることを期待しています。汎用的に使える自作関数を複数まとめてモジュールにして、それを使うことも紹介してゆきます。

ところで、Casio Basic で書いたものはプログラムと言い、そのファイルをプログラムファイルと言っています。一方、Python については、取扱説明書の表現によれば、スクリプトを書いて、スクリプトファイルに保存する、と言っているので、本連載でもこの表現を使ってゆきます。

Python は機械学習やDeep Learningで多用されている言語で、"1つのことをするには1つの書き方しかない" といったポリシーで、分かりやすく使いやすいことを目指して作られた言語です。様々な専門分野に特化したライブラリ(=モジュール)が様々な人によって多く提供されるようになり、初心者でも専門家でも使いやすい言語して大きな広がりをみせています。また、Webサービスでも使われています。Pythonを使えるようになるのは、今後を考えると有意義なことだと思います。

Casio Python は組込用Pythonである Micro Python Ver 1.94 がライセンスされ、カシオがカスタマイズして電卓に組み込んだものなので、一般にPCで使われている Python (CPython) に比べて Casio Python は機能が限定されています。

実際に管理人がCasio Pythonを使い始めて感じるのは、提供されている機能が限定されているので、覚えることが極端に少なく、理解が容易だということです。カシオ電卓に搭載されている Python は、一般のPythonのサブセット版ではありますが、最近のソフトウェア技術のパラダイムに沿った基本仕様は損なわれていませんので、Python の味見用のお試し版、学校での学習用の位置づけではないかと思います。

Casio Basic を使いこなして作ったプログラムを Casio Python に移植しようとすると、機能不足、関数不足のために、どうしても移植できないプログラムがあります。Casio Basic と同じ事を実現するために必要不可欠な機能が Casio Python にはまだ揃っていないのは残念なことですが、それなりに工夫すれば使いこなせると思います。
例えば、Casio BasicGetkeyLocate に相当するものが Casio Python にはありません。
Getkey ⇒ これは、どうにもなりません。カシオさんお願い...
Locate ⇒ さらに機能追加した locate() 関数を自分で作り使用中 ⇒ ユーザー関数 locate()

ところで、Python は本来オブジェクト指向言語です。ライブラリ(=モジュール)を呼び出して使う時に、オブジェクト指向特有のobject.method といった記述が必要なケースもありますが、オブジェクトとメソッドをあまり意識しないで、関数型言語のように使える側面もあると思います。特に Casio Python ではほとんど関数型言語のように使っても、とりあえずスクリプトが書けると思います。マニュアルの記述にもオブジェクトやメソッドという単語が殆ど出てこず、関数という表現が多用されています。



Casio Python の潜在能力 - fx-CG50 OS3.40以降

純正Casio Basic のプログラムから移植可能なものを 実際にCasio Python書き直してみました。すると、

Casio Basic に比べて実行速度が圧倒的に速い!

以下で紹介しますが、その爆速っぷりは、チョット驚きです!
俄然 Casio Python に興味を持ち、詳細を調べたくなりました。
先ずは、その速さを紹介します。


純正Casio BasicCasio Python の処理速度の比較 
Casio Python でグラフィックス描画、関数計算、加算計算をさせて、純正Casio Basic と比較してみます。

グラフィック描画の比較 - フラクタルでシダの葉

shida_py_s 

これは、以下の記事で取り上げたフラクタル図形の描画プログラムです。
fx-CG50 でフラクタル - シダの葉

そこで作成した 純正Casio Basic のプログラム SHIDACG.g3m を Python に移植しました。
Casio Basic プログラムファイルと Python スクリプトファイルのダウンロード

フラクタル図形描画時間の比較
純正Casio Basic40 分
Casio Python24.4 秒
アドイン版Casio Basic - C.Basic for CG Ver 1.45β43.0 秒
※ C.Basic for CG Ver 1.45βについては こちらを参照
※ C.Basic は、コードの最適化により、さらに10倍近く高速化可能です。
 比較のため下記の純正Casio Basic コードを走らせました。

Python は100倍速い!

Casio Basicオリジナルコードと Python スクリプトを以下に示します。できる限り忠実に再現していることが分かると思います。

shida_g3m_py_src3


関数計算の比較 - sin と cos の繰り返し計算

tstfinc2_g3m_py_src
TSTFUNC2.g3m と tstfunc2.py のダウンロード

N に1000 を入力した時 (1000回ループ) の関数計算の比較
純正Casio Basic6.5 秒
Casio Python0.6 秒
アドイン版Casio Basic - C.Basic for CG Ver 1.45β0.38 秒

Pythonでの関数計算が10倍速いことが分かります。


加算計算の比較 - 加算とループ

tstsum2_g3m_py_src2  
TSTSIM2.g3m と tstsum2.py のダウンロード

N に1000 を入力した時 (1000回ループ) の加算計算の比較
純正Casio Basic2.2 秒
Casio Python0.15 秒
アドイン版Casio Basic - C.Basic for CG Ver 1.45β0.02 秒

Pythonでの加算計算が15倍程度速い結果になりました。


Python が速い理由
加算計算や関数計算で10~15倍程度速いだけでなく、グラフィックス描画がさらに10倍程度速いことが分かりました。

ちなみに、シダの葉描画のPhytonスクリプトの最下行にある show_screen() を for ループの外に出すと、VRAM上で全てのデータが揃ってから一気に画面転送を行うことになり、この時は 20.8 秒とさらに速くなります。つまり、show_screen() による画面転送速度が極めて速いことが、シダの場描画の爆速の原因と言えます。

複雑な計算やグラフィックスやグラフ表示には極めて役立つという、Casio Pythonの潜在能力が明らかになったと思います。今後のグラフィックおよびグラフ描画のためのモジュール追加が期待されます。

なお、シダの葉描画で利用した clear_screen() メソッドは、OS3.40では隠し機能でマニュアルやカタログ機能にありませんでしたが、OS3.50 で両方に載るようになりました。但し、マニュアルやカタログ機能で見つからない隠し関数がまだ多くあり、これらを探索する方法があります。
Casio Python とは? の "2.2 Casio Python に組み込まれているオブジェクト名の一覧" を参照してください。



Casio Python 探索

Casio Python の潜在能力が意外に高いことが分ったので、これから色々と詳細を調べながら探索を進めたいと思っています。

fx-CG50 や fx-9750GIII のマニュアルを読んでも、Python モードでの操作法は書かれていますが、関数やメソッドの使い方は殆ど書かれていません。幸いなことに今ではネットで簡単に調べることができるのではありますが、それでもクラス、メソッドなどで構成される言語に慣れていないと分かりにくいと思います。

そこで、Casio Python を探索しながら、記事をいくつか書こうと思います。



 目 次

次の記事 - 1. はじめに:電卓で作る初めてのスクリプト




応援クリックをお願いします。励みになるので...
にほんブログ村 IT技術ブログ 開発言語へ


 


keywords: fx-CG50Pythonfx-9750GIIIfx-9860GIIIプログラム関数電卓

リンク集 | ブログ内マップ

テーマ : プログラム関数電卓
ジャンル : コンピュータ

Casio Pyhton - 目次

Python Casio Python
 Casioグラフ関数電卓の Python を使ってみる
     
- 目 次


初版:2020/06/13
更新:2021/03/23

目 次

リファレンス


連載記事

Python を基礎から勉強するためではなく、高級言語や Casio Basic でのプログラミング経験がある人 (管理人もその部類) がとりあえずスクリプトを書いて (C# や java の経験があるととりあえず書けてしまうのが Python の良いところです)、正しく動くものを作ってみる目的で記事を書いています。これは理解のためのポイントだと管理人が思う項目は、逐次解説を入れています。

カシオ電卓のPythonモード(Casio Python) はPCなどで使える普通の Python とは異なり、ベースとなる micro Python やカシオ独自の仕様が反映されているので、実際にスクリプトを書いて動作させながら、仕様を調べるのが、本連載のもう一つの目的です。

Python を基本から学びたい方は、ネットに良いものがありますので (例えば Python入門 など)、そちらで学習されることをお勧めします。


0. Casio Python のポテンシャル

1. はじめに: 電卓で作る初めてのスクリプト
  1.1 Casio Python での初めてのスクリプト
   1.1.1 fx-CG50 の準備
   1.1.2 Casio Python の起動
   1.1.3 スクリプトファイルの作成
   1.1.4 スクリプトの作成
   1.1.5 スクリプトの実行
  1.2 シェル画面
  1.3 スクリプト作成と編集の2つの方法

2. Casio Python とは?:MicroPython と CPython
  2.1 Casio Python の素性
   2.1.1 MicroPython 1.9.4
   2.1.2 CPython 3.5.4
  2.2 Casio Python に組み込まれているオプジェクト名の一覧
  2.3 Casio Python の学習方法

3. Casio Python の入出力:プログラムのデザイン
  3.1 入力 - input()
  3.2 出力 - print() とグラフィックス出力関数
  3.3 Casio Pyhton で作るプログラム

4. Casio Python への移植:モンテカルロ法(1)
  4.1 コメントの書き方
  4.2 モジュールの呼び出し - import
  4.3 メニュー機能の作成
  4.4 circle() 関数の作成
  4.5 グラフィックス出力
   4.5.1 自作関数 circle() で円を描画
   4.5.2 draw_string() 関数で文字列を出力
   4.5.3 while 文により繰り返し
   4.5.4 random() 関数により正方形内の座標(x, y)を取得
   4.5.5 set_pixel() 関数で (x, y) に点を描画
   4.5.6 点の座標から円周率を計算
   4.5.7 draw_string() 関数で回数と円周率を表示
   4.5.8 draw_string() 関数でシミュレーション終了の表示

5. 関数の作成と活用:grp_color() 関数の作成
  5.1 関数の引数について
  5.2 関数内の変数の有効範囲 - 変数のスコープ
  5.3 circle() 関数の拡張
  5.4 色指定関数を追加
   5.4.1 RGBによる色指定
   5.4.2 タプル型
   5.4.3 色指定方法の拡張
   5.4.4 グラフィック画面での色指定関数 grp_color() 関数の作成
  5.5 circle() 関数拡張の書き換え
  5.6 モンテカルロ法(2) - monteca2.py

6. グラフィックス出力関数の追加:line() とユーザーモジュールの作成
  6.1 line() 関数の作成
  6.2 グラフィックスユーザーモジュールの作成
  6.3 line() 関数の動作確認 - ckLine.py
  6.4 line() 関数の動作確認 - ckLine2.py

7. テキスト出力関数の追加:locate() をユーザーモジュールに追加
  7.1 locate() 関数の作成
   7.1.1 フォントピッチ dx, dy の設定
   7.1.2 column から x を、row から y を算出
   7.1.3 フォントサイズ指定の処理
   7.1.4 locate() の完成 
  7.2 locate() をユーザーモジュールに追加
  7.3 locate() の動作確認 - ckLocate.py
  7.4 モンカルロ法(2) - monteca3.py

8. シェル画面入力の工夫:高速素因数分解(1)
  8.1 高速素因数分解プログラム - FactorG のざっくりした分析
  8.2 入力と入力値のチェック
  8.3 変数の初期化
  8.4 素数が素因数かどうかをチェック
  8.5 "エラストテレスの篩い" で素因数の探索
  8.6 ckpwr() 関数の作成
  8.7 disp() 関数の作成 

9. Python らしい反復処理:高速素因数分解(2)
  9.1 frac() 関数をユーザーモジュール (u.py) に追加
  9.2 不要な処理を削除
  9.3 disp() 関数のスリム化
  9.4 素数が素因数かどうかをチェックする部分のスリム化
  9.5 素因数を探索する残りの部分のスリム化
  9.6 スリム化したスクリプト - FactorG2.py

10. 大きな数の計算:高速素因数分解(3)
  10.1 変更箇所の洗い出し
  10.2 増加する素因数に合わせて表示桁数を拡張する - disp() の変更
  10.3 入力桁数の制限を拡張し、それに合わせてエラー表示を変更する
  10.4 素因数探索回数の表示を追加する
  10.5 素因数が10桁を超えるときの表示の変更 - disp() の変更

11. 関数呼出オーバーヘッド:高速素因数分解(4)
  11.1 関数呼出あり (ver 3) となし (ver 4) の実行速度の違い
  11.2 Casio Python の関数呼出オーバーヘッド

12. 要素数の多いリスト:高速素因数分解(5)
  12.1 増分リストを拡張する
  12.2 要素数480個の探索数増分リスト
  12.3 探索数を拡大して15桁対応高速素因数分解のスクリプト
  12.4 探索数拡大による高速化
  12.5 増分リストをさらに拡張する

13. 10進数除算の出力と精度:高速素因数分解(7)
  13.1 Casio Python での除算結果のデータ型
  13.2 Casio Python での関数計算結果のデータ型
  13.3 Casio Python での除算計算の精度
  13.4 0 および負の整数入力の処理
  13.5 入力ルーチンの修正
   13.5.1 入力委が小数点以下 0 の小数の場合に処理を追加 
   13.5.2 入力値が負の数値の場合の処理
   13.5.3 入力が 0 の時の処理
  13.6 完成した FactorG7.py の動作確認

14. CGモデルとFXモデルのPythonモードの違い
  14.1 CGモデルからFXモデルへのスクリプトの移植
   14.1.1 CGモデルからFXモデルへのスクリプト移植のポイント
   14.1.2 CGモデルとFXモデルの判定方法
  14.2 FXモデルのPythonモードの制限 - バッファサイズ
   14.2.1 編集できるスクリプトファイルの行数制限
   14.2.2 print() 出力のスタック・バッファサイズ
   14.2.3 Shell 画面出力のスタック・バッファサイズ
  14.3 Casio Python自体の機能の違い
   14.3.1 CGモデル OS3.40 と OS3.50 の違い
   14.3.2 CGモデル OS3.50 と FXモデル OS3.40 の違い

15. RGB値による色設定
  15.1 RGB値の色を確認するプログラム
  15.2 fx-CG50 の高精細カラー液晶 - 実は 16bit カラー
  15.3 モノクロ液晶モデルでのRGB値による色設定
  15.4 RGB値でのピクセル色設定と読取りの実験 - CGモデルとFXモデル

16. circle() 関数のFXモデルへの拡張
  16.1 CGモデルとFXモデルの場合分け
  16.2 circle() 関数を使ってみる


17. shell 画面とグラフィック画面の活用:コラッツ問題
  17.1 コラッツ問題 (Collatz Problem) とは?
  17.2 コラッツ数列の計算
  17.3 出力の工夫ポイント
   17.3.1 input() での入力
   17.3.2 フォントサイズを変えて添え字の描画
   17.3.3 リアルタイムに変わる数値の表示
   17.3.4 グラフの描画
   17.3.5 数列の最大値を求めて表示
   17.3.6 shell 画面でのスクロール出力
  17.4 コラッツ問題で遊んでみる



リファレンス

Casio Python - オプジェクト一覧 (← をクリック)
Casio Python に実装されているものを網羅しています。
  • オブジェクト一覧には、個別に作成したリファレンスへのリンクがあります。
  • 個別のリファレンスには、Python 公式サイトでの解説へのリンクを掲載しています。
  • 個別のリファレンスには、Casio Python 独自の項目を優先して記述します。
  • 一旦作成したものは、随時追記。修正します。

e-Gadgetで作成したユーザーモジュール - オブジェクト一覧 (↓下記)
e-Gadget で作成したユーザーモジュールに含まれているもの

最新バージョン: ver 1.5 <ダウンロード>
サポートモデル / OS:
  ・CGモデル - fx-CG50 / OS3.40以降
  ・FXモデル - fx-9750GIII
/ OS3.40以降
  ・FXモデル - fx-9860GIII
/ OS3.40以降

ファイル名: u.py / ufx.py
  全サポートモデルに対応 - CGモデル/FXモデルの自動判定機能あり
  ・ufx.py は FXモデルで開くために、u.py から一部コメントを削除して150行以下にしたもの
実装オブジェクト
  isCG(): CGモデルかFXモデルかを判定
  grp_color(): グラフィック描画関数の様々な色設定をタプル型に変換
  circle(): グラフィックス画面に円を描画
  line(): グラフィック画面に線分を描画
  locate(): グラフィック画面にテキストを出力
  ・frac(): 浮動小数点型数値の小数部を取得




応援クリックをお願いします。励みになるので...
にほんブログ村 IT技術ブログ 開発言語へ


 


keywords: fx-CG50Pythonfx-9750GIIIfx-9860GIIIプログラム関数電卓

リンク集 | ブログ内マップ

テーマ : プログラム関数電卓
ジャンル : コンピュータ

Casio Python - Shell 画面とグラフィック画面の活用:コラッツ問題

Python Casio Python
 Casioグラフ関数電卓の Python を使ってみる
     - Shell 画面とグラフィック画面の活用:コラッツ問題 
目次


初版:2020/12/12


前の記事 - 16. circle() 関数のFXモデルへの拡張 |  次の記事 -


17. Shell 画面とグラフィック画面の活用:コラッツ問題

<fx-CG50 OS3.40以降にに対応>

Casio Python は、入出力に極めて大きな制限があります。そこで今回は、Casio Python で shell 画面とグラフィック画面 (描画画面) をうまく活用して、出力方法を工夫してみます。そのための題材として、コラッツ問題を取り上げます。


17.1 コラッツ問題 (Collatz Problem) とは?

数学の未解決問題の1つにコラッツ問題 (コラッツ予想) があります。

コラッツ問題
任意の正の整数 n に対して、関数 f(n) を以下のように定義する。
  f(n) = 3n+1  (n が奇数のとき)  
      n/2    (n が偶数のとき)
この関数を使って、繰り返し演算を続けて、以下のように数列 ai を作る。
 ・a0 = n
 ・ai = f(ai-1)
どのような 正の整数 n のときでも、数列は 1 に到達する。


要素が1になるまでの数列 ai (i≧0) をコラッツ数列と言います。小学生レベルの四則演算で計算できる一見単純そうな問題です。コンピュータを使って実際に計算して要素が1に到達しない例がまだ見つかっていません。5x260 まで計算されているそうです。しかし、数学的にまだ証明されていません。

数学的に証明されていない限り、もっと大きな整数でコラッツ予想が成り立たなくなる例が見つかるかも知れません。

実際に計算してみます。

(1) n=6 のときのコラッツ数列 (i≧0)
ai = 6, 3, 10, 5, 16, 8, 4, 2, 1
 最初に ai=1 になるとき i = 8 

(2) n=25 のときのコラッツ数列 (i≧0)
ai = 25, 76, 38, 19, 58, 29, 88, 44, 22, 11, 34, 17, 52, 26, 13, 40, 20, 10, 5, 16, 8, 4, 2, 1
 最初に ai=1 になるとき i = 23

(3) n=27 のときのコラッツ数列 (i≧0)
ai = 27, 82, 41, 124, 62, 31. 94, 47, 142, 71, 214, 71, 214, 107, 322, 161, 484, 242, 121, 364, 182, 91, 274, 137, 412, 206, 103, 310, 155, 466, 233, 700, 350, 175, 526, 263, 790, 395, 1186, 593, 1780, 890, 445, 1336, 668, 334, 167, 502, 251, 754, 377, 1132, 566, 283, 850, 425, 1276, 638, 319, 958, 479, 1438, 719, 2158, 1079, 3238, 1619, 4858, 2429, 7288, 3644, 1822, 991, 2734, 1367, 4102, 2051, 6154, 3077, 9232, 4616, 2308, 1154, 577, 1732, 866, 433, 1300, 650, 325, 976, 488, 244, 122, 61, 184, 92, 46, 23, 70, 35, 106, 53, 160, 80, 40, 20, 10, 5, 16, 8, 4, 2, 1
 最初に ai=1 になるとき i = 111 

n1 から 26 までは、手計算でも順調にできますが、n=27 になると途端に大変になり、9232 まで大きくなってしまうと、本当に 1 になるのか不安になりますが、結果として 1 に到達できました。

より詳しくは "コラッツ問題" や "コラッツ予想" で検索すれば、色々な情報が得られます。 

管理人としては、Casio Python で、
 ・計算しながら、数列の要素をグラフ表示させ、
 ・最後は、コラッツ数列の全ての要素を出力させる
そのようなスクリプトを作りたくなりました。

計算中は、数列の i の値と ai の値 (数列の要素の値) がどんどん変化している様子をリアルタイムで表示させながら、ai の値の推移をグラフ表示させ、計算が終われば ai の最大値も表示させます。最後に算出した全ての ai (数列の全ての要素) を出力させます。

今回は、以下のような動作をする Casio Python スクリプト - Collatz.py を作ります。




 Collatz.py  ダウンロード

"""Sample script

 Collatz Problem Simulation ver1.0;
 - Graphical and numerical output of
  Collatz sequence a0, a1, a2 ... an

 Collatz.py
  for fx-CG50 OS3.40 or later
   by Krtyski/e-Gadget
"""


from u import *

#Input integer >1
print('\n\n\n== Collatz Problem ==\n\n')
n=int(input('Input number:'))

#Initialize variables
offset=150
scale=50
a=[0]; a[0]=ni=0

#Draw initial screen
locate(0, 0, 'a :', size='l')
locate(3, 1, '0', size='s')
locate(3, 0, a[0], color='blue', size='l')
locate(0, 1, 'i :', size='l')
locate(0, 2, 'a :', size='l')
locate(3, 5, 'i', size='s')
locate(03, 'a :', size='l')
locate(2, 7, 'max', size='s')
line(offset, 191, 383, 191)

#Calculate and draw element of sequence
while a[i]!=1:
 #Increment index
 i+=1

 #Erase previous numbers
 locate(3, 1, i-1, color=0, size='l')
 locate(3, 2, a[i-1], color=0, size='l')

 #Calculate next element of sequence
 if a[i-1]%2:
  a.append(3*a[i-1]+1)
 else:
  a.append(a[i-1]//2)

 #Draw updated numbers
 locate(3, 1, i, size='l')
 locate(3, 2, a[i], size='l')

 #Draw (update) graph of sequence
 if a[i]//scale<192:
  line(offset+i-1, 191-a[i-1]//scale, offset+i, 191-a[i]//scale, (2550, 0))

#Draw finally maximum element
locate(3, 3, max(a), size='l')

#Output all elements of sequence
locate(0, 7, '[EXIT]', size='l')
for i in range(i+1):
 print('a' + str(i) + '=' + str(a[i]))


このスクリプトを見て、何をやっているのかスグに判る場合は、Collatz.py を実行して遊んでみて下さい。

これ以降は、Cllatz.py の解説を行い、少し遊んでみた結果を紹介します。



17.2 コラッツ数列の計算


今回は、コラッツ数列をリストで表現することにします。

つまり、
 数列:ai = a0, a1, a2, ... ,ai-1, ai

 リスト:a = [a0, a1, a2, ... ,ai-1, ai]

とすると、i 番目の要素は a[i] で得られます。

スクリプトの冒頭で、リスト a の定義 (a=[0]) を行います。この定義では、要素が 数値 0 が 1個しかないリストを定義しています。計算が終わらないと要素数が判らないので、次の要素を算出したときに、.append() 関数を使ってリストにその要素を追加し、i を1つづつ増やすようにします。

1つめの要素は、i=0 のときであり、その値は n なので、a[0]=n とします。
すると、リスト a[i] は以下のようにして求められます。

a=[0]; a[0]=ni=0

while a[i]!=1:
 #Increment index
 i+=1

 #Calculate next element of sequence
 if a[i-1]%2:
  a.append(3*a[i-1]+1)
 else:
  a.append(a[i-1]//2)


これを解説します;

a[i]
を求めるには、1つ前の a[i-1] を使って、

a[i-1] が奇数のとき、
 言い換えると、a[i-1] を 2で割ったときの余りが 0 でない時、
 つまり a[i-1]%2 が 0 でない時
  3 × a[i-1] +1 を リスト a の末尾に追加する、つまり
  a.append(3*a[i-1]+1)

・そうでない時、つまり a[i-1] が偶数のとき
  a[i-1] / 2 をリスト a の末尾に追加する、つまり
  a.append(a[i-1]//2)
  ※ ここで、リストの要素を整数にしたいので、// 演算子 (割り算の商 = 整数部を得る)を使う

この処理を a[i] が 1 でない限り繰り返せば良いので、上記の処理を whilr a[i] != 1: のループの中に記述しています。


17.3 出力の工夫ポイント

17.3.1 input() での出力
最初に1以上の整数を入力させるには、Casio Python に備わっている唯一の入力関数 input() を使います。
input() を使うと、キー入力した文字列や数値は、全て文字列として返ります。今回は、入力した数値を 数値変数 n に代入したいので、input() の戻り値を int() 関数で整数に変換します。

n=int(input())

入力時に表示する文字列は input()( ) 内に指定します。

n=int(input('Input number:'))

さて、今回は以下のような入力画面にしてみます。
Collatz_input 

チョット格好良いですね!
ちなみに、Collatz Problemコラッツ問題 のことです。

このような配置にするには、改行を利用します。文字列の中に \n を入れたら、そこで改行します。

そこで、入力画面は以下のようにします。

print('\n\n\n== Collatz Problem ==\n\n')
n=int(input('Input number:'
))



17.3.2 フォントサイズを変えて添え字の描画
コラッツ数列 ai は、添え字を使っているので、出力画面でも添え字が表現できると格好イイですよね!
そこで、以下のような表示にしてみました。

Collatz_27_1    

ここでは、a0aiamax を表示しています。

実は最初は以下のような表示にしましたが、添え字の 0i 、特に i が見づらいので却下。上のように見やすさを優先しました。

Collatz_27_2   


グラフィック画面では、大、中、小 の 3つのサイズのフォントが使えます。フォントピッチ (フォントの大きさ) を考えると、大フォントのピッタリ 1/2 が小フォントになっているので、大フォントの添え字に小フォントを使うとうまくいきそう。

ユーザー関数 locate() を作っているので、これを使います。スクリプトはチョット面倒に見えますが、コンセプトは簡単なので、地道に作ってゆきます。

 ユーザー関数 locate() の説明
 ※ locate() を使うには、ユーザーモジュール u.py をインポートします。
   from u import *


結果は、以下のようになります。

locate(0, 0, 'a :', size='l')
locate(3, 1, '0', size='s')
locate(3, 0, a[0], color='blue', size='l')
locate(0, 1, 'i :', size='l')
locate(0, 2, 'a :', size='l')
locate(3, 5, 'i', size='s')
locate(03, 'a :', size='l')
locate(2, 7, 'max', size='s')


ところで、コラッツ数列 ai を、今回のスクリプトでは リスト a で扱っていて、コラッツ数列の定義から、a0n なので、
 a[0] = n
となります。

上の locate() が並んだスクリプトの、上から3番目に a[0] で出てきていますが、これは n に置き換えても同じです。そして、このときだけ 青色で表示していることは、color='blue' がパラメータに入っていることで判ります (ユーザー関数 locate() の仕様)。


17.3.3 リアルタイムに変わる数値の表示
これを実現するには、グラフィック画面 (描画画面) に文字列を描画するしかありません。print() でshell 画面に出力するのでは、スクロールしてしまいます。

グラフィック画面で文字列を描画するには、casioplot モジュールに含まれる draw_string() 関数を使えば良いのですが、より使いやすくするために作ったユーザー関数 locate() を使います。

描画した文字列を新しい文字列で上書きするには、以下の手順が必要になります。
 (1) 色を指定して文字列を描画
 (2) 同じ位置に同じ文字列を白で再描画して、文字列を消去
 (3) 同じ位置に新たな文字列描画
手順 (2) が不可欠です。

locate() を使うには、ユーザーモジュール u.py をインポートします。
  from u import *


コラッツ数列を リスト a として得るための while ループの中に、i の値と a[i] の値をリアルタイムで描画するコードを追加すると、以下のようになります。

#Calculate and draw element of sequence
while a[i]!=1:
 #Increment index
 i+=1

 #Erase previous numbers ◀ この2行↓を追加
 locate(3, 1, i-1, color=0, size='l')
 locate(3, 2, a[i-1], color=0, size='l')

 #Calculate next element of sequence
 if a[i-1]%2:
  a.append(3*a[i-1]+1)
 else:
  a.append(a[i-1]//2)

 #Draw updated numbers ◀ この2行↓を追加
 locate(3, 1, i, size='l')
 locate(3, 2, a[i], size='l')



17.3.4 グラフの描画
横軸に i の値、縦軸に a[i] の値 をとって、ユーザー関数 line() を使ってグラフを描画します。

基本的には、点 (i-1, a[i-1]) と 点 (i, a[i]) を結んだ線分を描画することになりますが、画面全体の左上が原点になる座標系から、所定の長方形のエリアの中に収まるように変換して描画する必要があります。

今回は、点 (x, y) について以下のように変換します。

横座標 (x座標) は、
  左から 150 ドットの位置から右側を使うので、変数 offset = 150 として、
  offset だけ右にずらすので、xoffset+x に変換します。

縦座標 (y 座標) は、
  上下反転させ、変数 scale = 50 として 50分の1に圧縮して整数部を取得するので、
  y191-y//scale に変換します。 

ちなみに、// 演算は除算した時の商、つまり割り算結果の整数部を得るもので、これにより整数値が得られます。line() 関数に指定する座標パラメータは整数でなくてはならない(仕様です)ので、// 演算を用います。


以上から、
 点 (i-1, a[i-1]) は、点 (offset+i-1, 191-a[i-1]//scale)
 点 (i, a[i]) は、点 (offset+i, 191-a[i]//scale)
変換できます。

そして、今回はグラフを赤 (r, g,b) = (255, 0, 0) で描画することにしているので、line() 関数を使った線分描画は、以下のようになります。

 line(offset+i-1, 191-a[i-1]//scale, offset+i, 191-a[i]//scale, (255, 0, 0))


ところで、今回のグラフ描画エリアも画面に表示します。但し長方形のエリアではなく、エリアの下端のみを直線で示すことにします。そこで、line() 関数を使って、以下のようにします。
 line(offset, 191, 383, 191)


さて、今の line() 関数の仕様では、画面の外の座標を指定すると、エラーにならず見えない範囲に線分描画の動作を行うので、見えない線分描画に時間がかかってしまいます。そこで、今回指定したグラフ描画エリアの外の座標を設定する場合は描画をスキップするようにして、全体として高速化します。

具体的には、a[i]//scale の値が グラフ描画エリアの高さ(縦幅)に収まるときだけ、line() 関数で描画するようにしました。
 if a[i]//scale192:
  line(offset+i-1, 191-a[i-1]//scale, offset+1, 191-a[i]//scale, (255, 0, 0
))


Collatz_27_1  


ユーザー関数 line() の説明
  ※ line() を使うには、ユーザーモジュール u.py をインポートします。
    from u impot *  


offset=150            ◀ この行を追加
scale=50             ◀ この行を追加
a=[0]; a[0]=ni=0

line(offset, 191, 383, 191) ◀ この行を追加

#Calculate and draw element of sequence
while a[i]!=1:
 #Increment index
 i+=1

 #Erase previous numbers
 locate(3, 1, i-1, color=0, size='l')
 locate(3, 2, a[i-1], color=0, size='l')

 #Calculate next element of sequence
 if a[i-1]%2:
  a.append(3*a[i-1]+1)
 else:
  a.append(a[i-1]//2)

 #Draw updated numbers
 locate(3, 1, i, size='l')
 locate(3, 2, a[i], size='l')

 #Draw (update) graph of sequence ◀ この2行↓を追加
 if a[i]//scale<192:
  line(offset+i-1, 191-a[i-1]//scale, offset+i, 191-a[i]//scale, (2550, 0))



17.3.5 数列の値の最大値を求めて表示
冒頭で、整数 27 からコラッツ数列を求めたところ、最大 9232 とかなり大きな数になったのに、少し驚きました。そこで、求めた数列の最大の要素を表示したいと思いました。

そこで while ループ終了後に、a[i] の最大値を求めて表示します。
リスト a の要素から最大値を返す関数 max() を使えば、max(a) が最大値を返します。

ところで、max() 関数は、最大値が2つ以上有る時は、リストの中で一番左にあるものを返します。今回は最大値のみに感心があるので、max() 関数を使うだけで十分です。

この最大値を locate() 関数で、所定の位置に描画するには、以下のように記述します。
 locate(3, 3, max(a), size='l')

==========

#Initialize variables
offset=150
scale=50
a=[0]; a[0]=ni=0

#Calculate and draw element of sequence
while a[i]!=1:
 #Increment index
 i+=1

 #Erase previous numbers
 locate(3, 1, i-1, color=0, size='l')
 locate(3, 2, a[i-1], color=0, size='l')

 #Calculate next element of sequence
 if a[i-1]%2:
  a.append(3*a[i-1]+1)
 else:
  a.append(a[i-1]//2)

 #Draw updated numbers
 locate(3, 1, i, size='l')
 locate(3, 2, a[i], size='l')

 #Draw (update) graph of sequence
 if a[i]//scale<192:
  line(offset+i-1, 191-a[i-1]//scale, offset+i, 191-a[i]//scale, (2550, 0))

#Draw finally maximum element ◀ この1行を追加
locate(3, 3, max(a), size='l')



17.3.6 shell 画面でのスクロール出力
最後に、コラッツ数列の全要素を表示したいのですが、計算するまでは要素数が判らないので、グラフィック画面に出力するよりも、shell 画面でスクロールさせて出力するのが良いと思います。スクロールアップ/ダウンも自由にできるので、全要素を詳しく見るのに適しています。

冒頭の動画で判るように、グラフィック画面 と shell 画面を、適材適所でうまく使い分けると良いと思います。

shell 画面でリスト a の要素 a[0] から a[i] までを表示するには、

 for i in range(i+1):
  print(a[i])


とすればOKです。

ここでは、コラッツ数列の要素を表示するのが本来な目的なので、例えば i = 70 のときの ai を表示する場合は、

 a70=991

と表示しようと思います。そこで、スクリプトは以下のようになります。

 for i in range(i+1):
  print('a' + str(i) + '=' + str(a[i
]))



グラフィック画面 (描画画面) が表示されているとき、shell 画面に切り替える方法が2つあります。
[EXIT]キーを押す
input() を実行する

今回は、[EXIT] キーを押せば良いので、グラフィック画面のまま維持しておき、[EXIT] キーを押してもらってから、shell 画面に切り替わるようにします。そこで、操作ガイドの意味で、グラフィック画面の左下に [EXIT] と表示しようと思います。
 locate(0, 7, '[EXIT]', size='l')

これに続いて、上記の リスト a の要素表示のスクリプトを追加すれば、うまくゆきます。

最後に Casio Python の制限について、付記しておきます。fx-CG50 の Pythonモードでは、print() 関数の出力は最大 200行に制限されています。コラッツ数列 (= リスト a) の要素数が 200 個を超えると、全ての要素が表示されません。

色々と整数を変えて遊んでいる限りは、今のところ 200 行を超えるケースがないので、なんとかなりそうです。以下の動画では、整数が20000以下のときは、要素数 200 を超える場合が結構少ないことが判ります。





もし 200 行を超えるものが沢山でてくるならば、チョット面倒ですが工夫次第でなんとかできると思います。

==========

全要素表示の部分を最後に追加し、今回検討したスクリプト全体を以下に示します。

from u import *

#Input integer >1

print('\n\n\n== Collatz Problem ==\n\n')
n=int(input('Input number:'))

#Initialize variables
offset=150
scale=50
a=[0]; a[0]=ni=0

#Draw initial screen
locate(0, 0, 'a :', size='l')
locate(3, 1, '0', size='s')
locate(3, 0, a[0], color='blue', size='l')
locate(0, 1, 'i :', size='l')
locate(0, 2, 'a :', size='l')
locate(3, 5, 'i', size='s')
locate(03, 'a :', size='l')
locate(2, 7, 'max', size='s')
line(offset, 191, 383, 191)

#Calculate and draw element of sequence
while a[i]!=1:
 #Increment index
 i+=1

 #Erase previous numbers
 locate(3, 1, i-1, color=0, size='l')
 locate(3, 2, a[i-1], color=0, size='l')

 #Calculate next element of sequence
 if a[i-1]%2:
  a.append(3*a[i-1]+1)
 else:
  a.append(a[i-1]//2)

 #Draw updated numbers
 locate(3, 1, i, size='l')
 locate(3, 2, a[i], size='l')

 #Draw (update) graph of sequence
 if a[i]//scale<192:
  line(offset+i-1, 191-a[i-1]//scale, offset+i, 191-a[i]//scale, (2550, 0))

#Draw finally maximum element
locate(3, 3, max(a), size='l')

#Output all elements of sequence ◀ この3行↓を追加
locate(0, 7, '[EXIT]', size='l')
for i in range(i+1):
 print('a' + str(i) + '=' + str(a[i]))


17.4 コラッツ問題で遊んでみる

適当に入力する整数を選んで、コラッツ数列の要素数 i が大きい順に並べてみました。ここでは、170 ⇒ 160 ⇒ 151 を例に挙げています。当然これより大きな要素数になるケースは間違いなくあります。

Collatz_4321 Collatz_5433 Collatz_2345 

グラフの縦軸は、縮尺率 scale = 50 にしているから、191 x 50 = 9550 が表示できる amax の上限です。これを超える場合は表示されていません。上の3つの例では、いずれも amax は 9550 以上になっています。

これに続く要素数は、例えば以下のように 146 ⇒ 144 ⇒ 136 になりました。

Collatz_865 Collatz_654 Collatz_543 

要素数 i に関わらず要素の最大値 amax9232 になっています。最初に与える整数 a0 (=リストa[0]) が 27 の時も amax = 9232 で同じなので、数列要素の計算の終盤の計算パターンが同じになっているようですね。amax が 9550 以下なので、グラフは描画エリア内に収まっています。

さらに、これに続く要素数は、例えば以下のように 132 ⇒ 127 ⇒ 125 なりました。

Collatz_1234 Collatz_235 Collatz_345 

要素の最大値 amax9232 になる計算パターンは、a0 (=a[0]) が小さい時に多く見られるようです。

そして、amax9232 になるとき、a0 (=a[0]) の最小値は 27 なのは簡単に確認できます。その時の要素数 i111 です。

Collatz_27_3 

コラッツ問題 / コラッツ予想 で計算するコラッツ数列は、なんとも面白い動きをします。いずれ誰かが数学的に証明するのでしょう。その時が楽しみです。




目 次

前の記事 - 16. circle() 関数のFXモデルへの拡張

次の記事 -





応援クリックをお願いします。励みになるので...
にほんブログ村 IT技術ブログ 開発言語へ


 


keywords: fx-CG50Pythonfx-9750GIIIfx-9860GIIIプログラム関数電卓

リンク集 | ブログ内マップ

テーマ : プログラム関数電卓
ジャンル : コンピュータ

Casio Python - ユーザー関数 locate()

Python Casio Python
 Casioグラフ関数電卓の Python を使ってみる
     - リファレンス 
目次
locate()

ユーザーモジュール u.py に含まれる関数 u.py ver 1.5 ダウンロード

初版:2020/12/09

[対応モデル] - fx-CG50 OS3.20 以降、fx-9750GIII / fx-9860GIII OS3.21 以降

文字(列)  や 数値を、桁と行で位置を指定してグラフィック画面に描画します。
併せて、描画の色とフォントサイズを指定することもできます。
さらに、画面に転送するか、VRAMのみに転送するかの選択も可能です。

書 式locate(column, row, obj, color=1, size='m', show=1)

引 数
- 第1引数 - column: 描画する文字(列) / 数値 の座標 (桁, 行) の 桁座標、0 以上の整数
- 第2引数 - row: 描画する文字(列) / 数値の座標 (桁, 行) の 行座標、0 以上の整数
- 第3引数 - obj: 描画する文字(列) / 数値
- 第4引数 - color: 色を指定するパラメータ引数で、省略可能でデフォルトは 1 です。
       具体的な設定は grp_color() のリファレンスを参照。
- 第5引数 - size: フォントサイズを指定するパラメータ引数で、'l', 'm', 's' のいずれか。
       これは省略可能で、書略時は 'm' になります。
       'l' - 大フォント、'm' - 中フォント、's' - 小フォント
- 第6引数 - show: 1 のときは液晶に出力し、0 のときはVRAMへ転送するだけです。
       これは
パラメータ引数で、省略時は 1 になります。

引数を5つ、値だけを指定するとエラーになります ⇒ 関数の引数は こちら を参照。
引数を5つだけ設定する場合は、パラメータと共に設定します。
  例) line(5, 50, 30, 60, color=4) / line(5, 50, 30, 60, show=0)
フォントサイズの詳細 (CGモデル)
   ・'l' (大): 桁は 0~21、行は 0~7  | フォントピッチは 横16 dot、縦24 dot
   ・'m' (中): 桁は 0~31、行は 0~11 | フォントピッチは 横12 dot、縦16 dot
   ・'s' (小): 桁は 0~47、行は 0~17 | フォントピッチは 横 8 dot、縦12 dot
 フォントサイズの詳細 (FXモデル)
   ・'l' (大): 桁は 0~20、行は 0~7 | フォントピッチは 横 6 dot、縦 8 dot
   ・'m' (中): 桁は 0~31、行は 0~9 | フォントピッチは 横 4 dot、縦 6 dot
   ・'s' (小): 'm' と同じ (↑)


関数定義:

from u import *

def locate(column, row, obj, color=1, size='m', show=1):
 #set fontsize & pitch
 #from 4th argument(size)
 if isCG(): #CG model
  dx={'s':8, 'm':12, 'l':16}
  dy={'s':12, 'm':16, 'l':24}
 else: #FX model
  dx={'s':4'm':4, 'l':6}
  dy={'s':6, 'm':6, 'l':8}
 sz={'s':'small', 'm':'medium', 'l':'large'}

 #data transfer to VRAM
 draw_string(column*dx[size], row*dy[size], str(obj), grp_color(color), sz[size])

 if show: #data transfer to screen
  show_screen
()



スクリプトの解説:
  テキスト出力関数の追加:locate() をユーザーモジュールの追加 ⇒ こちら
  CGモデルとFXモデルの判定方法について ⇒ こちら
 ▶ grp_color() ユーザー関数 ⇒ こちら



応援クリックをお願いします。励みになるので...
にほんブログ村 IT技術ブログ 開発言語へ


 


keywords: fx-CG50Pythonfx-9750GIIIfx-9860GIIIプログラム関数電卓

リンク集 | ブログ内マップ

テーマ : プログラム関数電卓
ジャンル : コンピュータ

Casio Python - 要素数の多いリスト:高速素因数分解(5)

Python Casio Python
 Casioグラフ関数電卓の Python を使ってみる
     - 要素数の多いリスト:高速素因数分解(5) 
目次


初版:2020/08/03
大幅な追記修正:2020/08/06
Windowsアプリ修正:2020/08/13
追記修正:2020/11/14

前の記事 - 11. 関数のオーバーヘッド |  次の記事 - 13. 10進数除算の出力と精度


12. 要素数の多いリスト:高速素因数分解(5) <fx-CG50 OS3.4以降>

前回、関数呼出を無くして、15桁対応素因数分解を少し高速化しました。

高速素因数分解 15桁対応、僅かに高速化版 - FactorF4.py のダウンロード

以下のように素因数が15桁になるケースでは、実行時間が 23分30秒 もかかり、素因数探索のループを回る回数が 435万回程度になります (search:4349444)。
15digit_prime 

ここでは以下のロジックで探索数を決めています。今回は、同じロジックで、探索数を拡張する方法を検討します。

要素数48個の探索数増分リストを用いる方法
現在のロジックは、探索数の重複を効率よく減らせる点に特徴があり、
a) 素数 2, 3, 5, 7, 11 を探索数として除算を繰り返し、
b) 13222 の整数で 2, 3, 5, 7 のいずれかの倍数でない探索数で除算を繰り返し、
エラストテレスの篩いで素因数を調べる、という手法を採用しています。


12.1 増分リストを拡張する
同じロジックで、探索数を拡張してみます。それにより探索の無駄をさらに減らして、高速化を試みます。
要素数48個の増分リストを以下のように拡張します。

 探索数増分リストを拡張する方法
A) 素数 2, 3, 5, 7, 11, 13 を探索数として除算を繰り返し、
B) 172326 の整数で 2, 3, 5, 7, 11 のいずれかの倍数でない探索数で除算を繰り返し、
エラストテレスの篩いで素因数を調べる、


素因数探索で取りこぼしを防ぐために 172326 までの整数のうち 2, 3, 5, 7, 11 のいずれかの倍数でない整数を用い、その探索数は 480 個あります。 探索数増分リスト increment[ ] の要素は 480 個になります。スクリプトとして書き出す際は、1行8要素として、リストの定義に60行必要になります。

探索数を 17~2326 の範囲で探し、それが480個になる詳しい説明は、fx-5800P 素因数分解 - 高速化 を参照してください。
簡単にまとめると、
  • 2, 3, 5, 7, 11 の最小公倍数は 2310
  • 2, 3, 5, 7, 11 のいずれかの倍数でない探索数の増分は 2310 の周期
  • 探索数を17から調べるので、1周期の終わりは 2326
  • 2, 3, 5, 7, 11 のいずれの倍数でない整数の個数:480個
  • これの算出方法は以下の通り;
  •  2 の倍数の個数:2310/2 = 1155
  •  3 の倍数の個数:2310/3 = 770
  •  5 の倍数の個数:2310/5 = 462
  •  7 の倍数の個数:2310/7 = 330
  •  11 の倍数の個数:2310/11 = 210 ⇒ 以上の合計:2927個
  • .
  •  2と3の公倍数の個数:2310/6 = 385
  •  2と5の公倍数の個数:2310/10 = 231
  •  2と7の公倍数の個数:2310/14 = 165
  •  2と11の公倍数の個数:2310/11 = 105
  •  3と5の公倍数の個数:2310/15 = 154
  •  3と7の公倍数の個数:2310/21 = 110
  •  3と11の公倍数の個数:2310/33 = 70
  •  5と7の公倍数の個数:2310/35 = 66
  •  5と11の公倍数の個数:2310/55 = 42
  •  7と11の公倍数の個数:2310/77 = 30 ⇒ 以上の合計:1358個
  • .
  •  2と3と5の公倍数の個数:2310/2/3/5 = 77
  •  2と3と7の公倍数の個数:2310/2/3/7 = 55
  •  2と3と11の公倍数の個数:2310/2/3/11 = 35
  •  2と5と7の公倍数の個数:2310/2/5/7 = 33
  •  2と5と11の公倍数の個数:2310/2/5/11 = 21
  •  2と7と11の公倍数の個数:2310/2/7/11 = 15
  •  3と5と7の公倍数の個数:2310/3/5/7 = 22
  •  3と5と11の公倍数の個数:2310/3/5/11 = 14
  •  3と7と11の公倍数の個数:2310/3/7/11 = 10
  •  5と7と11の公倍数の個数:2310/5/7/11 = 6 ⇒ 以上の合計:288個
  • .
  •  2と3と5と7の公倍数の個数:2310/2/3/5 = 11
  •  2と3と5と11の公倍数の個数:2310/2/3/5/11 = 7
  •  2と3と7と11の公倍数の個数:2310/2/3/7/11 = 5
  •  2と5と7と11の公倍数の個数:2310/2/5/7/11 = 3
  •  3と5と7と11の公倍数の個数:2310/3/5/7/11 = 2 ⇒ 以上の合計:28個
  • .
  •  2と3と5と7と11の公倍数の個数:2310/2/3/5/7/11 = 1個
  • .
  •  以上から 2927 - 1358 + 288 - 28 + 1 = 480個
この480個の増分リスト increment[ ] の要素は、エクセルを使って半手動で求めることができますが、計算の正確さ、計算結果をスクリプトに転記する作業の正確さ(これが最大の目的)、そして時間の節約のために、アプリ SearchNum11.exe を Visual Studio 2019 Community の C# で作りました。C# はマウスクリックだけで殆どを自動作成してくれ、全部で30行程度のコーディングで出来ました。

要素数480の探索数増分リストの要素を出力するWindowsデスクトップアプリ
 - SearchNum11 のダウンロード


ダウンロードファイルには、アプリの実行ファイル (SearchNum11.exe) だけでなくソースコードも同梱していますので、ご自由にお使いください。

SearchNum11.exe を起動したところ:
SearchNum11_1

右上の [Calculate] ボタンを押すと、各要素とコンマ ,  がテキストボックスに出力されます。
SearchNum11_result
出力されると [Calculate] ボタンを押せなくしています。

1行あたり 8要素で改行するように出力しています。
テキストボックスでは、改行したり編集ができます。そこで、これをカットアンドペーストにより、エディタで開いた スクリプトに貼り付ければ、増分リストの要素を簡単に書き込めます。

例えば、上のテキストボックスの中でクリックしてから、[Ctrl]+[A] で全部を選択してから、[↑] キーを押して 最後の 12, のコンマの前まで (12 まで) 選択した状態で、[Ctrl]+[C] でコピーしてから、

エディタで、
increment=/
[ ]


] の中にカーソルを移動して、[Ctrl]+[V] とすればペーストできます。 


12.2 要素数480個の探索数増分リスト
先ずは、480個の探索数増分を要素数にするリスト increment[ ] を採用して、どの程度の実行速度改善が見られるのか、調べて見ます。

prime_list=[2,3,5,7,11,13]
for b in prime_list:
 search+=1
 d=a/b
 if frac(d)==0:
  e+=1
  z[e]=b
  while 1:
   a=int(d)
   z[e+21]+=1
   d=a/b
   if frac(d):
    break
  c=int(sqrt(a))
 if b>c:break


increment=\
[4,2,4,6,2,6,4,2,
 4,6,6,2,6,4,2,6,
 4,6,8,4,2,4,2,4,
 14,4,6,2,10,2,6,6,
 4,2,4,6,2,10,2,4,
 2,12,10,2,4,2,4,6,
 2,6,4,6,6,6,2,6,
 4,2,6,4,6,8,4,2,
 4,6,8,6,10,2,4,6,
 2,6,6,4,2,4,6,2,
 6,4,2,6,10,2,10,2,
 4,2,4,6,8,4,2,4,
 12,2,6,4,2,6,4,6,
 12,2,4,2,4,8,6,4,
 6,2,4,6,2,6,10,2,
 4,6,2,6,4,2,4,2,
 10,2,10,2,4,6,6,2,
 6,6,4,6,6,2,6,4,
 2,6,4,6,8,4,2,6,
 4,8,6,4,6,2,4,6,
 8,6,4,2,10,2,6,4,
 2,4,2,10,2,10,2,4
 2,4,8,6,4,2,4,6,
 6,2,6,4,8,4,6,8,
 4,2,4,2,4,8,6,4,
 6,6,6,2,6,6,4,2,
 4,6,2,6,4,2,4,2,
 10,2,10,2,6,4,6,2,
 6,4,2,4,6,6,8,4,
 2,6,10,8,4,2,4,2,
 4,8,10,6,2,4,8,6,
 6,4,2,4,6,2,6,4,
 6,2,10,2,10,2,4,2,
 4,6,2,6,4,2,4,6,
 6,2,6,6,6,4,6,8,
 4,2,4,2,4,8,6,4,
 8,4,6,2,6,6,4,2,
 4,6,8,4,2,4,2,10,
 2,10,2,4,2,4,6,2,
 10,2,4,6,8,6,4,2,
 6,4,6,8,4,6,2,4,
 8,6,4,6,2,4,6,2,
 6,6,4,6,6,2,6,6,
 4,2,10,2,10,2,4,2,
 4,6,2,6,4,2,10,6,
 2,6,4,2,6,4,6,8,
 4,2,4,2,12,6,4,6,
 2,4,6,2,12,4,2,4,
 8,6,4,2,4,2,10,2,
 10,6,2,4,6,2,6,4,
 2,4,6,6,2,6,4,2,
 10,6,8,6,4,2,4,8,
 6,4,6,2,4,6,2,6,
 6,6,4,6,2,6,4,2,
 4,2,10,12,2,4,2,10,
 2,6,4,2,4,6,6,2,
 10,2,6,4,14,4,2,4,
 2,4,8,6,4,6,2,4,
 6,2,6,6,4,2,4,6,
 2,6,4,2,4,12,2,12
]

while 1:
 for i in increment:
  search+=1
  b+=id=a/b
  if frac(d)==0:
   e+=1
   z[e]=b
   while 1:
    a=int(d)
    z[e+21]+=1
    d=a/b
    if frac(d):
     break
   c=int(sqrt(a))
  if b>cbreak
 if b>c
break


大きな要素数のリスト
ここで気になるのは、480個の整数要素のリストが正しく扱えるのかどうか?という点なので、チョット確認してみます。

というのも、以下のスクリプトを実行すると MemoryError が発生します。
lst = list(range(480))
print(lst)

ところで、要素数を 480 から 123 に変更して、
lst = list(range(123))
print(lst)
を実行すると、MemoryError が発生しなくなります。

つまり、この MemoryError は、print() が使用するメモリが不足して発生していることが分かり、これは Casio Python 特有の仕様だと思われます。その詳細はまだよく分かりません。


さて、以下のスクリプトを試します。

n=480; s=0
lst = list(range(n))
for i in range(n):
 s+=lst[i]
print(len(lst))
print(s)

これは、
・リスト lst = [0, 1, 2, ... , 477, 478, 479] を作り、
・各要素の和を計算して 変数 s に格納
・リスト lst の要素数と 変数 s の値を表示
するスクリプトです。

実行結果は、
480
114960

となり、正しい値さと分かります。正解は (0+480)×480/2 で簡単に検算できます。
従って、要素数 480 個の整数リスト自体は問題無く使えそうです。


12.3 探索数を拡大した15桁対応高速素因数分解のスクリプト

2, 3, 5, 7, 11 のいずれの倍数でもない数を探索数とし、探索数を 480 個に拡大したスクリプトは以下のようになります。

※ 高速素因数分解:15桁入力対応、探索数拡張版 - FactorG5.py のダウンロード

fx-CG50 Pythonモード:高速素因数分解 - FactorG5.py
"""Sample script

 Exercise;
 ported from Casio Basic
 "Prime Factor 15 digits"
  factorG.py
   ver 5.0

 by Krtyski/e-Gadget
"""

from u import *
digit=
15
search=0

def disp():
 global a,e,f,z
 clear_screen()
 locate(00f3'm'0)
 locate(160': ' str(len(inp)) + ' digits'3'm'0)
 locate(16, 11, 'search:'+str(search), 2, 'm', 0)
 line(0,15,383,15,1,0)
 for i in range(115):
  if i<=e:
   dx=int(i/12)*16
   dy=int(i/12)*11
   locate(0+dxi-dyz[i], 1'm'0)
   locate(10+dxi-dy'^('2'm'0)
   locate(12+dxi-dyz[i+21], 1'm'0)
   locate(15+dxi-dy')'2'm'0)
 locate
(00''0'm'1)


while
1:
 try:
  inp=str(eval(input('Number:')))
 except (SyntaxError,TypeError,NameError) as e:
  print(e)
  print('*must be number or\n expression')
  continue
 if '.' in inp:
  print('*must be integer')
  continue
 elif inp.isdigit():
  if len(inp)>digit:
   print('*must be '+str(digit)+' digit\n or less')
   continue
  else:
   f=int(inp)
   break
 else:
  continue

z=list(range(23))
for e in range(1,23):
 z[e]=0
e=0
a=f
c=int(sqrt(a))

prime_list=[2,3,5,7,11]
for b in prime_list:
 search+=1
 d=a/b
 if frac(d)==0:
  e+=1
  z[e]=b
  while 1:
   a=int(d)
   z[e+21]+=1
   d=a/b
   if frac(d):
    break
  c=int(sqrt(a))
 if b>c:break


increment=\
[4,2,4,6,2,6,4,2,
 4,6,6,2,6,4,2,6,
 4,6,8,4,2,4,2,4,
 14,4,6,2,10,2,6,6,
 4,2,4,6,2,10,2,4,
 2,12,10,2,4,2,4,6,
 2,6,4,6,6,6,2,6,
 4,2,6,4,6,8,4,2,
 4,6,8,6,10,2,4,6,
 2,6,6,4,2,4,6,2,
 6,4,2,6,10,2,10,2,
 4,2,4,6,8,4,2,4,
 12,2,6,4,2,6,4,6,
 12,2,4,2,4,8,6,4,
 6,2,4,6,2,6,10,2,
 4,6,2,6,4,2,4,2,
 10,2,10,2,4,6,6,2,
 6,6,4,6,6,2,6,4,
 2,6,4,6,8,4,2,6,
 4,8,6,4,6,2,4,6,
 8,6,4,2,10,2,6,4,
 2,4,2,10,2,10,2,4,
 2,4,8,6,4,2,4,6,
 6,2,6,4,8,4,6,8,
 4,2,4,2,4,8,6,4,
 6,6,6,2,6,6,4,2,
 4,6,2,6,4,2,4,2,
 10,2,10,2,6,4,6,2,
 6,4,2,4,6,6,8,4,
 2,6,10,8,4,2,4,2,
 4,8,10,6,2,4,8,6,
 6,4,2,4,6,2,6,4,
 6,2,10,2,10,2,4,2,
 4,6,2,6,4,2,4,6,
 6,2,6,6,6,4,6,8,
 4,2,4,2,4,8,6,4,
 8,4,6,2,6,6,4,2,
 4,6,8,4,2,4,2,10,
 2,10,2,4,2,4,6,2,
 10,2,4,6,8,6,4,2,
 6,4,6,8,4,6,2,4,
 8,6,4,6,2,4,6,2,
 6,6,4,6,6,2,6,6,
 4,2,10,2,10,2,4,2,
 4,6,2,6,4,2,10,6,
 2,6,4,2,6,4,6,8,
 4,2,4,2,12,6,4,6,
 2,4,6,2,12,4,2,4,
 8,6,4,2,4,2,10,2,
 10,6,2,4,6,2,6,4,
 2,4,6,6,2,6,4,2,
 10,6,8,6,4,2,4,8,
 6,4,6,2,4,6,2,6,
 6,6,4,6,2,6,4,2,
 4,2,10,12,2,4,2,10,
 2,6,4,2,4,6,6,2,
 10,2,6,4,14,4,2,4,
 2,4,8,6,4,6,2,4,
 6,2,6,6,4,2,4,6,
 2,6,4,2,4,12,2,12
]

while 1:
 for i in increment:
  search+=1
  b+=id=a/b
  if frac(d)==0:
   e+=1
   z[e]=b
   while 1:
    a=int(d)
    z[e+21]+=1
    d=a/b
    if frac(d):
     break
   c=int(sqrt(a))
  if b>cbreak
 if b>c
break

if a>1:
 e+=1
 z[e]=a
 a=1
 z[e+11]=
1

disp
()



12.4 探索数拡大による高速化
前回作成した Ver 4.0 と今回作成した Ver 5.0 の実行速度を比較しました。

比較対象は、前回と同じく以下の10個のケースを用います。

search_30445_f Search_39969_f seach_47824_f 
Case1_14digit search_90872_f search_189426_f.png 
search_263027_f.png search_343013_f.png search_2025804_f.png 
search_4349444_f.png 


表2- 探索数拡大による実行時間の違い - Ver 4.0 と Ver 5.0 の比較
Ver 4.0Ver 5.0
入力値serach実行時間
[秒]
search実行時間
[秒]
高速化
[秒]
search
減少率[%]
547,896,321,054,789 30,445 9.627,679 8.70.99.1
7,845,162,037,849 39,969 12.536,337 11.41.19.1
748,159,026,374,815 47,824 15.043,478 13.71.39.1
36,209,518,473,620 55,241 17.150,221 15.51.69.1
986,753,421,098,675 90,872 28.682,612 26.02.69.1
96,835,724,130,261 189,426 58.1172,207 52.75.49.1
748,159,026,396,835 263,027 80.5239,116 73.17.49.1
96,835,724,136,209 343,013 104.5311,832 95.39.29.1
748,159,026,336,209 2,025,804 644.0
(10分44秒)
1,841,642 585.0
(9分45秒)
59.09.1
362,095,184,736,209 4,349,444 1,361.0
(22分41秒)
3,954,042 1,235.6
(20分35.6秒)
125.49.1


search回数は、どの入力値でも一律に 9.1% 減少して、高速化していることが分かります。実行速度は 9~10%程度の高速化となりました。


12.5 増分リストをさらに拡張する 

ところで、さらに探索数を拡大して、2, 3, 5, 7, 11, 13 のいずれかの倍数でない探索数を得るための増分リストに拡張してみました。

要素数5760個の探索数増分リストを用いる方法
A-2) 素数 2, 3, 5, 7, 11, 13, 17 を探索数として除算を繰り返し、
B-2) 1930048 の整数で 2, 3, 5, 7, 11, 13 のいずえかの倍数でない探索数で除算を繰り返し、
エラストテレスの篩いで素因数を調べる。

この場合は、取りこぼしを具せずためには、2, 3, 5, 7, 11, 13 のいずれかの倍数でない整数は 5750 個あります。この個数の求め方は、上と同様に計算して得られます。

Casio Python のスクリプトファイルは、以前の記事で 300行が上限だと確認しています。
リスト increment[ ] の定義以外に103行使っているので、リスト定義を 196行以内に納める必要があります。
そこで、1行あたり30個の要素を並べることにすれば、要素の記述に192行必要となり、スクリプト全体で 295 行とぎりぎり収まります。

2, 3, 5, 7, 11, 13 のいずれかの倍数でない探索数の増分リストを出力するたツール SearchNum13.exe も作りました。
要素数 5760 の探索数増分リストの要素を出力するWindowsデスクトップアプリ
 - SearchNum13 のダウンロード


これを起動したところ;
SearchNum13_1

右上の [Calculate] ボタンをクリックすると、求める要素とコンマ , が出力されます。
SearchNum13_result

1行30要素で192行となるように出力されます。コピー&ペーストでエディタで追加すると楽です、というか唯一の現実的な方法ですね!

さて、5760個の整数リストはかなりメモリを使うので、Casio Python で使えるかどうかを調べます。

要素数 5760 個のリスト
n=5760; s=0
lst = list(range(n))
for i in range(n):
 s+=lst[i]
print(len(lst))
print(s)

を実行すると、
5760
16585920
と出力され、0 + 1 + 2 + ・・・ + 5758 + 5769 = 16585920 と、正しい計算結果になります。正解は、(0+5760)×5760/2 で簡単に検算できます。これにより、整数の要素数 5760 個に拡張しても、整数リストは問題無く動作することは確認できました。

5760個の整数リストが使えることが分かったので、実際にコピー&ペーストしてスクリプトを作りました。

高速素因数分解:15桁入力対応、さらに探索数拡張版 - FactorG6.py のダウンロード

これを実行すると、以下のようなエラーが表示されます。
Stack_Error 

エラーメッセージ:RuntimeError: pystack exhausted

Casio Python が動作する際に確保するスタックを使い果たした、という意味だと考えられます。

制御構造は、FactorG5.py と同じなので、構造制御のために使うスタック数は同じです。すると、スタックを使い果たすのは、長大なリストが原因になっている筈です。

スタックがどのように使われているのかは不明ですが、リスト increment[ ] の定義で使っている改行がスタックを増やしているかも知れないと考え、1行の文字数を長くして行数を減らし、何通りか試してみました。

すると、実行すると Invalid Data size エラー発生して実行出来ないケースと RuntimeError: pystack exhausted エラーが発生して実行できないケースの両方があることが分かりました。

Invalid Data size エラーが発生する時のスクリプトファイルサイズは、pystack exhausted エラーの時のファイルサイズよりも大きくなっています。Invalid Data size エラーが発生する条件が、必ずしもファイルサイズが大きいだけでもないことも分かりました。

今回は、アルゴリズムの要請から極端に大きな要素数のリストを使う事例となり、その結果エラーに阻まれて、探索数をさらに拡張できませんでした。Casio Python では、スクリプトファイルサイズの制限やスタックサイズの制限が厳しいことも判りました。


何かお気づきの方は、コメント欄に情報をお寄せください。



目 次

前の記事 - 11. 関数のオーバーヘッド

次の記事 - 13. 10進数除算の出力と精度





応援クリックをお願いします。励みになるので...
にほんブログ村 IT技術ブログ 開発言語へ


 


keywords: fx-CG50Pythonfx-9750GIIIfx-9860GIIIプログラム関数電卓

リンク集 | ブログ内マップ


テーマ : プログラム関数電卓
ジャンル : コンピュータ

Casio Python - 大きな数の計算:高速素因数分解(3)

Python Casio Python
 Casioグラフ関数電卓の Python を使ってみる
     - 大きな数の計算:高速素因数分解(3) 
目次


初版:2020/07/28
追記修正:202011/14

前の記事 - 9. Python らしい反復処理 |  次の記事 - 11. 関数呼出オーバーヘッド


10. 大きな数の計算:高速素因数分解(3)

今回は、Casio Python では、より大きな桁数の整数で精度を確保できることを紹介し、15桁入力に対応した高速素因数分解スクリプトに改造します。

前回までに作った高速素因数分解スクリプトは、オリジナルの Casio Basic プログラムから移植したものなので、入力値の桁数はオリジナルと同じように10桁が上限になっています。

Casio Basic は仕様上10進数10桁を超えると計算精度の保証がないので、入力値は最大10桁に制限されます。一方、Casio Python では10進数10桁を超えても精度が確保されます。


例えば、15桁の入力を2で除算してみます。
999,999,999,999,999,999 / 2 = 49,999,999,999,999.8
15digit 

入力15桁に対して、1桁多い16桁で出力されます。
正解は 499999999999999.5 なので、16桁の最下位に誤差が発生しています。


次に、16桁の入力を2で除算してみます。
9,999,999,999,999,999,999 / 2 = 50,000,000,000,000.0
16digit 

16桁の入力に対して、出力は1桁多い17桁になっています。
奇数を2で除算して整数になることから、少なくとも16桁の入力では、計算精度に問題があることが分かります。

続いて、1/3 の計算結果は、出力16 桁になります。
digit_one_third_shell 

円周率 pi も 16桁出力になります。
digit_pi_shell 

Casio Python での浮動小数点計算は、10進数15桁以内の入力なら、計算結果は15桁の精度が確保できるようです。

そこで、高速素因数分解の入力値の桁数を 10 桁から 15 桁に拡張することにします。


10.1 変更箇所の洗い出し

前回作った高速素因数分解 - 10桁対応版:FactorG2.py のダウンロード

15桁対応にするために、前回作った下記のスクリプトに対して、下記の2項目を変更します。

(1) 増加する素因数に合わせて表示行数を拡張する - disp() の変更
(2) 入力桁数の制限を拡張し、それに合わせてエラー表示を変更する


これら2つの項目に対応して変更箇所を以下に赤文字で示します。

fx-CG50 Pythonモード:高速素因数分解 - FactorG2.py
"""Sample script

 Exercise;
 ported from Casio Basic
 "Prime Factor 10 digits"
  FactorG.py
   ver 2.0

 by Krtyski/e-Gadget
"""

from u import *

def ckpwr():
 global a,b,c,d,e,z
 e+=1
 z[e]=b
 while 1:
  a=int(d)
  z[e+11]+=1
  d=a/b
  if frac(d):
   break
 c=int(sqrt(a))

def
disp():
 global a,e,f,z
 clear_screen()
 locate(00f3'm'0)
 locate(160': ' str(len(inp)) + ' digits'3'm'0)
 line(0,15,383,15,1,0)
 for i in range(112):  ◀ (1) 増加する素因数に対応して、表示行数を拡張する
  if i<=e:
   locate(,,z[i], 1'm'0)
   locate(10i'^('2'm'0)
   locate(12iz[i+21], 1'm'0)
   locate(15i')'2'm'0)
 locate
(00''0'm'1)


while
1:
 try:
  inp=str(eval(input('Number:')))
 except (SyntaxError,TypeError,NameError) as e:
  print(e)
  print('*must be number or\n expression')
  continue
 if '.' in inp:
  print('*must be integer')
  continue
 elif inp.isdigit():
  if len(inp)>10:    ◀ (2) 桁数制限を拡張する
   print('*must be 10 or less') ◀ (2) 拡張した桁数制限に合わせて変更
   continue
  else:
   f=int(inp)
   break
 else:
  continue

z=list(range(23))
for e in range(1,23):
 z[e]=0
e=0
a=f
c=int(sqrt(a))

prime_list=[2,3,5,7,11]
for b in prime_list:
 d=a/b
 if frac(d)==0ck_pwr()
 if b>c: break


increment=\
[2,4,2,4,6,2,6,4,
 2,4,6,6,2,6,4,2,
 6,4,6,8,4,2,4,2,
 4,8,6,4,6,2,4,6,
 2,6,6,4,2,4,6,2,
 6,4,2,4,2,10,2,10]

while 1:
 for i in increment:
  b+=id=a/b
  if frac(d)==0ck_pwr()
  if b>cbreak
 if b>c
break

if a>1:
 e+=1
 z[e]=a
 a=1
 z[e+11]=
1

disp
()



10.2 増加する素因数に合わせて表示行数を拡張する - disp() の変更

素数の積を調べてみます。

2から41までの素数の積;
2*3*5*7*11*13*17*19*23*29*31*37*41 = 304,250,263,527,210
15桁になります。

2から43までの素数の積;
2*3*5*7*11*17*19*23*29*31*37*41*43 = 21,306,211,311,576,906
17桁になります。

これらから、素数13個で15桁になることが確認でき、乗数を含めて出力には最大 13行必要になります。

前回作った disp() を示します。

def disp():
 global a,e,f,z
 clear_screen()
 locate(00f3'm'0)
 locate(160': ' str(len(inp)) + ' digits'3'm'0)
 line(0,15,383,15,1,0)
 for i in range(112):  ◀ (1) 増加する素因数に対応して、表示行数を拡張する
  if i<=e:
   locate(,,z[i], 1'm'0)
   locate(10i'^('2'm'0)
   locate(12iz[i+21], 1'm'0)
   locate(15i')'2'm'0)
 locate
(00''0'm'1)


Ver2_10lines1 

1画面あたり、素因数は11個まで表示できるので、for i in range(12): としました。
これを13行まで拡張するためには、以下のように2列にして13行まで表示できるように変更します。

13lines 

従って、反復処理の for 文は、for i in range(14): で良いですが、1行余裕をみて

for i in range(15):

としておきます。しかし、これだけでは、12行目 と 13行目 は、右側に表示されないので、4行ある locate() の引数を適切に変更してゆきます。
 
i の値が 12 未満では左側に表示、12 以上の場合は右側に表示すれば良いですね。
そこで、if i <12: / else: と場合分けをする方法もあるでしょうが、もう少しスッキリと簡潔に記述するために、locate() で指定する 桁数と行数の増分 dxdyi の値に応じて計算で求める方式を考えました。
 dx=int(i/12)*16  i<12 のとき dx=0i≧12 のとき dx=16
 dy=int(i/12)*11
  i<12 のとき dy=0i≧12 のとき dy=11

これを使って、スクリプトを以下のように変更しました。

def disp():
 global a,e,f,z
 clear_screen()
 locate(00f3'm'0)
 locate(160': ' str(len(inp)) + ' digits'3'm'0)
 line(0153831510)
 for i in range(115):  ◀ 変更部分
  if i<=e:
   dx=int(i/12)*16  ◀ 追加部分: i < 12 のとき dx=0、i ≧ 12 のとき dx=16
   dy=int(i/12)*11    ◀ 追加部分: i < 12 のとき dy=0、i ≧ 12 のとき dy=11
   locate(0+dxi-dyz[i], 1'm'0)    ◀ 第1引数に +dx を、第2引数に -dy を追加 
   locate(10+dxi-dy'^('2'm'0)    ◀ 第1引数に +dx を、第2引数に -dy を追加
   locate(12+dxi-dyz[i+21], 1'm'0) ◀ 第1引数に +dx を、第2引数に -dy を追加
   locate(15+dxi-dy')'2'm'0)    ◀ 第1引数に +dx を、第2引数に -dy を追加 
 locate
(00''0'm'1)


他にもスッキリと記述する方法があると思いますが、これで進めます。


10.3 入力桁数の制限を拡張し、それに合わせてエラー表示を変更する

前回までに作った "入力と入力チェックの記述" を以下に示し、今回変更する部分を赤文字で示します。

while 1:
 try:
  inp str(eval(input('Number:')))
 except (SyntaxErrorTypeErrorNameError) as e:
  print(e)
  print('*must be number or\n expression')
  continue
 if '.' in inp:
  print('*must be integer')
  continue
 elif inp.isdigit():
  if len(inp)>10:    ◀ (2) 桁数制限を拡張する
   print('*must be 10 or less') ◀ (2) 拡張した桁数制限に合わせて変更
   continue
  else:
   f int(inp)
   break
 else:
  continue

ここで、赤文字で追記している2箇所は、10 を 15 に変更すればOKです。

今回は、変数 digit を使って、これに 15 を代入して使うことにします。

digit=15

をスクリプトの冒頭に書いておけば、代入する数値を書き換えるだけで入力桁数の制限を簡単に変更できます。

そこで、上記2箇所は、以下のように変更します。

if ken(inp) > digit:
 print('*must be '+str(digit)+' digit\n or less')

 continue

出力文字列に エスケープシーケンス \n (バックスラッシn) による改行を入れているので、狭い画面の電卓で以下のように見やすくなります。
*must be 15
 or less



入力と入力チェック部分は、まとめると以下のようになります。

digit=15

・・・
・・・

while
1:
 try:
  inp str(eval(input('Number:')))
 except (SyntaxErrorTypeErrorNameError) as e:
  print(e)
  print('*must be number or\n expression')
  continue
 if '.' in inp:
  print('*must be integer')
  continue
 elif inp.isdigit():
  if len(inp)>digit:
   print('*must be '+str(digit)+' digit\n or less')
   continue
  else:
   f int(inp)
   break
 else:
  continue


10.4 素因数探索回数の表示を追加する

ここまで作ったスクリプトを実際に以下の2通りの入力値で走らせてみると、実行時間が大きく異なります。

NoSearch_ver3_fast_14digits NoSearch_ver3_slow_13digits 

左の14桁に比べて、右の13桁は、桁数が1つ少ないにもかかわらず、かなり時間がかかります。
素因数の桁数が大きいほど時間がかかると考えれば、素因数探索の while ループを回る回数が大きいほど時間がかかると理解できます。

そこで、実際に素因数探索の回数をカウントして、以下のように表示しようと思います。

Search_ver3_fast_14digits Search_ver3_slow_13digits 

search は探索回数です。

変数 search の初期化
スクリプト冒頭の digit=15 の直後に search=0 を追加し、変数 search を探索回数とします。

from u import *
digit=
15
search=0


 変数 search のカウントアップ
変数 search をカウントアップ (serch+=1) するのは、以下の赤文字で示した 2 箇所になります。

prime_list=[2,3,5,7,11]
for b in prime_list:
 search+=1    ◀ ここで search をカウントアップ
 d=a/b
 if frac(d)==0ck_pwr()
 if b>c: break


increment=\
[2,4,2,4,6,2,6,4,
 2,4,6,6,2,6,4,2,
 6,4,6,8,4,2,4,2,
 4,8,6,4,6,2,4,6,
 2,6,6,4,2,4,6,2,
 6,4,2,4,2,10,2,10]

while 1:
 for i in increment:
  search+=1   ◀ ここで search をカウントアップ
  b+=id=a/b
  if frac(d)==0ck_pwr()
  if b>cbreak
 if b>c
break



変数 search の表示
カウントアップした search の表示は、disp() 関数の内で、以下の赤文字で示したようにします。

def disp():
 global a,e,f,z
 clear_screen()
 locate(00f3'm'0)
 locate(160': ' str(len(inp)) + ' digits'3'm'0)
 locate(16, 11, 'search:'+str(search), 2, 'm', 0)   ◀ search の表示
 line(0,15,383,15,1,0)
 for i in range(115):
  if i<=e:
   dx=int(i/12)*16
   dy=int(i/12)*11
   locate(0+dxi-dyz[i], 1'm'0)
   locate(10+dxi-dy'^('2'm'0)
   locate(12+dxi-dyz[i+21], 1'm'0)
   locate(15+dxi-dy')'2'm'0)
 locate
(00''0'm'1)



10.5 素因数が10桁を超えるときの表示の変更 - disp() の変更

入力桁数を15桁まで拡張したので、見つかった素因数が10桁を超えることがあります。
例えば、以下のようなケースです。

11digit_prime 

4つめの素因数 58402449151 (11桁) が、乗数表示 ^(1 ) の表示と被ってしまっています。
これでは、ちょっと格好悪いので、素因数が10桁を超えた時は、超えた桁数に応じて乗数表示を右に移動して、被らないようにします。

disp() 関数内では、素因数は z[i] なので、素因数の桁数は、len(str(z[i])) で分かります。
この桁数が 10 を超えるとき、その超過分は len(str(z[i])) - 10 で求められます。これを変数 r に代入しておきます。

この r を使えば、乗数表示にかかわる locate() の第1引数に r を加えれば、被らなくなりなります。

この変更を行った disp() 関数は以下のようになります。変更部分は赤文字で示します。

def
disp():
 global a,e,f,z
 clear_screen()
 locate(00f3'm'0)
 locate(160': ' str(len(inp)) + ' digits'3'm'0)
 line(0,15,383,15,1,0)
 for i in range(115):
  if i<=e:
   r=0
   if len(str(z[i]))>10:
    r=len(str(z[i]))-10
   dx=int(i/12)*16
   dy=int(i/12)*11
   locate(0+dxi-dyz[i], 1'm'0)
   locate(10+dx+ri-dy'^('2'm'0)    ◀ 第1引数に r を加算する
   locate(12+dx+ri-dyz[i+21], 1'm'0)   ◀ 第1引数に r を加算する
   locate(15+dx+ri-dy')'2'm'0)     ◀ 第1引数に r を加算する
 locate
(00''0'm'1)



変更したスクリプトで、上記の "被ってしまった" ケースを入力すると、以下のように改善されたことが分かります。
Case1_14digit 


別のケースとして、素因数が15桁になるときは、以下のように乗数表示は17桁めに移動して、被りません。
15digit_prime 

なお、このケースでは、探索回数 search が、なんと434万回を超えています。

==========

以上で、15桁入力への拡張版が完成です。

※ 高速素因数分解 15桁入力対応版 - FactorG3.py のダウンロード

fx-CG50 Pythonモード:高速素因数分解 - FactorG3.py
"""Sample script

Exercise;
ported from Casio Basic
"Prime Factor 10 digits"
factorG.py
ver 3.0

by Krtyski/e-Gadget
"""

from u import *
digit=
15
search=0

def ckpwr():
 global a,b,c,d,e,z
 e+=1
 z[e]=b
 while 1:
  a=int(d)
  z[e+11]+=1
  d=a/b
  if frac(d):
   break
 c=int(sqrt(a))

def
disp():
 global a,e,f,z
 clear_screen()
 locate(00f3'm'0)
 locate(160': ' str(len(inp)) + ' digits'3'm'0)
 locate(16, 11, 'search:'+str(search), 2, 'm', 0)
 line(0,15,383,15,1,0)
 for i in range(115):
  if i<=e:
   dx=int(i/12)*16
   dy=int(i/12)*11
   locate(0+dxi-dyz[i], 1'm'0)
   locate(10+dxi-dy'^('2'm'0)
   locate(12+dxi-dyz[i+21], 1'm'0)
   locate(15+dxi-dy')'2'm'0)
 locate
(00''0'm'1)


while
1:
 try:
  inp=str(eval(input('Number:')))
 except (SyntaxError,TypeError,NameError) as e:
  print(e)
  print('*must be number or\n expression')
  continue
 if '.' in inp:
  print('*must be integer')
  continue
 elif inp.isdigit():
  if len(inp)>digit:
   print('*must be '+str(digit)+' digit\n or less')
   continue
  else:
   f=int(inp)
   break
 else:
  continue

z=list(range(23))
for e in range(1,23):
 z[e]=0
e=0
a=f
c=int(sqrt(a))

prime_list=[2,3,5,7,11]
for b in prime_list:
 search+=1
 d=a/b
 if frac(d)==0ck_pwr()
 if b>c: break


increment=\
[2,4,2,4,6,2,6,4,\
 2,4,6,6,2,6,4,2,\
 6,4,6,8,4,2,4,2,\
 4,8,6,4,6,2,4,6,\
 2,6,6,4,2,4,6,2,\
 6,4,2,4,2,10,2,10]

while 1:
 for i in increment:
  search+=1
  b+=id=a/b
  if frac(d)==0ck_pwr()
  if b>cbreak
 if b>c
break

if a>1:
 e+=1
 z[e]=a
 a=1
 z[e+11]=
1

disp
()






目 次

前の記事 - 9. Python らしい反復処理

次の記事 - 11. 関数呼出オーバーヘッド





応援クリックをお願いします。励みになるので...
にほんブログ村 IT技術ブログ 開発言語へ


 


keywords: fx-CG50Pythonプログラミング入門プログラミングプログラム関数電卓

リンク集 | ブログ内マップ

テーマ : プログラム関数電卓
ジャンル : コンピュータ

Casio Python - Python らしい反復処理:高速素因数分解(2)

Python Casio Python
 Casioグラフ関数電卓の Python を使ってみる
     - Python らしい反復処理:高速素因数分解(2) 
目次


初版:2020/07/25
修正:2020/07/26
追記修正:2020/11/14

前の記事 - 8. シェル画面入力の工夫 |  次の記事 - 10. 大きな数の計算


9. Python らしい反復処理:高速素因数分解(2) <fx-CG50 OS3.40以降>

前回は、Casio Basic プログラムを、とりあえず動かすことを優先させ、できるだけ忠実に Casio Python に移植しました。
291行もある長大なスクリプトです。似たような処理が何度も繰り返されており、無駄が満載です。

オリジナルの Casio Basic のコードは、その言語仕様から長大なソースになっていますが、今回は、Python らしい反復処理を用いて、大幅にスリム化します。主に、以下に示した5項目 の変更を行って、無駄を減らします。

(1) farc() 関数をユーザーモジュールに追加し、ここでの関数定義を削除する
(2) 不要な記述を削除する
(3) disp() 関数内を、反復処理に置換えて、記述を簡潔にする
(4) を反復処理に置き換えて、記述を簡潔にする
(5) を反復処理に書き換えて、記述を簡潔にする



前回作ったスクリプト - FactorG1.py のダウンロード

fx-CG50 Pythonモード:高速素因数分解 - FactorG1.py
"""Sample script

Exercise;
ported from Casio Basic
"Prime Factor 10 digits"
factorG.py
ver 1.0

by Krtyski/e-Gadget
"""

from u import *

def ckpwr():
 global a,b,c,d,e,z
 e+=1
 z[e]=b
 while 1:
  a=int(d)
  z[e+11]+=1
  d=a/b
  if frac(d):
   break
 c=int(sqrt(a))

def disp():
 global a,b,c,d,e,f,z
 if a>1:     ◀ (3) disp() は素因数探索が終わって最後に実行されるので、
  e+=1      ◀ (3) この if 文も 探索囚虜後に1度だけ実行される。
  z[e]=a    ◀ (3) この if 文は、最後に探索した素因数を z[ ] に格納する処理で、
  a=1       ◀ (3) 表示処理の一部ではないので、ここから削除し、
  z[e+11]=1   ◀ (3) disp() 関数実行の直前に記述することにする。
 d=int(e/6)   ◀ (2) 複数ページ切替のためのページ番号 d を計算⇒ 不要、削除する
 if e-6*d>0:   
 (2) 複数ページ切替のためのページ番号 d を計算⇒ 不要、削除する
  d+=1     ◀ (2) 複数ページ切替のためのページ番号 d を計算⇒ 不要、削除する
 c=1

 clear_screen()
 locate(00f3'm'0)
 b 6*(c-1)+1   (3) 複数ページ切替機能を使わないので、bは常に1 ⇒ b=1 と書き換える
 locate(0, 1, z[b], 1, 'm', 0) ◀ (3) ここから disp() 関数終わりまでは反復処理に置き換える
 locate(10, 1, '^(', 2, 'm', 0)                    ↓
 locate(12, 1, z[b+11], 1, 'm', 0)
 locate(15, 1, ')', 2, 'm', 0)
 if b+1<=e:
  locate(0, 2, z[b+1], 1, 'm', 0)
  locate(10, 2, '^(', 2, 'm', 0)
  locate(12, 2, z[b+12], 1, 'm', 0)
  locate(15, 2, ')', 2, 'm', 0)
 if b+2<=e:
  locate(0, 3, z[b+2], 1, 'm', 0)
  locate(10, 3, '^(', 2, 'm', 0)
  locate(12, 3, z[b+13], 1, 'm', 0)
  locate(15, 3, ')', 2, 'm', 0)
 if b+3<=e:
  locate(0, 4, z[b+3], 1, 'm', 0)
  locate(10, 4, '^(', 2, 'm', 0)
  locate(12, 4, z[b+14], 1, 'm', 0)
  locate(15, 4, ')', 2, 'm', 0)
 if b+4<=e:
  locate(0, 5, z[b+4], 1, 'm', 0)
  locate(10, 5, '^(', 2, 'm', 0)
  locate(12, 5, z[b+15], 1, 'm', 0)
  locate(15, 5, ')', 2, 'm', 0)
 if b+5<=e:
  locate(0, 6, z[b+5], 1, 'm', 0)
  locate(10, 6, '^(', 2, 'm', 0)
  locate(12, 6, z[b+16], 1, 'm', 0)
  locate(15, 6, ')', 2, 'm', 0)
 if b+6<=e:
  locate(0, 7, z[b+6], 1, 'm', 0)
  locate(10, 7, '^(', 2, 'm', 0)
  locate(12, 7, z[b+17], 1, 'm', 0)
  locate(15, 7, ')', 2, 'm', 0)
 if b+7<=e:
  locate(0, 8, z[b+7], 1, 'm', 0)
  locate(10, 8, '^(', 2, 'm', 0)
  locate(12, 8, z[b+18], 1, 'm', 0)
  locate(15, 8, ')', 2, 'm', 0)
 if b+8<=e:
  locate(0, 9, z[b+8], 1, 'm', 0)
  locate(10, 9, '^(', 2, 'm', 0)
  locate(12, 9, z[b+19], 1, 'm', 0)
  locate(15, 9, ')', 2, 'm', 0)
 if b+9<=e:
  locate(0, 10, z[b+9], 1, 'm', 0)
  locate(10, 10, '^(', 2, 'm', 0)
  locate(12, 10, z[b+20], 1, 'm', 0)
  locate(15, 10, ')', 2, 'm', 0)
 if b+10<=e:
  locate(0, 11, z[b+10], 1, 'm', 0)
  locate(10, 11, '^(', 2, 'm', 0)
  locate(12, 11, z[b+21], 1, 'm', 0)    
  locate(15, 11, ')', 2, 'm', 0)     ◀ (3) ここまでを反復処理で置き換える
 locate(0, 0, '', 0, 'm', 1)

def frac(x):          ◀ (1) この関数をユーザーモジュール u.py ver 1.4 に追加し
 retuen x-int(x)         これを削除する


while
1:
 try:
  inp=str(eval(input('Number:')))
 except (SyntaxError,TypeError,NameError) as e:
  print(e)
  print('*must be number or\n expression')
  continue
 if '.' in inp:
  print('*must be integer')
  continue
 elif inp.isdigit():
  if len(inp)>10:
   print('*must be 10 or less')
   continue
  else:
   f=int(inp)
   break
 else:
  continue

z=list(range(23))
for e in range(1,23):
 z[e]=0
z[1]=0    ◀ (2) 上の for文で 0 に初期化されているので不要、削除する
z[12]=0    
◀ (2) 上の for文で 0 に初期化されているので不要、削除する
e=0
a=f
c=int(sqrt(a))

b=2;d=a/b         ◀ (4) ここから whle 1: の上までを反復処理に置き換える
if frac(d)==0:ckpwr()        
if b>c:disp()          disp() 関数は、素因数探索が終わってから実行される
b=3;d=a/b           ので、ここで disp() を実行する必要が無い。
if frac(d)==0:ckpwr()     そこで、この部分の disp() は break に置換えても
if b>c:disp()          問題ない。disp() は最後に1度だけ実行すれば良い
b=5;d=a/b
if frac(d)==0:ckpwr()
if b>c:disp()
b=7;d=a/b
if frac(d)==0:ckpwr()
if b>c:disp()
b=11;d=a/b
if frac(d)==0:ckpwr()        
if b>c:disp()        ◀ (4) ここまでを反復処理で置き換える

while 1:
 b+=2;d=a/b       ◀ (5) ここから disp() の前までを反復処理に置き換える
 if frac(d)==0:ckpwr()       
 if b>c:break
 b+=4;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=2;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=4;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=6;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=2;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=6;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=4;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=2;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=4;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=6;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=6;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=2;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=6;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=4;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=2;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=6;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=4;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=6;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=8;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=4;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=2;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=4;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=2;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=4;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=8;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=6;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=4;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=6;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=2;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=4;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=6;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=2;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=6;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=6;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=4;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=2;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=4;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=6;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=2;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=6;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=4;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=2;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=4;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=2;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=10;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=2;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=10;d=a/b
 if frac(d)==0:ckpwr()       ↓
 if b>c:break        (5) ここまでを反復処理で書き換える

disp()



9.1 frac() 関数をユーザーモジュール (u.py) に追加

この関数は、浮動小数点の数値から小数部を取り出すもので、Casio Python には備わっていないません。便利な関数なのでユーザー関数としてユーザーモジュールに追加します。バージョンアップしたユーザーモジュールは下記からダウンロードできます。

ユーザーモジュール - u.py ver 1.4 のダウンロード


9.2 不要な処理を削除

上記の5箇所 (赤文字で表記した部分) を削除します。


9.3 disp() 関数のスリム化

disp() 関数の冒頭にある if 文
 if a>1:
  e+=1
  z[e]=a
  z[e+11]=1

を削除します。

というのは、この部分は最後に探索した素因数を リスト z[ ] に格納する処理なので、disp() から切り離し、disp() の直前に記述するように変更します。disp() は全ての素因数を見つけた後に、1度だけ実行されるので、この切り離しは機能に影響はなく、ソースの見やすさが改善されると思います。

さらに、
b = 6*(c-1)+1

b=1
に変更します。

というのも、変数 c は、結果表示のページを示し、c=1 の時は1ページ目、c=2 の時は2ページ目を示します。
現在の Casio Python ではインタラクティブなキー入力を取得する機能が実装されていないので、全ての素因数を1ページ内に表示する仕様に変更しています。従ってページ数のカウントは不要になり、常に c=1 となります。すると
b = 6*(c-1)+1 は、常に b=1 となるので、変数 b の初期化をこのように書き換えます。

これまでをまとめると、 disp() は下記になります。

def disp():
 global a,b,c,d,e,f,z
 clear_screen()
 locate(00f3'm'0)
 b=1                  (3) 複数ページ切替機能を使わないので、常に b=1、b=1 に書き換えた
 locate(0, 1, z[b], 1, 'm', 0) ◀ (3) ここから disp() 関数終わりまでは反復処理に置き換える
 locate(10, 1, '^(', 2, 'm', 0)                    ↓
 locate(12, 1, z[b+11], 1, 'm', 0)
 locate(15, 1, ')', 2, 'm', 0)
 if b+1<=e:
  locate(0, 2, z[b+1], 1, 'm', 0)
  locate(10, 2, '^(', 2, 'm', 0)
  locate(12, 2, z[b+12], 1, 'm', 0)
  locate(15, 2, ')', 2, 'm', 0)
 if b+2<=e:
  locate(0, 3, z[b+2], 1, 'm', 0)
  locate(10, 3, '^(', 2, 'm', 0)
  locate(12, 3, z[b+13], 1, 'm', 0)
  locate(15, 3, ')', 2, 'm', 0)
 if b+3<=e:
  locate(0, 4, z[b+3], 1, 'm', 0)
  locate(10, 4, '^(', 2, 'm', 0)
  locate(12, 4, z[b+14], 1, 'm', 0)
  locate(15, 4, ')', 2, 'm', 0)
 if b+4<=e:
  locate(0, 5, z[b+4], 1, 'm', 0)
  locate(10, 5, '^(', 2, 'm', 0)
  locate(12, 5, z[b+15], 1, 'm', 0)
  locate(15, 5, ')', 2, 'm', 0)
 if b+5<=e:
  locate(0, 6, z[b+5], 1, 'm', 0)
  locate(10, 6, '^(', 2, 'm', 0)
  locate(12, 6, z[b+16], 1, 'm', 0)
  locate(15, 6, ')', 2, 'm', 0)
 if b+6<=e:
  locate(0, 7, z[b+6], 1, 'm', 0)
  locate(10, 7, '^(', 2, 'm', 0)
  locate(12, 7, z[b+17], 1, 'm', 0)
  locate(15, 7, ')', 2, 'm', 0)
 if b+7<=e:
  locate(0, 8, z[b+7], 1, 'm', 0)
  locate(10, 8, '^(', 2, 'm', 0)
  locate(12, 8, z[b+18], 1, 'm', 0)
  locate(15, 8, ')', 2, 'm', 0)
 if b+8<=e:
  locate(0, 9, z[b+8], 1, 'm', 0)
  locate(10, 9, '^(', 2, 'm', 0)
  locate(12, 9, z[b+19], 1, 'm', 0)
  locate(15, 9, ')', 2, 'm', 0)
 if b+9<=e:
  locate(0, 10, z[b+9], 1, 'm', 0)
  locate(10, 10, '^(', 2, 'm', 0)
  locate(12, 10, z[b+20], 1, 'm', 0)
  locate(15, 10, ')', 2, 'm', 0)
 if b+10<=e:
  locate(0, 11, z[b+10], 1, 'm', 0)
  locate(10, 11, '^(', 2, 'm', 0)
  locate(12, 11, z[b+21], 1, 'm', 0)    
  locate(15, 11, ')', 2, 'm', 0)     ◀ (3) ここまでを反復処理で置き換える
 locate(0, 0, '', 0, 'm', 1)



変数 b に着目すると、b は リスト z[ ] のインデックスとして利用されていて、このインデックスは、1 から 11 まで 1 づつ増えています。そこで、それぞれの if 節にある locate() 関数4行分をワンセットとして、for 文による反復処理に書き換えます。

ついでに、入力値の桁数を表示する機能を追加します。

clear_screen() 以下を書き換えます。

 clear_screen()
 locate(00f3'm'0)
 locate(160': ' str(len(inp)) + ' digits'3'm'0) #入力値の桁数表示
 line(0,15,383,15,1,0)
 for i in range(112):
  if i<=e:
   locate(,,z[i], 1'm'0)
   locate(10i'^('2'm'0)
   locate(12iz[i+21], 1'm'0)
   locate(15i')'2'm'0)
 locate
(00''0'm'1)


さて、このように書き換えた結果、変数 b, c, d  はこの関数内では不要になりました。従ってグローバル変数の宣言において、不要になった変数 b, c, d を削除して、gobal 文を書き換えます。

range(1, 23) は、リスト [1, 2, 3, 4, 5, ... ,21, 22, 23] を生成して、その要素を最初から最後まで順に i に適用して反復を行います。range(23) は 要素が 0 から始まりますが、range(1, 23) とすれば、要素は 1 から始まります。

Note: range() の書式:
range(start=0, stop, step=1)
・引数が1つの場合は、range(stop) で、デフォルト start=0、step=1 が適用されます。
・引数が2つの場合は、range(start, stop) で、デフォルト step=1 が適用されます。
・引数が3つの場合は、range(start, stop, step) と全てを明示的に指定します。


完成した disp() 関数は下記のようになります。スッキリと無駄を省けました。

def disp():
 global a,e,f,z
 clear_screen()
 locate(00f3'm'0)
 locate(160': ' str(len(inp)) + ' digits'3'm'0)
 line(0,15,383,15,1,0)
 for i in range(112):
  if i<=e:
   locate(,,z[i], 1'm'0)
   locate(10i'^('2'm'0)
   locate(12iz[i+21], 1'm'0)
   locate(15i')'2'm'0)
 locate
(00''0'm'1)



9.4  素数が素因数かどうかをチェックする部分のスリム化

以下の部分を反復処理に置き換えます。

b=2;d=a/b
if frac(d)==0:ckpwr()
if b>c:disp()
b=3;d=a/b
if frac(d)==0:ckpwr()
if b>c:disp()
b=5;d=a/b
if frac(d)==0:ckpwr()
if b>c:disp()
b=7;d=a/b
if frac(d)==0:ckpwr()
if b>c:disp()
b=11;d=a/b
if frac(d)==0:ckpwr()
if b>c:disp
()


上の方針でも書きましたが、disp() は本来 (Casio Basic プログラムを見れば分かります) 全ての素因数探索が終わったあとで実行されます。従って、ここに disp() は不要で、break に置き換えることで、無駄な関数呼び出しを避けます。これにより実行時間も多少は短くなります。disp() は今回のスクリプトの最後で1度だけ実行されれば問題ありません。

なお、ここでは if 文 は1行で書いていますが、Python の公式サイトでは、この書き方を推奨していません。

変数 b に着目すると、b の値が 2, 4, 5, 7, 11 と異なるだけで、あとは同じ処理が列挙されています。これを反復処理に書き換えるには、リストを使った for 文が使えます。

この b の値は素数なので、これらの素数を要素としたリスト prime_list = [2, 3, 5, 7, 11] を定義して使うことにします。リストは、かぎ括弧を使って定義できます。

そして、以下のような for 文が書けるのは、Python のようなオブジェクト指向言語に 特有の便利なところです。

反復処理への書き換えは、以下のようになります。


prime_list=[2,3,5,7,11]
for b in prime_list:
 d=a/b
 if frac(d)==0ckpwr()
 if b>c: break


リスト prime_list の要素を最初から最後まで順に b に適用し、反復処理を行うことができます。リストの要素の最後が適用されれば反復が終わります。このように、反復回数を明示的に設定せずに、リストというシーケンス型オブジェクトに依存して for 文の反復が行えるのは、オブジェクト指向言語の大きな特徴です。Casio BasicC 言語のような関数型言語では出てこないと思います。

余談になりますが、for i in range(10): という記述は、i を 0 ~ 9 まで1つづつ適用するのですが、実はリストの要素を生成しながら、その要素を順次適用して反復を行うといった内部動作をしていて、リストの要素を最初おから最後まで順に適用しているわけです。

もう少し詳しく言えば、range(10) は、リスト [0, 1, 2, 3, 4, 5, 6, 7, 8, 9] の要素を順次生成し、その要素を 変数 i  に適用しながら反復処理を行いますが、for 文ではリスト全体を生成するのではなく、要素を生成しながら適用させるといった働き (ジェネレータの動作) をしています。実際にメモリ上にリストを生成しているのではないので、リストを格納するメモリが節約されます。リストが非常に大きな場合は、メモリ節約の効果が極めて大きくなり、このような時にはジェネレータの効果が絶大になります。

一方、
num = list(range(10))

と書くと、実際にリストを生成して、リスト num = [0, 1, 2, 3, 4, 5, 6, 7, 8, 9] が定義されます。


9.5 素因数を探索する残りの部分のスリム化

以下の部分 (while 1: 以下 144行の部分) を反復処理で書き換えます。
while 1:
 b+=2;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=4;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 ・・・
 ・・・
 ・・・
 b+=2;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=10;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break

変数 b に着目すると、b を増やしながら同じ処理を列挙しています。但し b の増分に簡単な規則性はなく、簡単な計算では求められません。この増分がどうやって決められているのかは、Casio Basic プログラムを公開している記事を参照してください。
fx-5800P 素因数分解 - 高速化

Python の for 文では、既に上で紹介したように、b の増分を要素としたリストを使えば簡単に反復処理ができます!

増分 (インクリメント) のリストなので、リストの名前を increment とし、以下のように定義します。

increment = [2, 4, 2, 4, 6, 2, 6, 4, 2, 4, 6, 6, 2, 6, 4, 2, 6, 4, 6, 8, 4, 2, 4, 2, 4, 8, 6, 4, 6, 2, 4, 6, 2, 6, 6, 4, 2, 4, 6, 2, 6, 4, 2, 4, 2, 10, 2, 10]

要素数は48個あります。

反復処理に書き換えると、以下のようになります。

increment=[2, 4, 2, 4, 6, 2, 6, 4, 2, 4, 6, 6, 2, 6, 4, 2, 6, 4, 6, 8, 4, 2, 4, 2, 4, 8, 6, 4, 6, 2, 4, 6, 2, 6, 6, 4, 2, 4, 6, 2, 6, 4, 2, 4, 2, 10, 2, 10]
while 1:
 for i in increment:
  b+=i; d=a/b
  if frac(d)==0: ck_pwr()
  if b>c: break
 if b>c: break

反復処理への書き換えに伴い、赤文字と青文字部分を追加しています。

赤文字は、for 文による反復のために追加した部分です。
リスト increment の要素を最初から最後まで順に 変数 i に適用しながら、反復を行います。
for 文の中にある if b>c: break の break は、b>c が成立するときに  for ループから抜けるだけで、while ループからの脱出ができません。ここでは一気に while ループからも脱出する動作にしたいところです (オリジナルの Casio Basic のソースでは Goto によるジャンプで、一気にループから脱出しています) 。

そこで、for ループから抜けたところで、もう一度青文字で追加した if b>c: break を実行することで、while ループからも一気に抜けるようにします。

=====

ところで、リスト increment の定義は、要素が48個もあるので、長々と記述します。電卓の狭い画面では見づらいので、改行して見やすくします。

スクリプトの動作に影響を与えずにソースを改行するためには、改行したいところにバックスラッシュ \ を入れます。

但し、括弧の中は改行してもスクリプトの動作に影響を与えないので、電卓の画面で見やすいように適宜改行を入れます。

increment=\
[2,4,2,4,6,2,6,4,
 2,4,6,6,2,6,4,2,
 6,4,6,8,4,2,4,2,
 4,8,6,4,6,2,4,6,
 2,6,6,4,2,4,6,2,
 6,4,2,4,2,10,2,10]

while 1:
 for i in increment:
  b+=id=a/b
  if frac(d)==0ckpwr()
  if b>cbreak
 if b>c
break


格段に無駄を省いて、スリム化できました!


9.6 スリム化したスクリプト - FactorG2.py

291行もあったスクリプトが、反復処理を使って書き換えることで 91行までスリム化できました。

Note: Casio Python スクリプトファイルのサイズ制限
今回、実際に長大なスクリプトを編集している際、Casio Python は、300 行を超えると電卓上ではこれ以上改行がができなくなることが分かりました。
さらに、PC上のエディタで 300 行を超えるスクリプトを書いてから電卓に転送して、このスクリプトファイルを開こうとすると、"Invalid Data Size" という通信エラーが発生し、ファイルを開けません。
300 行で改行ができなくてなっても、既存行内への入力はできます。おそらくファイルサイズの上限までは行内の入力は可能で、保存も可能です。ファイルサイズの上限はまだ確認できていません。
1ファイルあたり最大 300 行の制限は CGモデル (fx-CG50) の Casio Python の仕様
1ファイルあたり最大 150行の制限は FXモデル (fx-9750GIII, fx-9860GIII) の Casio Python の仕様


高速素因数分解 - 10桁対応版:FactorG2.py のダウンロード

fx-CG50 Pythonモード:高速素因数分解 - FactorG2.py
"""Sample script

Exercise;
ported from Casio Basic
"Prime Factor 10 digits"
factorG.py
ver 2.0

by Krtyski/e-Gadget
"""

from u import *

def ckpwr():
 global a,b,c,d,e,z
 e+=1
 z[e]=b
 while 1:
  a=int(d)
  z[e+11]+=1
  d=a/b
  if frac(d):
   break
 c=int(sqrt(a))

def
disp():
 global a,e,f,z
 clear_screen()
 locate(00f3'm'0)
 locate(160': ' str(len(inp)) + ' digits'3'm'0)
 line(0,15,383,15,1,0)
 for i in range(112):
  if i<=e:
   locate(,,z[i], 1'm'0)
   locate(10i'^('2'm'0)
   locate(12iz[i+21], 1'm'0)
   locate(15i')'2'm'0)
 locate
(00''0'm'1)


while
1:
 try:
  inp=str(eval(input('Number:')))
 except (SyntaxError,TypeError,NameError) as e:
  print(e)
  print('*must be number or\n expression')
  continue
 if '.' in inp:
  print('*must be integer')
  continue
 elif inp.isdigit():
  if len(inp)>10:
   print('*must be 10 or less')
   continue
  else:
   f=int(inp)
   break
 else:
  continue

z=list(range(23))
for e in range(1,23):
 z[e]=0
e=0
a=f
c=int(sqrt(a))

prime_list=[2,3,5,7,11]
for b in prime_list:
 d=a/b
 if frac(d)==0ck_pwr()
 if b>c: break


increment=\
[2,4,2,4,6,2,6,4,
 2,4,6,6,2,6,4,2,
 6,4,6,8,4,2,4,2,
 4,8,6,4,6,2,4,6,
 2,6,6,4,2,4,6,2,
 6,4,2,4,2,10,2,10]

while 1:
 for i in increment:
  b+=id=a/b
  if frac(d)==0ck_pwr()
  if b>cbreak
 if b>c
break

if a>1:
 e+=1
 z[e]=a
 a=1
 z[e+11]=1

disp
()




FactorG2_output 

この入力値での実行時間は 0.8秒後半の 0.9 秒で、前回の FactorG1.py より僅かに速くなっています。




目 次

前の記事 - 8. シェル画面入力の工夫

次の記事 - 10. 大きな数の計算





応援クリックをお願いします。励みになるので...
にほんブログ村 IT技術ブログ 開発言語へ


 


keywords: fx-CG50Pythonプログラム関数電卓

リンク集 | ブログ内マップ

テーマ : プログラム関数電卓
ジャンル : コンピュータ

Casio Python - シェル画面入力の工夫

Python Casio Python
 Casioグラフ関数電卓の Python を使ってみる
     - シェル画面入力の工夫:高速素因数分解(1) 
目次


初版:2020/07/23
追記修正:2020/11/14

前の記事 - 7. テキスト出力関数の追加 |  次の記事 - 9. Python らしい反復処理


8. シェル画面入力の工夫:高速素因数分解(1) <fx-CG50 OS3.40以降>

本連載は、Casio Basic で書いた フラクタル図形 - シダの葉描画プログラムを Casio Python に移植したことから始まっています。そして、その爆速っぷりから Casio Python に興味を持ち、本連載を始めたいきさつを書きました。

プログラミング習得は、実際に動く面白いプログラムを書きながら細かいところを覚えてゆくのが最善の方法だと思っています。そこで前回までは、以前 Casio Basic で作ったプログラム - モンテカルロ法シミュレーション を Casio Python に移植しながら、学習してきました。

今回の記事からしばらくは、以前 Casio Basic で作った高速素因数分解プログラム - FactorGCasio Python に移植してゆきます。

このプログラムの詳細は、fx-5800P 素因数分解 - 高速化 をご覧ください。高速素因数分解のロジックを詳しく説明しています。4倍の高速化が達成されています。

FactorG_CasioBasic FactorG4_comparer 
左の画面は、fx-CG50 で実行した Casio Basic プログラム FactorG の出力です。
右の画面は、Casio Python で作った FactorG スクリプト最終版での出力です。

オリジナルのCasio Basicプログラムは、1画面に全ての素因数が表示しきれない時は、特定キーを押して次のページ / 前のページ を切り替えて表示できるようにしています。一方、現バージョン (OS3.40) の Casio Python はキー操作でインタラクティブにキー入力を取得するために必要な関数が用意されていないので、全ての素因数を1画面に表示する必要があります。そのため中サイズのフォントでグラフィックス画面に出力しています。

FactorG4_10digit 
これは、7849516203 の素因数分解の結果ですが、fx-CG50 で走らせた Casio Basic のプログラムの実行時間は 12秒 で、Casio Python では、0.9 秒 となり、10倍以上高速に計算しています。やはり速いです!

FactorG4_15digit 
この画面は、入力値の桁数制限を15桁まで拡張したスクリプトでの結果です。

Casio Basic は10桁の精度がありますが、Casio Python は15桁の精度があるので、15桁まで拡張可能です。素因数の種類が12個以上の場合は、画面を分割して右側に表示するようにしています。

Casio Basic では入力時に数値だけでなく、計算式も入力でき、計算結果が入力値になります。そこで同じ機能を Casio Python スクリプトに実装しました。上の例では、素数の積の式 2*3*5*7*11*13*17*19*23*29*31*37*41 を入力した結果です。

さて、Casio Python は、シェル画面でしか入力ができません。今回は、とりあえず Casio Basic プログラムと同じように入力できるようにします。その上で、Casio Basic プログラムをできるだけ忠実に Casio Python へ移植して、とりあえず動くようにします。

そして、次回以降で Python らしいスクリプトに改造してゆくことにします。

今回作成するスクリプトは、以下からダウンロードできます;
fx-CG50 Casio Python 用 高速素因数分解スクリプト - FactorG1.py




8.1 高速素因数分解プログラム - FactorG のざっくりした分析

先ずは、グラフ関数電卓 (fx-9860G シリーズや fx-CGシリーズ) の Casio Basic で動作する FactorG のコードを眺めながら、Casio Python への移植作業の方針をざっくりと考えてみます。

FactorG は全く同じコードで fx-9860G シリーズと fx-CGシリーズで動作します。ここでは fx-CG用の g3m ファイルを具体的に眺めてみます。

このプログラムは結構な行数があるので、電卓でポチポチと入力するのは大変です。そこで、PC上でエディタを使って、編集&入力すべきでしょう。Casio Basic のコードをエディタで読み込み、それを Casio Python スクリプトに変更・編集しました。動作確認のためには、fx-CG50 に転送し、実際に走らせる必要があります。

以下からダウンロードした zip ファイルには グラフ関数電卓用の Casio Basic プログラムファイルに加えて、テキストファイルも含んでいるので、これをエディタで読み込み、Casio Python に書き換えてゆきます。

グラフ関数電卓用高速素因数分解プログラム - FactorG のダウンロード


FactorG_g3m_comment


WFSUB


以下のようなざっくりとした作戦で進めてゆこうと思います。
  1. 入力と入力値のチェックの部分は、Casio Basic と同様の動作になるように Casio Python 向けのロジックで書きます。 
  2. 変数の初期化は、サクッと移植します。
  3. 素数が素因数かどうかをチェックして素因数の乗数を計算するサブルーチン WFSUB を関数化し (ckpwr() の作成)、Goto 1 のジャンプ先の処理も関数化します (disp() の作成)。
  4. While ループのエラストテレスの篩い部分も、上で関数化したものを使って変更します。但し、While ループを Goto で脱出している部分は、Casio Python ではそのまま実装できません。Casio Python には goto が実装されておらず、while ループの脱出は break のみで行えるので、このあたりはロジックの変更が必要です。
  5. Lbl 1Lbl 2 の部分はまとめて関数化します (disp() の作成)。
  6. インタラクティブなキー入力は、Casio Python ではできないので、この部分は移植不可能です。そこでdisp() 関数の処理では、全ての素因数を1画面内で表示できるように、中サイズのフォントを使ってグラフィックス画面(描画画面) へ出力するようにします。
  7. 移植作業全体として、変数のスコープ、大域変数 (グローバル変数) の扱いと Goto の代わりの処理に留意します。

8.2 入力と入力値のチェック

この部分の Casio Basic のコードは以下です。

"NUMBER"?→F
If F<1 Or F≧1Exp10
Then "NUMBER MUST BE ≧1 And <1Exp10"
Stop
IfEnd
If F≠Int (F)
Then "NUMBER MUST BE AN INTEGER"
Stop
IfEnd

プログラムを強制終了させる Stop に相当する関数は、一般の Python では quit()exit() 関数が用意されていますが、 Casio Python には実装されていません。

そこで、ロジックを変更して、この入力部分を while 1:  無限ループの中に入れ、10桁以下の整数 の時だけ break でループから脱出し、10桁を超える入力や小数の入力の時はループを継続して、再入力させるようにします。つまり、正しく入力できるまでループが回り続け、先に進めないロジックに変更します。

さて、Casio Python では、シェル画面で使える input() が唯一の入力関数で、この関数は文字列を返します。
input()

さらに Casio Basic では ?→ 命令は、関数や式を入力した時の結果の数値を取得できますので、Casio Python でも同じような機能を盛り込みます。そこで、式 (文字列) を評価してその結果を返す eval() 関数を使います。

f=eval(input('Number:'))

入力値が10桁を超える場合、そして小数の場合、さらにテンキー以外のキー入力の場合は while ループを継続して、入力待ちの状態を継続し、それ以外の場合を整数入力になるように分岐方法を考え、その時は変数 f に整数を格納してループを抜けるロジックを考えました。

一例として、以下のようにしました。

while 1:
 inp=str(eval(input('Number:')))
 if '.' in inp:
  print('must be integer')
  continue
 elif inp.isdigit():
  if len(inp)>10:
   print('must be <= 10 digit')
   continue
  else:
   f=int(inp)
   break
 else:
  continue


これで、Casio Basic の入力とほぼ同じ機能が得られます。

小数の判定には、文字列中に小数点 . があるかどうかで調る
これは、Pythonらしい記述にしてみました。
if '.' in inp:
文字列 inp の中に、1文字だけの文字列 '.' があれば...という意味です。
inp に 小数点が含まれていれば、'.' in inpTrue (真) を返します。
小数点が含まれていれば、小数を入力したことになるので、continuewhile ループを継続し、入力待ちの状態を継続します。

入力が数値かどうかを判定後、桁数の判定を行う
入力が数値かどうかは、文字列が数値であるかどうかを判定する関数 isdigit() を使います。この関数は、文字列のまま、その内容が数値であるかどうかを調べられます。isdigit() は、文字列の帰属関数であり、今は文字列 inp について調べるので、
inp.isdigit()
と記述します。
isdigit() は数値であるなら、True (真) を返すので、if 文の条件として使いました。
print() での出力文字列は、電卓の画面に収まるように、簡潔に短くしました。

入力が小数でなく、11桁以上でない時、入力文字列を数値に変換して変数 f に格納
文字列 inp の長さを調べる関数 len() を使って、len(inp) が10を超えると print('must be <= 10') でエラー表示をして、continue でループを継続させ、入力待ちを継続するようにしました。 
整数の文字列の長さを10を超えない時は、10桁以下の整数であると判定し、文字列を整数に変換して 変数 f に代入し、このときだけ while 1: ループを抜けて、処理を先に進めるようにしました。

最後の else: 節は、上記以外の入力のとき while ループを継続することで、入力待ちを継続


さて、このスクリプトでは、テンキー以外を押した時に強制終了しますが、Casio Basic プログラムでは強制終了しません。

[EXE] を押しただけのときの挙動の違い - SyntaxError
Casio Basic では、数値を入力せず単に [EXE]を押すだけでは、入力待ちになります。一方、Casio Python の上のスクリプトでは SyntaxError となり、強制終了します。

 式として不完全な入力をした時の挙動の違い - SyntaxError
Casio Basic では、(、()、*、+、-、/、** など、式として不完全な入力後に [EXE] を押して入力を確定すると SyntaxError が発生しますが、プログラムは強制終了せず、修正して正しい式を入力して [EXE] を押すまで入力待ちが継続します。これは ?→  命令の仕様です。一方、Casio Python では、SyntaxError となり強制終了します。

引数に誤りのある関数を入力した時の挙動の違い - TypeError
Casio Basic では、例えば引数を指定しないで sin()in() などを入力して [EXE] を押すと SyntaxError が発生しますが、プログラムは強制終了せず、修正して正しく入力して [EXE] を押すまで入力待ちが継続します。これは ?→ 命令の仕様です。一方、Casio Python では、引数を正しく指定せずに入力すると TypeError が発生し、強制終了します。

アルファベット入力時の挙動の違い
Casio Basic  では、アルファベットを入力すると、アルファベットに何か数値が格納されているときはその数値を取得し、そうでない時は 0 を取得します。一方、Casio Python では、NameError が発生して強制終了します。 

Casio Python でのこれらの挙動は、上のスクリプトで利用している eval() で発生するエラーです。eval() は必要なのでこれを使いつつ、これらのエラーが発生しても強制終了させないために、スクリプトを改造してみます。赤文字が追加する部分です。

while 1:
 try:
  inp=str(eval(input('Number:')))
 except (SyntaxError, TypeError, NameError) as e:
  print(e)
  print('*must be number or\n expression')
  continue
 if '.' in inp:
  print('*must be integer')
  continue
 elif inp.isdigit():
  if len(inp)>10:
   print('*must be 10 digit\n or less')
   continue
  else:
   f=int(inp)
   break
 else:
  continue

エラーが発生すると、スクリプトの制御を失い、強制終了してしまいます。そこで、エラーを検知しても制御を取り戻すことで、強制終了を回避できます。そのような時には Try ステートメントが役立ちます。

エラーを検出し制御を取り戻す - Try / except
エラーの検出には、Try ステートメントと except ステートメントを利用します。
Try: 節except: 節 の間に、エラーを検出したい処理を挟みます。
今回、エラーを検出したい処理は
inp = str(eval(input('Number:')) 
で、インデントレベルを下げて記述します。

エラーは、eval() 関数で発生し、エラーの種類は上で調べたように、SyntaxErrorTypeErrorNameError の3種類です。これら3種類のエラーを検出するためには、
except (SyntaxError, TypeError, NameError): とタプル型で複数のエラーの種類を記述すれば、これら3種類のエラーのいずれかが発生した時、制御を取り戻せます。

さらに、
except (SyntaxError, TypeError, NameError) as e:
と記述すると、文字列変数 e にエラーの詳細が格納されます。
今回のスクリプトでは、
・SyntaxError 発生時は、e に "invalid syntax" が、
・TypeError 発生時は、e"function missing 1 required positional argument" が、
・NameError 発生時は、e"'x' is not defined" が、
格納されます。

そこで、print(e) により、エラーの詳細を出力します。
except 節は、Try 節と同じインデントレベルにします。print(e) もそのインデントが必要です。

続けて、エラーにならない入力のアドバイスを出力します。同じインデントレベルで下記を追加します。
print('*must be number or\n expression')
電卓の画面は狭いので、出力文字列の中にエスケープシーケンス \n を追加し、そこで改行しています。このエスケープシーケンスについては以前紹介しています。
出力は以下のようになります。
*must be number or
 expression


except 節の最後に continue を記述し、while へジャンプさせ、ループを継続させます。


上のスクリプトで、青色で示した
print('*must be 10 digit\n or less')
も、表現を変更し、さらに電卓画面に収まりやすいように、\n で改行します。出力は以下のようになります。
*must be 10 digit
 or less


以上で、入力と入力値のチェックの部分の移植が完了しました。


8.3 変数の初期化

上で検討した 入力と入力値のチェック では、Casio Basic の {1,22}→DIM Mat Z については無視しました。これは見つかった素因数とその乗数を格納するための行列 Z を定義して必要なメモリを確保するためのものです。メモリ確保は重要なのでプログラムの冒頭で実施してるわけです。

Casio Python への移植では、行列 Z の代わりに リストz (配列) を使います。リスト定義は、変数の初期化の部分で一緒に行うことにします。Casio Basic のコードを忠実に移植します。

z=list(range(23))
for e in range(1,23):
 z[e]=0
z[1]=0
z[12]=0
e=0
a=f
c=int(sqrt(a))



リスト Z の定義と初期化
先ず、要素数22個のリストz を以下のように定義します。
z=list(range(23))
この状態だと、リストの要素は [0, 1, 2, 3, 4,...20, 21, 22] となります。
次に、全ての要素を 0 にするため、for 文を使って、
for e in range(1,23):
 z[e]=0

続いて、z[1]=0z[12]=0 としています。z[1] は 1つめに見つかった素因数で、z[12] は1つめに見つかった素因数の乗数です。これらはすでに 0 で初期化されているので不要ですが、重要なので明示的に初期化を2回行っている Casio Basic に忠実に移植しました。

その他変数の初期化
リストz のインデックス e0 で初期化、入力値を素因数で割った値 a は 入力値 f で初期化、素因数探索範囲の最大値 c√a に近い整数で初期化しています。

8.1 と 8.2 の検討結果を、以下にまとめます。

while 1:
 try:
  inp=str(eval(input('Number:')))
 except (SyntaxError,TypeError,NameError) as e:
  print(e)
  print('*must be number or\n expression')
  continue
 if '.' in inp:
  print('*must be integer')
  continue
 elif inp.isdigit():
  if len(inp)>10:
   print('*must be 10 or less')
   continue
  else:
   f=int(inp)
   break
 else:
  continue

z=list(range(23))
for e in range(1,23):
 z[e]=0
z[1]=0
z[12]=0
e=0
a=f
c=int(sqrt(a))


8.4 素数が素因数かどうかをチェック

Casio Basic のコード;

2→B:A÷B→D:Frac (D)=0⇒Prog "WFSUB":B>C⇒Goto 1

Casio Python に移植するには、スクリプト構造を検討します。Casio Python には goto がないので、Goto 1 のジャンプ先である Lbl 1 とそれに続く Lbl 2 までを関数 disp() にまとめ、Goto 1disp() に置き換えます。

Prog "WFSUB" は、関数コールに置換え、WFSUB の内容を ckpwr() にまとめます。サブルーチン WFSUB は、見つかった素因数の乗数を求めるので、check power (乗数をチェック) の意味で、ckpwr() としました。

関数 dsip()ckpwr() は後で作るとして、先ずは上記の1行を Casio Python で書き換えます。

b=2;d=a/b
if frac(d)==0:ckpwr()
if b>c:disp()


なお、Casio Basic に実装されている関数 Frac() は、小数から小数部を得ます。ところが Casio Python には Frac() の相当する関数が無いので、新たに作成して使いました。

def frac(x):
 return x-int(x)


ところで、if frac(d)==0:ckpwr() は、本来下記のように記述することが Python 公式サイトで推奨されています。

if frac(d)==0:
 ckpwr()


しかし、1行に書いても正しく動作します。

b=2;d=a/b
という記述も、Python 公式サイトでは以下のように2行で書くことを推奨しています。
b=2
d=a/b



Note: Python の 単文 と 複文
b=2d=a/b は単文で、セミコロン ; で区切って、複数の単文を1行に記述できます。
しかし、if 文のような複文を セミコロン ; で区切って1行に記述しても正常に動作しません。
例えば、以下のように記述すると、エラーになります。
b=2;d=a/b;if frac(d):ckpwr();if b>c:disp() #エラー

これまでのところをまとめ、スクリプトの全体像は以下のようにします。

from u import *

def ckpwr():
 #あとで作ります

def disp():
 #あとで作ります

def frac(x):
 retuen x-int(x)

while
1:
 try:
  inp=str(eval(input('Number:')))
 except (SyntaxError,TypeError,NameError) as e:
  print(e)
  print('*must be number or\n expression')
  continue
 if '.' in inp:
  print('*must be integer')
  continue
 elif inp.isdigit():
  if len(inp)>10:
   print('*must be 10 or less')
   continue
  else:
   f=int(inp)
   break
 else:
  continue

z=list(range(23))
for e in range(1,23):
 z[e]=0
z[1]=0
z[12]=0
e=0
a=f
c=int(sqrt(a))

b=2;d=a/b
if frac(d)==0:ckpwr()
if b>c:disp()
b=3;d=a/b
if frac(d)==0:ckpwr()
if b>c:disp()
b=5;d=a/b
if frac(d)==0:ckpwr()
if b>c:disp()
b=7;d=a/b
if frac(d)==0:ckpwr()
if b>c:disp()
b=11;d=a/b
if frac(d)==0:ckpwr()
if b>c:disp()

#2の倍数、3の倍数、5の倍数、7の倍数以外の探索数で "エラストテレスの篩い" で素因数を探索
while 1:
 #この部分もあとで作ります


8.5 "エラストテレスの篩い" で素因数を探索

上記の「素数が素因数かどうかをチェックする」は、素数 2,3、5、7、11 探索数として、"エラストテレスの篩い" で素因数を探索しています。

それに続いて、2の倍数、3の倍数、5の倍数、7の倍数以外の探索数を用いて、素因数を探索します。Casio Basic のコードは以下のようになっています。

While 1
B+2→B:A/B→D:Frac (D)=0⇒Prog "WFSUB":B>C⇒Goto 1
B+4→B:A/B→D:Frac (D)=0⇒Prog "WFSUB":B>C⇒Goto 1
B+2→B:A/B→D:Frac (D)=0⇒Prog "WFSUB":B>C⇒Goto 1
B+4→B:A/B→D:Frac (D)=0⇒Prog "WFSUB":B>C⇒Goto 1
B+6→B:A/B→D:Frac (D)=0⇒Prog "WFSUB":B>C⇒Goto 1
B+2→B:A/B→D:Frac (D)=0⇒Prog "WFSUB":B>C⇒Goto 1
B+6→B:A/B→D:Frac (D)=0⇒Prog "WFSUB":B>C⇒Goto 1
・・・
・・・
(48行)
・・・
WhileEnd


Whileループの中のそれぞれ1行は、上記で作ったものと殆ど同じです。行頭で、探索数に 2, 4, 6, 8, 10 のいずれかを順に加算しているところだけが違っています。

この探索数の決め方についての数学的背景については、fx-5800P 素因数分解 - 高速化 に書いてあるので、気になる場合は参照してうださい。

While ループ内の一番上の行

B+2→B:A/B→D:Frac(D)=0⇒Prog "WFSUB":B>C⇒Goto 1

Casio Pyton に書き換えると、以下のようになります;

b+=2;d=a/b
if frac(d)==0:ckpwr()
if b>c:disp()


この書式は While ループの中に 48 行あります。そこで、Casio Python で while 1: の下に、上記の処理を48個記述すればそれで良いかと言えば、それがうまく動作しません

何が問題なのか?

Casio Basic で上記のコードが実行されると、Whle ループの中から Goto で対応する Lbl へジャンプすると同時にループを抜けます。しかし、Casio Python では、while ループの中で disp() を呼び出し、disp() の処理が終われば、呼出位置に戻り、ループを抜けられません。従って、うまく動作しません。

そこで、Goto 1 の代わりに break でループを抜け、抜けた直後に disp() を実行するように記述すれば、Casio Basic と同じ動作になります。

具体的には、以下のようになります。

while 1:
 b+=2;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=4;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=2;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 ・・・
 ・・・
 (48行)
 ・・・
disp()



これで、全体のロジックが決まったので、これまで作った部分 と これから作る部分 をまとめて、スクリプト全体の構成を以下に示します。


fx-CG50 Pythonモード:高速素因数分解 - FactorG1.py
from u import *

def ckpwr():
 #あとで作ります (8.5 参照)

def disp():
 #あとで作ります (8.6 参照)

def frac(x):
 retuen x-int(x)

while
1:
 try:
  inp=str(eval(input('Number:')))
 except (SyntaxError,TypeError,NameError) as e:
  print(e)
  print('*must be number or\n expression')
  continue
 if '.' in inp:
  print('*must be integer')
  continue
 elif inp.isdigit():
  if len(inp)>10:
   print('*must be 10 or less')
   continue
  else:
   f=int(inp)
   break
 else:
  continue

z=list(range(23))
for e in range(1,23):
 z[e]=0
z[1]=0
z[12]=0
e=0
a=f
c=int(sqrt(a))

b=2;d=a/b
if frac(d)==0:ckpwr()
if b>c:disp()
b=3;d=a/b
if frac(d)==0:ckpwr()
if b>c:disp()
b=5;d=a/b
if frac(d)==0:ckpwr()
if b>c:disp()
b=7;d=a/b
if frac(d)==0:ckpwr()
if b>c:disp()
b=11;d=a/b
if frac(d)==0:ckpwr()
if b>c:disp()

while 1:
 b+=2;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=4;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=2;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=4;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=6;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=2;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=6;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=4;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=2;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=4;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=6;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=6;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=2;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=6;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=4;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=2;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=6;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=4;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=6;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=8;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=4;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=2;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=4;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=2;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=4;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=8;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=6;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=4;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=6;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=2;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=4;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=6;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=2;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=6;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=6;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=4;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=2;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=4;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=6;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=2;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=6;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=4;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=2;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=4;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=2;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=10;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=2;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break
 b+=10;d=a/b
 if frac(d)==0:ckpwr()
 if b>c:break

disp()


8.6 ckpwr() 関数の作成

Casio Basic の WFSUB サブルーチンは下記です。

WFSUB

これを Casio Python の関数 ckpwr() にまとめます。この関数名は check power (乗数をチェック) からきています。

個々の処理を Casio Python のスクリプトに書き換えるのは、特に難しいところは無いと思います。

ここで非常に重要なポイントとして、変数のスコープについて考えなければなりません。というのも、Casio Basic では全ての変数は究極のグローバル変数 (大域変数) です。一般にグローバル変数とは、1つのプログラムの中の変数はどこからでもアクセスでき、どこからでもその値の変更が有効である、そのような変数です。対してローカル変数は関数やサブルーチンの中だけでアクセスや値の変更ができ、その関数やサブルーチンの外にある同じ名前の変数には影響を及ぼさない、そのような変数です。

Casio Basic では、変数は、メモリ内にある全てのプログラムやサブルーチンで共有され、どこからでも有効にアクセスでき、どこからでもその値を有効に変更できるもので、一般のグローバル変数よりもさらにグローバル、つまり究極のグローバル変数と言えます。メモリ内のいずれかのプログラムで変数の内容が変われば、全てのプログラムやサブルーチンで同じ名前の変数を参照すると、変化した値になっています。

ところが、Casio Python では、変数スコープについて、C や C# などとも少々事情が異なります。変数の代入が行われる場所が関数の中か外かでスコープが決まります。
  • 関数内で変数に代入が行われると、その変数はローカル変数 (局所変数) となり、その関数内だけで有効になります。その変数は関数外では不定になり、関数外で使えばエラーになります。
  • 関数外で変数に代入が行われると、その変数はグローバル変数となり、その変数は全ての関数内で有効になります。
変数が有効になる範囲を変数のスコープといいます。

Casio Python の変数スコープの仕様は一般的なもので、むしろ Casio Basic が特殊です。

但し、C, C++, C# に比べて Casio Python の変数スコープは少々曖昧に思えるところがあります。スクリプトを実行する際、まず最初にスクリプトの上から下までざっと解釈を行います。そして、関数を使う前(上)で関数定義を行う必要があります。

もし、関数定義の後(下)で関数の外で変数に代入すると、その変数はグローバル変数となり、本来は関数内でも有効になるはずですが、関数定義よりは後に変数に代入されているので、関数内で使おうとすると不定となり、エラーになります。この点で、Casio Python での変数スコープは曖昧で、
  • 変数への代入が関数内で行われるのか関数外で行われるかで、変数スコープ(グローバルかローカルか)が決まります。
  • 変数への代入が、どこで行われるかによって、有効になるか無効になるか (エラーになるか) が決まります。
そこで、Python では変数のスコープは原則ローカルだと思ってコーディングすることが重要です。もし関数内でグローバル変数を使いたい場合は、上の曖昧さを排除するために、関数内で global スレートメントを使って、明示的にグローバルだと宣言するのが良いと思います。

さて、関数 ckpwr() を定義する場合は、global ステートメントで 変数を列挙して、それらがグローバル変数であることを明示的に宣言します。

def ckpwr():
 global a,b,c,d,e,z
 e+=1
 z[e]=b
 while 1:
  a=int(d)
  z[e+11]+=1
  d=a/b
  if frac(d):
   break
 c=int(sqrt(a))


8.7 disp() 関数の作成

Casio Python スクリプトに書き換える対象の Casio Basic のコードは以下です。

Lbl 1
If A>1
Then Isz E
A→Mat Z[1,E]
1→A
1→Mat Z[1,E+11]
IfEnd
Int (E/6)->D
E-6*D>0⇒Isz D
1→C

Lbl 2
ClrText
Locate 1,1,F
Locate 12,1,C
Locate 13,1,":"
Locate 14,1,D
6*(C-1)+1⇒B
Locate 1,2,Mat Z[1,B]
Locate 11,2,"^("
Locate 13,2,Mat Z[1,B+11]
Locate 16,2,")"
If B+1≦E
Then Locate 1,3,Mat Z[1,B+1]
Locate 11,3,"^("
Locate 13,3,Mat Z[1,B+12]
Locate 16,3,")"
IfEnd
If B+2≦E
Then Locate 1,4,Mat Z[1,B+2]
Locate 11,4,"^("
Locate 13,4,Mat Z[1,B+13]
Locate 16,4,")"
IfEnd
If B+3≦E
Then Locate 1,5,Mat Z[1,B+3]
Locate 11,5,"^("
Locate 13,5,Mat Z[1,B+14]
Locate 16,5,")"
IfEnd
If B+4≦E
Then Locate 1,6,Mat Z[1,B+4]
Locate 11,6,"^("
Locate 13,6,Mat Z[1,B+15]
Locate 16,6,")"
IfEnd
If B+5≦E
Then Locate 1,7,Mat Z[1,B+5]
Locate 11,7,"^("
Locate 13,7,Mat Z[1,B+16]
Locate 16,7,")"
IfEnd

Casio Basic の Locate に対応する Casio Python の関数 locate() は既に作っており、ユーザーモジュール u.py (ver 1.3) から呼び出して使います。

locate() の書式:
locate(column, row, obj, color=1, size='m', show=1)

Casio Basic の Locate コマンドの第1引数、第2引数、第3引数は、そのまま column, row, obj に対応します。
予め移植の方針として、中サイズのフォントを使うことにしているので、第4引数には 'm' を指定します。
locate() を沢山使うので、第5引数には 0 を指定し、最後に locate(0, 0, '', 'm', 1) を追加して、第5引数で 1 を指定することで、全ての locate() の出力を一気に VRAM から画面へ転送することにします。

関数内で使う変数については、全て global ステートメントでグローバル変数の宣言を行います。

def disp():
 global a,b,c,d,e,f,z
 if a>1:
  e+=1
  z[e]=a
  a=1
  z[e+11]=1
 d=int(e/6)
 if e-6*d>0:
  d+=1
 c=1

 clear_screen()
 locate(00f3'm'0)
 b 6*(c-1)+1
 locate(0, 1, z[b], 1, 'm', 0)
 locate(10, 1, '^(', 2, 'm', 0)
 locate(12, 1, z[b+11], 1, 'm', 0)
 locate(15, 1, ')', 2, 'm', 0)
 if b+1<=e:
  locate(0, 2, z[b+1], 1, 'm', 0)
  locate(10, 2, '^(', 2, 'm', 0)
  locate(12, 2, z[b+12], 1, 'm', 0)
  locate(15, 2, ')', 2, 'm', 0)
 if b+2<=e:
  locate(0, 3, z[b+2], 1, 'm', 0)
  locate(10, 3, '^(', 2, 'm', 0)
  locate(12, 3, z[b+13], 1, 'm', 0)
  locate(15, 3, ')', 2, 'm', 0)
 if b+3<=e:
  locate(0, 4, z[b+3], 1, 'm', 0)
  locate(10, 4, '^(', 2, 'm', 0)
  locate(12, 4, z[b+14], 1, 'm', 0)
  locate(15, 4, ')', 2, 'm', 0)
 if b+4<=e:
  locate(0, 5, z[b+4], 1, 'm', 0)
  locate(10, 5, '^(', 2, 'm', 0)
  locate(12, 5, z[b+15], 1, 'm', 0)
  locate(15, 5, ')', 2, 'm', 0)
 if b+5<=e:
  locate(0, 6, z[b+5], 1, 'm', 0)
  locate(10, 6, '^(', 2, 'm', 0)
  locate(12, 6, z[b+16], 1, 'm', 0)
  locate(15, 6, ')', 2, 'm', 0)
 if b+6<=e:
  locate(0, 7, z[b+6], 1, 'm', 0)
  locate(10, 7, '^(', 2, 'm', 0)
  locate(12, 7, z[b+17], 1, 'm', 0)
  locate(15, 7, ')', 2, 'm', 0)
 if b+7<=e:
  locate(0, 8, z[b+7], 1, 'm', 0)
  locate(10, 8, '^(', 2, 'm', 0)
  locate(12, 8, z[b+18], 1, 'm', 0)
  locate(15, 8, ')', 2, 'm', 0)
 if b+8<=e:
  locate(0, 9, z[b+8], 1, 'm', 0)
  locate(10, 9, '^(', 2, 'm', 0)
  locate(12, 9, z[b+19], 1, 'm', 0)
  locate(15, 9, ')', 2, 'm', 0)
 if b+9<=e:
  locate(0, 10, z[b+9], 1, 'm', 0)
  locate(10, 10, '^(', 2, 'm', 0)
  locate(12, 10, z[b+20], 1, 'm', 0)
  locate(15, 10, ')', 2, 'm', 0)
 if b+10<=e:
  locate(0, 11, z[b+10], 1, 'm', 0)
  locate(10, 11, '^(', 2, 'm', 0)
  locate(12, 11, z[b+21], 1, 'm', 0)
  locate(15, 11, ')', 2, 'm', 0)
 locate(0, 0, '', 0, 'm', 1)


このスクリプトで、グラフィックス画面へ中サイズのフォントで出力すると最大11個の素因数まで出力表示できます。これで十分かどうかを検証します。

1)
素因数として、2、3,5,7,11,13,17,19,23,29 の 10個の場合、これを得る数値はこれらの積で、
2 x 3 x 5 x 7 x 11 x 13 x 17 x 19 x 23 x 29 = 6,469,693,230 は10桁の数です。つまり結果表示は10桁あれば十分で、このスクリプトの表示行数の範囲内です。

2)
素因数として、2、3、5、7、11、13、17、19、23、29、31 の 11個の場合、これを得る数値はこれらの積で、
2 x 3 x 5 x 7 x 11 x 13 x 17 x 19 x 23 x 29 x 31 = 200,560,490,130 は12桁の数です。つまり結果表示は11桁あれば十分で、このスクリプトの表示行数の範囲内です。

3)
なお、今回作ったスクリプトは、入力値が最大10桁に制限される Casio Basic プログラムと同等なので、入力値が10桁に制限される場合の素因数の種類は最大数は 10 個 です。

以上から、10桁制限の今回のスクリプトでは、全ての素因数を表示できることが確認できました。

FactorG1_output1 



目 次

前の記事 - 7. テキスト出力関数の追加

次の記事 - 9. Python らしい反復処理





応援クリックをお願いします。励みになるので...
にほんブログ村 IT技術ブログ 開発言語へ


 


keywords: fx-CG50Pythonfx-9750GIIIfx-9860GIIIプログラム関数電卓


テーマ : プログラム関数電卓
ジャンル : コンピュータ

Casio Python - ユーザー関数 line()

Python Casio Python
 Casioグラフ関数電卓の Python を使ってみる
     - リファレンス 
目次
line()

ユーザーモジュール u.py に含まれる関数 u.py ver 1.5 ダウンロード

初版:2020/11/08

[対応モデル] - fx-CG50 OS3.20 以降、fx-9750GIII / fx-9860GIII OS3.21 以降

線分を描画します。

書 式line(x1, y1, x2, y2, color=1, show=1)

引 数
- 第1引数 - x1: 描画する線分の端の座標 (x1, y1) の x座標、0 以上の整数
- 第2引数 - y1: 描画する線分の端の座標 (x1, y1) の y座標、0 以上の整数
- 第3引数 - x2: 描画する線分のもう一方の端の座標 (x2, y2) の x座標、0 以上の整数
- 第4引数 - y2: 描画する線分のもう一方の端の座標 (x2, y2) の y左表、0 以上の整数
- 第5引数 -  color: パラメータ引数で、省略可能でデフォルトは 1 です。
       具体的な設定は grp_color() のリファレンスを参照。
- 第6引数 - show: パラメータ引数で、省略可能でデフォルトは 1 です。

引数を5つ、値だけを指定するとエラーになります ⇒ 関数の引数は こちら を参照。
引数を5つだけ設定する場合は、パラメータと共に設定します。
  例) line(5, 50, 30, 60, color=4) / line(5, 50, 30, 60, show=0)


関数定義:

from casioplot import *

def
line(x1y1, x2, y2, color=1, show=1):
 rgb = grp_color(color)
 dx x2 x1; dy = y2 - y1
 if dx==0 and dy==0: #avoid division by 0 error
  set_pixel(x1, y1, rgb)
  return
 if abs(dx)>abs(dy):
  if dx: #when dx is not 0
   k = int(dx/abs(dx)) #k=1 or -1
   slope = dy/dx
   for x in range(0, dx, k):
    set_pixel(x1+x, y1+int(x*slope), rgb)
 else:
  if dy: #when dy is not 0
   k = int(dy/abs(dy)) #k=1 or -1
   slope = dx/dy
   for y in range(0, dy, k):
    set_pixel(x1+int(y*slope), y1+y, rgb)
 if show:
  show_screen
()



スクリプトの解説:
  グラフィック出力関数の追加:line() 関数とユーザーモジュールの作成 ⇒ こちら
  grp_color() ユーザー関数について ⇒ こちら




応援クリックをお願いします。励みになるので...
にほんブログ村 IT技術ブログ 開発言語へ


 


keywords: fx-CG50Pythonfx-9750GIIIfx-9860GIIIプログラム関数電卓

リンク集 | ブログ内マップ

テーマ : プログラム関数電卓
ジャンル : コンピュータ

Casio Python - グラフィックス出力関数の追加

Python Casio Python
 Casioグラフ関数電卓の Python を使ってみる
     - グラフィックス出力関数の追加:line() とユーザーモジュールの作成 
目次


初版:2020/07/11
修正:2020/07/17
追記修正:2020/11/08

前の記事 - 5. 関数の作成と活用 |  次の記事 - 7. テキスト出力関数の追加


<fx-CG50 OS3.40 以降、fx-9750GIII, fx-9860GIII OS3.40以降>


6. グラフィックス出力関数の追加:line() とユーザーモジュール の作成

前回は、グラフィックス描画関数で使うための 色指定関数 grp_color() 関数を作り、それを使って 円描画関数 circle() を拡張しました。その上で、モンテカルロ法シミュレーションをカラー化した monteca2.py を作りました。

Casio Python の現在の言語仕様では、Casio Basic で作った実用プログラムを完全に移植できないレベルの貧弱な仕様で、発展途上と言えます。Casio Basic で実用プログラムが作れる最大の要因は、望む位置に数値や文字列を出力できる Locate コマンドと電卓のほぼ全てのキー入力を取得可能な Getkey コマンドの存在です。一方 Casio Python には LocateGetkey に相当するものがありません。

従って、プログラム=スクリプトを実行する際、シェル画面において 唯一の入力関数 input() を用いて必要な入力を済ませ、その後グラフィックス画面(描画画面)で出力を行うといった構造にするのが、現状では実用的なスクリプトを作成する最善の方法だと考えています。

Getkey 相当の関数については、カシオによる今後のOSアップデートに期待するしかありません。Locate に相当する locate() 関数は自作できるので、これについては次回紹介する予定です。

ところで、既に circle() 関数を作っていますので、今回は線分描画の line() 関数を作って、さらにユーザーモジュールに追加します。


6.1 line() 関数の作成 [2020/07/12 大幅修正]

[2020/07/12 ロジック見直し] 読者のK様からのご指摘により、垂直に近い傾きの大きな線分は 反復回数が少なくなるため飛び飛びの点描画になり、きちんと線分が描画されないことが判明しました。そこで、K様のご提案に従って、線分の傾きが常に45度を超えない処理2つを組み合わせることで、十分な反復回数を確保し、確実に線分を描画できるように修正しました。K様、ありがとうございます。

line_schematic 

点描画を反復して線分を描く際、x 座標 を1だけ増減させた時 y 座標は 傾きの分を増減させます。
最初に与えられる (x1, y1)(x2, y2) の2つの座標から、
dx = x2 - x1
dy = y2 - y1
とすると、傾き slope = dy/dx となります。

x 座標が x だけ変化したとき、y 座標は x*slope だけ変化します。(x1, y1) 座標を基準にして、x 座標 が x だけ変化したときの座標は、(x1+x, y1+x*slope) になります。この座標への点の描画を反復して線分を描くのですが、その反復回数を range(dx) で与えます。

線分の傾きが大きいと dx は小さな値になり、range(dx) で十分な反復回数が得られません。画面の左端から右端にわたる線分を描画する場合の反復回数は 191回未満になると、理屈の上では線分がまばらになります。傾きが垂直に近いと range(dx) は数回以下になってしまい、適切に線分を描けなくなります。

反復回数の一番厳しい条件として 最低でも 191 回確保するためには、傾きの絶対値が45度未満が条件で、別の表現をすれば、|dx| > |dy| が条件になります。| | は絶対値です。Python での表現では、適切な線分描画の条件は、abs(dx) > abs(dy) となります。

では、傾きの絶対値が 45 度以上になる場合は、x 座標と y 座標を入れ替えると、傾きの絶対値が 45 度未満になります。直交座標で2つの座標を入れ替えても数学的に等価であることは保証されているので、スクリプトも x 座標と y 座標を入れ替えて、数学的に等価になるようにすれば良いことになります。このときの傾きは、x と y を入れ替えて sllope = dx/dy となります。y 座標を y だけ変化させるとき、x 座標は y*slope だけ変化します。この変化は、(x1, y1) を基準にすると、座標 (x1+y*slope, y1+y) になります。

この問題は、以下の簡単なスクリプトで確認できます。

from u import *
for x in range(383):
 line(0, 0, x, 191)


三角形の領域が塗りつぶされるべきところ、そうなっていなかったことから、K様のご指摘を検証できました。今回の修正の結果、三角形の領域が正しく塗りつぶされることが確認できました。

なお、show が1の場合にVRAMから画面にデータ転送を行う処理は、関数定義の一番最後に一括して記述します。
----------

(x1, y1) と (y1, y2) の間に線分を描きます。返値はなしです。すると関数定義の1行目は次のようになります。

def line(x1, y1, x2, y2, color=1, show=1):

ここで、第5引数の color はデフォルトで 1、前回作った grp_color() の引数仕様をそのまま引き継ぎます。
第6引数の show は、VRAMから画面への転送を行うかどうかを指定し、デフォルトで 1、つまり転送します。

但しここで注意しなければならないのは、dx = 0 の時、つまり垂直な線分を描画しようとするとき、0で除算できないのでエラーになります。これについては、とりあえずスクリプトを書いた後に検討します。

abs(dx) > abs(dy)、つまり線分の傾きが45度未満の場合

def line(x1, y1, x2, y2, color=1, show=1):
 rgb = grp_color(color)
 dx = x2 - x1
 dy = y2 - y1
 if abs(dx)>abs(dy):
  if dx: #when dx is not 0
   slope = dy/dx
   for x in range(dx):
    set_pixel(x1+x, y1+x*slope, rgb)
 else:
   (slope≧1 の時の処理)

 if show:
  show_screen()


最初は、line() の第5引数 color から、タプル型の rgb 値を得るために、grp_color() を使っています。というのも set_pixel() の第3引数にタプル型の rgb 値を使って色指定するのが仕様になっているからです。
 
さて、点の描画を反復させるために for 文を使い、range(dx) により、x を 0 から dx-1 まで1づつ増やしながら set_pixel() で点を描画させます。

この反復処理では、最初の x = 0 のときに (x1, y1) に点を描画します。次の描画、つまり隣の点がどうなるかと言えば、x が 1 増え、点描画の x 座標は x1 + x   なので実際は x1 + 1 となります。x 座標の増分が x なので、そのときの y座標の増分は、x に傾き slope を掛けた値になります。  つまり y1 + x*slope です。

但し、ピクセル単位で点を描画するので、set_pixel() の引数は整数でなければなりません。 しかし slope は整数とは限りませんので、x*slope は小数点を切り捨てて整数にしてやる必要があります。

追加した部分は赤文字にしています。
def line(x1, y1, x2, y2, color=1, show=1):
 rgb = grp_color()
 dx = x2 - x1
 dy = y2 - y1
 if abs(dx)>abs(dy):
  if dx: #when dx is not 0
   slope = dy/dx
   for x in range(dx):
    set_pixel(x1+x, y1+int(x*slope))
 else:
  (slope≧1 の時の処理)

0 での除算エラーは dx=0 で発生します。このとき、abs(dx)=0 なので、if の条件 abs(dx)>abs(dy)0 > abs(dy) となります。絶対値は常に 0 以上なので、この条件は False (偽) となります。つまり除算エラーが 発生する dx=0 の時は、この if 文以下は実行されず、else 節にジャンプします。除算エラーについては、else 節を書く時に改めて検討します。


実は、range() の使い方で、考えなければならない問題があります。

for x in range(dx): の動作を細かくみてゆきます。
仮に dx = 10 だとすると、range(dx) はリスト [0, 1, 2, 3, 4, 5, 6, 7, 8, 9] を作り、左から順に x に適用しながら反復処理を行います。実は、dx が正か負かで動作が変わってしまう問題があります。

dx が正の時は、問題ありません。 しかし dx = x2 - x1 なので、dx は常に正である保証はありません。

dx が負の場合、例えば -10 のときは、点描画の最初の座標 (x1, y1) は (x2, y2) よりもにあり、点描画は右から左へ進んでゆきます。つまり、range(dx) は、リスト [0, -1, -2, -3, -4, -5, -6,  -7, -8, -9] を作り、左から順に x に適用しながら反復処理を行う必要があります。

従って、range(dx) は、dx が正の場合は range(0, dx, 1) で、dx が負の場合は range(0, dx, -1) でないといけないことが分かります。

ちなみに、range(stop) は、range(start, stop, step) と書き換えることができ、start のデフォルトは 0 、step のデフォルトは 1 です。
range() のリファレンス参照

そこで、dx が正でも負での同じ動作にするためには、range(dx) と記述する代わりに range(0, dx, k) とし、dxが正の場合は k=1dxが負の場合は k=-1 となるように k を決めれば、問題を解消できます。

k の決め方は色々な方法があると思いますが、ここでは dxdxの絶対値で割り算すると 1 か -1 かになるので、この方法で k を決定します。但し range() の仕様上 k は整数でなければならないので、int() 関数で小数点以下を切り捨て整数にします。

k = int(dx/abs(dx))

abs() 関数は引数の絶対値を返す関数なので、これを使いました。

追加した部分は赤文字にしています。
def line(x1, y1, x2, y2, color=1, show=1):
 rgb = grp_color()
 dx = x2 - x1
 dy = y2 - y1
 if abs(dx)>abs(dy):
  if dx: #when dx is not 0
   k = int(dx/abs(dx))
   slope = dy/dx
   for x in range(0, dx, k):
    set_pixel(x1+x, y1+int(x*slope))
 else:
  (slope≧1 の時の処理)

 if show:
  show_screen()


これで、abs(dx)>abs(dy) の時の処理は完成です。

from casioplot import*
def
line(x1, y1, x2, y2, color=1, show=1):
 rgb = grp_color()
 dx = x2 - x1
 dy = y2 - y1
 if abs(x)>abs(y):
  if dx: #when dx is not 0
   k = int(dx/abs(dx))
   slope = dy/dx
   for x in range(0, dx, k):
    set_pixel(x1+x, y1+int(x*slope))
 else:
  #abs(dx)>abs(dy) でない時の処理

 if show: 
#data transfer to screen
  show_screen()


abs(dx)>abs(dy) でない時、つまり線分の傾きが 45度以上の場合
  
 [2020/07/12 追記]

直交座標系では、2つの座標を入れ替えても数学的に等価であることは保証されています。そして、x 座標と y 座標を入れ替えると、傾きは 45度未満になって、適切に線分を描けます。

そこで、x と y を入れ替え、dx と dy を入れ替えて作ったスクリプトを書きます。但し、set_pixel() の引数については、上で検討したように、単なる文字の入れ替えではなく、数学的な意味で等価になるように変更します。

from casioplot import*
def
line(x1, y1, x2, y2, color=1, show=1):
 rgb = grp_color()
 dx = x2 - x1
 dy = y2 - y1
 if abs(dx)>abs(dy):
  if dx: #when dx is not 0
   k = int(dx/abs(dx))
   slope = dy/dx
   for x in range(0, dx, k):
    set_pixel(x1+x, y1+int(x*slope))
 else:
  if dy: #when dy is not 0
   k = int(dy/abs(dy))
   slope = dx/dy
   for y in range(0, dy, k):
    set_pixel(x1+int(y*slope), y1+y)

 if show: 
#data transfer to screen
  show_screen()

今追加した else 節について、0 除算のエラーについて検討します。上で検討したように、dx=0 の場合は else 節の処理にジャンプします。dx=0 の時は slope=0 となるので、問題ないことが分かります。

では、dy=0 になる場合を検討します。最初の if の条件 abs(dx)>abs(dy) は、abs(dx)>0 となるので、これは常に True (真) です。つまり dy=0 の場合は、最初の if 文で処理され、slope=0 となり問題ありません。

除算エラーになる最後のケースは、dx=0 かつ dy=0 の場合です。この条件が成り立つ時は、x1 = x2 かつ y1 = y2 です。このとき描画するのは線分でなくて、点になることが分かります。この条件が成り立つ時は、if 文の条件には合わないので else 節へジャンプし、そこで slope の計算で除算エラーが発生します。

そこで、if 文の上に以下の処理を追加して、dxdy がともに 0 の時、座標 (x1, y1) に点を描画した後、return で終了し、それ以下を実行しないようにします。

from casioplot import*
def
line(x1, y1, x2, y2, color=1, show=1):
 rgb = grp_color()
 dx = x2 - x1
 dy = y2 - y1
 if dx==0 and dy==0: # avoid division by 0 error
  set_pixel(x1, y1, rgb)
  return
 if abs(dx)>abs(dy):
  if dx: #when dx is not 0
   k = int(dx/abs(dx))
   slope = dy/dx
   for x in range(0, dx, k):
    set_pixel(x1+x, y1+int(x*slope))
 else:
  if dy: #when dy is not 0
   k = int(dy/abs(dy))
   slope = dx/dy
   for y in range(0, dy, k):
    set_pixel(x1+int(y*slope), y1+y)

 if show: 
#data transfer to screen
  show_screen()


6.2 グラフィックスユーザーモジュールの作成

前回作成した grp_color()circle() そして今回作った line() は、グラフィックス画面に出力する汎用関数なので、いちいち関数定義をコピーして使うのは煩わしいです。そこで、calioplotmath, random モジュールのように呼び出して簡単に使うようにします。

そのためには、自分で作ったユーザー関数の定義を1つのスクリプトファイルに収め、そのスクリプトをモジュールとして呼びだして使います。今回は、このスクリプトファイル名を u.py とします。このモジュール(ユーザーモジュール)を使うには、スクリプトの冒頭に

from u import *

と書くだけです。

u.py (ユーザーモジュール ver 1.5) のダウンロード [2020/11/08 修正]

このモジュールには、grp_color()circle()line() が含まれます。これらは、モノクロ液晶モデル (FXモデル) にも対応しています。但し、line() 関数は FXモデル対応のために変更する必要がなく、そのままです。

これまで作った関数のスクリプトは、電卓内に保存されている筈です。電卓をPCとリンクし、エディタを使ってPC上でカット&ペーストするのが、間違いなく、簡単に u.py を作れます。
エディタについては、Casio Python - はじめに:電卓で作る初めてのスクリプト の 1.3 スクリプと作成と編集の2つの方法 を参照してください。


6.3 line() 関数の動作確認 - ckLine.py [fx-CG50 OS3.4以降専用]

line() 関数を実際に使ってみます。

 ckLine.py のダウンロード  - このスクリプトはCGモデル専用です。

from u import *

line(0,0,383,0,'black')
line(383,0,383,191,'blue')
line(383,191,0,191,'red')
line(0,191,0,0,'magenta')
line(10,10,373,10,'green')
line(373,10,373,181,'cyan')
line(373,181,10,181,'yellow')
line(10,181,10,10,1)
line(0,0,383,191,2)
line(0,191,383,0,3)

実行結果は以下のようになります;

lines_output 

FXモデル (fx-9750GIII, fx-9860GIII OS3.40以降) では、各座標値を3で割った整数部に置き換える必要があります。
例えば、383383//3 に変更します。[2020/11/08 追記]


6.4 line() 関数の動作確認 - ckLine2.py [fx-CG50 OS3.40 以降専用]

line() 関数を使った別のスクリプトです。
左上の座標 (0, 0) を基点にして、扇状に線分を -90度から0度まで反復描画してみます。反復する際にカラーコードを 1 ~ 7 まで繰り返し変更してみると、なんとも不思議な模様が得られました。

ckLine2.py のダウンロード - このスクリプトはCGモデル専用です

from u import *

for x in range(383):
 line(0, 0, x, 191, x%7)
for y in range(191, 0, -1):
 line
(0, 0, 383, y, y%7)



Note: 除算の余りを求める - %
整数 x7 で除算した時の余りは x%7 で得られます。この計算により 0 ~ 7 の整数を算出できます。


ckLine2_3 

ckLine2_5 

ckLine2_6 


FXモデル (fx-9750GIII, fx-9860GIII OS3.40以降) では、各座標値を3で割った整数部に置き換える必要があります。
例えば 383 を 383//3 に変更します。[2020/11/08 追記]



目 次

前の記事 - 5. 関数の作成と活用

次の記事 - 7. テキスト出力関数の追加





応援クリックをお願いします。励みになるので...
にほんブログ村 IT技術ブログ 開発言語へ


 


keywords: fx-CG50Pythonプログラム関数電卓

リンク集 | ブログ内マップ


テーマ : プログラム関数電卓
ジャンル : コンピュータ

最新記事
検索フォーム
最新コメント
カテゴリ
C# (3)
Visitors
Online Counter
現在の閲覧者数:
プロフィール

やす (Krtyski)

Author:やす (Krtyski)
since Oct 30, 2013


プログラム電卓は、プログラムを作って、使ってナンボ!

プログラム電卓を実際に使って気づいたこと、自作プログラム、電卓での Casio Basic や Casio Python プログラミングについて書いています。

なお管理人はカシオ計算機の関係者ではなく、Casio Basicが面白いと感じる1ユーザーです。


写真: 「4駆で泥んこ遊び@オックスフォード郊外」

リンク
月別アーカイブ
Sitemap

全ての記事を表示する

ブロとも申請フォーム

この人とブロともになる

QRコード
QR