穴あき薄板の応力集中

2.5.1 Step を定義する

次の Create Step（Stepを作成） アイコン をクリックして Create Step（Stepを作成） ダイアログを開いて新しい Step を作成します。名前を Apply-Roark-Conditionsとし、 General procedure category が選ばれていることを確認して、 Static, General を手続き一覧から選択します。次図を参照してください。

次に Continue ボタンをクリックすると Edit Step（Stepを編集） ダイアログが表示されます。既定設定のままで構いません。主な設定は、（擬似）時間 period が 1.0、小変形仮定（Nlgeom=off）を使用し、 Other tab では direct equation solver を使用する設定です。説明文は任意で、入力しなくても構いません。

Assign Section（Sectionを割り当て） OK 確定ボタンを押して Step の作成を完了します。 Step Manager アイコン を、 Create Step（Stepを作成） アイコン右側から開くと、新しい Step の横に緑のチェックマーク（）が表示されます。次図を参照してください。

2.5.2 Field output を定義する

現在、Coreform IGA for Abaqusがフィールド出力としてサポートするのは変位（U）、応力（S）、相当塑性ひずみ（PEEQ）です。また現時点では、対応するすべての解析でこれらを常に出力するよう固定されています。将来はより柔軟に設定できるよう改善予定ですが、現状では field output を変更しないことを推奨します。

2.5.3 History output を定義する

この例題では Abaqus の history-output request は不要です。次節で定義する Coreform IGA line Probe が、後処理に使用する応力分布を記録します。

2.5.4 Coreform IGA Probe を定義する

Coreform IGA for Abaqus の強力な機能の1つが Coreform IGA Probe です。通常の Probe や history output は node や quadrature point の位置（疎な場合があります）でしか結果を得られなかったり、出力ファイルに保存済みの結果に対してしか使えなかったりします。これに対して Coreform IGA Probe は、任意位置でさまざまな結果を照会できます。² これらの値は解析中に高精度で計算されます。³現在 Coreform IGA for Abaqus には3種類の Probe があります。

• Point Probe: 1つの (X,Y,Z) 位置に配置する Probe。
• Multipoint Probe: 複数の (X,Y,Z) 位置に配置する Probe。
• Line Probe: 2つの (X,Y,Z) 位置で指定した line 上に、等間隔で並ぶ複数の点を配置する Probe。

この例題では Line Probe Line Probe を板厚中央面に作成します。始点は穴の半径位置（Y軸上）とし、Y軸方向に板外表面まで延ばします。Coreform IGA Probe を作成するには Create Coreform IGA Probes アイコン をクリックし、下表と図に示す parameter-value pair を設定します。

Assign Section（Sectionを割り当て） OK 確定ボタンを押して Coreform IGA Probe の作成を完了します。

以上で Step モジュールの作業は完了です。

2.7.1 Pinned displacement boundary condition

次の Create Boundary Condition（境界条件を作成） アイコン をクリックして Create Boundary Condition（境界条件を作成） ダイアログを開き、名前 Hold-XMinで新しい boundary condition を作成します。 Step 選択欄では、作成済みの Apply-Roark-Conditions Step が選択されていることを確認し、category に Mechanical を、type に Symmetry/Antisymmetry/Encastre を選択します。

Assign Section（Sectionを割り当て） Continue ボタンをクリックします。

viewport で boundary condition の region を選択するよう求められます。前述の X-minimum surface を選択します。set 名として XMin.

Assign Section（Sectionを割り当て） Done（完了） ボタンをクリックします。

を推奨します。続く Edit Boundary Condition（境界条件を編集） ダイアログで、 ダイアログで PINNED

Assign Section（Sectionを割り当て） OK を選択し、surface の displacement を0に拘束します。 確定ボタンを押して boundary condition の作成を完了します。 アイコン を、 Create Boundary Condition（境界条件を作成） Boundary Condition Manager アイコンを開くと、新しい boundary condition の横に緑のチェックマーク（）が表示されます。次図を参照してください。

）が表示されます。また viewport には、拘束された surface と degree of freedom（DOF）を示す boundary condition の表示が現れます。

2.7.2 Pressure load condition

次の Create Load Condition（荷重条件を作成） アイコン をクリックして Create Load Condition（荷重条件を作成） ダイアログを開き、名前 Pull-XMaxで新しい boundary condition を作成します。 Step 選択欄では、作成済みの Apply-Roark-Conditions Step が選択されていることを確認し、category に Mechanical を、type に Pressure を選択します。

Assign Section（Sectionを割り当て） Continue ボタンをクリックします。

で新しい load condition を作成します。viewport で region を選択するよう求められるので、前述の X-maximum surface を選択します。set 名として XMax.

Assign Section（Sectionを割り当て） Done（完了） ボタンをクリックします。

を推奨します。続く を推奨します。続く Edit Load（荷重を編集） ダイアログで、下表と図に示す parameter-value pair を設定します。

Assign Section（Sectionを割り当て） OK 確定ボタンを押して load condition の作成を完了します。

viewport には、先ほどの boundary condition に加えて pressure load condition が表示されます。圧力を表す vector が surface から外向き（+X方向）であることを確認してください。surface 内向き（-X方向）なら、Magnitude 欄の - 符号を付け忘れている可能性があります。

以上で Load モジュールの作業は完了です。

2.10.1 Abaqus ODB file を扱う

Abaqus ODB file を表示すると、 XMin と XMax set を定義した位置に2つの triangle mesh と、多数の ✕ user element marker が表示されます。

次に Create Display Group アイコン をクリックして Create Display Group ダイアログを開き、 RNODE3D 対象の element type を選択して Remove（削除） アイコン をクリックし、user-element marker を表示対象から除外します。

その後 Dismiss（閉じる） をクリックしてダイアログを閉じます。

