張興剛 戴丹

1)(貴州大學物理學院,貴陽 550025)

2)(貴州大學計算機科學與技術(shù)學院,貴陽 550025)

二維顆粒堆積中壓力問題的格點系統(tǒng)模型?

張興剛1)?戴丹2)

1)(貴州大學物理學院,貴陽 550025)

2)(貴州大學計算機科學與技術(shù)學院,貴陽 550025)

糧倉效應,應力分布,格點系統(tǒng),三對角矩陣

為了便于從理論上探究糧倉效應產(chǎn)生的機理,處理筒倉中顆粒介質(zhì)的應力分布等問題,將二維顆粒堆積簡化為格點系統(tǒng),并且以隨機堆積為物理背景提出了一個由吸收系數(shù)p及側(cè)向傳遞系數(shù)q決定的力傳遞模型.給出了矩陣形式的力傳遞方程,提出基于二階差分方程的方法同時求解傳遞系數(shù)矩陣的特征值和特征向量,從理論上導出了一種典型情況下容器底部壓力分布與頂部壓力分布的關(guān)系式.對有效質(zhì)量隨總質(zhì)量變化關(guān)系的理論分析表明,該模型可以給出與Janssen模型類似的結(jié)果.對無負載情況下的底部應力分布進行了理論計算,結(jié)果表明容器底部中央應力最大,離中央越遠應力越小.運用數(shù)值計算討論了p與q對容器底部壓力隨堆積高度變化曲線的影響.

1 引 言

顆粒物質(zhì)靜力學主要研究靜態(tài)顆粒體系的堆積結(jié)構(gòu)、顆粒間的相互作用、顆粒介質(zhì)中的應力分布等問題,在糧倉的設(shè)計、堆石壩的穩(wěn)定性、地基的承載能力等工程問題中有直接應用,山崩、滑坡、塌陷等災害的理解也與其密切相關(guān).有許多有趣的顆粒物質(zhì)靜力學現(xiàn)象,例如糧倉效應和沙堆應力凹陷[1].人們很早就發(fā)現(xiàn)了糧倉效應,即糧倉底部所受壓力隨糧食堆積高度的增加趨于飽和的現(xiàn)象.早在1895年,德國工程師Janssen就提出了一個基于某些假設(shè)的連續(xù)介質(zhì)模型,給出了糧倉效應的簡單解釋[1,2].在糧倉效應的相關(guān)實驗中,通常關(guān)注容器底部測得的有效質(zhì)量Me(即容器底部壓力與重力加速度g的比值)與容器中顆?？傎|(zhì)量Mt的關(guān)系.根據(jù)Janssen模型可建立顆粒介質(zhì)的力學平衡方程,在容器頂部無負載的情況下,可以求得有效質(zhì)量為

式中Ms為飽和質(zhì)量.(1)式對于二維和三維的情況都適用[3].在一定實驗條件下,所得到的實驗結(jié)果可以用Janssen模型予以解釋.

然而,顆粒堆中的應力問題并不簡單.一般情況下需要結(jié)合應力平衡方程?·σ+ρg=0(其中σ為顆粒介質(zhì)的應力張量場,ρ為介質(zhì)的密度分布)、介質(zhì)的本構(gòu)關(guān)系、邊界條件等因素,才能從理論上確定介質(zhì)中的應力分布.但是,目前還沒有普遍適用的關(guān)于靜態(tài)顆粒介質(zhì)的本構(gòu)關(guān)系,實際顆粒堆中的邊界條件也往往難以確定;更為困難的是,顆粒介質(zhì)的宏觀性質(zhì)往往與其制備方式、加載條件等因素密切相關(guān).事實上,顆粒物質(zhì)靜力學中有許多基本的物理問題需要深入研究[4].直到最近,仍有不少關(guān)于糧倉效應、Janssen模型等經(jīng)典問題的研究報道,人們從實驗[3,5?10]、計算機模擬[11]、理論模型等[3,12?14]方面對筒倉中顆粒堆的壓力問題進行了更深入的探討.在大部分研究工作中,容器底部所受壓力隨顆粒堆積高度的變化(可間接用Me(Mt)函數(shù)表征)仍是一個基本問題.不過,人們用不同的方式、從不同的角度加深了對這一問題的認識.Vanel等[3]發(fā)現(xiàn)頂部負載較大時,Me(Mt)會呈現(xiàn)非單調(diào)變化,這與Janssen模型預言的單調(diào)變化不相符,用OSL(oriented stress linearity)模型可比較好地理解這一現(xiàn)象.Wambaugh等[7]用光彈性實驗對不同頂部負載下的Me(Mt)關(guān)系進行研究,觀察到了非單調(diào)變化、Me沒有趨向飽和、Janssen模型中的特征長度并非常數(shù)等現(xiàn)象.另外,研究者通過勻速移動筒倉側(cè)壁[5]、上提或下壓側(cè)壁等[6,14]方式改變實驗條件,研究了Me(Mt)關(guān)系、顆粒對側(cè)壁的最大靜摩擦力等問題.最近,有人提出新的實驗方法測量顆粒對側(cè)壁的作用力[9,10],例如文獻[9]中通過圖像關(guān)聯(lián)技術(shù)測量彈性容器壁的應變場,從而間接測量側(cè)壁所受的應力分布.

