PFCD中文講解例題

1、通過例子學(xué)習(xí)幾種常見命令例二：PFC3D目錄下的GuideStartfooting.dat，程序如下（注意：本程序與原程序不同，特加了plot set rotation (24.0,352.0,340.0)這句以使得模型轉(zhuǎn)過一定角度，讓用戶看的更清楚，另外還將原程序一分為二，主要是為了更好地學(xué)習(xí)每個(gè)命令的作用），以下是程序代碼。;fname: footing.DAT (tutorial example for PFC3D)new ; clear program state to begin new problemset random ; reset random-number generato

2、rtitle 'Tutorial Example'wall id=1 face ( 0, 0, 0) (10, 0, 0) (10, 0 -5) ( 0, 0,-5)wall id=2 face ( 0, 0, 0) ( 0, 0,-5) ( 0, 2,-5) ( 0, 2, 0)wall id=3 face (10, 0,-5) (10, 0, 0) (10, 2, 0) (10, 2,-5)wall id=4 face (10, 2,-5) (10, 2, 0) ( 0, 2, 0) ( 0, 2,-5)wall id=5 face ( 0, 0,-5) (10, 0,-5

3、) (10, 2,-5) ( 0, 2,-5)wall id=6 face ( 0, 0, 0) ( 0, 2, 0) (10, 2, 0) (10, 0, 0)gen id=1,750 rad 0.15,0.20 x=0,10 y=0,2 z=-5,0;plot create Footingplot set title text 'Basic collection of particles in a box'plot set rotation (24.0,352.0,340.0)plot add ball yellowplot add wall whiteplot add a

4、xes brownplot showpause;wall id=1 kn=1e8 ks=1e8wall id=2 kn=1e8 ks=1e8wall id=3 kn=1e8 ks=1e8wall id=4 kn=1e8 ks=1e8wall id=5 kn=1e8 ks=1e8wall id=6 kn=1e8 ks=1e8;prop density 2000 kn 1e8 ks 1e8ini rad mul 1.47;set hist_rep=5hist ball zvel 3 1 0hist diagnostic mufset dt dscalecycle 3000以下是程序中逐條命令的解釋

5、：1、我們首先定義了模型的邊界，使用WALL命令，如下：wall id=1 face ( 0, 0, 0) (10, 0, 0) (10, 0 -5) ( 0, 0,-5)wall id=2 face ( 0, 0, 0) ( 0, 0,-5) ( 0, 2,-5) ( 0, 2, 0)wall id=3 face (10, 0,-5) (10, 0, 0) (10, 2, 0) (10, 2,-5)wall id=4 face (10, 2,-5) (10, 2, 0) ( 0, 2, 0) ( 0, 2,-5)wall id=5 face ( 0, 0,-5) (10, 0,-5) (10

6、, 2,-5) ( 0, 2,-5)wall id=6 face ( 0, 0, 0) ( 0, 2, 0) (10, 2, 0) (10, 0, 0)這些命令創(chuàng)造了一個(gè)由六個(gè)墻面圍成的箱子，每個(gè)墻面都賦予了一個(gè)獨(dú)一無二的ID號(hào)。每個(gè)墻面四邊形通過給出四邊形的四個(gè)頂點(diǎn)坐標(biāo)來定義，即四個(gè)頂點(diǎn)連成四條線，然后再圍成一個(gè)四邊形墻面。注意：頂點(diǎn)的順序決定了墻面是否被激活，激活的墻面可以自動(dòng)識(shí)別與小球的接觸，不激活的話就會(huì)被小球穿透，2、然后在箱子里填充小球，用Generate命令，也可簡(jiǎn)寫，如下：gen id=1,750 rad 0.15,0.20 x=0,10 y=0,2 z=-5,0這個(gè)命令規(guī)定了

