PFC3D中文講解例題

1、通過例子學習幾種常見命令例二：PFC3D目錄下的GuideStartfooting.dat，程序如下（注意：本程序與原程序不同，特加了plot set rotation (24.0,352.0,340.0)這句以使得模型轉過一定角度，讓用戶看的更清楚，另外還將原程序一分為二，主要是為了更好地學習每個命令的作用），以下是程序代碼。;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)造了一個由六個墻面圍成的箱子，每個墻面都賦予了一個獨一無二的ID號。每個墻面四邊形通過給出四邊形的四個頂點坐標來定義，即四個頂點連成四條線，然后再圍成一個四邊形墻面。注意：頂點的順序決定了墻面是否被激活，激活的墻面可以自動識別與小球的接觸，不激活的話就會被小球穿透，2、然后在箱子里填充小球，用Generate命令，也可簡寫，如下：gen id=1,750 rad 0.15,0.20 x=0,10 y=0,2 z=-5,0這個命令規(guī)定了

7、小球填充的空間范圍，小球的大小以及填充的總個數(shù)。Id=1,750的含義是：ID號從1到750其實就生成了750個小球，并給每個小球指定了各自的ID號。Rad 0.15,0.20的含義是：小球半徑從0.15m到0.20m隨機選擇x=0,10 y=0,2 z=-5,0的含義是：小球填充的空間是0<x<10，0<y<2，-5<z<0的一個立方體空間。如果填充的空間太小或者小球的半徑設置的太大，以至于這個空間容不下這么多的小球，那系統(tǒng)會自動減少生成的球的個數(shù)，并且系統(tǒng)會提示實際生成的球的個數(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顯示相關的命令的前面都用plot引出，然后后面再加相關的命令，如create Footing的意思就是創(chuàng)建一個名字為Footing的可視化窗口（這個在菜單欄點擊“window”，下拉菜單可發(fā)現(xiàn)多了一項vi

9、ew 1Footing，用戶可自己看下）；set title text是設定標題， set rotation 是設定模型的角度，以方便觀察；add ball yellow是設定球的顏色為黃色；add wall white是設定墻面（也就是邊界面）為白色；add axes brown是設定坐標軸為棕色； 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，可簡寫為prop，定義密度density為2000kg/m3，法向剛度和切向剛度均為108N/m6、改變小球半徑ini rad mul 1.47ini應該是initial的簡寫，用于初始化小球半徑（有待考證）rad表示半徑mul是multiply的簡寫，是倍數(shù)的意思。1.47是倍數(shù)因子，表示現(xiàn)在半徑是原來半徑的1.

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

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

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

14、ll wall，如運行完畢后想顯示別的結果，如接觸力則plot cforce（contact force的簡寫）?？偠灾琾lot后面加你想顯示的東東即可，用戶需靈活使用。此外，程序中的注釋用英文的分號; ，后面加注釋的內容即可。第二天的講解中小球都是自由狀態(tài)，無外力的作用，因此有必要討論施加外力的情況。（1） 對小球施加重力作用;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第一句是設定重力加速度為9.81m/s2，方向為Z軸負方向，grav是重力gravity的簡寫；第二句是設定球與球之間的摩擦系數(shù)，fric是friction的簡寫，在該程序中設定摩擦系數(shù)為1.0。運行結果如下圖所示：（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運行結果如圖所示與只施加重力的情況相比，這個程序多加了如下幾句代碼：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是刪除編號為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是建立一個ID為10的墻面，設定面的法向剛度、切向剛度、摩擦系數(shù)和Z方向的速度，這個速度的單位是m/step，因為開啟了density scaling（程序set dt dscale的作用），如果是完全的動力學分析的話，單位將會是m/s（這個地方可解釋為

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

25、，這句話的意思是重置隨機數(shù)，因為generate生成小球的時候，半徑經常是個范圍，需要隨機生成，如果不加這句的話，隨機得到的半徑大小可能一直是某一組隨機數(shù)，故需重置一下，這與vb編程中的隨機數(shù)重置設定是一樣的。（2）plot wall id=on 顯示所有wall的ID（3）定義字符串變量命令：macro（手冊上的例子均是將某字符串賦給某個變量，因此我取名字符串變量，如果是賦予具體數(shù)值，比如a=5之類的功能，應該有別的命令可實現(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以上程序需注意兩點：1、 macro命令的結構形式：macro+變量名+字符串，字符串要加單引號；2、 倒數(shù)第二句是重新定義了變量size，需要注意的是size重新定義的時候需加單引號，第一次定義不需要加。（4）一般分析流程這個流程應該是先建立模型，設定接觸和材料屬性，施加載荷和約束，然后再運行程序使達到穩(wěn)態(tài)，最后根據(jù)需要改變條件繼續(xù)求解，從而得出結果。（5）生成球體的兩個命令命令generate和命令ballGenerate（簡寫為gen）在第二天的講解中

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