顆粒體系的宏觀力學性質(zhì)應該由其微觀結(jié)構(gòu)、相互作用的統(tǒng)計特征決定.兩個宏觀上看起來相同的顆粒堆,它們的微觀狀態(tài)一般不同,其中重要微觀因素的不同有可能導致宏觀上力學性質(zhì)的重大差異.前人的研究表明,對于同樣的問題,不同的實驗過程往往會得到不同的實驗結(jié)果[5?10];另外,已有宏觀理論模型只在一定條件下可以用于理解部分實驗結(jié)果,目前不是很清楚這些模型適用(或者不適用)的物理機理.可以看到,這些問題主要來源于顆粒體系在微觀結(jié)構(gòu)和相互作用方面的多樣性、不確定性、易變性,以及尚不明確的體系的微觀特征與宏觀力學性質(zhì)之間的關(guān)系.因此,加強對微觀結(jié)構(gòu)、力網(wǎng)、力傳遞等問題的研究,建立微觀特征與宏觀性質(zhì)之間的聯(lián)系,對于顆粒靜力學來說非常重要[4,15].文獻[13]建立了圓盤的晶格堆積模型,基于力學平衡條件得到了相鄰兩層顆粒間力傳遞的規(guī)則,并且利用力傳遞規(guī)則給出了容器底部總壓力會隨著堆積層數(shù)增加趨向飽和的結(jié)果.實際的顆粒堆一般是非晶格結(jié)構(gòu),難以直接從微觀上對其進行理論研究,且更細致的研究應該討論應力分布,而不只是總的壓力.標量q模型[16,17]將顆粒堆積抽象為格點系統(tǒng),通過上下兩層格點間重量的隨機傳遞簡化實際堆積中顆粒間的相互作用;結(jié)合馬爾科夫過程理論,q模型可以比較好地用于理解隨機堆積中力的概率分布等問題.

受到q模型的啟發(fā),為了便于從理論上探究頂部壓力分布與底部壓力分布之間的關(guān)系,本文將二維顆粒體系簡化、抽象為帶有邊界的格點系統(tǒng),以隨機堆積為物理背景建立了一個由p和q這兩個參數(shù)決定的力傳遞模型,其中吸收系數(shù)p反映了邊界由于摩擦而對顆粒產(chǎn)生的支撐作用,側(cè)向傳遞系數(shù)q描述了顆粒向其正下方左右兩側(cè)傳遞力的能力.該模型的解析計算涉及到一種三對角矩陣的相似對角化,我們提出基于二階差分方程的方法同時求解三對角矩陣的特征值和特征向量,對p=1/2,q=1這一典型的情況進行了解析計算.最后,通過模型的求解討論了Me(Mt)關(guān)系、底部壓力分布等問題.該研究加深了對微觀力傳遞規(guī)則與宏觀壓力分布之間關(guān)系的認識.

2 格點系統(tǒng)模型