7、小球填充的空間范圍，小球的大小以及填充的總個(gè)數(shù)。Id=1,750的含義是：ID號(hào)從1到750其實(shí)就生成了750個(gè)小球，并給每個(gè)小球指定了各自的ID號(hào)。Rad 0.15,0.20的含義是：小球半徑從0.15m到0.20m隨機(jī)選擇x=0,10 y=0,2 z=-5,0的含義是：小球填充的空間是0<x<10，0<y<2，-5<z<0的一個(gè)立方體空間。如果填充的空間太小或者小球的半徑設(shè)置的太大，以至于這個(gè)空間容不下這么多的小球，那系統(tǒng)會(huì)自動(dòng)減少生成的球的個(gè)數(shù)，并且系統(tǒng)會(huì)提示實(shí)際生成的球的個(gè)數(shù)。3、模型的顯示用命令plot，程序如下：plot create Footi

8、ngplot set title text 'Basic collection of particles in a box'plot set rotation (24.0,352.0,340.0)plot add ball yellowplot add wall whiteplot add axes brownplot show從上面命令可以看出plot的用法規(guī)律，即與plot顯示相關(guān)的命令的前面都用plot引出，然后后面再加相關(guān)的命令，如create Footing的意思就是創(chuàng)建一個(gè)名字為Footing的可視化窗口（這個(gè)在菜單欄點(diǎn)擊“window”，下拉菜單可發(fā)現(xiàn)多了一項(xiàng)vi

9、ew 1Footing，用戶可自己看下）；set title text是設(shè)定標(biāo)題， set rotation 是設(shè)定模型的角度，以方便觀察；add ball yellow是設(shè)定球的顏色為黃色；add wall white是設(shè)定墻面（也就是邊界面）為白色；add axes brown是設(shè)定坐標(biāo)軸為棕色； show是顯示模型的命令。4、定義墻面的法向剛度（normal stiffness）和切向剛度（shear stiffness）wall id=1 kn=1e8 ks=1e8wall id=2 kn=1e8 ks=1e8wall id=3 kn=1e8 ks=1e8wall id=4 kn=1e

10、8 ks=1e8wall id=5 kn=1e8 ks=1e8wall id=6 kn=1e8 ks=1e8其中kn表示法向剛度，ks表示切向剛度，單位為N/m5、定義球的密度，法向剛度和切向剛度prop density 2000 kn 1e8 ks 1e8使用命令property，可簡(jiǎn)寫為prop，定義密度density為2000kg/m3，法向剛度和切向剛度均為108N/m6、改變小球半徑ini rad mul 1.47ini應(yīng)該是initial的簡(jiǎn)寫，用于初始化小球半徑（有待考證）rad表示半徑mul是multiply的簡(jiǎn)寫，是倍數(shù)的意思。1.47是倍數(shù)因子，表示現(xiàn)在半徑是原來半徑的1.

11、47倍。7、設(shè)定監(jiān)測(cè)變量和采樣間隔時(shí)間set hist_rep=5hist ball zvel 3 1 0hist diagnostic mufhist是history歷史記錄是簡(jiǎn)寫，set hist_rep=5表示每隔5個(gè)運(yùn)算步記錄一次變量的值（默認(rèn)的是10個(gè)運(yùn)算步記錄一次）；本程序設(shè)定了兩個(gè)監(jiān)測(cè)對(duì)象，一個(gè)是hist ball zvel 3 1 0，表示的是中心坐標(biāo)最接近x=3，y=1，z=0的球的z方向的速度被監(jiān)測(cè)和記錄（在這個(gè)例子中，282號(hào)球最接近這個(gè)坐標(biāo)，因此作為監(jiān)測(cè)球）；另一個(gè)是hist diagnostic muf，是用來監(jiān)測(cè)系統(tǒng)的不平衡力，主要用于靜力學(xué)分析，當(dāng)這個(gè)值趨于很小時(shí)

