胡昌成，趙 耀，袁 華

(華中科技大學(xué) 船舶與海洋工程學(xué)院，湖北 武漢 430074)

一種雙向曲率板成型應(yīng)變分布的力學(xué)計(jì)算方法

胡昌成，趙 耀，袁 華

(華中科技大學(xué) 船舶與海洋工程學(xué)院，湖北 武漢 430074)

在雙向曲率板線加熱成型及滾壓成型的工藝設(shè)計(jì)研究中，一項(xiàng)十分重要的工作就是依據(jù)所需目標(biāo)形狀計(jì)算應(yīng)變分布并以此確定加工路徑和工藝參數(shù)，因此應(yīng)變分布計(jì)算的精度直接影響到加工路徑和工藝參數(shù)的準(zhǔn)確性。本文介紹了一種基于力學(xué)方法的應(yīng)變分布計(jì)算方法，首先對(duì)平板施加節(jié)點(diǎn)位移場(chǎng)進(jìn)行彈塑性計(jì)算，由于存在回彈現(xiàn)象，需對(duì)節(jié)點(diǎn)位移場(chǎng)進(jìn)行反變形修正，本文將得到的形狀與目標(biāo)進(jìn)行對(duì)比后，再利用兩者偏差對(duì)位移場(chǎng)進(jìn)行迭代修正，直到兩者偏差滿足精度要求，最后將得到的應(yīng)變場(chǎng)以初應(yīng)變的形式輸入到平板模型上，用初應(yīng)變法的結(jié)果對(duì)應(yīng)變的計(jì)算方法進(jìn)行驗(yàn)證，證明該方法的精度，為后續(xù)加工路徑和工藝參數(shù)的準(zhǔn)確確定提供基礎(chǔ)。

雙向曲率板；應(yīng)變分布；力學(xué)方法；面內(nèi)應(yīng)變；彎曲應(yīng)變

0 研究背景

船舶建造等行業(yè)中常用的曲面板一般是通過機(jī)械加工或者線加熱成型的方式得到，從力學(xué)成型機(jī)理上看，曲面板的成型是將彎曲應(yīng)變和（或者）面內(nèi)應(yīng)變施加到平板上得到的，對(duì)于單向曲率板只需要彎曲應(yīng)變即可得到，對(duì)于雙向曲率板則需要彎曲應(yīng)變和面內(nèi)應(yīng)變的共同作用。在實(shí)際加工中，單向曲率板可以用三輥彎板機(jī)等設(shè)備直接加工得到，對(duì)于雙向曲率板則需要借助于線加熱成型、線滾壓成型、模壓成型等方式實(shí)現(xiàn)目標(biāo)形狀板的加工。

線加熱成型和線滾壓成型是船舶外板等雙向曲率板建造中十分常用的成型方法，其原理是通過熱源或滾輪的移動(dòng)在被加工板上逐漸產(chǎn)生塑性變形，通過將局部范圍內(nèi)的變形擴(kuò)展至整條加工路徑上最終產(chǎn)生整體的變形，其整體形狀是通過對(duì)局部變形的控制來實(shí)現(xiàn)的，由于其整體變形與加工路徑、工藝參數(shù)之間沒有明顯的關(guān)系，很難直接通過所需的目標(biāo)形狀確定加載路徑和工藝參數(shù)。為了實(shí)現(xiàn)雙向曲率板成型的自動(dòng)化，需要根據(jù)所需的目標(biāo)形狀確定對(duì)應(yīng)的加工路徑和工藝參數(shù)，大量科研工作者的研究[1 – 3]表明，線成型加工產(chǎn)生的應(yīng)變分布與加工路徑之間具有明顯的關(guān)系：面內(nèi)應(yīng)變及彎曲應(yīng)變的主應(yīng)變方向與加工路徑垂直。根據(jù)這一特性，提出了一種自動(dòng)化線成型的方案[3 – 5]：首先根據(jù)目標(biāo)形狀計(jì)算應(yīng)變分布，再根據(jù)應(yīng)變確定加工路徑的位置，最后確定相應(yīng)的工藝參數(shù)。

