基于ANSYS的自適應網格劃分

何為網格自適應劃分？

ANSYS程序提供了近似的技術自動估計特定分析類型中因為網格劃分帶來的誤差。（誤差估計在ANSYS Basic Analysis Procedures Guide第五章中討論。）通過這種誤差估計，程序可以確定網格是否足夠細。如果不夠的話，程序將自動細化網格以減少誤差。這一自動估計網格劃分誤差并細化網格的過程就叫做自適應網格劃分，然后通過一系列的求解過程使得誤差低于用戶指定的數值（或直到用戶指定的最大求解次數）。

 

自適應網格劃分的先決條件

ANSYS軟件中包含一個預先寫好的宏，ADAPT.MAC，完成自適應網格劃分的功能。用戶的模型在使用這個宏之前必須滿足一些特定的條件。（在一些情況下，不滿足要求的模型也可以用修正的過程完成自適應網格劃分，下面還要討論。）這些要求包括：

(1)    標準的ADAPT過程只適用于單次求解的線性靜力結構分析和線性穩態熱分析。

(2)    模型最好應該使用一種材料類型，因為誤差計算是根據平均結點應力進行的，在不同材料過渡位置往往不能進行計算。而且單元的能量誤差是受材料彈性模量影響的。因此，在兩個相鄰單元應力連續的情況下，其能量誤差也可能由于材料特性不同而不一樣。在模型中同樣應該避免殼厚突變，這也可能造成在應力平均是發生問題。

(3)    模型必須使用支持誤差計算的單元類型。（見表1）

(4)    模型必須是可以劃分網格的：即模型中不能有引起網格劃分出錯的部分。

 

表1  自適應網格劃分可用單元

 

2-D Structural Solids

PLANE2  2-D 6-Node Triangular Solid

PLANE25 Axisymmetric Harmonic Solid

PLANE42 2-D 4-Node Isoparametric Solid

PLANE82 2-D 8-Node Solid

PLANE83 Axisymmetric Harmonic 8-Node Solid

 

3-D Structural Solids

SOLID45 3-D 8-Node Isoparametric Solid

SOLID64 3-D Anisotropic Solid

SOLID73 3-D 8-Node Solid with Rotational DOF

SOLID92 3-D 10-Node Tetrahedral Solid

SOLID95 3-D 20-Node Isoparametric Solid

 

3-D Structural Shells

SHELL43 Plastic quadrilateral Shell

SHELL63 Elastic Quadrilateral Shell

SHELL93 8-Node Isoparametric Shell

 

2-D Thermal Solids

PLANE35 2-D 6-Node Triangular Solid

PLANE75 Axisymmetric Harmonic Solid

PLANE55 2-D 4-Node Isoparametric Solid

PLANE77 2-D 8-Node Solid

PLANE78 Axisymmetric Harmonic 8-Node Solid

 

3-D Thermal Solids

SOLID70 3-D 8-Node Isoparametric Solid

SOLID87 3-D 10-Node Tetrahedral Solid

SOLID90 3-D 20-Node Isoparametric Solid

 

3-D Thermal Shells

SHELL57 Plastic Quadrilateral Shell

 

如何使用自適應網格劃分：基本過程

進行自適應網格劃分的基本過程包括如下步驟：

1.       象其他線性靜力分析或穩態熱分析一樣，先進入前處理器（/PREP7或Main Menu>Preprocessor）。然后指定單元類型，實參和材料特性，要滿足上面提到的條件。

2. 用實體建模過程建立模型，用可以劃分網格的面或體建模。用戶不需指定單元大小也不用劃分網格，ADAPT宏會自動劃分網格。（如果要同時劃分面和體網格，生成ADAPTMSH.MAC用戶子程序－見后。）

3. 在PREP7中或在SOLUTION(/SOLU或Main Menu>Solution)中指定分析類型，分析選項，載荷和載荷步選項。在一個載荷步中僅施加實體模型荷載和慣性荷載（加速度，角加速度和角速度）。（通過ADAPTBC.MAC用戶子程序       可以施加有限單元載荷，固連和約束方程。通過ADAPTSOL.MAC用戶子程序       可以加入多個載荷步。這些子程序在后面還要討論。）

4. 如果在PREP7中，退出前處理器[FINISH]。（可以在SOLUTION或在初始狀態下激活ADAPT宏）。

5. 用下列方法激活自適應求解。

   Command: ADAPT

   GUI: Main Menu>Solution>Adaptive Mesh