12、，表示系統(tǒng)達(dá)到穩(wěn)定狀態(tài)。8、設(shè)定求解選項(xiàng)set dt dscale因?yàn)檫@是一個(gè)靜態(tài)分析，因此我們開啟density-scaling選項(xiàng)（簡(jiǎn)寫為dscale）來優(yōu)化求解效率，dt應(yīng)該是這個(gè)選項(xiàng)中的某個(gè)模式（有待考證）。9、設(shè)定求解步數(shù)cycle 3000當(dāng)寫入這句程序后，計(jì)算就開始了，計(jì)算的步數(shù)，時(shí)間步，總時(shí)間，不平衡力平均值和最大不平衡力這些信息均顯示在屏幕上，并且每隔5個(gè)計(jì)算步更新一次（這是set hist_rep=5的作用），圖像動(dòng)畫每隔20個(gè)計(jì)算步更新一次（默認(rèn)），這個(gè)更新時(shí)間可用命令set pinterval來改變。另外，阻尼已經(jīng)在運(yùn)動(dòng)方程中被默認(rèn)設(shè)定，用于吸收振動(dòng)能量，更快地達(dá)到平衡

13、狀態(tài)。10、查看結(jié)果求解完畢后，用戶可輸入命令plot his 1來查看歷史記錄中第一個(gè)變量的數(shù)值曲線。同理可輸入命令plot his 2來查看歷史記錄中第二個(gè)變量的數(shù)值曲線。由曲線可知系統(tǒng)在cycle=1000時(shí)已趨于穩(wěn)定，故可以將之前的3000步改為1000步，以減少計(jì)算時(shí)間（個(gè)人愚見，僅供參考）。先對(duì)plot補(bǔ)充點(diǎn)說明：1、plot create *命令可用于創(chuàng)建多個(gè)圖像顯示窗口，但每次只能看一種結(jié)果，可以通過菜單欄里的window選項(xiàng)卡來切換窗口，查看想看的結(jié)果2、如只想顯示小球，則使用命令plot ball，如只想顯示墻面，則使用命令plot wall，如二者都想顯示則plot ba

14、ll wall，如運(yùn)行完畢后想顯示別的結(jié)果，如接觸力則plot cforce（contact force的簡(jiǎn)寫）?？偠灾?，plot后面加你想顯示的東東即可，用戶需靈活使用。此外，程序中的注釋用英文的分號(hào); ，后面加注釋的內(nèi)容即可。第二天的講解中小球都是自由狀態(tài)，無外力的作用，因此有必要討論施加外力的情況。（1） 對(duì)小球施加重力作用;fname: footing.DAT (tutorial example for PFC3D)new ; clear program state to begin new problemset random ; reset random-number genera

15、tortitle 'Tutorial Example'wall id=1 face ( 0, 0, 0) (10, 0, 0) (10, 0 -5) ( 0, 0,-5)wall id=2 face ( 0, 0, 0) ( 0, 0,-5) ( 0, 2,-5) ( 0, 2, 0)wall id=3 face (10, 0,-5) (10, 0, 0) (10, 2, 0) (10, 2,-5)wall id=4 face (10, 2,-5) (10, 2, 0) ( 0, 2, 0) ( 0, 2,-5)wall id=5 face ( 0, 0,-5) (10, 0,

16、-5) (10, 2,-5) ( 0, 2,-5)wall id=6 face ( 0, 0, 0) ( 0, 2, 0) (10, 2, 0) (10, 0, 0)gen id=1,750 rad 0.15,0.20 x=0,10 y=0,2 z=-5,0;plot create Contact_Forceplot set title text 'Model at equilibrium'plot add ball yellowplot add wall whiteplot add cforce blackplot add axes brownplot show;wall i

17、d=1 kn=1e8 ks=1e8wall id=2 kn=1e8 ks=1e8wall id=3 kn=1e8 ks=1e8wall id=4 kn=1e8 ks=1e8wall id=5 kn=1e8 ks=1e8wall id=6 kn=1e8 ks=1e8;prop density 2000 kn 1e8 ks 1e8ini rad mul 1.47;set hist_rep=5hist ball zvel 3 1 0hist diagnostic mufset dt dscale;set grav 0 0 -9.81prop fric 1.0;cyc 1000此程序與第二天的程序幾乎