在上述方案中，第1步就是根據(jù)所需目標(biāo)形狀計(jì)算應(yīng)變分布，因此計(jì)算得到的應(yīng)變分布精度直接影響到后續(xù)確定加工路徑和工藝參數(shù)的準(zhǔn)確性。根據(jù)目標(biāo)形狀計(jì)算應(yīng)變分布的方法主要有力學(xué)展開法和微分幾何方法，Letcher[6]通過求解普遍的含高斯曲率的possion方程得到應(yīng)變場(chǎng)，Shin[7]強(qiáng)調(diào)了面內(nèi)應(yīng)變的作用，提出了面內(nèi)應(yīng)變的計(jì)算方法；Yu[8]通過求解一個(gè)約束非線性規(guī)劃方程得到最小應(yīng)變場(chǎng)的分布；張雪彪[9]利用參數(shù)曲面表達(dá)的方法進(jìn)行了帆形板的展開。然而Cheng[1]指出利用微分幾何方法通常只能得到面內(nèi)應(yīng)變分布，對(duì)于彎曲應(yīng)變則很難得到，而基于力學(xué)的展開方法則可以同時(shí)準(zhǔn)確得到這2種應(yīng)變的分布情況。Ueda[3 – 5]利用大變形彈性有限元法計(jì)算從初始平板到最終形狀所需的應(yīng)變，可是并沒有給出具體的求解過程；Liu[10]將目標(biāo)形狀板放置于2塊剛性體平板之間，目標(biāo)形狀板與平板之間無摩擦，兩平板之間距離不斷減小直至等于目標(biāo)形狀板的板厚，此方法得到的應(yīng)變分布與所需的應(yīng)變分布符號(hào)相反，并且沒有考慮回彈對(duì)應(yīng)變結(jié)果的影響。

本文介紹了一種基于迭代修正的力學(xué)方法進(jìn)行應(yīng)變分布的求解，通過彈塑性有限元方法計(jì)算得到應(yīng)變的分布，同時(shí)考慮回彈的影響，利用迭代修正法對(duì)目標(biāo)形狀施加反變形以減少回彈的影響，能夠快速準(zhǔn)確地計(jì)算得到成型雙向曲率板所需的應(yīng)變分布，最后以初應(yīng)變法計(jì)算得到的應(yīng)變分布所對(duì)應(yīng)的變形形狀，用初應(yīng)變法的結(jié)果對(duì)應(yīng)變分布的計(jì)算方法進(jìn)行驗(yàn)證，證明了應(yīng)變計(jì)算方法的準(zhǔn)確性。

1 總體方案

本文介紹的計(jì)算應(yīng)變分布力學(xué)方法按照?qǐng)D1所示的流程實(shí)現(xiàn)：

1）首先根據(jù)所需的目標(biāo)形狀，對(duì)目標(biāo)形狀曲面進(jìn)行三維曲面擬合，并根據(jù)擬合的結(jié)果計(jì)算有限元模型中每個(gè)節(jié)點(diǎn)的位移場(chǎng)；

2）將得到的位移場(chǎng)以節(jié)點(diǎn)位移載荷的形式施加到彈塑性模型上，由于平板模型在發(fā)生大撓度變形的同時(shí)會(huì)伴隨著板長和板寬方向的位移，因此在施加節(jié)點(diǎn)位移時(shí)只施加板厚方向的位移，在板長和板寬方向則處于自由狀態(tài)，同時(shí)對(duì)整塊板進(jìn)行剛體位移約束；

3）利用Abaqus有限元軟件進(jìn)行彈塑性計(jì)算，進(jìn)行平衡迭代計(jì)算出平板的變形，在此步驟中會(huì)計(jì)算回彈對(duì)變形的影響，得到的應(yīng)變場(chǎng)和變形場(chǎng)均為回彈后的結(jié)果；

4）將得到的變形結(jié)果與所需形狀進(jìn)行對(duì)比，由于兩模型存在著坐標(biāo)不一致等問題，本文利用基于最小二乘的原理對(duì)模型進(jìn)行空間坐標(biāo)變換，實(shí)現(xiàn)兩模型坐標(biāo)的統(tǒng)一，并計(jì)算兩者之間的總體偏差；