注意，可以在熱或結構分析中使用ADAPT宏，但不能在一次自適應分析中同時進行這兩種不同類型的計算。在自適應網格劃分的迭代過程中，單元的大小將作調整（在FACMN和FACMX指定的范圍內）以減小或增加單元能量誤差，直到誤差滿足指定的數值（或指定的最大求解次數）為止。

 6. 當自適應網格計算收斂時，程序自動將單元形狀檢查打開[SHPP,ON]。然后返回SOLUTION或初始狀態，這取決于激活ADAPT的狀態。接下來可以進       入POST1用標準操作進行后處理。

修改基本過程

 

選擇自適應性

如果用戶清楚某個部分網格劃分的誤差相對影響較小時（如應力水平較低且變化較小），可以將這些區域從自適應網格劃分中排除以加快分析速度。同樣，用戶也許想將接近應力奇異點的部分（如集中載荷）排除掉。選擇邏輯操作可以解決這類問題。

 

圖 選擇自適應能改進有應力集中的模型

 

如果用戶選擇了一個關鍵點集，ADAPT宏仍將包含進所有的關鍵點（在選擇的和未選擇的關鍵點都作網格改動），除非將ADAPT命令（Main Menu>Solution>Adaptive Mesh）中KYKPS設為1。

如果用戶選擇了一個面或體集，ADAPT宏將只在選擇的區域調整網格大小。此時必須在激活ADAPT宏之前在PREP7中對整個模型進行網格劃分。

 

用用戶子程序定制ADAPT宏

標準的ADAPT宏并不能滿足特定的分析需要。例如，用戶可能想同時對面和體進行網格劃分，這在標準宏當中是不可以的。對于這種或其他一些類似情況，可以對ADAPT宏進行修改使之適用于特定的分析。ANSYS程序用宏這種方式完成自適應網格劃分，本身就使得用戶可以對其進行相應的修改以適應不同的要求。方便的是，用戶不用總是通過修改ADAPT代碼的方式來定制宏。宏的三個部分可以用用戶子程序的方法來修改，這個方法將ADAPT宏和用戶文件分開，用戶可以生成子程序由ADAPT宏來調用。這三個部分是：網格劃分命令序列，邊界條件命令序列和求解命令序列。相應的用戶子程序名為ADAPTMSH.MAC,ADAPTBC.MAC和ADAPTSOL.MAC。下面看一下這三個子程序的功能：

 

生成用戶網格劃分子程序（ADAPTMSH.MAC）

缺省情況下，如果模型中包含一個或多個體時，ADAPT宏將只對體劃分網格而不對面進行劃分。如果當前選擇集中沒有體，宏才對面進行劃分。如果要同時對面和體進行劃分的話，就要生成用戶子程序ADAPTMSH.MAC來提供相應的操作。在重新劃分網格之前要清楚所有這些實體中劃分過網格的實體。子程序大致如下：

C*** Subroutine ADATMSH.MAC-Your name-Job Name-Data Created

TYPE,1             !指定劃分網格的單元類型屬性

ACLEAR,3,5,2       !清除本程序中要重新劃分網格的面和體的網格

VCLEAR,ALL

AMESH,3,5,2        !對面3和5劃分網格（ADAPT不對其他面劃分網格）

TYPE,2             !改變單元類型，劃分體

VMESH,ALL          !對所有體劃分網格

請查閱TYPE,ACLEAR,VCLEAR,AMESH和VMESH命令得到更詳細的解釋。

     我們建議將C***行包含在文件中以區別不同的宏。這一行將在任務的輸出中出現，就可以確認ADAPT宏正確地調用了用戶子程序。

 

生成用戶邊界條件子程序（ADAPTBC.MAC）

ADAPT宏在每次求解循環時都要清除并重新劃分網格，因此模型的結點和單元也要不斷的改變。這就使得所有的有限單元載荷，自由度固連，約束方程等施加在結點和單元的邊界條件都不能使用。如果要包含這些有限單元約束時，就要使用用戶子程序，ADAPTBC.MAC。在這個子程序中，可以選擇結點然后定義有限單元載荷，自由度固連和約束方程。下面是一個ADAPTBC.MAC的例子：

C*** Subroutine ADAPTBC.MAC-Your name-Hob Name-Data CReated

NSEL,S,LOC,X,0      !選擇X坐標為0的結點

D,ALL,UX,0          !令選擇的所有結點UX為0

NSEL,S,LOC,Y,0      !選擇Y坐標為0的結點

