Pythonの平面トラス解析プログラムの入力データ作成方法について説明します。慣れないうちはトラスの形状や荷重を配列で表現するのが難しいかもしれませんが、図を見ながら理解してください。各変数の内容を覚えなければプログラム全体が見えてきません。
なお、この記事は以下のリンク先で紹介した「平面トラス解析プログラム」における入力データ作成部分を解説したものです。
０． FEMなど数値解析シリーズ
この記事は「FEMなど数値解析シリーズ」の記事です。一連の記事は、以下のリンク集を参照してください。
１．入力データ（サンプル1）
それでは、具体的な解析モデルから入力データを作成します。
(1) 解析モデル
下図は解析モデルの節点番号、節点座標、部材番号、部材断面積、部材のヤング係数を示します。
なお、節点0はピン支持、節点2はローラー支持です。これらの節点の拘束条件はトラスの節点変位を固定するものであり、境界条件ともいいます。
荷重条件は下図の通りです。節点1に鉛直方向下向きに大きさ3の荷重を加えます。
(2) 入力データの作成
コード01は、平面トラス解析プログラムの入力データを作成するdef関数です。平面トラス解析プログラムのコード全体は以下のリンクで紹介されています。
なお、このリンク先の「平面トラス解析プログラムのサンプルコード全体」をこの記事では「トラスサンプルコード」と呼ぶことにします。
#コード01
def input_data_1():
'''入力データ例2（3節点、3部材トラスの例） ※適用範囲は32,767節点以下、65,535部材以下。
input_data: [INPUT][DATA]
(1) 入力データ
siz_nod(=2):節点自由度。x、y方向変位の2種類
var_max(=65535):部材数などで指定した変数の型が取り得る最大値
uint16→0～65535、unit36→0～4294967295、unit64→0～18446744073709551615
総節点数(nod_tn)は（var_max-1)/siz_nod以下、総部材数(nod_tn)は(var_max)以下。
nod_tn:総節点数
mem_tn:総部材数
bc_tn:境界条件を与える節点の数
nl_tn:外力が作用する節点の数
mem_ym :ヤング係数
nod_coo[i, j]:i節点の座標(j=0:x方向, j=1:y方向)
mem_are[i] :i部材の断面積
mem_nn[i, j]:i部材を構成する節点番号(j=0, 1)
bc_nn[i]:境界条件を与えられている節点番号
bc_ff[i, j]:fixの方向:0, freeの方向:1(j=0:x方向, J=1:y方向)
nl_nn[i]:外力が作用する節点番号
nl_mag[i, j]:外力の大きさ(j=0:x方向, J=1:y方向)
'''
#関数からの出力
global siz_nod, var_max
global nod_tn, mem_tn, bc_tn, nl_tn
global mem_ym
global nod_coo, mem_are, mem_nn, bc_nn, bc_ff, nl_nn, nl_mag
siz_nod = 2
var_max = 65535
nod_tn = 3
mem_tn = 3
bc_tn = 2
nl_tn = 1
mem_ym = 1.
#適用範囲の確認
if (nod_tn * siz_nod) > (var_max - 1):
print('error:節点数が適用範囲を超えている。nod_tn * siz_nod > var_max -1 ')
if mem_tn > var_max:
print('error:部材数が適用範囲を超えている。mem_tn > var_max')
#'uint16' 0~65535
#numpy配列の宣言
nod_coo = np.zeros((nod_tn, siz_nod), dtype = 'float64')
mem_are = np.zeros(mem_tn, dtype = 'float64')
mem_nn = np.zeros((mem_tn, 2), dtype = 'uint16')
bc_nn = np.zeros(bc_tn, dtype = 'uint16')
bc_ff = np.zeros((bc_tn, siz_nod), dtype = 'bool')
nl_nn = np.zeros(nl_tn, dtype = 'uint16')
nl_mag = np.zeros((nl_tn, siz_nod), dtype = 'float64')
nod_coo[:] = [[-100, 0],
[0, 100],
[100, 0]]
mem_are[:] = 100.
mem_nn[:] = [[0,1],
[1,2],
[0,2]]
bc_nn[:] = [2,0]
bc_ff[:] = [[1, 0],
[0, 0]]
nl_nn[0] = 1
nl_mag[0, :] = [0., -3.]データー数が多いためreturnで変数を返すのではなく、27～30行目でグローバル変数として扱っています。
また、40～44行目では、想定以上の節点数や部材数のモデルを入力したときにエラーを出力します。
48～54行目ではNumPy配列を定義しています。Pythonのfloat型は64ビットの浮動小数点数表現（倍精度）しか指定できませんが、NumPyではfloat16, float32, float64, float128の4種類を指定できます。精度が高いほど計算精度はあがりますが、必要とするメモリ容量が多くなり、計算速度は低下します。必要に応じて調整します。
なお、Python3のint型は上限、下限のないint型しか指定できませんが、NumPyではint8, int16(short), int32(long), int64の符号付整数と、uint8, uint16, uint32, uint64の符号なし整数を指定できます。いずれも上限、加減が決まっており、表現できる数の範囲が大きいほど必要とするメモリ容量が多くなります。
mem_nn, bc_nn, nl_nnは、負の値とならないためuint16としました。uint16は0～65535の範囲の符号なしの整数です。int32(long)ではなくuint16を指定することによりメモリを節約します。
bc_ffは、境界条件を表現する配列であり、固定（fix）の場合は0(False)、自由（free）の場合は1(True)の2種類しかないため、メモリを節約するためにint32(long)ではなくbool型を適用しています。
解析モデルが大きくなると大きなメモリが必要となりますので、解析精度とのバランスを考え、できるだけメモリを節約した方がよいと思います。
32行目のsiz_nodは節点自由度の総数であり、平面トラスではsiz_nod = 2、立体トラスではsiz_nod = 3です。自由度を示す値をsiz_nodと明示することで立体トラスへの改良を容易にします。
33行目のvar_maxは、「トラスサンプルコード」のdef bc_index()において節点数が制限値である(655535 - 1) / 2を超えていないか判定したり、ind_fix[]やind_free[]の算定で使用します。
２．節点座標などを配列で表現する方法
それでは、節点座標、部材断面、荷重条件、境界条件を、配列で表現する方法を説明します。
(1) 総節点数、節点番号、節点座標
・nod_tn = 総節点数
・nod_coo[節点番号, (x座標:0, y座標:1)] = 節点座標
つまり、「nod_coo[0, 0] = -100.」とは、節点番号0の節点のx座標が-100.であることを示します。
(2) 総部材数、部材番号
・mem_tn = 総部材数
・mem_nn[部材番号, (部材0端:0, 部材1端:1)] = 節点番号
つまり、「mem_nn[0, 1] = 1」とは、部材番号0の1端の節点番号が1であることを示します。
なお、「部材0端」「部材1端」とは部材の両側に0、1の番号をふり、部材のどちら側の端部かを番号により区別します。
(3) 部材断面積、ヤング係数
・mem_ym = ヤング係数（共通）
・mem_are[部材番号] = 断面積
つまり、mem_are[1] = 100.とは部材番号1の部材の断面積が100.であることを示します。
なお、「トラスサンプルコード」では、トラス部材は同様の材質である場合を想定してヤング係数は各部材共通の値を用いています。
(4) 境界条件
・bc_tn = 境界条件を与える節点の数
・bc_nn[境界条件番号] = 節点番号
・bc_ff[境界条件番号, (x方向:0, y方向:1)] = fix:0, free:1
つまり、bc_nn[0] = 2とは、境界条件が与えられる0番目の節点のは節点番号が2であることを示します。
また、bc_ff[0, 0] = 1とは、境界条件が与えられる0番目の節点のx方向(0)はfree(1)であることを示します。
(5) 荷重条件
・nl_tn = 荷重が作用する節点の数
・nl_nn[荷重番号] = 節点番号
・nl_mag[荷重番号, (x方向:0, y方向:1)] = 荷重
つまり、nl_nn[0] = 1とは、荷重が作用する0番目の節点のは節点番号が1であることを示します。
また、nl_mag[0, 1] = -3とは、荷重が作用する0番目の節点のy方向(1)荷重は-3であることを示します。
３．入力データ（サンプル２）
それでは練習です。サンプル１を参考に下の解析モデルの入力データを作成してください。
(1) 解析モデル
以下、総節点数5、総部材数7の解析モデルです。
(2) 入力データの作成
上記、解析モデル（サンプル２）の入力データはコード02です。いかがでしょうか。
#コード02
def input_data_2():
'''入力データ例2（5節点、7部材トラスの例） ※適用範囲は32,767節点以下、65,535部材以下。
適用範囲は32,767節点以下、65,535部材以下。
'''
#関数からの出力
global siz_nod, var_max
global nod_tn, mem_tn, bc_tn, nl_tn
global mem_ym
global nod_coo, mem_are, mem_nn, bc_nn, bc_ff, nl_nn, nl_mag
siz_nod = 2
var_max = 65535
nod_tn = 5
mem_tn = 7
bc_tn = 2
nl_tn = 2
mem_ym = 1.
#適用範囲の確認
if (nod_tn * siz_nod) > (var_max - 1):
print('error:節点数が適用範囲を超えている。nod_tn * siz_nod > var_max -1 ')
if mem_tn > var_max:
print('error:部材数が適用範囲を超えている。mem_tn > var_max')
#'uint16' 0~65535
#numpy配列の宣言
nod_coo = np.zeros((nod_tn, siz_nod), dtype = 'float64')
mem_are = np.zeros(mem_tn, dtype = 'float64')
mem_nn = np.zeros((mem_tn, siz_nod), dtype = 'uint16')
bc_nn = np.zeros(bc_tn, dtype = 'uint16')
bc_ff = np.zeros((bc_tn, siz_nod), dtype = 'bool')
nl_nn = np.zeros(nl_tn, dtype = 'uint16')
nl_mag = np.zeros((nl_tn, siz_nod), dtype = 'float64')
nod_coo[:] = [[0., 0.],
[100., 100.],
[200., 0.],
[300., 100.],
[400., 0.]]
mem_are[:] = 100.
mem_nn[:] = [[0, 1],
[1, 2],
[2, 3],
[3, 4],
[0, 2],
[2, 4],
[1, 3]]
bc_nn = [0, 4]
bc_ff[:] = [[0, 0],
[1, 0]]
nl_nn[:] = [1, 3]
nl_mag[:] = [[0., -3.],
[0., -3.]]４．ファイルによる入力
「トラスサンプルコード」は学習用なので、小回りが効くように関数の中で直接入力データを代入しています。しかし、大きな解析モデルではファイルで入出力しなければ大変です。
ファイルで入出力する方法は以下のリンク先を参考にしてください。
私が実際に購入した教材のご紹介
以下、私が実際に購入したPythonの教材をまとめてみました。 Pythonを学習する上で、少しでもお役に立つことができればうれしいです。
・Python♪私が購入したPythonの書籍のレビュー
・UdemyのPythonの動画講座を書籍を買う感覚で購入してみた
その他
液晶ペンタブレットを買いました
(1) モバイルディスプレイを買うつもりだったのに激安ペンタブレット購入
以下、私が光回線を導入した時の記事一覧です。
(1) 2020年「光回線は値段で選ぶ」では後悔する ←宅内工事の状況も説明しています。
(2) NURO光の開通までWiFiルーターを格安レンタルできる
(3) NURO光の屋外工事の状況をご紹介。その日に開通！
(4) 光回線開通！実測するとNURO光はやっぱり速かった
(5) ネット上のNURO光紹介特典は個人情報がもれないの？
