LLMで情報収集中
Pynastranで結果の取得をしたい。用件を箇条書きで示します。▪︎計算はノーマルモード解析▪︎プロパティに含まれる要素ごとに合計値を出力▪︎出力したい結果はひずみエネルギー▪︎最終的にはモデル全体のひずみエネルギーの合計値とプロパティ単位の合計値の比を求め、プロパティ単位での全体ひずみエネルギーに対する割合(寄与度)を取得したい。 必要情報があれば先に聞いてください。 ネット検索必要なら利用してください。
了解です。pyNastranで SOL103（ノーマルモード） の 要素ひずみエネルギー（ESE） を読み、Property（PID）単位に合計 → 全体合計に対する寄与度（比率） を出す、という流れで組めます。
まず前提：NASTRAN側で「ESE」を出力しておく OP2やPCHに ESE が出ていないと、pyNastran側では取得できません。一般に Case Control に例えば下記のような指定を入れます（GUIだと「Strain Energy」をONにするのと同義）。 ESE(PUNCH)=ALL（PCH出力） （OP2で欲しいなら）解析のI/O設定で OP2 を有効化（環境によって PARAM,POST,-1 等が必要な場合あり）
※ OP2内では「モードの要素ひずみエネルギー」系テーブル（例：ONRGY1/ONRGY2）として扱われます。
pyNastranでやること（設計）
BDFを読む：eid -> pid の対応表を作る（要素からPIDを取る）
OP2を読む：strain_energy（ESE）結果を取り出すpyNastranは op2.op2_results.strain_energy.cquad4_strain_energy など「要素タイプ別」に保持します。
モードごとに 要素ひずみエネルギーを eid で取り出し pid に集約して合計 全体合計で割って寄与度（割合）算出
CSV等に出力
サンプルコード（OP2 + BDFから、PID別寄与度を出す）
なるべく「要素タイプが混在しても動く」ように、strain_energy配下を introspection して総なめします。また、ESEの列は「energy / percent / density」のような構成になり得るので（古いWikiにもその記述があります）、基本は get_headers() が取れればそれに従い、取れない場合は先頭列を energy とみなす、で堅めにしています。
from future import annotations from collections import defaultdictfrom pathlib import Path import numpy as npimport pandas as pd from pyNastran.bdf.bdf import BDFfrom pyNastran.op2.op2 import OP2
def build_eid_to_pid_map(bdf: BDF) -> dict[int, int]: """要素EID -> Property ID(PID) の対応を作る。PIDを持たない要素は除外。""" eid_to_pid: dict[int, int] = {} for eid, elem in bdf.elements.items(): # RBE2等 pid を持たない要素もある pid = getattr(elem, "pid", None) if pid is None: continue try: pid_int = int(pid) except Exception: continue eid_to_pid[int(eid)] = pid_int return eid_to_pid
def iter_strain_energy_result_objects(op2: OP2): """ op2.op2_results.strain_energy 配下のテーブル(dict)を列挙する。 例: cquad4_strain_energy, ctetra_strain_energy ... など。4 """ se = op2.op2_results.strain_energy for name in dir(se): if not name.endswith("_strain_energy"): continue obj = getattr(se, name, None) if isinstance(obj, dict) and obj: # subcase -> result_object のdict yield name, obj
def pick_energy_column(result_obj) -> int: """ result_obj.data のどの列が strain energy かを推定する。 headersが取れるなら、それを優先。取れないなら 0 列目を採用。 """ if hasattr(result_obj, "get_headers"): headers = result_obj.get_headers() # 例: ['strain_energy', 'percent', 'strain_energy_density'] のようなケースを想定 for i, h in enumerate(headers): h_low = str(h).lower() if "strain" in h_low and "energy" in h_low and "percent" not in h_low and "density" not in h_low: return i if h_low in ("energy", "strain_energy"): return i # それでも見つからなければ、percent/density以外の最初の列を採用 for i, h in enumerate(headers): h_low = str(h).lower() if "percent" in h_low or "density" in h_low: continue return i return 0
def compute_property_contribution( bdf_path: str | Path, op2_path: str | Path, subcase: int | None = None, modes: list[int] | None = None, # 0-based index（pyNastranのtime/mode index想定）) -> pd.DataFrame: # --- read BDF (for eid->pid) --- bdf = BDF(debug=False) bdf.read_bdf(str(bdf_path), xref=False) eid_to_pid = build_eid_to_pid_map(bdf) # --- read OP2 (for strain energy) --- op2 = OP2(debug=False) op2.read_op2(str(op2_path)) rows = [] for table_name, subcase_dict in iter_strain_energy_result_objects(op2): for key, result in subcase_dict.items(): # subcase指定がある場合：keyが int とは限らない（タプルキーになる例あり）5 # ここは「intなら一致」「タプルなら先頭のsubcaseを使う」方針で吸収 if subcase is not None: k_subcase = None if isinstance(key, int): k_subcase = key elif isinstance(key, tuple) and len(key) >= 1 and isinstance(key[0], int): k_subcase = key[0] if k_subcase != subcase: continue # result.element : (nelements,) # result.data : (nmodes, nelements, ncols) を想定 eids = np.asarray(getattr(result, "element", []), dtype=int) if eids.size == 0 or not hasattr(result, "data"): continue data = np.asarray(result.data) if data.ndim != 3: continue energy_col = pick_energy_column(result) nmodes = data.shape[0] mode_indices = range(nmodes) if modes is None else [m for m in modes if 0 <= m < nmodes] for mi in mode_indices: energy = data[mi, :, energy_col].astype(float) pid_sum = defaultdict(float) total = 0.0 for eid, e in zip(eids, energy): pid = eid_to_pid.get(int(eid)) if pid is None: continue pid_sum[pid] += float(e) total += float(e) if total == 0.0: continue for pid, se_val in pid_sum.items(): rows.append( { "table": table_name, "subcase_key": str(key), "mode_index": mi, "pid": pid, "strain_energy": se_val, "contribution": se_val / total, "total_strain_energy_mode": total, } ) df = pd.DataFrame(rows) if df.empty: return df # 見やすく：寄与度[%] df["contribution_pct"] = 100.0 * df["contribution"] # 典型的には table が要素タイプごとに分かれてるので、 # 「PID単位に最終集約（要素タイプを跨いで合算）」した集計表も出しておく agg = ( df.groupby(["subcase_key", "mode_index", "pid"], as_index=False) .agg( strain_energy=("strain_energy", "sum"), total_strain_energy_mode=("total_strain_energy_mode", "max"), ) ) agg["contribution"] = agg["strain_energy"] / agg["total_strain_energy_mode"] agg["contribution_pct"] = 100.0 * agg["contribution"] agg = agg.sort_values(["subcase_key", "mode_index", "contribution_pct"], ascending=[True, True, False]) return agg
if name == "main": # 入力 bdf_path = "model.bdf" op2_path = "model.op2" # 例：subcase=1 の全モードで実行（modesは 0-based index） df = compute_property_contribution(bdf_path, op2_path, subcase=1, modes=None) print(df.head(30)) df.to_csv("strain_energy_contribution_by_pid.csv", index=False)
ここまでで満たせる要件 ✅ ノーマルモード解析（SOL103）のESE結果を取得（OP2に入っている前提） ✅ PID（プロパティ）単位で合計値 ✅ モデル全体の合計に対する PIDごとの比（寄与度） ✅ CSV出力まで
先に確認したい「必要情報」（最小限） あなたの環境に合わせてコードの“ツメ”を調整したいので、次だけ教えてください（1行でOKです）：
OP2は出ていますか？（ESEを入れて再計算済みか）
使用ソルバは MSC / NX / OptiStruct のどれですか？（OP2の癖対策用。ESEでキーがタプル化する例もあります）
寄与度は 「各モードごと」 が欲しいですか？それとも 「指定モード範囲で総和」（例：1〜6次モード合算）ですか？
この3点が分かれば、あなたの出したい最終形式（例：モード別シート、合算1枚、PID名も併記…など）に合わせて完成形に寄せます。
OP2にひずみエネルギーは出力済み。ソルバはNX。寄与度は各モードごとで欲しいです。
了解です（NX／OP2にESE出力済み／寄与度は各モードごと）。 NX Nastran の ESE（要素ひずみエネルギー） を pyNastran で読むときの “ハマりどころ” は主にこれです： strain_energy の辞書キーが int（=subcase番号）ではなく、tuple になることがある（例：(1, 2, 1, 0, 0, '', '') みたいなキー）→ 「subcase=1だけ欲しい」のに単純に dict[subcase] ができないことがあるので、キーから subcase を抽出してフィルタします。
以下に、NXのtupleキーも吸収して、(モード×PID) の寄与度をCSVにするコードを完成形として置きます。
from __future__ import annotations
from collections import defaultdict
from pathlib import Path
import numpy as np
import pandas as pd
from pyNastran.bdf.bdf import BDF
from pyNastran.op2.op2 import OP2
def build_eid_to_pid_map(bdf: BDF) -> dict[int, int]:
"""要素EID -> PID。pidを持たない要素(RBE等)は除外。"""
out: dict[int, int] = {}
for eid, elem in bdf.elements.items():
pid = getattr(elem, "pid", None)
if pid is None:
continue
try:
out[int(eid)] = int(pid)
except Exception:
continue
return out
def get_subcase_from_key(key) -> int | None:
"""
NXでありがちな strain_energy のキー（tuple）から subcase を取り出す。
例: (1, 2, 1, 0, 0, '', '') の先頭が isubcase。1
"""
if isinstance(key, int):
return key
if isinstance(key, tuple) and len(key) >= 1 and isinstance(key[0], int):
return key[0]
return None
def iter_strain_energy_tables(op2: OP2):
"""
op2.op2_results.strain_energy 配下の *_strain_energy dict を列挙。
要素タイプ別（cquad4/ctetra等）に分かれている。2
"""
se = op2.op2_results.strain_energy
for name in dir(se):
if not name.endswith("_strain_energy"):
continue
dct = getattr(se, name, None)
if isinstance(dct, dict) and dct:
yield name, dct
def pick_energy_column(result_obj) -> int:
"""dataのどの列が strain energy か推定（headersがあれば優先）。"""
if hasattr(result_obj, "get_headers"):
headers = [str(h).lower() for h in result_obj.get_headers()]
# まず "strain" & "energy" を含み、percent/density ではない列を探す
for i, h in enumerate(headers):
if ("strain" in h and "energy" in h) and ("percent" not in h) and ("density" not in h):
return i
# 次点：energy/strain_energy
for i, h in enumerate(headers):
if h in ("energy", "strain_energy"):
return i
# それでもダメなら percent/density以外の最初
for i, h in enumerate(headers):
if ("percent" in h) or ("density" in h):
continue
return i
return 0
def get_mode_numbers(result_obj, nmodes: int) -> list[int]:
"""
pyNastran側の「モード番号（1始まり）」を可能なら取得し、
取れなければ 1..nmodes を採用。
"""
# 多くの結果オブジェクトは _times に固有値/周波数/モード番号が入ることがある
modes = getattr(result_obj, "modes", None)
if modes is not None and len(modes) == nmodes:
try:
return [int(x) for x in modes]
except Exception:
pass
times = getattr(result_obj, "_times", None)
if times is not None and len(times) == nmodes:
# SOL103だと "mode number" が入る場合と固有値が入る場合があるので、
# ここは「整数っぽければそれを採用」くらいに留める
try:
as_int = [int(x) for x in times]
# 1..nmodesに近いならモード番号扱い
if set(as_int) == set(range(1, nmodes + 1)):
return as_int
except Exception:
pass
return list(range(1, nmodes + 1))
def compute_modal_pid_contribution(
bdf_path: str | Path,
op2_path: str | Path,
subcase: int = 1,
) -> pd.DataFrame:
# --- read BDF ---
bdf = BDF(debug=False)
bdf.read_bdf(str(bdf_path), xref=False)
eid_to_pid = build_eid_to_pid_map(bdf)
# --- read OP2 ---
op2 = OP2(debug=False)
op2.read_op2(str(op2_path))
# --- accumulate (mode, pid) across all element types ---
# NXでは element-type別テーブルに分かれるので、最後に跨って合算する
accum = defaultdict(float) # (mode_no, pid) -> sum_energy
totals = defaultdict(float) # mode_no -> total_energy
for table_name, dct in iter_strain_energy_tables(op2):
for key, res in dct.items():
if get_subcase_from_key(key) != subcase:
continue
eids = np.asarray(getattr(res, "element", []), dtype=int)
if eids.size == 0 or not hasattr(res, "data"):
continue
data = np.asarray(res.data)
if data.ndim != 3:
continue
nmodes = data.shape[0]
energy_col = pick_energy_column(res)
mode_nos = get_mode_numbers(res, nmodes)
# data: (nmodes, nelems, ncols)
for mi, mode_no in enumerate(mode_nos):
energy = data[mi, :, energy_col].astype(float)
# PID集約
mode_total = 0.0
for eid, e in zip(eids, energy):
pid = eid_to_pid.get(int(eid))
if pid is None:
continue
accum[(mode_no, pid)] += float(e)
mode_total += float(e)
totals[mode_no] += mode_total
# --- build dataframe ---
rows = []
for (mode_no, pid), se in accum.items():
total = totals.get(mode_no, 0.0)
if total == 0.0:
continue
rows.append({
"subcase": subcase,
"mode_no": mode_no, # 1始まりのモード番号（推定）
"pid": pid,
"strain_energy": se,
"total_strain_energy": total,
"contribution": se / total,
"contribution_pct": 100.0 * (se / total),
})
df = pd.DataFrame(rows)
if df.empty:
return df
df = df.sort_values(["mode_no", "contribution_pct"], ascending=[True, False]).reset_index(drop=True)
return df
if __name__ == "__main__":
bdf_path = "model.bdf"
op2_path = "model.op2"
subcase = 1
df = compute_modal_pid_contribution(bdf_path, op2_path, subcase=subcase)
print(df.head(30))
df.to_csv("modal_strain_energy_contribution_by_pid.csv", index=False)
NX向け：SOL103 ESE → PID合計 → 寄与度（各モード） from future import annotations from collections import defaultdictfrom pathlib import Pathimport numpy as npimport pandas as pd from pyNastran.bdf.bdf import BDFfrom pyNastran.op2.op2 import OP2
def build_eid_to_pid_map(bdf: BDF) -> dict[int, int]: """要素EID -> PID。pidを持たない要素(RBE等)は除外。""" out: dict[int, int] = {} for eid, elem in bdf.elements.items(): pid = getattr(elem, "pid", None) if pid is None: continue try: out[int(eid)] = int(pid) except Exception: continue return out
def get_subcase_from_key(key) -> int | None: """ NXでありがちな strain_energy のキー（tuple）から subcase を取り出す。 例: (1, 2, 1, 0, 0, '', '') の先頭が isubcase。1 """ if isinstance(key, int): return key if isinstance(key, tuple) and len(key) >= 1 and isinstance(key[0], int): return key[0] return None
def iter_strain_energy_tables(op2: OP2): """ op2.op2_results.strain_energy 配下の *_strain_energy dict を列挙。 要素タイプ別（cquad4/ctetra等）に分かれている。2 """ se = op2.op2_results.strain_energy for name in dir(se): if not name.endswith("_strain_energy"): continue dct = getattr(se, name, None) if isinstance(dct, dict) and dct: yield name, dct
def pick_energy_column(result_obj) -> int: """dataのどの列が strain energy か推定（headersがあれば優先）。""" if hasattr(result_obj, "get_headers"): headers = [str(h).lower() for h in result_obj.get_headers()] # まず "strain" & "energy" を含み、percent/density ではない列を探す for i, h in enumerate(headers): if ("strain" in h and "energy" in h) and ("percent" not in h) and ("density" not in h): return i # 次点：energy/strain_energy for i, h in enumerate(headers): if h in ("energy", "strain_energy"): return i # それでもダメなら percent/density以外の最初 for i, h in enumerate(headers): if ("percent" in h) or ("density" in h): continue return i return 0
def get_mode_numbers(result_obj, nmodes: int) -> list[int]: """ pyNastran側の「モード番号（1始まり）」を可能なら取得し、 取れなければ 1..nmodes を採用。 """ # 多くの結果オブジェクトは _times に固有値/周波数/モード番号が入ることがある modes = getattr(result_obj, "modes", None) if modes is not None and len(modes) == nmodes: try: return [int(x) for x in modes] except Exception: pass times = getattr(result_obj, "_times", None) if times is not None and len(times) == nmodes: # SOL103だと "mode number" が入る場合と固有値が入る場合があるので、 # ここは「整数っぽければそれを採用」くらいに留める try: as_int = [int(x) for x in times] # 1..nmodesに近いならモード番号扱い if set(as_int) == set(range(1, nmodes + 1)): return as_int except Exception: pass return list(range(1, nmodes + 1))
def compute_modal_pid_contribution( bdf_path: str | Path, op2_path: str | Path, subcase: int = 1,) -> pd.DataFrame: # --- read BDF --- bdf = BDF(debug=False) bdf.read_bdf(str(bdf_path), xref=False) eid_to_pid = build_eid_to_pid_map(bdf) # --- read OP2 --- op2 = OP2(debug=False) op2.read_op2(str(op2_path)) # --- accumulate (mode, pid) across all element types --- # NXでは element-type別テーブルに分かれるので、最後に跨って合算する accum = defaultdict(float) # (mode_no, pid) -> sum_energy totals = defaultdict(float) # mode_no -> total_energy for table_name, dct in iter_strain_energy_tables(op2): for key, res in dct.items(): if get_subcase_from_key(key) != subcase: continue eids = np.asarray(getattr(res, "element", []), dtype=int) if eids.size == 0 or not hasattr(res, "data"): continue data = np.asarray(res.data) if data.ndim != 3: continue nmodes = data.shape[0] energy_col = pick_energy_column(res) mode_nos = get_mode_numbers(res, nmodes) # data: (nmodes, nelems, ncols) for mi, mode_no in enumerate(mode_nos): energy = data[mi, :, energy_col].astype(float) # PID集約 mode_total = 0.0 for eid, e in zip(eids, energy): pid = eid_to_pid.get(int(eid)) if pid is None: continue accum[(mode_no, pid)] += float(e) mode_total += float(e) totals[mode_no] += mode_total # --- build dataframe --- rows = [] for (mode_no, pid), se in accum.items(): total = totals.get(mode_no, 0.0) if total == 0.0: continue rows.append({ "subcase": subcase, "mode_no": mode_no, # 1始まりのモード番号（推定） "pid": pid, "strain_energy": se, "total_strain_energy": total, "contribution": se / total, 出力されるもの CSVに、モードごとに total_strain_energy（モデル全体合計） pidごとの strain_energy 合計 contribution（0〜1） contribution_pct（%）
が入ります。
NXでよくある注意（重要） tupleキー問題：上の get_subcase_from_key() で吸収しています（NXのESEでよく起きる挙動として報告あり）。 「モード番号」表示：pyNastranの結果オブジェクトによっては modes や _times に “モード番号そのもの” が入らず、固有値など別の値が入る場合があります。その場合このコードは 安全側で 1..nmodes を mode_no として採番します（寄与度計算自体は問題なし）。
次に、精度を上げるために1点だけ確認（質問というより“チェック項目”） もし「mode_no が 1,2,3… と綺麗に出てない」場合は、あなたのOP2内で モード番号/固有値がどの配列に入っているかを合わせ込みます。その時は、あなたの環境で次を1回だけ実行して出力を貼ってください（10行程度でOK）： from pyNastran.op2.op2 import OP2op2 = OP2()op2.read_op2("model.op2") se = op2.op2_results.strain_energy# 何か1個（例: cquad4_strain_energy）を選んで最初のresultを覗くdct = se.cquad4_strain_energykey0 = list(dct.keys())[0]res = dct[key0]print("key0 =", key0)print("get_headers =", res.get_headers() if hasattr(res,"get_headers") else None)print("has modes =", hasattr(res,"modes"))print("modes =", getattr(res,"modes", None))print("times =", getattr(res,"_times", None))print("data shape =", res.data.shape) （この情報があると、mode_no を確実に OP2 のモード番号に一致させられます。） ここまでで、要件どおり NX SOL103の“各モードごとのPID寄与度” を取れるはずです。