D,ALL,UX,0          !令選擇的所有結點UY為0

NSEL,ALL            !選擇所有結點

 

生成用戶求解子程序（ADAPTSOL.MAC）

ADAPT宏中的求解命令序列很簡單：

/SOLU

SOLVE

FINISH

這個缺省的命令序列只能求解單個載荷步。用戶可以將其他的命令序列加入用戶子程序ADAPTSOL.MAC中。

 

關于用戶子程序的一些其他說明

用戶可以象生成其他文件一樣生成這些子程序。就是說可以用APDL語言的*CREATE命令（Utility Menu>Macro>Create Macro）和APDL的*END命令，也可以用外部的文件編輯器。當ADAPT宏調用這些子程序時，軟件先搜索ANSYS根目錄，再搜索用戶根目錄，最后是當前目錄。因此，要確保其他目錄中不包括與所用文件同名的文件。解釋行（C***）會在輸出文件中出現，可以通過它來檢查是否使用了正確的文件。另外，通過在運行ADAPT宏之前用/PSEARCH,OFF(Utility Menu>Macro>Macro Search Path)指定軟件只搜索ANSYS根目錄和用戶當前目錄可以從某種程度上減少文件混用的可能。不管這些子程序在什么位置，它們都能被找到，除非將ADAPT命令的KYMAC選項設為1。

 

定制ADAPT宏(UADAPT.MAC)

有些情況下用戶需要修改ADAPT宏但不能通過單獨的用戶子程序的方式，那么就需要直接修改ADAPT宏的主體。但是，因為某些原因，我們不推薦直接對ADAPT宏進行修改。（例如，別的用戶和你同時使用一個軟件，在調用ADAPT宏時會發現宏被修改了！）因此，在ANSYS安裝中支持一個宏的拷貝文件UADAPT.MAC，便于用戶修改。

如果對UADAPT.MAC文件進行了修改，我們建議對修改后的文件取一個新的文件名。然后在調用時輸入這個文件名。要知道的是，如果新文件名是一個“unknown command”，ANSYS將搜索上級目錄，然后是登錄的目錄，最后是工作目錄，直到找到這個宏為止。如果修改的宏只能為一個用戶使用，那么存儲的位置應在用戶登錄目錄的層次之下（不能等于或高于這個目錄層次）。這樣，存儲的低層次的文件可以通過*USE命令（Utility Menu>macro>Execute Data Block）來調用。

 

自適應網格劃分的一些說明

下面的建議可能有助于自適應網格劃分的使用：

不需指定初始網格大小，但指定大小可能有利于自適應收斂。如果用戶指定了關鍵點網格大小，ADAPT宏在第一次循環時使用這個值，然后在隨后的循環中進行調整。用下列命令指定單元大小：

Command: KESIZE

GUI: Main Menu>Preprocessor>-Meshing-Size Cntrls>-Keypoints-All KPs

    Main Menu>Preprocessor>-Meshing-Size Cntrls>-Keypoints-Picked KPs

如果定義了線分段數或大小比例，ADAPT宏將在每次循環中都使用這個數值而不作改變。如果沒有定義任何形式的網格份數，在初始網格劃分時將使用缺省的網格大小[參見SMRTSIZE和DESIZE命令]。用下列方法指定線分段數或大小比例：

Command: LESIZE

GUI: Main Menu>Preprocessor>-Meshing-Size Cntrls>-Lines-All Lines

    Main Menu>Preprocessor>-Meshing-Size Cntrls>-Lines-Picked Lines

?         映射網格劃分適用于2-D實體和3-D殼單元。但面的映射劃分效果不明顯。

?         映射網格劃分適用于3-D實體。對體進行映射劃分比自由劃分效果要好的多。

?         總體上說，在自適應網格中有中間結點的單元比線性單元要好。

?         不要用集中載荷或尖角等引起奇異性的結構，因為此時ADAPT在這些奇異點處能量值將不收斂。如果模型中有集中載荷時，將其用施加在一個小面上的壓力等效。（或通過選擇將奇異部分排除在自適應網格劃分之外。）

?         在許多情況下，用一系列相對小的區域替代少數幾個大的區域將得到更好的網格劃分。

?         如果最大響應位置已知或事先可以推測，就在附近放置一個關鍵點。

?         如果是在交互方式下運行ADAPT，而ANSYS在沒有提示出錯信息時突然退出，可以在Jobmame.ADPT文件中查看自適應網格劃分部分以確定出錯原因。    同樣，在批處理方式下運行ADAPT時，可以看Jobname.ADPT確定出錯原因。