18、相同，唯一的不同是多加了兩句程序：set grav 0 0 -9.81prop fric 1.0第一句是設(shè)定重力加速度為9.81m/s2，方向?yàn)閆軸負(fù)方向，grav是重力gravity的簡(jiǎn)寫；第二句是設(shè)定球與球之間的摩擦系數(shù)，fric是friction的簡(jiǎn)寫，在該程序中設(shè)定摩擦系數(shù)為1.0。運(yùn)行結(jié)果如下圖所示：（2） 將上頂面換做一小塊下降的板，壓縮箱子里的小球程序如下：;fname: footing.DAT (tutorial example for PFC3D)new ; clear program state to begin new problemset random ; reset

19、random-number generatortitle 'Tutorial Example'wall id=1 face ( 0, 0, 0) (10, 0, 0) (10, 0 -5) ( 0, 0,-5)wall id=2 face ( 0, 0, 0) ( 0, 0,-5) ( 0, 2,-5) ( 0, 2, 0)wall id=3 face (10, 0,-5) (10, 0, 0) (10, 2, 0) (10, 2,-5)wall id=4 face (10, 2,-5) (10, 2, 0) ( 0, 2, 0) ( 0, 2,-5)wall id=5 fac

20、e ( 0, 0,-5) (10, 0,-5) (10, 2,-5) ( 0, 2,-5)wall id=6 face ( 0, 0, 0) ( 0, 2, 0) (10, 2, 0) (10, 0, 0)gen id=1,750 rad 0.15,0.20 x=0,10 y=0,2 z=-5,0;plot create Contact_Forceplot set title text 'Model at equilibrium'plot add ball yellowplot add wall blueplot show;wall id=1 kn=1e8 ks=1e8wall

21、 id=2 kn=1e8 ks=1e8wall id=3 kn=1e8 ks=1e8wall id=4 kn=1e8 ks=1e8wall id=5 kn=1e8 ks=1e8wall id=6 kn=1e8 ks=1e8;prop density 2000 kn 1e8 ks 1e8ini rad mul 1.47;set hist_rep=5hist ball zvel 3 1 0hist diagnostic mufset dt dscale;set grav 0 0 -9.81prop fric 1.0;delete wall 6wall id 10 face (0,0,0) (0,2

22、,0) (2,2,0) (2,0,0)wall id 10 ks 1e8 kn 1e8 fric 1.0wall id 10 zvel=-0.5e-3property n_bond=5e5 s_bond=5e5hist wall zforce id 10cyc 6000運(yùn)行結(jié)果如圖所示與只施加重力的情況相比，這個(gè)程序多加了如下幾句代碼：delete wall 6wall id 10 face (0,0,0) (0,2,0) (2,2,0) (2,0,0)wall id 10 ks 1e8 kn 1e8 fric 1.0wall id 10 zvel=-0.5e-3property n_bond

23、=5e5 s_bond=5e5hist wall zforce id 10其中，delete wall 6是刪除編號(hào)為6的面，即箱子的頂面。wall id 10 face (0,0,0) (0,2,0) (2,2,0) (2,0,0)wall id 10 ks 1e8 kn 1e8 fric 1.0wall id 10 zvel=-0.5e-3是建立一個(gè)ID為10的墻面，設(shè)定面的法向剛度、切向剛度、摩擦系數(shù)和Z方向的速度，這個(gè)速度的單位是m/step，因?yàn)殚_啟了density scaling（程序set dt dscale的作用），如果是完全的動(dòng)力學(xué)分析的話，單位將會(huì)是m/s（這個(gè)地方可解釋為