為了研究筒倉中顆粒物質(zhì)對容器底部產(chǎn)生的總壓力以及壓力的分布,可以考慮圖1所示二維顆粒體系.如果顆粒的尺寸遠小于容器的尺寸,則可以將容器中的顆粒體系看作連續(xù)的顆粒介質(zhì).對于二維的情況,顆粒介質(zhì)的密度可以用ρ=Mt/(LxLy)表征,其中Mt為顆粒介質(zhì)的總質(zhì)量;邊界處壓力的分布可以用應力(即單位長度邊界所受的壓力)的分布表征.于是顆粒介質(zhì)頂部所受的負載可用應力函數(shù)σ(x,0)描述,容器底部的應力分布可用函數(shù)σ(x,Ly)描述.二維顆粒堆積中關(guān)于壓力分布的一個基本問題,是根據(jù)σ(x,0)以及顆粒體系的一些特征求出σ(x,Ly).實際的靜態(tài)顆粒體系往往是非晶格的堆積,很難從微觀上對其進行細致的幾何和力學分析.為了避免直接分析堆積結(jié)構(gòu)帶來的困難,將二維顆粒堆積簡化為圖2所示的規(guī)則排列格點系統(tǒng).在格點系統(tǒng)中,一個格點等效于一小團顆粒,第m行第n列的格點可用符號P(m,n)表示,設(shè)每個格點的質(zhì)量均為m0,則每個格點都會受到重力m0g的作用.記豎直方向上相鄰格點的間距為Δy,水平方向上的間距為Δx,則

一般情況下,不能簡單地使用力傳遞的概念理解物理問題,除非這一概念有力學規(guī)律的支撐.已有的研究表明[12,13,16,17],在一定條件下可以采用力傳遞這種簡單、直觀的方式處理靜態(tài)顆粒體系中的許多力學問題,但應該注意的是,不同的微觀結(jié)構(gòu)往往會導致不同的力傳遞規(guī)則.類似于標量q模型,以隨機堆積為物理背景,定義格點系統(tǒng)中一種典型的力傳遞規(guī)則,其大致情況如圖2所示.規(guī)定同一行的格點間沒有豎直方向的分力,不相鄰的兩行格點間也沒有相互作用.對于任意格點P(m,n)

而言,稱其所受的重力、頂部對其產(chǎn)生的載荷以及第m?1行格點(如果存在的話)對其產(chǎn)生的作用力為輸入力,將P(m,n)對第m+1行格點及其對容器的作用力稱為輸出力.對于每一個格點,只有上述5種豎直方向分力的存在,為了保證格點在豎直方向的力平衡條件,總輸入力必須等于總輸出力.

1)當m=M時,格點只向容器底部輸出力,有fM=wM;

2)當m/=M,n/∈{1,N}時,P(m,n)將其總輸入力的q/2輸出給P(m+1,n?1),將總輸入力的q/2輸出給P(m+1,n+1),剩下的(1?q)wmn輸出給P(m+1,n),稱q為側(cè)向傳遞系數(shù),滿足0＜q≤1;

3)當m/=M,n∈{1,N}時,P(m,1)將其總輸入力的(1?p)q輸出給P(m+1,2),將總輸入力的(1?p)(1?q)輸出給P(m+1,1),剩下的pwm1輸出給容器左側(cè),稱p為吸收系數(shù),滿足0≤p＜1,P(m,N)的力傳遞情況與P(m,1)類似.

根據(jù)上述定義的力傳遞規(guī)則,當1≤m≤M?1時有力傳遞方程

其中N階方陣

稱為傳遞系數(shù)矩陣,這是一個三對角矩陣[18].矩陣中的吸收系數(shù)p反映了容器側(cè)壁與顆粒之間的摩擦所導致的格點總輸入力被器壁吸收的程度,若容器側(cè)壁與顆粒無摩擦,則p=0,有摩擦時0＜p＜1.矩陣中的側(cè)向傳遞系數(shù)q反映了顆粒微團向其左右兩側(cè)傳遞力的程度,主要由微觀的堆積結(jié)構(gòu)決定,若顆粒微團只會向其正下方傳遞力,則q=0,但是這種情況在隨機堆積中一般不存在.為表述方便,將該模型稱為pq模型.為了由(4)式和(5)式推出wM,作變換可以得到

式中E為N階單位矩陣.假設(shè)

則容易推出

假設(shè)A相似于對角矩陣Λ,即存在可逆矩陣Q使A=QΛQ?1,并且假設(shè)矩陣A?E可逆,那么結(jié)合(4)式、(8)式和(9)式可以得到