?         如果模型中有些區域有過度的扭曲時，在網格劃分中就會出錯。在這種情況下，用KESIZE命令(Main Menu>Preprocessor>-Meshing-Size Cntrls>-Keypoints-Picked KPs)中SIZE域指定扭曲區域附近關鍵點的最大單元長度。同時，ADAPT命令中的FACMX將設為1，阻止過度扭曲部分單元大小增加。

?         應當存儲結果文件（Jobname.RST或Jobname.RTH）。在ADAPT運行過程中程序如果發生中斷，結果文件中將保存ADAPT過程已完成求解的內容。

?         在自適應網格運行之前應輸入SAVE命令（Utility Menu>File>Save as Jobname.db）。在程序出錯中斷時，可以用Jobname.db重新啟動計算。

《ANSYS Verification Manual》中的實例：

VM193               二維熱對流的自適應分析

VM205               均布載荷下橢圓薄膜的自適應分析

下面以帶圓孔矩形平板的構應力集中分析來說明自適應網格方法。

如圖1所示，一個帶圓孔的矩形薄板左右兩邊受均布拉力。

幾何尺寸及材料屬性如下：w=h=10mm,R=0.5mm，E=2e5MPa，μ=0.3，q=1MPa。

由于模型和載荷具有對稱性，因此只需要考慮1/4模型。

下面分別是利用菜單和命令流方式進行有限元分析的方法。

圖1  帶圓孔矩形平板的拉伸示意圖

方法一、通過加密網格實現應力集中的精確模擬

1．GUI菜單建模分析過程

第一步，清除內存準備分析。

（1）清除內存。

GUI：Utility Menu>File>Clear& Start New，單擊OK按鈕。

（2）更換工作文件名。

GUI：Utility Menu>File>Change Jobname，輸入the plate under tension。

（3）定義標題。

GUI：Utility Menu>File>Change Title，輸入the plate under tension analysis。

第二步，創建有限元模型。

（1）進入前處理器。

GUI：Main Menu>Preprocessor

（2）定義單元類型。

GUI：Main Menu>Preprocessor>Element Type>Add/Edit/Delete，選擇plane42單元。

（3）定義材料。

GUI：Main Menu>Preprocessor>Material Props>Material Models，彈出對話框，連續雙擊Structural>linear>Elastic>Isotropic，彈出對話框，在值域EX處輸入2e5，在值域PRXY處輸入0.3。

（4）創建矩形。

GUI：Main Menu>Preprocessor>Modeling>Create>Areas>Rectangle，選擇By Dimensions，彈出對話框，在x1、x2值域中分別輸入0和5，y1、y2值域中分別輸入0和5。

（5）創建實心圓。

GUI：Main Menu>Preprocessor>Modeling>Create>Areas>Circle，選擇Solid Circle，彈出對話框，在WP x和WP y域中分別輸入0和0，在Radius域中輸入0.5。

（6）布爾運算。

GUI：Main Menu>Preprocessor>Modeling>Operate>Booleans>Subtract，選擇Areas，彈出對話框，在Subtract Areas對話框中輸入1，單擊OK按鈕，接著彈出對話框，在對話框中輸入2。

（7）設置單元屬性。

GUI：Main Menu>Preprocessor>Meshing>>Mesh Attributes>All Areas，彈出對話框，MAT對應項設置為1，TYPE對應項設置為1 PLANE42。

第三步，靜力學分析。

（1）進入求解器。

GUI：Main Menu>Solution

（2）選擇分析類型。

GUI：Main Menu>Solution-Analysis Type- New analysis，選中Static。

（3）施加對稱邊界條件。

GUI：Main Menu>Solution>Define Loads>Apply> Structural>Displacement>Symmetry B.C.>

On Lines，彈出對話框，選擇9和10兩條線。

（4）施加均布拉力。

GUI：Main Menu>Solution>Define Loads>Apply>Structural>Pressure>On Lines，彈出對話框，選擇線2，然后單擊OK按鈕，彈出Apply PRES on lines對話框，在value欄里輸入-1。

（5）自適應執行求解。

