2D数値計算の境界条件を画像から作るPythonコード例

  • このエントリーをはてなブックマークに追加

2Dの流体解析をする時、計算領域内部に置いた構造物の影響を考慮したくなります。ここでは画像の輪郭から流体解析に使用できる境界条件を作る1つの例を、Pythonの流体解析コードと画像処理コードをコラボしながら考察してみたいと思います。

こんにちは。wat(@watlablog)です。本日はネタ回!特に役に立たない画像処理と流体解析コードのコラボをしてみた結果をメモしました

この記事のモチベーションと目標

2D数値計算の例

当ブログでは「Pythonで2次元ラプラス方程式を数値計算する方法」の記事で2次元(2D)空間の数値計算を行いました。

紹介した記事では学習目的であることから計算領域の中に何も条件を設定せず、4辺への境界条件を与えるというシンプルなシミュレーションを行っていました(下図)。

2D数値計算の例

今回はこのシミュレーションをアレンジする内容として、計算領域内部に構造物を置く方法について検討してみます。

計算領域数値計算境界条件…といった用語が曖昧な方は是非先ほど紹介した記事をご覧になってください。

計算領域内に配置する構造物のイメージ

ここでいう構造物とは下図のように、計算領域内に空間とは異なる格子点を配置することを意味しています。

この図では星形の構造物を粗い点(赤点)で表現しているので不格好ですが、十分な解像度を持てば複雑な形状も計算領域内に表現可能です。

構造物のイメージ

矩形や円形ならまだしも、複雑な形状を置きたい場合はプログラミング時に苦労します。例えば上記直交格子の場合は、\(j,k\)番目の反復計算時に何らかの分岐を作る必要があります。

目標:画像から境界条件を設定する格子点を自動抽出する

そこで本題です。この記事では2Dの数値計算において、画像から境界条件を設定する格子点を自動抽出することを目標とします。

下図のように、画像の複雑な輪郭を自動抽出し座標情報に変換してしまえば、数値計算時の条件分岐にそのまま使うことができます。

画像内オブジェクトから輪郭線を抽出するイメージ図

本ページの目標を以下の動画に示します。複雑な境界条件を画像から自動設定して2Dシミュレーションに活用するところとします。複雑な形状…という意味で漢字にしてみました。

目標動画

今回はノリだけで書いています。正式な数値解析のメッシュ分割法にはもっと良いやり方があると思いますが、それらを無視してできそうだったのでやってみたというだけ…。

画像内オブジェクト情報を使って2D数値計算するPythonコード

動作環境

今回は以下の環境でコーディングしました。

Python環境

Python Python 3.9.6
PyCharm (IDE) PyCharm CE 2020.1
Numpy 1.21.1
matplotlib 3.4.3
pillow 7.1.2

PC環境

Mac OS macOS Catalina 10.15.7
CPU 1.4[GHz]
メモリ 8[GB]

サンプル画像

境界条件を設定する「境界」の文字が書かれた(描かれた?)画像ファイルimg-kyoukai.pngをこちらに置きます。Pythonプログラムの直下に置いて是非検証等にお使いください。ちなみにパワポで適当に文字を書いてスクショしただけです。

①画像サイズを格子に合わせる

コピペ用の全コードは最後に置きますが、まずは要所の手順説明から始めます。

先ほどの画像は適当なサイズでスクショしたものなので、計算領域との対応をとるためにまずは格子サイズの情報を使ってリサイズ(ついでにグレースケール化)を行います。

OpenCVを用いたリサイズは以下の記事でも扱いました。

Python/OpenCVで画像をリサイズする方法

②二値化する

次に画像処理ではおなじみ、二値化処理を行います。ここでは文字部分を白く(輝度=255)したかったのでcv2.THRESH_BINARY_INVを使っています。

ちなみに、二値化にも色々種類があります。画像の状態に合わせて以下の記事も参考にしてみてください。

Python/OpenCVで画像の二値化をする方法
Python/OpenCVの適応的閾値処理で綺麗な二値化!

二値化するとこのようになります。

二値化画像

③輪郭抽出をする

輪郭抽出は以下のコードで行います。抽出した輪郭は階層構造hierarchy)になっているので、それを一度平坦化して\(x,y\)座標のリストにしています。画像の座標と一般的なデカルト座標はy軸の方向が逆なので絶対値をかけています。

こちらが抽出結果(緑線が輪郭)です。

二値化画像の輪郭を抽出

輪郭抽出は以下の記事で詳しく説明しています。

Python/OpenCVで画像内オブジェクトの輪郭抽出をする

④境界条件設定に構造の条件を追記

Pythonで2次元ラプラス方程式を数値計算する方法」の記事で作成した境界条件を設定する関数に、今回輪郭情報を抽出した座標の条件を追記します。

ここでは辺の境界条件を0、構造の境界条件を1にしています(構造が一定温度で発熱しているイメージ)。

全コード

以下に全コードを示します。ほとんどが既存記事のコードのままです。

しかし今回のコードはSOR法最適な緩和係数計算が通用しません。そのため\(\omega\)は1.5を手入力しています。計算結果がおかしいと思った時は緩和係数を1にしてガウス・ザイデル法で解いてみることをお勧めします。

反復計算時、img_binaryの輝度値が255の位置における計算はしないようにしています。

計算には関係していないですが、post_img_rendering関数で計算された分布に対して元画像をレンダリングして境界をはっきりと見せています。

上記コードを実行すると、冒頭の動画が作成されます。

境界条件を設定した解析結果

以下のような離れた場所にオブジェクトが存在する場合は、それぞれから色が発生するという結果も得られます。ラプラス方程式はポテンシャル計算に使われるため、全体を滑らかにしようとする作用があります。

異なる画像の解析結果

星の画像img-star.pngはこちら。コードを実行する時に「img-laplace2d」フォルダを空にしてください(動画作成に前の画像が入ってしまうため)。

まとめ

このページでは過去に学んだ2次元ラプラス方程式の数値計算コードに対して、構造物の境界を画像からつくる方法を紹介しました。

正式なメッシュ分割技法ではないと思うのであくまで参考…遊び的な内容ですが、一応思ったとおりの動作はしていると思います。

ただし、反復計算時の境界付近の取り扱いはかなり適当になっているため、計算が不安定になっているかも知れません。それらは今後の学習次第で修正を入れようと思います。

ちょっと思いついてやってみました!
Twitterでも関連情報をつぶやいているので、wat(@watlablog)のフォローお待ちしています!

SNSでもご購読できます。

コメントを残す

*