大規(guī)模水閘非線性有限元并行計(jì)算研究

張向陽，秦忠國(河海大學(xué)工程力學(xué)系，江蘇南京210098)*

大規(guī)模水閘非線性有限元并行計(jì)算研究

張向陽，秦忠國
(河海大學(xué)工程力學(xué)系，江蘇南京210098)*

借助開源的PETSc并行數(shù)值算法包，使用C ++編寫了一個基于消息傳遞模式的非線性彈性有限元的并行計(jì)算程序.程序通過多次調(diào)用PETSc中的求解器組件KSP完成了Duncan-Chang雙曲線模型的并行有限元分析.最后通過一個實(shí)際工程水閘模型進(jìn)行測試，結(jié)果表明多線程的計(jì)算效率比單線程的計(jì)算效率要高，四個進(jìn)程的加速比達(dá)到了1.91.同時也表明，通過借助于成熟的數(shù)學(xué)工具箱進(jìn)行非線性有限元并行程序開發(fā)是可行的，有效的.

非線性有限元;增量法;并行程序; PETSc

0 引言

由于工程技術(shù)和工業(yè)制造領(lǐng)域的不斷進(jìn)步和發(fā)展，在預(yù)測模型的構(gòu)造和模擬［1］、大型實(shí)際工程的計(jì)算與數(shù)值仿真、計(jì)算流體動力學(xué)、全球氣象預(yù)報［2］、圖像處理、非線性的動力學(xué)系統(tǒng)等領(lǐng)域中都對計(jì)算的精度和速度提出了極高的要求.對這些重大的科學(xué)計(jì)算問題進(jìn)行非線性有限元數(shù)值仿真時，模型網(wǎng)格經(jīng)常達(dá)到幾百萬甚至更多，計(jì)算量大且時間較長，再加上計(jì)算機(jī)內(nèi)存容量的限制，工程分析人員不得不縮小計(jì)算規(guī)模，導(dǎo)致計(jì)算結(jié)果的可靠性降低.當(dāng)前，采用并行有限單元法［3］可以很好的解決這些問題，即提高大型工程的求解規(guī)模和加快求解速度.

隨著并行計(jì)算的發(fā)展，一些優(yōu)秀的并行算法包引起了科研人員的關(guān)注.這些并行算法包提供并行的數(shù)據(jù)結(jié)構(gòu)和成熟的算法及標(biāo)準(zhǔn)的調(diào)用接口.基于所編譯得到的可移植可擴(kuò)展科學(xué)計(jì)算工具箱PESTc( Portable Extensible Toolkits for Scientific Computing)［4］，使用C ++語言編寫了面向?qū)ο蟮姆蔷€性彈性有限元分級加載的并行求解程序.程序采用Duncan-Chuang雙曲線本構(gòu)模型，其求解方法為增量法中的初始剛度法，并用實(shí)際工程算例進(jìn)行測試.

1 并行有限元的設(shè)計(jì)及實(shí)現(xiàn)

1.1科學(xué)計(jì)算工具箱PESTc的介紹

PETSc是由美國Argonne國家實(shí)驗(yàn)室Staish Balay等人開發(fā)的科學(xué)實(shí)驗(yàn)計(jì)算工具箱.PETSc的數(shù)據(jù)結(jié)構(gòu)主要包括向量Vec、矩陣Mat、結(jié)構(gòu)網(wǎng)格DA等［5］.PETSc工具箱中的方程組求解器組件主要有三種:線性方程求解器( KSP)、非線性求解器( SNES)、時間步進(jìn)求解器( TS).其中KSP是關(guān)鍵，SNES和TS都是以此為基礎(chǔ)的.基于這三個求解器，用戶可以開發(fā)自己的應(yīng)用程序.

1.2有限元并行化的策略

并行計(jì)算要解決計(jì)算時存儲空間的不足和問題求解時間長的兩大難題［11］.目前分布式環(huán)境的并行策略有基于系統(tǒng)方程組求解的并行有限元、并行EBE法、并行子結(jié)構(gòu)法、結(jié)構(gòu)非線性分析與動力分析的并行算法等.其中有限元并行子結(jié)構(gòu)法的基本思想是:將大規(guī)模的問題分解為若干個小問題，由各個處理器負(fù)責(zé)一個或多個小問題的計(jì)算［6］，計(jì)算過程中處理器之間進(jìn)行必要的數(shù)據(jù)交換，最后由各個處理器協(xié)調(diào)工作得出整體解［7］.