于是計算fM(即wM)相當于考慮傳遞系數(shù)矩陣A的相似對角化以及線性方程組(8)式的求解(或者矩陣A?E的可逆性).

3 傳遞系數(shù)矩陣的相似對角化

著重考慮N階傳遞系數(shù)矩陣A的相似對角化問題.在p和q取一般數(shù)值的情況下,通過繁雜的數(shù)學討論可以看到:矩陣A(p,q)的相似對角化本質(zhì)上涉及到非特殊的高次代數(shù)方程的求解,無法找到其解析的結(jié)果.因此,以特殊但在物理上屬于典型的p=1/2,q=1的情形進行理論分析,這相當于筒倉系統(tǒng)中顆粒與器壁的摩擦不是很小也不是很大,并且顆粒團由于成拱效應明顯而幾乎不向其正下方傳遞力的情況.此時矩陣A(1/2,1)是一個對稱的三對角矩陣,根據(jù)線性代數(shù)中對稱矩陣的相似對角化理論易知,存在N階正交矩陣以及實對角矩陣Λ=diag(λ1,λ2,···,λN)使得

式中λk為A(1/2,1)的特征值,pNk)T為與λk相對應的單位特征向量.

為了給出矩陣P與Λ,需要考慮A(1/2,1)的特征值與特征向量.設(shè)xT=(x1,x2,···,xN)T是A(1/2,1)的特征向量,λ是相應的特征值,則

傳統(tǒng)的線性代數(shù)理論是先由行列式給出矩陣的特征多項式,利用其對應的特征方程求解特征值,然后將特征值代入(12)式計算出相應的特征向量.這里采用直接分析齊次線性方程組(12)式的方法,同時進行特征值與特征向量的計算從而簡化理論分析過程.將A(1/2,1)的矩陣形式代入(12)式可得帶有邊值條件的差分方程為

由于A(1/2,1)是實對稱矩陣,其特征值都是實數(shù),因此只需在λ∈R的范圍內(nèi)討論(13)式.首先,根據(jù)二階常系數(shù)差分方程理論求(13)式中差分方程的通解.該差分方程對應的特征方程為

以x1和x2為邊值,則(15)式中的待定系數(shù)c1與c2應該滿足

由θ/=k0π易知該方程組中系數(shù)矩陣的行列式不為零,容易得到向量(c1,c2)T的惟一解,將其解代入(15)式有

利用(17)式可得xN?1,xN關(guān)于x1,x2的表達式,將它們代入(13)式中的邊值條件有

該方程組的系數(shù)矩陣行列式可化簡為

由上述分析過程可知,|λ|＜1時矩陣A(1/2,1)的特征值由方程det(θ,N)=0中θ的解的情況決定.顯然,當且僅當時det(θ,N)=0,其中結(jié)合θ/=k0π,λ=cosθ以及cosθ函數(shù)的基本性質(zhì)可以得到A(1/2,1)的N個不同的特征值

這也是A(1/2,1)的所有特征值.因此p=1/2,q=1的情況下不用討論λ=±1,|λ|＞1的情形.結(jié)合(17)式和(18)式可得與λk相對應的一個特征向量

利用歐拉公式及等比數(shù)列求和不難求出的長度的平方

于是得到N個滿足正交歸一條件的特征向量

利用(20)式和(23)式可以給出(11)式中的矩陣P與Λ.至此,矩陣A(1/2,1)的相似對角化完成.

4 結(jié)果與討論

由(20)式可知,A(1/2,1)的所有特征值都滿足?1＜λk＜1,由于其特征值不為1,因此矩陣A(1/2,1)?E是可逆的.只要令Q=P,就可以利用(10)式給出容器底部壓力分布fM與頂部壓力分布f0的關(guān)系,即

由于

因此

顯然(27)式中的求和部分只有當k取奇數(shù)時才有貢獻.下面采用近似的方法討論(27)式,一方面是因為不能嚴格地計算(27)式,另一方面是因為在實際問題中更關(guān)心M和N遠大于1時Me隨Mt變化的總趨勢.記

如果M的取值與N+1的平方有相同的量級,可設(shè)M=h(N+1)2,其中h可以反映豎直方向堆積的高度.則(29)式在N?1時可簡化為

將其代入(27)式,用求和部分的第1項代替整個求和部分,可以得到