5）若總體偏差滿足精度要求，則直接從有限元計(jì)算結(jié)果中提取應(yīng)變結(jié)果，若不滿足精度要求，則根據(jù)步驟3中得到的變形與所需形狀之間的偏差對(duì)位移場(chǎng)施加反向變形，重新進(jìn)行步驟1～步驟4中的過程，直至偏差滿足精度要求。

2 計(jì)算實(shí)例

2.1 模型建立

以目標(biāo)形狀為鞍形的平板成型為例，對(duì)力學(xué)方法計(jì)算應(yīng)變分布的過程進(jìn)行詳細(xì)說明。初始平板的尺寸為2 000 mm×1 000 mm×20 mm，楊氏模量為2.1×1011Pa，泊松比為0.3，彈塑性屬性參照1010 steel鋼，為了完整地得到整塊板的應(yīng)變分布，利用實(shí)體單元對(duì)模型進(jìn)行網(wǎng)格劃分，板厚方向劃分4層網(wǎng)格，網(wǎng)格大小為40 mm×40 mm×5 mm，建立的有限元模型及網(wǎng)格劃分如圖2所示。

2.2 位移場(chǎng)施加

本文以鞍形板的成型為例，定義鞍形板的目標(biāo)形狀為z=x2/4–y2/9，如圖3所示。

2.3 彈塑性計(jì)算

本文使用Abaqus有限元軟件進(jìn)行彈塑性的數(shù)值計(jì)算，Abaqus軟件可以根據(jù)輸入的節(jié)點(diǎn)位移載荷自動(dòng)進(jìn)行平衡迭代求解變形結(jié)果，設(shè)置2個(gè)分析步：第1步施加位移載荷進(jìn)行彈塑性計(jì)算；第2步將位移載荷卸載，此時(shí)被加工板會(huì)由于殘余應(yīng)力的作用而發(fā)生回彈作用，在回彈后所得到的形狀與目標(biāo)形狀不同。在回彈前后板的塑性最大主應(yīng)變分布如圖5所示，從圖中可以看出，回彈前后被加工板的塑性應(yīng)變分布幾乎不變，即回彈前后2種狀態(tài)的變形場(chǎng)對(duì)應(yīng)著同一種塑性應(yīng)變狀態(tài)，分析其原因，是由于在回彈前，板中存在殘余應(yīng)變，節(jié)點(diǎn)受外部載荷作用，此時(shí)板中的應(yīng)變分布不協(xié)調(diào)，在此狀態(tài)下的應(yīng)變分布所對(duì)應(yīng)的形狀是回彈后的板形狀，因此若根據(jù)此塑性應(yīng)變進(jìn)行加工路徑和工藝參數(shù)的確定，則與目標(biāo)形狀之間存在著偏差。

2.4 結(jié)果對(duì)比

圖6為彈塑性計(jì)算得到的回彈后形狀與目標(biāo)形狀對(duì)比結(jié)果圖，虛線表示回彈后的變形形狀，實(shí)線表示目標(biāo)形狀，從圖中可以看出由于回彈的影響，導(dǎo)致所得到的形狀與目標(biāo)形狀存在著較大的差別，同時(shí)彈塑性模型的坐標(biāo)與目標(biāo)形狀的坐標(biāo)不在同一個(gè)坐標(biāo)系中，所以在進(jìn)行結(jié)果對(duì)比之前必須進(jìn)行坐標(biāo)匹配，本文使用最鄰近點(diǎn)迭代（Iterative Closet Point，ICP）方法[11]進(jìn)行坐標(biāo)匹配，在完成坐標(biāo)匹配之后，使用式（1）進(jìn)行總體偏差的計(jì)算：

其中：σ為總體偏差，n為模型中性面上節(jié)點(diǎn)的數(shù)量，Di為第i個(gè)節(jié)點(diǎn)距離目標(biāo)形狀的最小距離。

2.5 反變形迭代修正

