一、問題描述

有半徑為a中心孔的均勻薄板受到單軸壓力，應力為1000MPa，中心孔半徑a = 0.5 in., 薄板高2h，寬2w，h = 3 in., w = 6 in., 彈性模量E = 2(10)⁶ psi，泊松比v=0.3，解決平面應力問題，并將有限元的近似解與基于彈性力學理論的精確解進行對比。

二、理論分析

考慮這類中心開孔方板，受到單軸拉力，圓孔圓心和矩形中心重合，處于平面應力狀態，使用有限元求解此結構的變形圖。

首先對此結構進行單元剖分，確定單元按照有限元的分析流程，要先對此結構進行單元剖分，確定單元與結點編號、以及單元的自由度編號。因為這里是平面應力問題，所以可以采用常應變三角形單元進行網格劃分，并且采用的是非結構化的網格。

在該結構中采用的是常應變三角形單元，在整體坐標系中單元剛度矩陣均用下式進行計算：

對于三個節點，有

還可以得到

代入可以得到：

其中A為三角形單元面積，給定三角形節點坐標，可以通過下式得到：

三角形單元應變能為

在三角形單元中，外力做功為

單元總勢能為

單元處于平衡狀態，總勢能最小，有

代入上式得到

其中k即為單元剛度矩陣。

三、圓孔的孔邊應力集中理論

 

五、網格劃分

六、應力云圖

七、對比分析

有限元解(數值解)，最終輸出的應力極值為3096MPa；彈性力學書上的理論解為3100MPa，原因是有限元網格劃分所存在的誤差，導致計算結果存在一定的誤差，但由于誤差不超過數值的5%，證明有限元仿真結果的準確性。

八、總結

有限元分析的最大特點就是標準化和規范化，這種特點時使大規模分析和計算成為可能。實現有限元分析標準化和規范化的載體就是單元，通過構造具有代表性的單元裝配成復雜的結構。在此例中構造三角形單元，先列出小剛度矩陣，再進行裝配，最后帶入邊界條件和外力得出位移、應力、應變的解，并且畫出云圖。云圖中體現了應力集中，即在小孔周圍網格密集（即應力大），在遠離小孔的地方應力網格稀疏，符合了彈性力學的理論。彈性力學中，孔邊的邊界條件是極坐標下的正應力與切應力均為零。而通過有限元方法求出的孔邊應力大致符合彈性力學理論

附錄代碼

#%%

import meshpy.triangle as triangle

import numpy as np

import numpy.linalg as la

import pandas as pd

import matplotlib.pyplot as pt

pt.figure(figsize=(12,6))#畫布的大小

# Quadratic element demo, by Aravind Alwan

def round_trip_connect(start, end):##劃分網格

   return [(i, i + 1) for i in range(start, end)] + [(end, start)]

mesh_points=[]

def main(R):

   JUZHEN = []

   # R = 0.5

   points = [(6, 0), (6, 3), (-6, 3), (-6, -3), (6, -3), (6, 0)]

   facets = round_trip_connect(0, len(points) - 1)

   circ_start = len(points)

   points.extend((R * np.cos(angle), R * np.sin(angle)) for angle in np.linspace(0, 2 * np.pi, 50, endpoint=False))

   facets.extend(round_trip_connect(circ_start, len(points) - 1))