這與Janssen模型得到的(1)式類似,隨著堆積高度的增加,有效質(zhì)量逐漸增大,最后達到飽和.飽和質(zhì)量滿足Ms≈16m0N3/π4,總質(zhì)量Mt≈m0hN3,(31)式可寫為

這是格點系統(tǒng)模型在M和N遠大于1并且取一級近似(即關(guān)于Me的求和式中只取k=1這一項)時得到的結(jié)果,與Janssen模型得到的(1)式在指數(shù)函數(shù)的系數(shù)部分有微小的差異,這說明pq模型中的力傳遞機理可以用于理解Janssen行為產(chǎn)生的原因.一方面,成拱效應使得顆粒能向其正下方的左右兩側(cè)傳遞力;另一方面,摩擦使得傳遞到邊界的作用力可以被吸收.這兩方面的因素相結(jié)合,產(chǎn)生了壓力飽和.

傳統(tǒng)的Janssen模型主要考慮容器底部總的壓力[1,2].利用pq模型,可以從理論上考慮容器底部的壓力分布問題.為簡單起見,仍討論頂部無負載的情形.由(24)式可得

類似于有效質(zhì)量的討論,這里也采用近似的方法計算(33)式,在N遠大于1且取一級近似的情況下,得到

可以通過數(shù)值計算檢驗(34)式的近似程度,為此,設(shè)格點所受重力m0g為單位1,分別利用(33)式和(34)式計算壓力分布fM(n)隨M的變化.圖3所示為精確解和近似解的數(shù)值計算結(jié)果,顯然,只需取一級近似就能比較好地解決壓力分布問題.從圖3還可以看到,容器底部中央所受壓力最大,從中央向兩側(cè)移動時壓力逐漸減小.另外,在容器底部的某一固定位置處,隨著堆積高度的增加,該位置處所受壓力逐漸增大,最終會達到飽和.pq模型中力的傳遞類似于一維隨機行走.在圖2所示格點系統(tǒng)中,y軸相當于隨機行走問題中的時間軸.隨機行走問題在宏觀上會導致擴散方程.實際上,將pq模型連續(xù)化后會得到帶有吸收邊界且含源的擴散方程.由于邊界的吸收作用會隨著堆積高度的增加逐漸由邊界擴散到中央,因此產(chǎn)生了圖3所示壓力分布.如果將顆粒體系當作連續(xù)介質(zhì),應該采用應力分布函數(shù)σ(x,Ly)描述容器底部的壓力分布情況.類似于有效質(zhì)量的討論,取M=h(N+1)2,利用(34)式、(2)式和(3)式,可以導出

圖3 壓力分布隨豎直方向?qū)訑?shù)M的變化Fig.3.Pressure distribution varies with the vertical layer number M.

式中β=Δy/Δx2是一個涉及格點系統(tǒng)模型微觀特征的量,需要從微觀上對顆粒堆積結(jié)構(gòu)及相互作用進行分析才能確定.

由于一般情況下不能求出矩陣A(p,q)相似對角化后的解析結(jié)果,因此為了研究p和q對壓力的影響,需要采用數(shù)值計算方法.由(4)式和(5)式可得

設(shè)m0g為單位1,在頂部無負載的情況下,利用(36)式可得底部總壓力為

在Δx,Δy確定的情況下,可以用M和N分別代表顆粒堆積在豎直方向的高度和水平方向的長度.

下面用數(shù)值計算討論p和q對壓力-高度曲線F(M)的影響.如圖4所示,當吸收系數(shù)或側(cè)向傳遞系數(shù)為0時,壓力只會隨高度線性增加,不會出現(xiàn)壓力飽和,這進一步說明了成拱效應與邊界摩擦是產(chǎn)生壓力飽和的主要原因.從圖中還可以看到,隨著吸收系數(shù)、側(cè)向傳遞系數(shù)的增加,飽和壓力的數(shù)值會減小.吸收系數(shù)越大代表邊界因摩擦而產(chǎn)生的支撐作用越強,側(cè)向傳遞系數(shù)越大代表內(nèi)部顆粒的重力能更容易地傳遞到邊界從而被支撐,因此產(chǎn)生了上述結(jié)果.

