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