2022年8月30日火曜日

MEMS

先週、MEMSセンサーを購入。

今まで省電力、小型、現場で長く持つ、などと聞いて来ましたが、実物を見たのは初めてでした。

本当に、小さい!

加速度センサーだったのですが、5~6㎜四方。これ、はんだ付けする自信はありません。基盤についた商品を購入してよかった。手に負えないところでした。スマホを分解したことはありませんが、このような小型センサーが所狭しと並んでいるのでしょう。モノづくり技術者にリスペクトです。

制御ボード、インターフェースボードなどを加えて独立した観測機器とするのでしょうね。以前よりラズパイでの制御もネットで見ますが、それではセンサーの小型化といった特徴を損なっています。小型パーツを組み合わせ、工作されるのでしょう。ノウハウを学びたい。

初心者ですので、無理はしません。PCにつなげて初期不良のないことは確認できましたので、料理法をゆっくり考えましょう。

 


2022年8月20日土曜日

位相 unwrapping

位相の unwrapping について調べていました。

SAR ではなく、単なる信号の再現(逆フーリエ)でこの話題が出てきました。

調べてみると、ありますね。

位相アンラッピング(Phase unwrapping)
http://retrofocus28.blogspot.com/2013/12/phase-unwrapping_26.html

def mywrap(y): for i in range(1, len(y)): div = y[i] - y[i-1] if div < -pi: y[i:] += 2*pi elif div > pi: y[i:] -= 2*pi else: pass return y

numpy なら1行です。

np.unwrap() 

簡単ですね。
理想的な条件なら再現できるようですが、実際のデータでは難しいようです。


2022年8月14日日曜日

Cython

Pandas の公式に速度改善に関するチュートリアルがあります。
https://pandas.pydata.org/pandas-docs/stable/user_guide/enhancingperf.html

CPython ではなく、Cython の利用です。
C で書いたライブラリを Python で読むのではなく、Python で書いた関数を C に渡してコンパイルする形。渡す際に型を定義したり、型チェックを外すことで速くすることが可能になっているようです。これは、使い勝手が良さそう。

プロファイルは動かなかったのですが、それ以外は動きました。遅い処理はこのような簡易な処理方法があることも覚えておきましょう。


C++ & Python

C++ で Pythonのモジュールを作りたい。

初めは Python と C++ の exe をバッチファイルで繋ごうかと思いましたが、ライブラリにすればPython から直接読めます。

チュートリアルがあり、GitHub にテンプレートが公開されていました。それをトレースしましょう。https://docs.microsoft.com/ja-jp/visualstudio/python/working-with-c-cpp-python-in-visual-studio?view=vs-2022

まずは環境設定。
VS2019 は Python3.7 までのサポートのようです。これは古い。VS2019 を削除し、2022に入れなおし。これでデフォが3.9?でしょうか。デフォを3.10にしたかったのですが、入れ替え方がわかりません。
空のPythonプロジェクトを作り、新たに conda 環境を Python=3.10 で作成。これを新しいプロジェクトに対する既定の環境に設定。あらためて空の Python プロジェクトを作り、テンプレートのコードをコピぺしました。が、これではダメ。再度開くと3.9に戻ります。
Google先生に伺ったところ、[全般][インタープリター]で変更できました。
Python のみで動くように初期コードに書き換えた後、3.10で動作確認。OK でした。

次に C++。
既存のプロジェクトとして ’superfastcode’ を追加。VS2017 製でしたので[プロジェクト][ソリューションの再ターゲット]にてアップデート。’superfastcode’ というライブラリとして Python から呼び出せるとのこと。ソースファイル module.cpp は既に追加されている状態です。
プロジェクトのプロパティから、いくつかの設定を変更します。
・[全般][プロジェクトの既定値][構成の種類] dllを指定。
・[C/C++][全般][追加のインクルード ディレクトリ]新たに作成したenv下の include フォルダーを追加。
・[C/C++][コード生成][ランタイムライブラリ]マルチスレッドDLL (/MD
・[リンカー][全般][追加のライブラリディレクトリ] ]新たに作成したenv下の libs フォルダーを追加。

これで一度ビルド。OKです。
’superfastcode2’ も同様に取り込んで設定。採用していた conda 環境に pybind11 をインストールしておきます。念のため、Win11 を再起動してからビルド。
うーん、引っかかります。Python3.6が必要と出ます。どこでしょうか。

結局、原因がわからなかったので’superfastcode2’ はアンロード。
リリース環境でも同じ設定でビルド。で、試算。