24、：靜力學(xué)只考慮系統(tǒng)達(dá)到平衡穩(wěn)態(tài)的情況，故單位與時(shí)間無關(guān)；動(dòng)力學(xué)關(guān)注的是每一時(shí)刻的狀態(tài)，故單位與時(shí)間有關(guān)）property n_bond=5e5 s_bond=5e5是為所有的接觸面定義法向和切向強(qiáng)度，以此來模擬初始的接觸力（應(yīng)該是類似于預(yù)緊力的作用，值得商榷），單位為Nhist wall zforce id 10用于監(jiān)測(cè)下降的板的Z方向的力的變化情況，計(jì)算完畢后如想查看，只需輸入plot hist 3即可。結(jié)果如圖所示這個(gè)力的曲線圖與英文幫助手冊(cè)上的圖明顯不同，但不知道是什么原因出的錯(cuò)，有待思考。下面附上英文幫助手冊(cè)上的圖。（1） 補(bǔ)充一句程序說明在每個(gè)程序的開頭，經(jīng)常有句set random

25、，這句話的意思是重置隨機(jī)數(shù)，因?yàn)間enerate生成小球的時(shí)候，半徑經(jīng)常是個(gè)范圍，需要隨機(jī)生成，如果不加這句的話，隨機(jī)得到的半徑大小可能一直是某一組隨機(jī)數(shù)，故需重置一下，這與vb編程中的隨機(jī)數(shù)重置設(shè)定是一樣的。（2）plot wall id=on 顯示所有wall的ID（3）定義字符串變量命令：macro（手冊(cè)上的例子均是將某字符串賦給某個(gè)變量，因此我取名字符串變量，如果是賦予具體數(shù)值，比如a=5之類的功能，應(yīng)該有別的命令可實(shí)現(xiàn)，待后續(xù)）Macro xx x 0 10Macro yy y 0 10Macro zz z 0 10Macro size id 1 100 rad 0.5 0.5Mac

26、ro ball xx yy zz sizeMacro size id 1001 1100 rad 0.1 0.1 ;注意！Gen ball以上程序需注意兩點(diǎn)：1、 macro命令的結(jié)構(gòu)形式：macro+變量名+字符串，字符串要加單引號(hào)；2、 倒數(shù)第二句是重新定義了變量size，需要注意的是size重新定義的時(shí)候需加單引號(hào)，第一次定義不需要加。（4）一般分析流程這個(gè)流程應(yīng)該是先建立模型，設(shè)定接觸和材料屬性，施加載荷和約束，然后再運(yùn)行程序使達(dá)到穩(wěn)態(tài)，最后根據(jù)需要改變條件繼續(xù)求解，從而得出結(jié)果。（5）生成球體的兩個(gè)命令命令generate和命令ballGenerate（簡(jiǎn)寫為gen）在第二天的講解中

27、已經(jīng)說過，它可一次生成大量的球填充到某個(gè)空間內(nèi)，ball命令是一次只生成一個(gè)球，這也是與generate的區(qū)別。Ball命令的格式為：ball rad * id * x * y * z *（*為各參數(shù)對(duì)應(yīng)的值，需用戶自己輸入）（6）施加邊界條件1、 給墻面施加速度和角速度Wall id 1 xvel 1e-3給ID號(hào)為1的墻面施加x方向速度，大小為1e-3Wall id 2 yspin 1e-3 x=2 y=2 z=0給ID號(hào)為2的墻面施加繞y軸轉(zhuǎn)動(dòng)的角速度，旋轉(zhuǎn)中心坐標(biāo)為（2,2,0）2、 給球施加恒定不變的速度或角速度為球施加速度或角速度，可用命令initial，change，property，但這種定義只能對(duì)球的速度進(jìn)行初始化，不能保證以后球一直以恒定的速度運(yùn)動(dòng)下去，因此需要用到fix命令。程序如下：Ini xvel 5e-3 zspin 1e-3 range id 40Fix x zspin range id 40上面程序的意思就是先給ID為40的球賦予x方向的速度和z方向的角速度，之后用fix命令保證這兩個(gè)速度恒定不變。3、 給球初始化運(yùn)動(dòng)參數(shù)Initial命令可初始化小球的力，