メインコンテンツへスキップ
見出し画像

書記が物理やるだけ#428 塑性力学における降伏条件と硬化則

    塑性力学の基本事項を見ていく。


    問題

    画像


    画像


    説明

    無負荷の状態から荷重をかけ始めると、特定の応力値(弾性限度)に達するまで、応力とひずみは比例関係を保つ。しかし、ひずみがさらに大きくなると材料が降伏し、この比例関係が崩れて応力とひずみの関係は非線形となる。この領域では、除荷しても変形が完全には戻らず、残留ひずみが生じる。このような永久的な変形を塑性変形と呼ぶ。

    画像


    バウシンガー効果(Bauschinger effect)は、一度ある方向に塑性変形を与えたのち、逆方向の荷重を加えると、再び同方向に荷重を加えたときより塑性変形が低い応力でおこる効果である。

    画像


    降伏条件には、主に「ミーゼスの降伏条件」と「トレスカの降伏条件」がある。

    ミーゼスの降伏条件(せん断ひずみエネルギー説)は、材料内部に蓄えられたひずみエネルギーが一定の限界値に達したときに降伏が生じるとする考え方である。鋼などの延性材料における実験結果と良好に一致するため、実務設計で最も広く用いられている。主応力空間における降伏曲面は円筒形で表される。

    一方、トレスカの降伏条件(最大せん断応力説)は、材料内の最大せん断応力が限界値に達したときに降伏するとする考え方である。主応力空間では、ミーゼスの円筒に内接する正六角柱として表される。これはミーゼスよりも降伏を厳しめ(安全側)に判定する傾向がある。

    画像


    画像


    解答

    降伏条件の計算には,主偏差応力を用いる。

    画像


    ミーゼスの降伏条件は、等価応力が偏差応力の不変量(2乗和に相当)で定義されるため、主応力空間における降伏曲面は静水圧軸(中心軸)からの距離が一定の円筒形となる。一方、トレスカの降伏条件は、最大主応力と最小主応力の差(最大せん断応力)を制限するものである。この条件は各主応力の大小関係に応じて6つの線形な境界式で定義されるため、形状はそれらが組み合わさった正六角柱となる。

    画像


    σ1-σ2平面に射影した降伏境界を以下に可視化する。

    import numpy as np
    import matplotlib.pyplot as plt
    
    def plot_yield_surfaces_2d(sigma_y=150):
        """
        sigma_1 - sigma_2 平面 (sigma_3 = 0) 上にミーゼスとトレスカの降伏曲線をプロットする。
        """
        # 応力範囲の設定(降伏応力の±1.5倍程度)
        limit = sigma_y * 1.5
        s_max = np.linspace(-limit, limit, 500)
        s1, s2 = np.meshgrid(s_max, s_max)
    
        # --- 1. ミーゼスの降伏条件 (von Mises) ---
        # 式: sigma_1^2 - sigma_1*sigma_2 + sigma_2^2 = sigma_Y^2
        # sigma_3 = 0 の場合の等価応力を計算
        sigma_eq_mises = np.sqrt(s1**2 - s1*s2 + s2**2)
    
        # --- 2. トレスカの降伏条件 (Tresca) ---
        # 式: max(|sigma_1 - sigma_2|, |sigma_2 - sigma_3|, |sigma_3 - sigma_1|) = sigma_Y
        # sigma_3 = 0 の場合: max(|sigma_1 - sigma_2|, |sigma_2|, |sigma_1|) = sigma_Y
        # sigma_2 と sigma_3=0 の差の絶対値
        # sigma_3 と sigma_1 の差の絶対値
        tresca_12 = np.abs(s1 - s2)
        tresca_23 = np.abs(s2 - 0)
        tresca_31 = np.abs(0 - s1)
        # 3つの最大値をとる
        sigma_eq_tresca = np.maximum(tresca_12, np.maximum(tresca_23, tresca_31))
    
        # --- 作図 ---
        plt.figure(figsize=(8, 8))
    
        # ミーゼスの楕円 (等高線として描画)
        # sigma_eq が sigma_y に等しい線のみを描く
        contour_mises = plt.contour(s1, s2, sigma_eq_mises, levels=[sigma_y], colors='b', linewidths=2)
        # plt.clabel(contour_mises, inline=True, fontsize=10, fmt='von Mises') # ラベルを追加したい場合
    
        # トレスカの六角形 (等高線として描画)
        contour_tresca = plt.contour(s1, s2, sigma_eq_tresca, levels=[sigma_y], colors='r', linewidths=2, linestyles='--')
        # plt.clabel(contour_tresca, inline=True, fontsize=10, fmt='Tresca') # ラベルを追加したい場合
    
        # --- グラフの装飾 ---
        plt.title(f'Yield Surfaces in $\sigma_1$-$\sigma_2$ Plane ($\sigma_3=0$)\n$\sigma_Y = {sigma_y}$ MPa', fontsize=14)
        plt.xlabel('Principal Stress $\sigma_1$ (MPa)', fontsize=12)
        plt.ylabel('Principal Stress $\sigma_2$ (MPa)', fontsize=12)
    
        # 軸線 (x=0, y=0)
        plt.axhline(0, color='black', linewidth=1)
        plt.axvline(0, color='black', linewidth=1)
    
        # グリッド
        plt.grid(True, linestyle=':', alpha=0.6)
    
        # 軸の範囲とアスペクト比を正方形に
        plt.xlim(-limit, limit)
        plt.ylim(-limit, limit)
        plt.gca().set_aspect('equal', adjustable='box')
    
        # 凡例の作成 (contourは直接凡例にできないため、ダミーの線を作成)
        plt.plot([], [], color='b', label='von Mises (Ellipse)')
        plt.plot([], [], color='r', linestyle='--', label='Tresca (Hexagon)')
        plt.legend(loc='upper right')
    
        plt.show()
    
    # --- 実行 ---
    if __name__ == "__main__":
        # これまでの演習と同じ降伏応力 150 MPa でプロット
        plot_yield_surfaces_2d(sigma_y=150)


    画像


    π 平面は、3次元の主応力空間において、静水圧軸 n=(1,1,1)^T に垂直な平面である。

    画像


    import numpy as np
    import matplotlib.pyplot as plt
    from matplotlib.patches import Circle
    
    def to_pi_plane_2d(s1, s2, s3):
        """
        3次元の主応力(または偏差応力)をpi平面上の2次元座標(x, y)に変換する。
        """
        x = (2*s1 - s2 - s3) / np.sqrt(6)
        y = (s2 - s3) / np.sqrt(2)
        return x, y
    
    def plot_pi_plane_and_point(sigma_y=150, target_sigma=(160, 60, 60)):
        """
        pi平面上にvon Mises(円)とTresca(六角形)の降伏曲面を描き、
        特定の応力状態(点P)をプロットする。
        """
        s1_t, s2_t, s3_t = target_sigma
    
        # --- 1. ミーゼスの円 (von Mises) ---
        R = np.sqrt(2/3) * sigma_y
    
        # --- 2. トレスカの六角形 (Tresca) ---
        tx_bound = (2*sigma_y - 0 - 0) / np.sqrt(6)
        ty_bound = 0
        angles = np.linspace(0, 2*np.pi, 7)
        tx = tx_bound * np.cos(angles)
        ty = tx_bound * np.sin(angles)
    
        # --- 3. 特定の応力状態(点P)の変換 ---
        sigma_m = sum(target_sigma) / 3
        s1, s2, s3 = s1_t - sigma_m, s2_t - sigma_m, s3_t - sigma_m
        px, py = to_pi_plane_2d(s1, s2, s3)
        r_p = np.sqrt(px**2 + py**2)
    
        # --- 作図 ---
        fig, ax = plt.subplots(figsize=(8, 8))
    
        # von Mises 円
        circle_mises = Circle((0, 0), R, color='blue', fill=False, linestyle='-', linewidth=2, label='von Mises')
        ax.add_patch(circle_mises)
    
        # Tresca 六角形
        ax.plot(tx, ty, 'red', linewidth=2, label='Tresca')
    
        # 特定の応力状態(点P)
        ax.plot(px, py, 'ko', markersize=10)
        
        # テキスト表示 (Raw文字列 fr"" を使用して LaTeX エラーを回避)
        label_p = fr"P $(160, 60, 60)$" + "\n" + fr"dist $r \approx {r_p:.1f}$"
        ax.text(px, py + 8, label_p, fontsize=10, ha='center', fontweight='bold')
    
        # --- グラフの装飾 ---
        ax.set_title(fr"$\pi$-Plane Yield Surfaces ($\sigma_Y = {sigma_y}$ MPa)", fontsize=14)
        ax.set_xlabel(r"$\pi$-axis $x$ ($\sigma_1$ direction)", fontsize=12)
        ax.set_ylabel(r"$\pi$-axis $y$ ($\sigma_2 - \sigma_3$ direction)", fontsize=12)
    
        ax.axhline(0, color='black', linewidth=0.8)
        ax.axvline(0, color='black', linewidth=0.8)
        ax.grid(True, linestyle=':', alpha=0.6)
    
        # 描画範囲の設定
        limit = max(R, r_p) * 1.3
        ax.set_xlim(-limit, limit)
        ax.set_ylim(-limit, limit)
        ax.set_aspect('equal')
        ax.legend(loc='upper right')
    
        # 判定結果の表示
        if r_p < R:
            res = "Result: Inside (Elastic)"
            res_col = "green"
        else:
            res = "Result: Outside (Yield)"
            res_col = "red"
        
        ax.text(0.05, 0.05, res, transform=ax.transAxes, fontsize=12, 
                color=res_col, bbox=dict(facecolor='white', alpha=0.8))
    
        plt.tight_layout()
        plt.show()
    
    if __name__ == "__main__":
        plot_pi_plane_and_point(sigma_y=150, target_sigma=(160, 60, 60))


    画像


    等方硬化では、降伏曲面が原点を中心に相似形を保って拡大する。π 平面におけるミーゼスの降伏円の半径について、降伏応力 σY  の増加に比例して降伏圏が拡大する。例えば、σY  が 150MPa から 200MPa に上昇した場合、半径は約 1.33 倍となる。このとき、負荷応力が 180MPa であれば、更新された降伏応力 200MPa を下回るため、塑性変形は進行せず弾性域内として扱われる。
    一方、移動硬化では、降伏曲面はその形状と大きさを維持したまま、応力空間内をシフトする。例えば、引張方向に 50MPa シフトした場合、引張側の降伏点は 150+50=200MPa に上昇するが、逆方向の圧縮降伏点は −150+50=−100MPa に低下する。このように、「一方向への塑性変形が逆方向の降伏応力の低下を招く現象」がバウシンガー効果である。

    import numpy as np
    import matplotlib.pyplot as plt
    
    def plot_mises_ellipse(ax, sigma_y, alpha1=0, alpha2=0, color='b', label='', linestyle='-'):
        """
        sigma_1 - sigma_2 平面 (sigma_3 = 0) 上に、
        背応力 (alpha1, alpha2) を考慮したミーゼスの降伏楕円をプロットする。
        式: (s1-a1)^2 - (s1-a1)(s2-a2) + (s2-a2)^2 = sigma_Y^2
        """
        # 描画用データ生成(角度をパラメータとする)
        theta = np.linspace(0, 2*np.pi, 200)
        
        # 標準的なミーゼス楕円(原点中心)の座標を生成
        # 幾何学的な変換を用いて楕円の長軸・短軸方向から座標を求める
        r_major = np.sqrt(2) * sigma_y # 長軸方向の半径
        r_minor = np.sqrt(2/3) * sigma_y # 短軸方向の半径
        
        # 45度回転した座標系での楕円
        x_rot = r_major * np.cos(theta)
        y_rot = r_minor * np.sin(theta)
        
        # 元の sigma_1, sigma_2 座標系に戻す(45度回転)
        s1_standard = (x_rot - y_rot) / np.sqrt(2)
        s2_standard = (x_rot + y_rot) / np.sqrt(2)
        
        # 背応力分だけ平行移動
        s1 = s1_standard + alpha1
        s2 = s2_standard + alpha2
        
        ax.plot(s1, s2, color=color, linestyle=linestyle, linewidth=2.5, label=label)
    
    # --- パラメータ設定 ---
    sigma_y0 = 150  # 初期降伏応力 (MPa)
    
    # ==================================================================
    # 可視化 [問題1:等方硬化 (Isotropic Hardening)]
    # ==================================================================
    def visualize_isotropic_hardening():
        sigma_y_new = 200  # 硬化後の降伏応力 (MPa)
        
        fig, ax = plt.subplots(figsize=(8, 8))
        
        # 1. 初期降伏曲面(原点中心、半径 小)
        plot_mises_ellipse(ax, sigma_y0, color='blue', label=f'Initial ($\sigma_Y={sigma_y0}$)')
        
        # 2. 等方硬化後の降伏曲面(原点中心、半径 大)
        plot_mises_ellipse(ax, sigma_y_new, color='cyan', linestyle='--', label=f'Isotropic Hardening ($\sigma_Y={sigma_y_new}$)')
        
        # 演習問題の点 P(180, 0) をプロット
        point_p = (180, 0)
        ax.plot(point_p[0], point_p[1], 'ro', markersize=10)
        ax.text(point_p[0]+5, point_p[1]+5, f'P {point_p}', color='red', fontsize=12, fontweight='bold')
    
        # グラフの装飾
        ax.set_title(f'Isotropic Hardening ($\sigma_1$-$\sigma_2$ Plane)\nEllipse Expands Uniformly', fontsize=14)
        ax.set_xlabel('$\sigma_1$ (MPa)', fontsize=12)
        ax.set_ylabel('$\sigma_2$ (MPa)', fontsize=12)
        
        limit = sigma_y_new * 1.3
        ax.set_xlim(-limit, limit)
        ax.set_ylim(-limit, limit)
        ax.axhline(0, color='black', linewidth=1)
        ax.axvline(0, color='black', linewidth=1)
        ax.grid(True, linestyle=':', alpha=0.6)
        ax.set_aspect('equal', adjustable='box')
        ax.legend(loc='upper right')
        
        plt.tight_layout()
        plt.show()
    
    # ==================================================================
    # 可視化 [問題2:移動硬化 (Kinematic Hardening)]
    # ==================================================================
    def visualize_kinematic_hardening():
        # 背応力テンソル (alpha1, alpha2, alpha3=0)
        alpha = (50, 0) # sigma_1 方向に 50 MPa シフト
        
        fig, ax = plt.subplots(figsize=(8, 8))
        
        # 1. 初期降伏曲面(原点中心)
        plot_mises_ellipse(ax, sigma_y0, color='blue', label=f'Initial (Center=(0,0))')
        
        # 2. 移動硬化後の降伏曲面(中心が移動、サイズは同じ)
        plot_mises_ellipse(ax, sigma_y0, alpha1=alpha[0], alpha2=alpha[1], 
                           color='magenta', linestyle='-.', label=f'Kinematic Hardening (Center={alpha})')
        
        # 中心位置(背応力)をプロット
        ax.plot(alpha[0], alpha[1], 'k+', markersize=15, markeredgewidth=2)
        ax.text(alpha[0]+5, alpha[1]+5, f'Back Stress {alpha}', color='black', fontsize=11)
    
        # バウシンガー効果の引張・圧縮降伏点を強調
        ax.plot(sigma_y0 + alpha[0], 0, 'go', markersize=8) # 引張側: 200
        ax.text(sigma_y0 + alpha[0]+5, 10, 'Tensile Yield (200)', color='green')
        
        ax.plot(-sigma_y0 + alpha[0], 0, 'go', markersize=8) # 圧縮側: -100
        ax.text(-sigma_y0 + alpha[0]+5, 10, 'Comp. Yield (-100)', color='green')
    
        # グラフの装飾
        ax.set_title(f'Kinematic Hardening ($\sigma_1$-$\sigma_2$ Plane)\nEllipse Translates (Bauschinger Effect)', fontsize=14)
        ax.set_xlabel('$\sigma_1$ (MPa)', fontsize=12)
        ax.set_ylabel('$\sigma_2$ (MPa)', fontsize=12)
        
        limit = (sigma_y0 + alpha[0]) * 1.3
        ax.set_xlim(-limit, limit)
        ax.set_ylim(-limit, limit)
        ax.axhline(0, color='black', linewidth=1)
        ax.axvline(0, color='black', linewidth=1)
        ax.grid(True, linestyle=':', alpha=0.6)
        ax.set_aspect('equal', adjustable='box')
        ax.legend(loc='upper right')
        
        plt.tight_layout()
        plt.show()
    
    # --- 実行 ---
    if __name__ == "__main__":
        visualize_isotropic_hardening() # 問題1 の可視化
        visualize_kinematic_hardening() # 問題2 の可視化


    画像


    画像


    本記事のもくじはこちら:


    あなたへのおすすめ