Running benchmarks with COUNT = 500000
[tanh(x) for x in d] (Python implementation) took 0.512 seconds
[fast_tanh(x) for x in d] (CPython C++ extension) took 0.123 seconds
やはりCの方が速い。
チュートリアルを見ると pybind11 を使うと ’superfastcode’ より遅くなっていますので、ここまでにしておきましょう。

C++ でライブラリを作り、Python で読み込んで計算するチュートリアルでしたが、分かり難いですね。これ、gcc でも設定できるでしょうか。自信ないですね。

gcc と エンコード

C ++の初心者です。

C、C#とともに何度か入門しましたが、実務ではそれほど用がなかったため、覚える前に辞めてしまう状態が続いています。此度、実務(Windows10)で使うことになり、あらためて入門です。

Visual Studio 2019 を使ってコンパイルしていると、scanf でエラー。

error C4996: 'scanf': This function or variable may be unsafe. Consider using scanf_s instead. To disable deprecation, use _CRT_SECURE_NO_WARNINGS.
そうでした。https://phreeqc.blogspot.com/2021/05/c.html
scanf_s を使うと通りましたが、デバッグでエラー。
test.exe (プロセス *****) は、コード -****** で終了しました。
忘れていますね。scanf_sで文字列を入力する場合には第三引数で最大配列数を指定する必要がありました(引数不足のwarning が出ていましたが、無視していました)。 

引数を指定し、できたexe を実行!
動きました。

表示がすぐに消えてしまうので、「Enterkey を押して終了させたい。」をコピペ。
https://teratail.com/questions/140273
これでOK。

次はgcc
tdm64-gcc-10.3.0-2 を Win11 に入れて、ターミナル起動。gcc -v で動作確認。