表示されている surface は、boundary condition と load condition で選択した surface に対応します。

2.10.2 IGA ODB file を扱う

IGA ODB file に切り替えると、Abaqus/CAE が user element の直接可視化をサポートしていないことへの代替として、解形状を近似する output mesh が表示されます。trim されていない element は C3D8 element で、trim された element は C3D4 element でメッシュ化されています。

応力集中の評価で着目する応力不変量である最大主応力の contour plot を作成できます。task pane の Plot Contours on Deformed Shape アイコン をクリックし、続いて toolbar の Field Output Dialog を開きます。 Field Output ダイアログで S を Output Variable と として選択し、Invariant に を Max. Principal を指定します。次図を参照してください。

Assign Section（Sectionを割り当て） OK ボタンをクリックします。

viewport には output mesh とともに、最大主応力で色付けされた解析結果が表示されます。またモデル内の最大主応力値を示す legend も表示されます。

Abaqus/CAE の Probe 機能も利用できます。応力集中が予想される位置を見やすい方向へ view を調整し、 Probe values アイコン を task pane から開きます。

を推奨します。続く をクリックします。 Probe Values ダイアログで Nodes を Probe 欄から選択します。次に viewport で、予想される最大応力位置に最も近い node を選択します。

すると をクリックします。 ダイアログ内に、その node に対応する最大主応力を示す行が追加されます。

最後に、 probe_data.json file を Python から Kernel Command Line Interface で開き、応力成分を取得して最大主応力を計算することもできます。現在の Probe schema では、sample は寄与する IGA part または leaf ごとにグループ化されています。 source_location_indices の値は、それぞれの group 内 sample を、要求した line-Probe location へ対応付けます。以下の script は final frame の全 sample を集め、その順序を復元し、境界で複数の leaf が寄与する場合は複数の値を保持します。

import json
import os

import numpy as np

job_dir = r"C:\Users\Owner\Documents\AbaqusTemp\CoreformIGA\ExampleProblems\PWH-Kt\JobDir"
probe_filename = os.path.join(job_dir, "probe_data.json")
with open(probe_filename) as probe_file:
    probe_data = json.load(probe_file)

stress_probe = probe_data[
    "Apply-Roark-Conditions"
]["history"]["Probe-Stress-Profile"]
components = ("XX", "YY", "ZZ", "XY", "XZ", "YZ")
samples = []

for group in stress_probe["part_groups"]:
    source_indices = group["source_location_indices"]
    locations = group["locations"]
    stress = group["Stress"]
    component_values = [stress[name][-1] for name in components]
    expected_count = len(source_indices)
    if len(locations) != expected_count or any(
        len(values) != expected_count for values in component_values
    ):
        raise ValueError("probe locations and final stress values are misaligned")

    for slot, (source_index, location) in enumerate(
        zip(source_indices, locations)
    ):
        sxx, syy, szz, sxy, sxz, syz = (
            values[slot] for values in component_values
        )
        stress_tensor = np.array([
            [sxx, sxy, sxz],
            [sxy, syy, syz],
            [sxz, syz, szz],
        ])
        max_principal = np.linalg.eigvalsh(stress_tensor)[-1]
        samples.append((
            source_index,
            group.get("part_id", -1),
            location,
            max_principal,
        ))

if not samples:
    raise ValueError("probe contains no final-frame stress samples")

samples.sort(key=lambda sample: (sample[0], sample[1]))
source_indices = np.asarray([sample[0] for sample in samples])
xyz = np.asarray([sample[2] for sample in samples], dtype=float)
s_max = np.asarray([sample[3] for sample in samples], dtype=float)
hole_edge_values = s_max[source_indices == 0]
if not hole_edge_values.size:
    raise ValueError("probe has no sample at requested source location 0")
hole_edge_stress = np.max(hole_edge_values)
print("Maximum principal stress at the hole edge:", hole_edge_stress)

受け入れ値となる最大主応力は、Roark Table 17.1, case 7a の net-section stress-concentration factor を用いて計算します。 [1]ここで \(D\) は板幅、 \(r\) は穴半径、 \(t\) は板厚、 \(P\) は作用圧力の大きさ、 \(A = Dt\) は荷重端面積です。

\[ K_t = 3.00 - 3.13\left(\frac{2r}{D}\right) + 3.66\left(\frac{2r}{D}\right)^2 - 1.53\left(\frac{2r}{D}\right)^3 \]

\[ \sigma_\text{nom} = \frac{F}{t\left(D - 2r\right)} = \frac{PA}{t\left(D - 2r\right)} \]

\[ \sigma_\text{max} = K_t \, \sigma_\text{nom} \]

この例題の寸法と荷重を代入すると、 \(2r/D = 0.4\), \(K_t = 2.23568\), \(F = 5.0000 \times 10^{-7}\,\text{lbf}\), \(\sigma_\text{nom} = 1.6667 \times 10^{-6}\,\text{lbf}/\text{in}^2\)、 \(\sigma_\text{max} = 3.7261 \times 10^{-6}\,\text{lbf}/\text{in}^2\).

以下では、この解析解による計算値と、本手順で得られた数値解析の観測値を比較します。

Probe に沿った最大主応力は、 XYData object を次のように作成して plot することもできます。

session.XYData(
    name="S_max_vs_Y",
    data=tuple(map(tuple, np.column_stack((xyz[:, 1], s_max)))),
    xValuesLabel="Y position",
    yValuesLabel="Maximum principal stress",
)

次に XY Data Manager アイコンをクリックして XY Data Manager ダイアログを開き、作成した S_max_vs_Y data を選択して Plot ボタンをクリックします。

をクリックします。viewport に、次図と同様の Probe の line plot が表示されます。

以上で Visualization モジュールの作業は完了です。