GUI：Main Menu>Solution>Solve>Adaptive Mesh，如圖2所示，發現自適應網格求解選項并未出現，這時候點擊Unabridged Menu如圖3所示，自適應網格求解選項就出現了。點擊Adaptive Mesh，彈出Adaptive Meshing and Solution對話框如圖4所示，在NSOLN即迭代求解次數輸入1，STARGT即結構誤差輸入2,其中2代表2%，FACMN即最小單元尺寸輸入0.05，FACMX即最大單元尺寸輸入1，點擊OK按鈕按鈕自適應網格劃分求解。

            圖2  未點擊Unabridged Menu之前              圖3  Adaptive Mesh選項出現

圖4  Adaptive Mesh自適應網格劃分求解設置

（6）退出求解器。

GUI：Main Menu>Finish

第四步，執行后處理。

（1）進入后處理器。

GUI：Main Menu>Solution>General Postproc

（2）切換結果坐標系到柱坐標下。

GUI：Main Menu>Solution>General Postproc>Options for Outp，彈出對話框，選擇全局柱坐標Global cylindric。

（3）柱坐標下觀察結果。

GUI：Main Menu>Solution>General Postproc>Plot Results>Contour Plot-Nodal Solu，彈出對話框，選擇Y-Component of stress，應力云圖如圖5所示。

圖5  柱坐標下的σ_y分布

（4）顯示單元應變能。

GUI： Main Menu>General Postproc>Plot Results>Contour Plot>Element Solu，彈出 Contour Element Solution Data對話框，選擇Energy-Strain energy，應變能云圖如圖6所示。

圖6  單元應變能云圖

（5）顯示單元能量偏差。

GUI： Main Menu>General Postproc>Plot Results>Contour Plot>Element Solu，彈出 Contour Element Solution Data對話框，選擇Error Estimation-Structural Error Energy，單元能量偏差云圖如圖7所示。

                                          圖7  單元能量偏差云圖

（6）顯示單元應力偏差。

GUI： Main Menu>General Postproc>Plot Results>Contour Plot>Element Solu，彈出 Contour Element Solution Data對話框，選擇Error Estimation-Absolute Maximum Stress Variation，單元應力偏差云圖如圖8所示。

                                            圖8  單元應力偏差云圖

結果討論：

 

                                           圖9  迭代次數為1時的單元  

                                    圖10  迭代次數為2時的單元

                                   圖11  迭代次數為3時的單元

                                   圖12  迭代次數為4時的單元

                                    圖13  迭代次數為5時的單元

圖14  迭代次數為6時的單元

表1 自適應網格不同迭代次數結果比較

迭代次數	節點數	單元數	symax	symin
1	49	34	2.827	-0.865265
2	394	347	3.094	-1.06
3	361	314	3.074	-1.062
4	365	318	3.087	-1.064
5	360	313	3.102	-1.063
6	354	307	3.101	-1.065

表2各種方法結果比較

方法	一般方法(esize=0.1)	一般方法(esize=0.05)	局部加密法	p單元法	自適應網格法
節點數	2418	8894	944	542	354
單元數	2318	8695	873	163	307
symax	2.975	3.077	3.107	3.092	3.101
symin	-.966125	-1.046	-1.063	-1.069	-1.065

從表1可以看出隨著迭代次數的增加結果趨于穩定，應力集中系數穩定在3.10左右，單元數及節點數隨著迭代次數的增加呈現波動趨勢，原因在于設定的結構誤差目標1%沒辦法滿足，這點可以從output window看出。針對這種情況可以通過減小最小單元最小尺寸來實現。

從表2可以看出p單元法節點數最少，自適應網格法單元數最小，它們均在代價最小的情況下得到了非常好的結果；一般網格加密的方法必須要整體加密才能得到較好的結果，原因在于應力梯度不大的區域加密網格效果也不明顯，只有在應力梯度較大的區域加密網格才能達到較好的效果，因此局部加密網格法應運而生，在應力梯度較大的位置加密網格，既能達到分析的精度又能避免浪費計算資源，可以說是針對應力集中不錯的一種方法。當然局部加密法的前提是要能知道哪些區域應力較大，因此在實際操作中可以先采用一般網格法，得到一個粗糙的結果，根據結果再來有目的性的局部加密網格，從而達到較好的分析效果。

綜上，在實際分析應力集中的時候可以采用p單元法、自適應網格法和局部網格加密法，這里需要說明的是針對三維實體單元p單元法以及自適應網格法計算規模可能會比較大，實體單元采用局部加密計算規模相對比較容易控制。

附件為自適應網格對于的命令流以及完整的文檔，象征性收費1元。

視頻鏈接：http://www.yqgqt.org.cn/college/video/c10819，有興趣的可以觀看，謝謝