該并行程序設(shè)計(jì)的方法是:從讀入數(shù)據(jù)后，到得到非線性問題的位移解的過程，都完成了并行化處理，使用并行子結(jié)構(gòu)法的思想［12］，但是不形成子結(jié)構(gòu)，而是把結(jié)構(gòu)的總體數(shù)據(jù)信息劃分為大概相同的子塊，分配到各個進(jìn)程當(dāng)中去，各個進(jìn)程預(yù)先估算所需要的存儲空間，并創(chuàng)建存儲空間.組集總體剛度矩陣時，每個進(jìn)程當(dāng)中都分發(fā)一份單元結(jié)點(diǎn)總體信息，這樣就消除了在形成總體剛度矩陣時，各個進(jìn)程之間的通信壓力.程序的并行處理主要針對初始應(yīng)力場計(jì)算和分級加載計(jì)算過程，這兩個過程有明顯的界限，但是都被并行化處理，最后得到位移、應(yīng)力解在進(jìn)程中保存結(jié)果.

1.3并行程序的分析及設(shè)計(jì)

非線性彈性是由本構(gòu)關(guān)系的非線性引起，即應(yīng)力σ與應(yīng)變ε矢量之間不再是線性關(guān)系.用增量形式寫出其本構(gòu)關(guān)系［8］:

物理方程

式中: DT表示為切線彈性矩陣.

圖1 程序的總體流程圖

再考慮到非線性彈性有限元的幾何方程和平衡方程可得:

式中: KT( u)是切線剛度矩陣，它是節(jié)點(diǎn)位移u的函數(shù); ce是選擇矩陣.

程序采用Duncan-Chang雙曲線本構(gòu)模型，要考慮初始應(yīng)力場的計(jì)算和荷載的分級加載求解，其程序計(jì)算流程如圖1所示［9］.程序從讀入數(shù)據(jù)到得到位移解的過程中均實(shí)現(xiàn)并行化處理，特別是在初始應(yīng)力場計(jì)算、分布加載求解步驟上予以高度并行化處理.因?yàn)檫@兩個部分所耗時間占整個程序的耗時比重比較大，通常能占到整個程序運(yùn)行時間的百分之八十［10］，所以該程序的并行效率是比較高的.

2 算例分析

2.1計(jì)算模型及運(yùn)行環(huán)境

對某水閘系統(tǒng)進(jìn)行三維非線性有限元的靜力分析，工況為水閘系統(tǒng)的完建期，所受荷載為重力荷載.通過此算例主要測試該求解程序的并行加速效果，測試環(huán)境為一臺普通的32位臺式機(jī)，操作系統(tǒng)為Windows 7，處理器為四核Intel( R) Core( TM) i3 －2120CPU，主頻3.3 GHz，可使用內(nèi)存3.03 GB.

水閘系統(tǒng)由鋼筋混凝土閘身和軟土地基兩部分組成.其中閘身的材料特性如下:彈性模量為2.7e10 Pa，泊松比為0.167，密度為2 500 kg/m3.地基材料參數(shù)采用Duncan-Chang雙曲線模型來確定，其中K =3.5e3，表示初始彈性模量的基數(shù); n取0.5，反應(yīng)E隨σ增長的關(guān)系;φ=30°; c = 0.3 kPa; Rf為破壞比，值為0.8; G為初始泊松比的基數(shù)，值為0.3; F為泊松比隨圍壓增長降低的速率，值為0.01; D為泊松比隨應(yīng)變的增長關(guān)系，值為10.0; Ku為土體卸載的參數(shù)，值為4.2e3.將水閘系統(tǒng)模型離散為101 748個單元(八節(jié)點(diǎn)六面體單元)，總節(jié)點(diǎn)數(shù)為112 670個.

2.2計(jì)算結(jié)果及分析

水閘系統(tǒng)的有限元模型如圖2所示，對上述有限單元網(wǎng)格模型進(jìn)行1～4個進(jìn)程的非線性有限元靜力計(jì)算.首先計(jì)算初始應(yīng)力場，然后將荷載平均等分為10級來施加，其計(jì)算結(jié)果如表1所示，取得了明顯的并行加速效果，其四個進(jìn)程計(jì)算的加速比和效率達(dá)到1.91和0.48.

圖2 水閘系統(tǒng)的有限元模型

表1 1－4進(jìn)程并行計(jì)算結(jié)果

計(jì)算過程中，程序的兩個主要并行部分初始應(yīng)力場計(jì)算和分級加載計(jì)算過程的運(yùn)行時間，加速比和效率如表2所示，并繪成折線圖，如圖3、4所示，其加速效果也較為明顯.

表2 兩個主要部分的并行計(jì)算時間及加速效果

圖3 進(jìn)程個數(shù)隨時間變化的關(guān)系圖

圖4 初始應(yīng)力場及分級加速計(jì)算的關(guān)系圖

對于用有限單元法分析一個實(shí)際問題時，方程組的求解耗時通常會占到百分之八十的比重.在上述計(jì)算過程中，每次求解都包含了11次的線性方程組的求解，即初始應(yīng)力場的計(jì)算和10次加載計(jì)算.得到求解方程的相關(guān)信息如表3所示，可見解方程的耗時能達(dá)到整個求解過程時間的86%以上，該程序采用PETSc工具箱中的并行求解器，取得了較好的效率和明顯的加速效果.

