有限元專場 | 有限元模擬分析誤差來源之剪力自鎖+沙漏效應
固有頻率怎麼算出來的
導讀：2020模擬知識周
即將收官，受到使用者高度好評！7月16日（今晚）20時有限元分析專場（第9場直播），我們邀請
模擬秀優秀講師
，
有限元培訓超過5000+
學員的
張老師
直播，關注
Ansys2020R2新版功能解讀
和
案例實操-面向工程應用的結構工程師思維方法，
感興趣的朋友一定要抓緊報名（抽取20個10000核時平行計算資源+10000元實物精品+免費領取58套模擬秀影片等）。以下是正文：
自我上一篇有限元模擬分析誤差來源剖析有限元模擬分析誤差來源之網格劃分，已經過去30多天了
，不少朋友在後臺留言，是否還會有續集？今天正式告訴大家，我還會繼續寫下去，
歡迎大家繼續關注我在模擬秀官網和APP釋出的原創文章，真心希望能引發工程師朋友們的共鳴。
溫馨提示：些新來的朋友可以透過以下文章連結檢視我以往的文章，老朋友也可以溫故而知新哦，任何疑問，歡迎大家留言交流。
有限元模擬分析誤差來源之網格劃分
有限元模擬分析誤差來源之固有頻率
有限元模擬分析誤差來源之邊界條件設定-動載荷
有限元模擬分析誤差來源之邊界條件，約束和point mass
有限元模擬分析誤差來源之材料引數設定，小心為妙！
而今天，我將繼續給大家帶來最新撰寫的原創文章
《有限元模擬分析誤差來源之剪力自鎖+沙漏效應》
。眾所周知，有限元分析中單元設定不當會造成計算的誤差，其中最常見的就是剪力鎖死和沙漏效應。下面我就這剪力自鎖和沙漏效應來討論下有限元模擬分析誤差的來源。如有不當，歡迎大家批評指正。
一、高斯積分
要搞清楚沙漏效應和剪力自鎖，先從單元平衡方程說起。咱們知道有限元的單元平衡方程是透過最小勢能原理推出來的。最小勢能原理長這個樣子。
看著好簡單，但簡單的公式往往蘊含著極深刻的道理。這個公式在說平衡狀態的位移場，使得總勢能取極值。總勢能包括應變能和外力勢能，這樣公式就變成這樣。
從表示式來看，應變能：
體積力勢能：
面力勢能：
發現這些都需要求積分。對計算機來說，求積分是困難的，但求和就很方便。於是工程師就開始找數學工具，就找到了高斯積分。高斯在幾百年前，就發現一個定積分可以近似等於多項式求和，
高斯積分公式如下所示。
其中：
叫做高斯積分權；
叫高斯積分點。
二、剪力自鎖
下一個概念完全積分，就是高斯積分點的個數可以對單元的剛度矩陣可以精確積分，其中一次單元每個方向兩個，二次單元每個方向三個。示意圖如下。
對於受彎載荷單元，單元變形應該如下所示。
完全積分的一次單元由於一次單元無法模擬彎曲，導致單元的剛度變大，示意圖如下。
如圖所示，在積分點處水平虛線和豎直虛線的夾角不再為90度，這樣剪應變就產生了。可能會導致計算結果不可信。
1、
算例演示
給出一個懸臂樑例的算例，如下所示。
透過材料力學，我們知道懸臂樑端部的位移為：
利用workbench來模擬，
看看剪力自鎖到底會帶來什麼樣的誤差？
workbench中單元屬性一般是自動賦予的，大家可能不知道怎麼設定單元屬性，首先設定完全積分，要在Gemetry中將Element control設定成manual，如下所示。
這個時候才可以對模型進行單元屬性設定，找到相關模型中Brick Integration scheme設定成Full，如下圖所示。
這是對完全積分設定，下面對一次單元進行設定，在mesh中將Element Midside Nodes改成Drop就是一次單元了，如下所示。
計算得到結果如下
端部位移為0。93mm，誤差達到了70%，是不可接受的。在將網格細化後，再計算，結果如下。
端部位移為0。0015mm，誤差為50%。
將一次單元變為二次單元，用同樣的網格條件，看看計算結果，如下。
計算結果為3mm，誤差3%，可以接受。
將應變在不同的單元提取出來，一次單元的應變情況如下。
看剪應變的絕對值，即0。0007。
二次單元的應變情況如下。
看剪應變的絕對值，即
。
看到在完全積分時，一次單元產生了較大的剪應變，即發生了剪力自鎖現象。當然剪力自鎖不一定對計算結果產生那麼大影響，我們細化網格再計算，結果如下。
計算結果為2。85mm，誤差8%，還是可以接受的。
2、剪力自鎖結論
1）受彎曲載荷時，採用完全積分的一次單元，會產生剪力自鎖，對計算精度產生影響；採用完全積分的二次單元，則不會產生剪力自鎖，計算精度高。
2）網格細化，降低網格長細比，可以提高計算精度，降低剪力自鎖帶來的影響。
三、沙漏效應
以上我們討論了剪力自鎖，接下來，我們討論下有限元模擬過程中經常會遇到的另一個問題-沙漏效應。
1、模態計算
不知道大家在做模態分析時，會不會遇到這樣的情況。
還是以前文中懸臂樑的例子，從圖中可以看出一階固有頻率只有6。8Hz，這個結果可信嗎？我們驗證下，固有頻率的計算公式如下。
K
是模態剛度
， M
是模態質量。
這兩個引數好像都不會計算，不著急，先計算剛度。對圖中振型進行觀察，發現第一階固有頻率是整體表現的起伏。那我們就可以透過撓度來計算，懸臂樑的撓度計算公式，大家沒忘吧。
把這個公式變形下，公式如下所示。
那麼模態剛度 ：
計算得到剛度：
還有質量呢？我們還是觀察振型，發現這個振型幾乎是整個模型參與的，那麼就認為這階模態質量就是整個模型的質量：
代入公式1中，得到一階固有頻率為 ：
這個和模擬的數值差太遠，為什麼會這樣呢？
從單元設定上，我將單元設定為減縮積分的一次單元，會不會跟這個有關呢？
2、減縮積分
上篇文章中討論了完全積分，這次說下減縮積分。減縮積分比完全積分在每個方向上少用一個積分點，那麼減速積分狀態下的一次單元和二次單元積分點的分佈，如下所示。
這樣當單元受彎時，就會產生這樣的情況，如下所示。
這樣中間的橫線，實際上既沒有伸長，也沒有縮短，所以單元沒有產生正應力。此外，橫線和豎線的夾角也沒有發生變化，所以單元也沒有產生剪應力。單元什麼應力也沒有，所以單元不能承受任何載荷，即剛度為零。這就是沙漏效應。
3、案例演示
把一次單元換成二次單元，在計算下模態。
模擬得到結果為91。2Hz。同理論計算值相比，只有0。6%的誤差。
減縮積分的二次單元是不是就不存在沙漏現象呢？如果把網格化的粗一些，計算結果變成如圖所示。
一階固有頻率變成0Hz，不光是一階固有頻率，其他階的也都變成0或者很小的數值，如下所示。
為什麼會這樣呢？實際上跟單元的長細比有關，如果單元特別細長，積分點就會非常接近，極限情況兩條積分線重合在一起。如下所示。
我們看變形後的單元形狀，二次單元可以彎曲，那麼豎線和橫線的夾角是90度，所以單元沒有剪應力。單元特別細長，橫線基本位於中性層，那麼它也不會伸長或縮短，所以沒有正應力。那麼沙漏效應又產生了。
4、沙漏效應總結
在減縮積分情況下，不管是一次單元還是二次單元，都有可能產生沙漏效應。要消除沙漏效應，就要提高網格質量，減小網格的長細比。
減縮積分可以有效的排除剪力自鎖的影響，不會產生剛度過大的情況。
全文完！
來自2020模擬周的邀請函
《2020模擬知識周》已於7月6日準時開幕！
此次活動，我們邀請到了
10位資深名企工程師
做客模擬秀直播間，分別從理論知識，軟體使用，實操講解，案例分析等方面幫助CAE模擬技術人員高效學習模擬知識，攻克模擬技術應用的難點。為使用者帶來各大行業真正接地氣的行業應用經驗。點選檢視詳情《2020模擬知識周：聚焦6大模擬行業》
7月16日（今晚）20時有限元分析專場專場，我們邀請
模擬秀優秀講師
，
有限元培訓超過5000+
學員的
張老師
直播，關注
Ansys2020R2新版功能解讀
和
案例實操-面向工程應用的結構工程師思維方法，講課大綱
如何
：