>gcc -v
gンン spec gpオト「ワキB
COLLECT_GCC=gcc
COLLECT_LTO_WRAPPER=C:/TDM-GCC-64/bin/../libexec/gcc/x86_64-w64-mingw32/10.3.0/lto-wrapper.exe
・・・
T|[gウト「 LTO ウkASY€: zlib zstd
gcc o[W 10.3.0 (tdm64-1)

文字化けします。
コマンドプロンプト表示をUTF-8に変更

>chcp 65001
Active code page: 65001

>gcc -v
Using built-in specs.
COLLECT_GCC=gcc
COLLECT_LTO_WRAPPER=C:/TDM-GCC-64/bin/../libexec/gcc/x86_64-w64-mingw32/10.3.0/lto-wrapper.exe
・・・
algorithms: zlib zstd
gcc version 10.3.0 (tdm64-1)

OKです。では、コンパイル。コードはUTF-8。

>g++ -o test test.c

できたexeを起動すると、文字化けしています。Win の実行ファイル側では shift-JIS 固定なのでしょう。
実行ファイルでのエンコードをCP932として指定します。

>g++ -o test test.c -fexec-charset=CP932

うまく変換できているように見えます。もう少しテストが必要かな。
ひとまず、環境構築は完了です。

2022年8月7日日曜日

PC Building Simulator

PC Building Simulator というゲーム。

購入したもののずっとインストールすらしておらず、昨晩、初めて触りました。
日々、PC のトラブルが持ち込まれ、それに対応していくモードが主でしょうか。ウイルス駆除だの、電源交換だの、メモリ増設だの。これ、やっていて辛くなってきます。実世界のシステム担当者の仕事を代わりにやっているようで、遊んでいる気になれません。ええ、あらためて担当者に感謝です。

フリーで作成できるモードは遊んでいる感じがします。高級な部品を思う存分使用できます。ケースが狭くてリザーバーとラジエーターが干渉したり、電源が足りなかったり。実機で組む前に試す、なんて使い方もできそうです。
3DMarkのベンチマークが動く、エラーで落ちるなんてところまでシミュレートしてくれます。CPU, GPU ともに水冷で組んでみましたが、GPUの温度はどこに出るのでしょう。クロック変更はどこかな?

これなら壊す心配がありません。不慣れな人の練習用に使っても良いかも。

2022年8月6日土曜日

テレワーク その2

第7波が全国的に猛威を振るっています。

第1波から現在まで、勤務先の私を含む周辺では一度もテレワークが許可されていません。executives はイロイロと対応できない?しない?まま、基本、全員出社が平常運転となっています。いえ、体調不良者や濃厚接触者、妊婦のいる家庭等一部の方は execcutives の許可によってテレワークが認められています。ただし在宅手当なし、と通達されています。
一部の支店ではテレワークを導入しているようですので、リクルート時には「在宅勤務導入」を表に出しているようです。ま、嘘ではないですね。

今日、知人から電話がありコロナ対応状況について話していました。
同業他社でも「小学生以下の子供がいる人はテレワーク」「普通にテレワーク」など。7波の特徴に応じた対応を取っている会社は、さすがと思います。就活では本質を見抜き、そういった会社を選んでいただきたいものです。


CIMへの対応? その2

他支店で、地質CIMに対応せざるを得ない状況に陥った方々が現れました。

受注者希望だったそうですが、発注者要望を受けているという少し変わった状況。で、地質の executives 経由ではなく、CIM の executives 経由で相談が届きました。
地質調査業務であれば、詳細度はさて置き3次元可視化する程度しか要望されないと思います。であれば、特に難しい点はありません。が、皆さん3次元化を未だに避けていらっしゃるので、できる人を探している状なのでしょう。
https://phreeqc.blogspot.com/2021/08/cim.html

ちょうど GEORAMA の更新時期になりました。今年から地質部門は新しい executive に変わったのですが、保守更新については判断保留になりました。私が地質を離れてからどなたも使われていませんし、地質担当の方々は未だに2次元CAD 主体ですので、即決できなかったのでしょう。国が CIM を進めている中、何らかの理由でこれまで3次元化を避けてきた、今後も積極的にはかかわらないというのが、新しい executives が引き継いだ経営方針なのかもしれません。

設計部署は数年前から Autodesk 製品を使って人を育てていますし、CIM に対応しています。現実的に地質CIMでは GEORAMA が納品への近道の一つですので、私なら更新は即決ですが、この会社でその判断は誤っているのかもしれません。

さて、相談は来ても、ソフトはない、人は育っていない、誰も手を動かそうとしないでは、弱りました。いえ、executives は予測していた結果でしょうから、それほど困っていないのでしょう。
本当に困っているのは担当者のみ。あれ?昔から変わっていないですね、この悪習。
さて、どうしましょうね。

弾性論とFEM

FrontISTR の ハンズオンに参加しました。

FrontISTR のビルドやインストール部分は時間の都合で割愛。メッシング は商用ソフトからのインポートが主流になりつつあるとのことで、これも割愛。ジョブを投げて計算を回すのに慣れるといった5時間でした。うーん。東大のスパコン?を使わないのでそこのジョブ管理を覚える必要はなかったのですが。ま、スパコンの使い方は似たようなものだと確認はできました。

機能紹介や理論の部分に METIS による領域分割、FEMにおける弾性論など懐かしい話題が含まれていました。FEMでの展開を忘れかけており、帰宅してから一通り復習です。


**************************************
線形変位場を仮定
→要素変位ー節点変位関係・・・補間関数[N]を使用
{ue}=[N]{u}
→物理空間座標系と計算空間座標系(正規化座標系)の関係・・・写像関数
アイソパラメトリック要素:補間関数と写像関数が同一(形状関数)

ひずみ-変位関係(支配方程式1)
→FEMでは要素ひずみー節点変位関係
{ε}=[B]{u}・・・(要素内の)ひずみは(節点)変位の空間勾配
Bマトリックス(ひずみ-変位マトリックス):[B]=∇[N]・・・補間関数の空間微分

応力-ひずみ関係(フックの法則、支配方程式2)
→FEMでは要素応力-節点変位関係
σ=E{ε}・・・1次元
{σ}=[D]{ε}=[D][B]{u}
Dマトリックス(応力-ひずみマトリックス):ばね。ラメ定数でも表現可。

仮想仕事の原理(支配方程式3)
外力の仮想仕事(力x変位)と
δW=1/2{u}T{f}:節点変位ベクトルと節点外力ベクトルの積
仮想ひずみエネルギー(ひずみx応力)の
δU=1/2∫{ε}T{σ}dV=∫[B]T{u}T[D][B]{u}dV
つり合い(要素剛性方程式)
{u}T{f}={u}T∫[B]T[D][B]{u}dV
{f}=[k]{u}
kマトリックス:k=∫[B]T[D][B]dxdydz
→FEMではヤコビアン|J|+正規化領域でガウスの数値積分:k=∫[B]T[D][B]|J|dξdηdζ

全体剛性マトリクス
{F}=[K]{U}

変位の計算
{U}=[K]-1{F}

節点変位→要素ひずみ→要素応力
応力は積分点で計算。
節点変位は要素間で連続(適合)するが、ひずみと応力の連続性は保障されない。
応力・ひずみは変位の勾配に比例する物理量であるため、変位より精度が低くなる。
実用的には積分点の値を節点に外挿し平均した節点平均応力、節点平均ひずみが用いられる。

※アイソパラメトリック1次要素→せん断ロッキング→低減積分で回避→アワーグラスモード→非適合要素(要素間の変位の連続性が保たれないものの、実用的には問題ない程度)

※四辺形1次要素の完全積分:4点(2次積分)、2次要素9点(3次積分)→低減積分:1次要素1点、2次要素4点

Pyrhon + GDAL

Raster を複数使用した安定計算。

解像度の違いでエラーを吐きます。
Arc で ポリゴンからラスターへ変換していたのですが、解像度が完全には一致しないようでした。仕方がないので、リサンプリングすることに。
以前、GDAL を使って Raster を読み込んでいたので、そこに追記。

# Open Raster
def open_src(file,band):
    src = gdal.Open(file, gdal.GA_ReadOnly)
    if src is None:
        print('file not found')
    else:
        print(src.GetDescription() )
        xy=src.GetGeoTransform()#
        print(xy)
        print(src.GetProjection())
        print(src.GetMetadata() )
        arr_b1=src.GetRasterBand(band).ReadAsArray(0)
        print(arr_b1.shape)
        nodata=src.GetRasterBand(band).GetNoDataValue())
        print(nodata)
        return src,xy,arr_b1,nodata


#図化
def plot_raster(arr_b1, nodata):
    arr_b1 = np.where(arr_b1==nodata, np.nan, arr_b1)
    print('max:',np.nanmax(arr_b1))
    print('min:',np.nanmin(arr_b1))
    print(arr_b1.shape)   
    plt.imshow(arr_b1)
    plt.colorbar()
    return arr_b1

    
#測地系変換
#EPSG:4612 から 4326に変換して保存
ds = gdal.Warp(file_out, file_in, srcSRS="EPSG:4612", dstSRS="EPSG:4326")
ds=None


#2つのラスターの領域、解像度を(ほぼ)揃える
#切り出し領域を決める
minx=136
miny=35
maxx=137
maxy=36
#領域 をforium で確認
fmap = folium.Map(location=[(miny+maxy)/2,(minx+maxx)/2], zoom_start=8)
line = [[miny,minx],[miny,maxx],[maxy,maxx],[maxy,minx],[miny,minx]]
folium.PolyLine(locations=line).add_to(fmap)
fmap
#ピクセル数を決める
width=18000
height=11100
#resampling
cut_file0="./output/cut0.tif"
cut_file1="./output/cut1.tif"
ds=gdal.Translate(cut_file0, file_in0, projWin=[minx,maxy,maxx,miny], width=width, height=height, noData=nodata)
ds=gdal.Translate(cut_file1, file_in1, projWin=[minx,maxy,maxx,miny], width=width, height=height, noData=nodata)
ds=None


#保存
def save_dst(rasterfile, arr_new, nodata, src):
    driver = gdal.GetDriverByName("GTiff") 
    dst = driver.Create(rasterfile,
                       xsize=arr_new.shape[1],
                       ysize=arr_new.shape[0],
                       bands=1,
                       eType = gdal.GDT_Float32)
    dst.SetGeoTransform(src.GetGeoTransform())
    dst.SetProjection(src.GetProjection())
    band1=dst_ds.GetRasterBand(1)
    band1.WriteArray(arr_new)
    band1.SetNoDataValue(nodata)
    dst.FlushCache()
 

Arc でラスターを2枚作成・書き出していましたので、マルチバンドにしてから1枚出力すれば解像度が完全に一致したのですが。
https://pro.arcgis.com/ja/pro-app/2.8/tool-reference/data-management/composite-bands.htm

ま、今後使えるかもしれませんので残しておきましょう。

 

2022年8月2日火曜日

float, double

C#で作成されたソフトと、Pythonで書かれたスクリプト。
同じテキストファイルを読んで計算したのですが、答えが異なりました。

原因は変数の型。

C# は float(32bit) で読み込むように作られていました。Python の float は C の double (64bit)で作られているそうですので、そのまま読み込むと前者の桁が小さくなります。https://docs.python.org/ja/3/library/stdtypes.html

後者を np.float32 で変更すれば良いのでしょうが、そもそも、input ファイルの段階で有効桁数以上の数字を消していなかったのがマズイ。不要な桁を消して読み込めば、一致しました。

手を抜くとダメですね。すぐに帰ってきます。