表3 解方程占整個求解時間的信息

對該水閘系統(tǒng)有限元網(wǎng)格模型進(jìn)行不同數(shù)目的分級加載計(jì)算，計(jì)算進(jìn)程數(shù)都為4，得到的計(jì)算結(jié)果如表4所示.由每次計(jì)算得到的位移最大值可知，計(jì)算的結(jié)果是符合初始剛度法對于Duncan-Chang模型的結(jié)果規(guī)律，即分級數(shù)目越多，計(jì)算結(jié)果從上方趨于精確結(jié)果.但要得到比較精確的結(jié)果，還需要更大的計(jì)算量，也從側(cè)面說明并行計(jì)算對于這樣非線性問題的精確求解的重要性.

表4 不同荷載步時的計(jì)算總耗時和位移最大值

3 結(jié)論

通過借助于PETSc提供的高層應(yīng)用開發(fā)平臺，使用面向?qū)ο笳Z言C ++，編寫了一個非線性彈性有限元的并行計(jì)算程序，并對該程序進(jìn)行了測試.結(jié)果顯示，該程序取得了明顯的并行加速效果，其四個進(jìn)程計(jì)算的加速比和效率分別達(dá)到了1.91和0.48.同時也表明借助于PETSc并行數(shù)值算法包來開發(fā)非線性有限元并行程序，可以不用考慮并行數(shù)據(jù)的分布存儲結(jié)構(gòu)和復(fù)雜的通信問題，降低了因程序并行化帶來的編程困難，縮短了并行程序的開發(fā)周期，提高了程序的實(shí)用性.

［1］KARTAM N，F(xiàn)LOOD I.Construction simulation using parallel computing environments［J］.Automation in Construction，2000，10( 1) : 69-78.

［2］XIE J，YANG C，ZHOU B，et al.High-performance computing for the simulation of dust storms［J］.Computers，Environment and Urban Systems，2010，34 ( 4) : 278-290.

［3］BABU?KA I，ELMAN H C.Some aspects of parallel implementation of the finite-element method on message passing architectures［J］.Journal of Computational and Applied Mathematics，1989，27( 1-2) : 157-187.

［4］張健飛.采用PETSc的有限元并行計(jì)算實(shí)現(xiàn)與優(yōu)化［J］.計(jì)算機(jī)工程與應(yīng)用，2010，46( 10) : 57-59.

［5］高懷玉，張峰，秦忠國.基于PETSc的有限元高性能求解方法［J］.河海大學(xué)學(xué)報(自然科學(xué)版)，2013，41 ( 4) : 365-370.

［6］張健飛，姜弘道.有限元子結(jié)構(gòu)并行算法的性能分析［J］.力學(xué)與實(shí)踐，2002，24( 5) : 35-37.

［7］胡寧，張汝清.有限元子結(jié)構(gòu)并行算法的效能分析［J］.應(yīng)用力學(xué)學(xué)報，1993，10( 1) : 124-128.

［8］任青文.非線性有限單元法［M］.南京:河海大學(xué)出版社，2000: 7-30.

［9］姜弘道，張健飛，秦忠國，等.工程科學(xué)中的高性能計(jì)算［M］.北京:科學(xué)出版社，2013.

［10］MICHAEL J Q.MPI與OpenMP并行程序設(shè)計(jì)［M］.北京:清華大學(xué)出版社，2004: 34-40.

［11］姜弘道.水利高性能計(jì)算的進(jìn)展［J］.水利水電科技進(jìn)展，2006，26( 2) : 70-76.

［12］謝旻.高可用MPI并行編程環(huán)境及并行程序開發(fā)方法的研究與實(shí)現(xiàn)［D］.長沙:國防科學(xué)技術(shù)大學(xué)，2007.

Study of Nonlinear Finite Element Parallel Computing based on Large-Scale Lock

ZHANG Xiangyang，QIN Zhongguo
( Department of Engineering Mechanics，Hehai University，Nanjing 210098，China)

With the help of the open source PETSc parallel numerical algorithm package，C ++ is used to write a nonlinear elastic finite element parallel computing program based on the model of message passing.The procedure is completed in parallel finite element analysis of Duncan-Chang hyperbola model by the KSP solver component repeatedly calls in PETSc.Through a practical engineering of sluice model to test，the result show that the calculation efficiency of multi thread is higher than the calculation efficiency of single thread，and the accelerate ratio of the four processes reach 1.91.It also show that using the mature mathematical toolbox to development nonlinear finite element parallel program is feasible and effective.

nonlinear finite element; incremental method; parallel program; PETSc

A

1673-9590( 2016) 01-0040-04

2015-03-16

張向陽( 1991－)，男，碩士研究生;

秦忠國( 1961－)，男，教授，博士，主要從事計(jì)算力學(xué)及工程仿真、水工軟件開發(fā)、高性能計(jì)算的研究

E-mail: 339681123@ qq.com.