圖4 (a)吸收系數(shù)p與(b)側(cè)向傳遞系數(shù)q對底部總壓力隨堆積高度變化曲線的影響Fig.4.E ff ects of parameters(a)p and(b)q on bottom pressure variation with packing height.

5 結(jié) 論

研究筒倉中顆粒物質(zhì)的微觀力傳遞規(guī)則與宏觀壓力分布之間的關(guān)系,是深入了解顆粒介質(zhì)力學性質(zhì)的重要途徑.傳統(tǒng)的q模型可以比較好地解釋隨機堆積中力的概率分布、力的漲落等統(tǒng)計問題,這說明q模型至少在一定程度上正確反映了隨機堆積中力傳遞的特征.為了研究筒倉中顆粒物質(zhì)的壓力、應力分布問題,本文針對二維隨機堆積建立了一個類似于q模型的格點系統(tǒng)模型——pq模型.通過對p=1/2,q=1這一典型情況的解析計算結(jié)果可以看到:對于Me(Mt)關(guān)系,pq模型能給出與Janssen模型類似的結(jié)果;另外,pq模型還可用于討論容器底部應力分布、點載荷下的響應等問題.不過應該注意到,盡管晶格堆積可以表現(xiàn)出與隨機堆積不一樣的力傳遞特征,但是晶格堆積模型也能導致與Janssen模型類似的Me(Mt)關(guān)系.如果統(tǒng)一采用格點模型理解問題,導致類似的結(jié)果來源于兩種模型的共同性:格點能向其正下方左右兩側(cè)的格點傳遞力,邊界處格點的總輸入力會被吸收.可以看到,格點系統(tǒng)模型用于理解顆粒靜力學問題時有其簡單、方便之處,同時也可以擁有豐富的內(nèi)容.本文討論的pq模型只是格點系統(tǒng)模型中一種可能的力傳遞形式,對于不同的堆積結(jié)構(gòu)和環(huán)境條件,可能具有不同的力傳遞規(guī)則,在數(shù)學上一般表現(xiàn)為不同的傳遞系數(shù)矩陣.這自然引入了一個值得深入研究的問題:力傳遞規(guī)則與堆積結(jié)構(gòu)等因素之間有何種關(guān)系？另外,在p與q取一般數(shù)值的情況下,只要對A(p,q)進行適當?shù)南嗨谱儞Q,就可將其轉(zhuǎn)化為對稱矩陣,并可利用文中提出的方法對其進行相似對角化,從理論上分析p,q等因素對有關(guān)壓力問題的影響.但是這一過程的嚴格數(shù)學討論十分繁瑣,而且由于涉及非特殊性高次代數(shù)方程的求解,一般不能給出解析結(jié)果.如果將離散的格點系統(tǒng)模型連續(xù)化,利用其中的力傳遞規(guī)則給出應力函數(shù)σ(x,y)所服從的微分方程(相當于力傳遞方程連續(xù)化后的結(jié)果),pq模型會給出拋物線方程,利用偏微分方程理論可以比較容易地對問題進行求解.不過,離散的格點系統(tǒng)模型也有其優(yōu)勢,例如可以很方便地進行數(shù)值計算.

[1]Lu K Q,Liu J X 2004Physics33 629(in Chinese)[陸坤權(quán),劉寄星2004物理33 629]

[2]de Gennes P G 1999Rev.Mod.Phys.71 374

[3]Vanel L,Claudin P,Bouchaud J P,Cates M E,Wittmer J P 2000Phys.Rev.Lett.84 1439

[4]Sun Q C,Hou M Y,Jin F 2011The Physics and Mechanics of Granular Matter(Beijing:Science Press)(in Chinese)[孫其誠,厚美瑛,金峰2011顆粒物質(zhì)物理與力學(北京:科學出版社)]

[5]Bertho Y,Frédérique G D,Hulin J P 2003Phys.Rev.Lett.90 144301

[6]Peng Z,Li X Q,Jiang L,Fu L P,Jiang Y M 2009Acta Phys.Sin.58 2090(in Chinese)[彭政,李湘群,蔣禮,符力平,蔣亦民2009物理學報58 2090]

[7]Wambaugh J F,Hartley R R,Behringer R P 2010Eur.Phys.J.E32 135

[8]Li X Q,Jiang Y M,Peng Z 2010J.Shandong Univ.(Natural Science)45 101(in Chinese)[李湘群,蔣亦民,彭政2010山東大學學報45 101]