由于回彈的原因造成變形結(jié)果與目標(biāo)形狀的不一致，使用反變形法進(jìn)行節(jié)點(diǎn)位移場(chǎng)的修正，具體方法是利用回彈后變形結(jié)果與目標(biāo)形狀之間的差別預(yù)估反變形量，將反變形量施加到節(jié)點(diǎn)位移場(chǎng)上，以此與回彈量相抵消，減小回彈后形狀與目標(biāo)形狀之間的差別。

根據(jù)式（1）計(jì)算得到總體偏差為σi，若σi大于最大允許偏差σmax，則需要進(jìn)行下一步迭代，則在第（i+1）次迭代中，反變形形狀為：

其中k為修正系數(shù)，其目的是保證每一次迭代的收斂性，利用反變形對(duì)節(jié)點(diǎn)位移場(chǎng)進(jìn)行修正，修正后的節(jié)點(diǎn)位移場(chǎng)為：

利用以上迭代算法對(duì)圖3所示的鞍形目標(biāo)形狀進(jìn)行計(jì)算，設(shè)置最大允許偏差σmax=2×10–4，利用程序進(jìn)行迭代計(jì)算，迭代4次后滿足精度要求，每次迭代所得到的總體偏差如表1所示。

表 1 迭代次數(shù)及總體偏差大小Tab. 1 The number of iterations and the total deviation

將每次迭代計(jì)算所得到的形狀與目標(biāo)形狀進(jìn)行對(duì)比，以圖2中縱向中心線Line1和橫向中心線Line2的變形結(jié)果進(jìn)行表示，如圖8所示，可以看出隨著迭代次數(shù)的增加，成型所得到的形狀逐漸逼近目標(biāo)形狀，在通過4次迭代計(jì)算之后，所得到的形狀與所需形狀之間的差別已經(jīng)很小，可以滿足精度要求。完成4次迭代后得到的塑性應(yīng)變分布如圖9所示。

3 計(jì)算結(jié)果驗(yàn)證

根據(jù)應(yīng)變與變形之間的關(guān)系，若已知應(yīng)變分布，就可以計(jì)算出該應(yīng)變分布所對(duì)應(yīng)的變形形狀，初應(yīng)變方法是將應(yīng)變分布輸入到模型上計(jì)算該應(yīng)變分布所對(duì)應(yīng)變形的方法，其可以通過簡單彈性有限元計(jì)算得到結(jié)構(gòu)的應(yīng)力和變形，利用這一方法，可以將迭代計(jì)算得到的塑性應(yīng)變以初應(yīng)變的形式輸入到平板模型上，進(jìn)行彈性計(jì)算得到該應(yīng)變所對(duì)應(yīng)的形狀，通過將得到的形狀與所需形狀進(jìn)行對(duì)比，進(jìn)行精度的驗(yàn)證。

將圖9中塑性應(yīng)變分布輸入到平板模型上，進(jìn)行彈性有限元計(jì)算，得到的形狀與所需形狀對(duì)比如圖10（a）所示，兩者撓度之差如圖10（b）所示，在圖10（a）中由于兩形狀區(qū)別非常小，所以得到的2個(gè)三維圖已經(jīng)重合，從圖（b）中可以看出可以兩形狀的最大偏差僅0.6 mm，說明經(jīng)過4次迭代計(jì)算得到的應(yīng)變分布所對(duì)應(yīng)的形狀與所需的形狀一致，因此可以證明本文介紹的應(yīng)變分布的計(jì)算方法精確。

4 結(jié) 語

本文利用彈塑性有限元方法計(jì)算成型雙向曲率目標(biāo)形狀板所需的應(yīng)變分布，并利用迭代法對(duì)結(jié)果進(jìn)行修正從而提高計(jì)算精度，通過本文研究得出以下結(jié)論：

1）以節(jié)點(diǎn)位移的形式可以將目標(biāo)形狀作為位移載荷輸入到有限元模型中，進(jìn)行彈塑性計(jì)算可以得到應(yīng)變和變形結(jié)果；

2）使用反變形迭代修正法可以減少回彈對(duì)變形的影響，提高應(yīng)變分布的精度，隨著迭代次數(shù)的增加，回彈后的形狀逐步逼近于所需的形狀，計(jì)算得到的應(yīng)變場(chǎng)精度也隨之提高；

3）通過初應(yīng)變方法對(duì)計(jì)算得到的應(yīng)變分布進(jìn)行彈性計(jì)算可以得到應(yīng)變分布對(duì)應(yīng)的變形，結(jié)果表明本文介紹的方法可以精確地計(jì)算得到所需形狀的應(yīng)變場(chǎng)，為后續(xù)加工路徑和工藝參數(shù)的確定提供了基礎(chǔ)。

[ 1 ]JIN C, Y L YAO. Process design of laser forming for threedimensional thin plates [J]. Transactions of the ASME. Journal of Manufacturing Science and Engineering, 2004. 126(2):217–25.

[ 2 ]XU Zhao-kang. Discuss on forming method for complex curved plate and its applicability [J]. Journal of Wuhan Institute of Shipbuilding Technology, 2006(03): 9–12.

[ 3 ]UEDA Y. Development of computer-aided process planning system for plate bending by line heating (report 1): relation between the final form of plate and the inherent strain(machanics, strength & structural design) [J].Transactions of JWRI, 1991. 20(2): 275–285.

[ 4 ]UEDA Y. Development of computer-aided process planning system for plate bending by line heating (report 2): practice for plate bending in shipyard viewed from aspect of inherent strain[J]. Journal of Ship Production, 1994. 10(4): 239–239.

[ 5 ]UEDA Y. Development of computer-aided process planning system for plate bending by line heating (report 3): relation between heating condition and deformation [J]. Journal of Ship Production, 1994. 10(4): 248–248.

[ 6 ]LETCHER J S. Lofting and fabrication of compound-curved plates [J]. Journal of Ship Research, 1993. 37(2): 166–175.

[ 7 ]SHIN J G, C H RYU. Nonlinear kinematic analysis of the deformation of plates for ship hull fabrication [J]. Journal of Ship Research, 2000. 44(3): 270–277.

[ 8 ]YU G, N M PATRIKALAKIS, T MAEKAWA. Optimal development of doubly curved surfaces [J]. Computer Aided Geometric Design, 2000. 17(6): 545–577.

[ 9 ]ZHANG Xue-biao. Research on ship-hull curved plate forming by pure heating [D]. Dalian: Dalian University of Trchnology,2006.

[10]LIU C, Y L YAO. FEM-based process design for laser forming of doubly curved shapes [J]. Journal of Manufacturing Processes, 2005. 7(2): 109–121.

[11]ZHANG Xue-chang, XI Jun-tong, YAN Jun-qi. Study of inspection technologies based on point cloud for complex surfaces and its system development [D]. Shanghai: Computer Integrated Manufacturing System, 2005(5): 770–774, 780.

Calculating strain field to forming a double-curved shape based on mechanic theory

HU Chang-cheng, ZHAO Yao, YUAN Hua
(School of Naval Architecture and Ocean Engineering Huazhong University of Science and Technology, Wuhan 430074, China)

The determination of process paths and parameters according to the strain field that required to forming a double-curved desired shape is an important work during the process design of line-heating and rolling forming in the process design of line heating and rolling forming. Therefore, the precision of Calculated strain field will affect the accuracy of process paths and parameters. This paper presents a method to calculate the strain filed required to develop a desired shape to a planar shape depend on mechanic theory. Firstly, nodal displacements is applied to the flat plate via elastic-plastic finite element analysis. Owing to the spring back, the shape calculated by elastic-plastic finite element analysis is different between the required shape, then the deviation between calculated shape and required shape is calculated, and reverse deformation is added to nodal displacements, this process continues until the deviation can meet the requirement of precision to terminate the loop. Finally, the obtained strain field is applied to flat model as initial strain to verify the accuracy of this method, and results shows that it has very high calculation precision, this method provide the foundation of the determination of process paths and parameters.

double-curved plate；strain field；mechanics methods；in-plane strain；bending strain

U671.3

A

1672 – 7649(2017)09 – 0006 – 06

10.3404/j.issn.1672 – 7649.2017.09.002

2016 – 11 – 03；

2016 – 12 – 15

胡昌成(1988 – )，男，博士研究生，主要從事板成型和非線性有限元方法方面的研究工作。