[9]Cambau T,Hure J,Marthelot J 2013Phys.Rev.E88 022204

[10]Li Z F,Peng Z,Jiang Y M 2014Acta Phys.Sin.63 104503(in Chinese)[李智峰,彭政,蔣亦民2014物理學報63 104503]

[11]Landry J W,Grest G S,Silbert L E,Plimpton S J 2003Phys.Rev.E67 274

[12]Marconi U M B,Petri A,Vulpiani A 2000Physica A280 279

[13]Zhang X G,Hu L,Long Z W 2006Chin.J.Comput.Phys.23 642(in Chinese)[張興剛,胡林,隆正文 2006計算物理23 642]

[14]Jiang Y M,Zheng H P 2008Acta Phys.Sin.57 7360(in Chinese)[蔣亦民,鄭鶴鵬 2008物理學報 57 7360]

[15]Gendelman O,Pollack Y G,Procaccia I 2016Phys.Rev.Lett.116 078001

[16]Liu C,Nagel S R,Schecter D A,Coppersmith S N,Majumdar S,Narayan O,Witten T A 1995Science269 513

[17]Coppersmith S N,Liu C,Majumdar S,Narayan O,Witten T 1996Phys.Rev.E53 4673

[18]Yang S L 2010Math in Practice and Theory40 155(in Chinese)[楊勝良2010數(shù)學的實踐與認識40 155]

Lattice model for pressure problems in two-dimensional granular columns?

Zhang Xing-Gang1)?Dai Dan2)

1)(College of Physics,Guizhou University,Guiyang 550025,China)
2)(College of Computer Science and Technology,Guizhou University,Guiyang 550025,China)

20 January 2017;revised manuscript

29 May 2017)

In order to make it easier to investigate some problems such as the mechanism of Janssen e ff ect and the stress distribution in granular medium,we simplify a granular column into a lattice system,in which a lattice point represents a small lump of granular medium and only neighbor interactions are considered.To study the disordered granular columns,a force propagation lattice model determined by the absorption coefficientpand the lateral transfer coefficientqis proposed,and this model is analyzed from the theoretical view.Firstly,the equation of force propagation in the matrix form is given,and this equation is determined by a tridiagonal matrixA(p,q)that is called transfer coefficient matrix.Based on the force transfer equation,the bottom force distribution varying with the top force distribution and the layer of lattice system is deduced,and its analytical solution refers to the similarity diagonalization of matrixA(p,q).Then,a method based on the second order di ff erence equations is proposed to obtain the eigenvalues and eigenvectors of the transfer coefficient matrix.The eigenvalues and eigenvectors ofA(p,q)can be rigorously deduced for a typical case,and with these results the pressure distribution relationship between the top and the bottom of the container is given.Based on these theoretical expressions,the relationship between the e ff ective mass and the total mass of granular medium is deduced,and it means that the force propagation model and the Janssen model can lead to similar results.Moreover,the bottom stress distribution is calculated without the top load.Calculations show that the stress distribution reaches a maximum at the center bottom and drops down to either side.Finally,numerical calculations are performed to investigate the e ff ects of parameterspandqon the relation between bottom pressure and packing height.Numerical results show that the saturated value of pressure decreases while parameterporqincreases.

Janssen e ff ect,stress distribution,lattice systems,tridiagonal matrix

(2017年1月20日收到;2017年5月29日收到修改稿)

10.7498/aps.66.204501

?貴州省科技合作計劃(批準號:20157641)、貴州大學引進人才科研基金(批準號:201334)和貴州大學教育教學改革研究項目(批準號:JGYB201517)資助的課題.

?通信作者.E-mail:xgzhang@gzu.edu.cn

?2017中國物理學會Chinese Physical Society

http://wulixb.iphy.ac.cn

PACS:45.70.–n,83.80.Fg,02.10.UdDOI:10.7498/aps.66.204501

*Project supported by the Science and Technology Cooperation Project of Guizhou Province,China(Grant No.20157641),the Scienti fi c Research Foundation for Talents of Guizhou University,China(Grant No.201334),and the Education and Teaching Reform Project of Guizhou University,China(Grant No.JGYB201517).

?Corresponding author.E-mail:xgzhang@gzu.edu.cn