可降解高聚物的降解跨尺度建模研究

高曉浩, 張?zhí)壹t, 楊智勇, 孫樂(lè)冉, 樊春美, 張德政

(1. 北京科技大學(xué) 計(jì)算機(jī)與通信工程學(xué)院, 北京 100083；

2. 材料領(lǐng)域知識(shí)工程北京市重點(diǎn)實(shí)驗(yàn)室, 北京 100083)

可降解高聚物的降解跨尺度建模研究

高曉浩1,2, 張?zhí)壹t1,2, 楊智勇1,2, 孫樂(lè)冉1,2, 樊春美1,2, 張德政1,2

(1. 北京科技大學(xué) 計(jì)算機(jī)與通信工程學(xué)院, 北京 100083；

2. 材料領(lǐng)域知識(shí)工程北京市重點(diǎn)實(shí)驗(yàn)室, 北京 100083)

可降解聚合物由于其良好的應(yīng)用前景吸引了廣大科研工作者的興趣。但降解過(guò)程復(fù)雜，使得建立降解過(guò)程的模型研究總是局限于某些片段。如在尺度上或微觀、或宏觀；研究方法上或離散、或連續(xù)；模型也或隨機(jī)或確定。本文在蒙特卡洛方法、元胞自動(dòng)機(jī)方法基礎(chǔ)上提出新的建模思路，即將不同尺度上離散化與連續(xù)性相耦合、隨機(jī)性與確定性相耦合的元胞蒙特卡洛自動(dòng)機(jī)方法(CMCA)。將連續(xù)問(wèn)題以元胞的方法進(jìn)行離散，在不同元胞上以概率實(shí)現(xiàn)隨機(jī)事件的發(fā)生，由元胞的更新迭代進(jìn)行演化，揭示了從微介觀到宏觀的降解變化過(guò)程。本文對(duì)具體的步驟和執(zhí)行算法進(jìn)行了詳細(xì)的解釋?zhuān)⒁跃廴樗岜P(pán)進(jìn)行實(shí)例計(jì)算，計(jì)算結(jié)果表明建模方法正確。

蒙特卡洛方法；元胞自動(dòng)機(jī)；生物材料；可降解高聚物；建模；元胞蒙特卡洛自動(dòng)機(jī)方法

1 前 言

生物可降解高分子(聚乳酸，聚丙交酯等)由于在生物體吸收性縫合材料、骨科固定、藥物緩釋系統(tǒng)及組織工程支架中的顯著優(yōu)點(diǎn)而備受廣泛關(guān)注。這些應(yīng)用能否成功在很大程度上取決于高分子的降解速率能否得到控制。高聚物的降解過(guò)程一般認(rèn)為是經(jīng)歷兩個(gè)過(guò)程：首先是在水分子的攻擊下，高分子長(zhǎng)鏈水解斷裂成短鏈，然后再在酶的作用下進(jìn)一步降解，最終生成無(wú)害的水和二氧化碳[1,2]。其中第一階段水解反應(yīng)將直接影響高聚物的形狀及強(qiáng)度變化，所以高聚物降解率的控制主要取決于第一階段即水解階段。用實(shí)驗(yàn)方法研究高聚物的降解過(guò)程周期太長(zhǎng)，且由于水解過(guò)程中的自催化反應(yīng)[3]使得不同尺寸的高聚物降解速率都不一樣，從而增加了工作量與不確定性。計(jì)算機(jī)建模的方法因可以使周期縮短，且可以提供更豐富的降解過(guò)程數(shù)據(jù)而引起了很多科研工作者的關(guān)注。建模方法也從剛開(kāi)始單一的擴(kuò)散、化學(xué)反應(yīng)、質(zhì)量傳遞等現(xiàn)象出發(fā)的宏觀機(jī)理模型[4～7]，微(或介)觀的蒙特卡洛模型(MC)[8～11]、元胞自動(dòng)機(jī)模型(CA)[12～14]等發(fā)展到宏觀與微(或介)觀相結(jié)合[15]或隨機(jī)與質(zhì)量傳遞相結(jié)合的模型[16]。降解機(jī)理復(fù)雜，過(guò)程涉及化學(xué)反應(yīng)、物理擴(kuò)散與力學(xué)性能的改變，體現(xiàn)在微觀尺度上的鏈隨機(jī)斷裂，介觀尺度上微孔的形成，到宏觀尺度上高聚物質(zhì)量的減少、結(jié)構(gòu)的坍塌等。目前的研究大多涉及一個(gè)至兩個(gè)尺度的建模，從研究的方法上，也是分別采用離散建模方法和連續(xù)建模方法。本文在蒙特卡洛方法與元胞自動(dòng)機(jī)方法的基礎(chǔ)上提出元胞蒙特卡洛自動(dòng)機(jī)法(CMCA)、并將其用于模擬微觀上隨機(jī)發(fā)生的高分子鏈斷裂、介觀尺度上孔隙的形成及宏觀上質(zhì)量的減少，模擬結(jié)果與實(shí)驗(yàn)比較，以證明本文方法的有效性。由此CMCA方法將微介觀離散的方法、隨機(jī)的方法相統(tǒng)一并與宏觀的擴(kuò)散方程相耦合，從不同尺度揭示高聚物的鏈斷裂、低聚體擴(kuò)散、孔洞形成等降解過(guò)程，并相統(tǒng)一并與宏觀的擴(kuò)散方程相耦合尺度、多方法相耦合的解決思路與算法。

2 建模理論與思路

2.1 模擬化學(xué)反應(yīng)過(guò)程的蒙特卡洛模型(MC)

蒙特卡洛方法(Monte Carlo method，簡(jiǎn)稱(chēng)MC)，也稱(chēng)統(tǒng)計(jì)模擬方法，是二十世紀(jì)四十年代中期由于科學(xué)技術(shù)的發(fā)展和電子計(jì)算機(jī)的發(fā)明，而被提出的一種以概率統(tǒng)計(jì)理論為指導(dǎo)的一類(lèi)非常重要的數(shù)值計(jì)算方法。它是指使用隨機(jī)數(shù)(或更常見(jiàn)的偽隨機(jī)數(shù))來(lái)解決很多計(jì)算問(wèn)題的方法。對(duì)于化學(xué)反應(yīng)的動(dòng)態(tài)過(guò)程，有兩個(gè)主要問(wèn)題：1)下一個(gè)反應(yīng)將在何時(shí)發(fā)生；2)下一個(gè)反應(yīng)是哪一種類(lèi)型。

蒙特卡洛隨機(jī)算法可以解決這兩個(gè)問(wèn)題，如下：

(1) 輸入M個(gè)反應(yīng)的反應(yīng)常數(shù)π1, π2, …… , πM和N種反應(yīng)分子的初始分子數(shù)x1, x2, …… , xN，同時(shí)將時(shí)間t設(shè)為零；

(2) 根據(jù)現(xiàn)時(shí)各反應(yīng)分子的分子數(shù)，計(jì)算α1, α2, …… ,αv,……, αM，對(duì)于二級(jí)反應(yīng)如：A+B→C，計(jì)算公式為

其中πv表示第v個(gè)反應(yīng)的反應(yīng)常數(shù)，xA，xB為參加第v個(gè)反應(yīng)的反應(yīng)物分子數(shù)；vα為該反應(yīng)的概率；(3) 產(chǎn)生兩個(gè)單位區(qū)間均勻分布的隨機(jī)數(shù)r1、r2，分別由公式

計(jì)算得到 ?t、μ，其中 ?t為兩個(gè)反應(yīng)的時(shí)間間隔，μ為下一時(shí)刻該發(fā)生的第μ個(gè)反應(yīng)；

(4) 使反應(yīng)時(shí)間增加 ?t并按第μ個(gè)反應(yīng)，調(diào)整有關(guān)分子的個(gè)數(shù)(反應(yīng)物分子數(shù)減一，生成物分子數(shù)加一)，然后回到(2)。

2.2 元胞自動(dòng)機(jī)算法

元胞自動(dòng)機(jī)(Cellular Automata，簡(jiǎn)寫(xiě)為CA)，是由馮·諾依曼在1948年提出的，它是一種動(dòng)態(tài)的網(wǎng)格模型，是定義在一個(gè)由具有離散、有限狀態(tài)的元胞空間上，并按照一定的局部規(guī)則，在離散的時(shí)間維度上演化的動(dòng)力系統(tǒng)。散布在規(guī)則網(wǎng)格(Lattice Grid)中的每一元胞(Cell)取有限的離散狀態(tài)，遵循同樣的作用規(guī)則，依據(jù)確定的局部規(guī)則作同步更新。大量元胞通過(guò)簡(jiǎn)單的相互作用而構(gòu)成動(dòng)態(tài)系統(tǒng)的演化。

元胞自動(dòng)機(jī)具有如下幾個(gè)要素：

(1) 元胞

元胞是元胞自動(dòng)機(jī)的基本單元，每個(gè)元胞都有自己的狀態(tài)，所有元胞可能出現(xiàn)的狀態(tài)構(gòu)成的集合稱(chēng)為元胞狀態(tài)空間。元胞離散地分布在線性，平面或空間維度上，根據(jù)鄰居元胞的狀態(tài)和轉(zhuǎn)化規(guī)則在離散的時(shí)間距離上進(jìn)行狀態(tài)變化。

(2) 元胞狀態(tài)

一般情況下，在真實(shí)的計(jì)算模擬應(yīng)用中，元胞狀態(tài)被定義為一個(gè)整數(shù)形式的離散集{k1,k2,…,kn-1,kn}，其中n是所有可能出現(xiàn)的狀態(tài)數(shù)。每一個(gè)整數(shù)都有一個(gè)被映射的現(xiàn)實(shí)意義，如狀態(tài)“1”代表“未降解”，狀態(tài)“0”表示“已降解”。

(3) 元胞空間

所謂元胞空間是指將元胞以網(wǎng)點(diǎn)的形式排列在線性，平面或三維歐幾里德空間之上所形成的集合。如何確立合適的元胞空間，對(duì)應(yīng)用元胞自動(dòng)機(jī)模型模擬復(fù)雜系統(tǒng)問(wèn)題來(lái)說(shuō)是非常重要的環(huán)節(jié)。

(4) 轉(zhuǎn)換規(guī)則

轉(zhuǎn)換規(guī)則是指根據(jù)此時(shí)刻當(dāng)前元胞的狀態(tài)和其鄰居元胞的狀態(tài)來(lái)推算下一時(shí)刻當(dāng)前元胞狀態(tài)的函數(shù)。該函數(shù)的輸入是當(dāng)前元胞及其鄰居的狀態(tài)，輸出是當(dāng)前元胞下一時(shí)刻應(yīng)該轉(zhuǎn)變成的狀態(tài)。簡(jiǎn)單地講，就是一個(gè)狀態(tài)轉(zhuǎn)移函數(shù)。

元胞自動(dòng)機(jī)的算法如圖1所示。

圖1 元胞自動(dòng)機(jī)方法示意圖Fig.1 Scheme of the Cellular Automata

2.3 基于MC與CA相結(jié)合的高聚物降解建模方法—元胞蒙特卡洛自動(dòng)機(jī)法(CMCA)

2.3.1 CMCA方法的思想與步驟

本文在蒙特卡洛方法與元胞自動(dòng)機(jī)方法的基礎(chǔ)上，提出元胞蒙特卡洛自動(dòng)機(jī)法(Cellular Monte Carlo Automatamethod，簡(jiǎn)稱(chēng) CMCA)，即將高聚物在幾何上進(jìn)行離散，在離散的元胞上模擬隨機(jī)鏈斷裂，實(shí)現(xiàn)微觀與介觀的降解過(guò)程。方法思想為：將要求解的物體離散為等份的元胞，元胞的尺寸越小，就越將問(wèn)題引入相應(yīng)的介觀或微觀尺度；在元胞上引入蒙特卡洛方法，用概率的方法隨機(jī)決定元胞內(nèi)的過(guò)程變化，元胞的演化除了自身的反應(yīng)，還要參考鄰居對(duì)這個(gè)元胞的影響，這樣所有元胞并行迭代演化，揭示介觀乃至微觀的動(dòng)態(tài)演化過(guò)程。CMCA算法的過(guò)程如圖2所示，具體方法如下。

圖2 元胞蒙特卡洛自動(dòng)機(jī)法 (CMCA) 執(zhí)行過(guò)程示意圖Fig.2 Scheme of the Cellular Monte Carlo Automata Method

以高聚物的降解為例，將高聚物離散為n×n等份的元胞(本文設(shè)n = 1000)，在每個(gè)元胞內(nèi)假設(shè)有一條高分子鏈，所有元胞的高分子鏈分子量呈正態(tài)分布；元胞的分子鏈按照蒙特卡洛方法的概率隨機(jī)選出進(jìn)行水解斷裂，一條分子鏈斷裂成兩條新的分子鏈，當(dāng)形成單體或低聚體(一般認(rèn)為聚合度在8以下)時(shí)，水解產(chǎn)物參與到水解反應(yīng)，出現(xiàn)了自催化現(xiàn)象，即鄰居元胞內(nèi)的低聚體會(huì)影響元胞的鏈斷裂方式及下一步水解反應(yīng)的概率，且低聚體會(huì)向外擴(kuò)散，形成高聚物的內(nèi)部孔洞，將此定為元胞的自動(dòng)演化規(guī)則，逐次迭代出高聚物降解的動(dòng)態(tài)演化過(guò)程。具體步驟為：

(1) 高聚物離散為n×n的元胞，元胞(i，j)在t時(shí)刻的狀態(tài)x(i，j，t)為：

元胞的鄰居選擇馮·諾依曼的4鄰居型，將元胞的上、下、左、右緊鄰的四個(gè)元胞作為當(dāng)前元胞的鄰居，如圖1所示，中間的黑色網(wǎng)格代表當(dāng)前元胞，周?chē)膫€(gè)灰色的網(wǎng)格表示其鄰居。

(2) 遍歷所有元胞，統(tǒng)計(jì)含有不同重復(fù)單元的分子鏈的條數(shù)xm(m = 1,2，…，M，M表示最長(zhǎng)鏈的聚合度)。

(3) 由式(2)、 (3)和隨機(jī)數(shù)r1, r2確定下一反應(yīng)的時(shí)刻 Δt及何種反應(yīng)μ，其中水解反應(yīng)是二級(jí)反應(yīng)，反應(yīng)式為：

故公式(1)可變?yōu)棣羦= πvxvxw，πv表示第v個(gè)反應(yīng)的反應(yīng)常數(shù)，xv是第v個(gè)反應(yīng)的聚合物分子數(shù)，xw代表水，假設(shè)水是充足的。用隨機(jī)的方法選取一條鏈將反應(yīng)μ對(duì)應(yīng)到元胞[i,j]。

(4) 元胞[i,j]發(fā)生鏈斷裂反應(yīng)，分子鏈數(shù)改變(新產(chǎn)生的鏈以鏈表表示)，分子量分布也隨之改變，元胞根據(jù)演變規(guī)則發(fā)生狀態(tài)改變。狀態(tài)演變規(guī)則制定如下。

分子量的變化將依據(jù)鏈的斷裂情況，本模型在鏈的斷裂模擬上采取末端斷裂與隨機(jī)斷裂相結(jié)合的方式，根據(jù)元胞自身和鄰居元胞有無(wú)低聚體產(chǎn)生作為判斷依據(jù)：有則認(rèn)為產(chǎn)生催化作用，鏈的斷裂采用隨機(jī)斷裂的方式；無(wú)則認(rèn)為無(wú)催化作用，鏈斷裂采用末端斷裂的方式。鏈斷裂后重新統(tǒng)計(jì)不同分子鏈的個(gè)數(shù)，重新計(jì)算各個(gè)元胞的分子量，如果元胞沒(méi)有降解，則狀態(tài)仍為固態(tài)；一旦元胞內(nèi)發(fā)生鏈的斷裂，則狀態(tài)由固態(tài)變?yōu)榻到鈶B(tài)。隨著降解進(jìn)一步進(jìn)行，如果元胞降解到一定程度，分子量低到某臨界值，分子鏈發(fā)生溶解，則元胞變成液體，由此認(rèn)為產(chǎn)生了孔隙。

統(tǒng)計(jì)液態(tài)元胞的數(shù)量Z(t)：

則可得到重量的損失率

(5) 元胞內(nèi)低聚體的擴(kuò)散

由第4步元胞內(nèi)高分子鏈斷裂，產(chǎn)生了新的低聚體，低聚體會(huì)向外擴(kuò)散，擴(kuò)散遵從Fick’s第二擴(kuò)散定律，采用Wang等[6]擴(kuò)散方程：

ce是聚合物酯鍵濃度(不能擴(kuò)散)，對(duì)應(yīng)于高分子材料的分子量(由鏈數(shù)可計(jì)算)；col是降解產(chǎn)生的低聚體(可以擴(kuò)散)的濃度；k1為無(wú)催化作用下水解反應(yīng)的速率常數(shù)，k2為有催化作用下水解反應(yīng)的速率常數(shù)，D為擴(kuò)散系數(shù)。D的計(jì)算是：D = D0+(1.3ε2-0.3ε3)*(D1-D0)，D0為低聚體在聚合物內(nèi)的擴(kuò)散系數(shù)，D1為低聚體在孔洞中的擴(kuò)散系數(shù)，ε為孔隙率(可由元胞的狀態(tài)統(tǒng)計(jì)計(jì)算)：

將n×n的元胞劃分粗粒度的塊，如本文n取1000時(shí)，劃分10*10個(gè)元胞為一塊，在塊與塊之間進(jìn)行擴(kuò)散。將微分方程離散為塊內(nèi)的差分計(jì)算，式(6)改寫(xiě)為：

即

?t為時(shí)間步長(zhǎng)。

(6) 更新時(shí)間t = t + ?t，回到第2步。

2.3.2 CMCA的算法描述

算法描述為：

3 計(jì)算與討論

3.1 計(jì)算

按以上程序執(zhí)行步驟，計(jì)算出數(shù)均分子量與初始數(shù)均分子量之比、重量損失、分子量分布隨時(shí)間變化如圖3所示。根據(jù)元胞狀態(tài)的演變得到孔隙的形成過(guò)程，如圖4所示。將模擬結(jié)果與實(shí)驗(yàn)值對(duì)比可見(jiàn)，本文提出的元胞蒙特卡洛自動(dòng)機(jī)方法(CMCA)模擬高聚物降解是正確可行的。

圖 3 CMCA 模型計(jì)算值與文獻(xiàn)[3]的實(shí)驗(yàn)值比對(duì)Fig.3 Comparison between the calculated values of the CMCA model and experimental values in the literature [3] (a) normalized molecular weight (b) weight loss (c) molecular weight distribution

圖4 CMCA方法模擬高聚物降解出現(xiàn)孔洞的過(guò)程Fig.4 Pore forming processes during the simulation of polymer degradation by CMCA

3.2 討論

3.2.1 分子鏈的斷裂方式

由模擬結(jié)果(圖3)可見(jiàn)，本文對(duì)于斷裂方式的選擇是比較合理的，即在鏈的斷裂模擬上采取末端斷裂與隨機(jī)斷裂相結(jié)合的方式。若完全采用末端斷裂的方式，將無(wú)法解釋分子量減小的現(xiàn)象，這是因?yàn)槟┒藬嗔逊绞礁赢a(chǎn)生短鏈，而且短鏈可能擴(kuò)散，使得平均分子量的減少十分緩慢。另一方面，若完全采用隨機(jī)斷裂的方式，將無(wú)法解釋質(zhì)量損失的現(xiàn)象，這是因?yàn)殡S機(jī)斷裂的方式產(chǎn)生可擴(kuò)散的短鏈太少。

3.2.2 對(duì)自催化反應(yīng)的看法

對(duì)于步驟(2)確定下一反應(yīng)的時(shí)間 Δt及何種反應(yīng)μ中的主導(dǎo)方程，選擇水攻擊高分子鏈的二級(jí)反應(yīng)，即在計(jì)算αv時(shí)，有αv= πvxvxw,，沒(méi)有選用文獻(xiàn)[6]中的反應(yīng)方程：

這項(xiàng)方程右邊第一項(xiàng)表示高聚物與水反應(yīng)，第二項(xiàng)表示產(chǎn)生的低聚體促進(jìn)了高聚物的降解，即自催化反應(yīng)，在元胞的分子鏈這樣微介觀尺度上，并不是每個(gè)元胞都能斷裂出低聚體，因此把部分元胞可能存在的自催化反應(yīng)納入下一反應(yīng)時(shí)刻的概率是不合理的。本文在所有高分子鏈都會(huì)受到水攻擊的均等機(jī)會(huì)上計(jì)算概率，確定下一反應(yīng)時(shí)刻，將元胞的鄰居元胞是否存在低聚體作為是否產(chǎn)生自催化反應(yīng)的依據(jù)，使低聚體攻擊高聚物產(chǎn)生隨機(jī)斷裂。

3.2.3 同種材料不同尺寸的對(duì)比

將CMCA模型用于預(yù)測(cè)與上例同種材料但不同尺寸的降解過(guò)程，結(jié)果如圖5所示。由圖可見(jiàn)，尺寸大的盤(pán)降解反而快于尺寸小的盤(pán)，這是由于在降解過(guò)程中，逐漸產(chǎn)生低聚體，呈酸性，相較于尺寸小的設(shè)備來(lái)說(shuō)，尺寸大的設(shè)備內(nèi)部的低聚體更加不易擴(kuò)散到外部，導(dǎo)致內(nèi)部的酸性較高，而本文又考慮酸的自催化作用，所以尺寸大的盤(pán)自催化作用更強(qiáng)，故而降解速度更快。

圖5 與文獻(xiàn)[3]同種材料的 2 mm 和 6 mm 盤(pán)降解過(guò)程模擬Fig.5 Degradation simulation of 2 mm and 6 mm discs with the same material in the literature [3] (a) normalized molecular weight (b) weight loss

3.2.4 模擬方法的效率

本文依據(jù)CMCA模型進(jìn)行高分子降解過(guò)程的仿真，采用計(jì)算機(jī)高級(jí)程序設(shè)計(jì)語(yǔ)言Java編寫(xiě)仿真程序，在配置為Intel(R) Core(TM) i7-4790 CPU，20.0 GB內(nèi)存的臺(tái)式PC機(jī)上，用時(shí)13 h完成了圖3所示的30周的降解模擬過(guò)程(包括元胞狀態(tài)圖圖4的繪制與計(jì)算數(shù)據(jù)文件的存儲(chǔ))。相較于傳統(tǒng)實(shí)驗(yàn)室方法，本研究方法充分利用了計(jì)算機(jī)速度快、成本低、可視化的特點(diǎn)，直觀的呈現(xiàn)實(shí)驗(yàn)結(jié)果并快速正確地揭示降解機(jī)理。

4 結(jié) 論

本文針對(duì)高聚物的降解建模中尺度單一、方法各異的現(xiàn)狀，提出了跨尺度的、離散與連續(xù)相結(jié)合的、隨機(jī)與確定相匹配的元胞蒙特卡洛自動(dòng)機(jī)方法(CMCA)。闡述了CMCA方法的建模思路與實(shí)現(xiàn)步驟，對(duì)高聚物DLPLA實(shí)例進(jìn)行了降解計(jì)算。從計(jì)算結(jié)果可見(jiàn)，在微(介)觀模擬高分子鏈斷裂、孔隙的形成及宏觀考慮低聚體擴(kuò)散等方面與實(shí)驗(yàn)結(jié)果吻合得很好，表明本文方法有效，并可以用于類(lèi)似問(wèn)題的解算。

[1] Kasuya K, Takagi K, Ishiwatari S, et al. Biodegradabilities of various aliphatic polyesters in natural waters [J]. Polymer Degradation and Stability, 1998, 59(1): 327-332.

[2] Siparsky G L, Voorhees K J, Miao F. Hydrolysis of polylactic acid (PLA) and polycaprolactone (PCL) in aqueous acetonitrile solutions: autocatalysis [J]. Journal of Environmental Polymer Degradation, 1998, 6(1): 31-41.

[3] Grizzi I, Garreau H, Li S, et al. Hydrolytic degradation of devices based on poly(DL-lactic acid) size-dependence [J]. Biomaterials, 1995, 16(4): 305-311.

[4] Rothstein S N, Federspiel W J, Little S R. A unified mathematical model for the prediction of controlled release from surface and bulk eroding polymer matrices [J]. Biomaterials, 2009, 30(8): 1657-1664.

[5] Li S, Garreau H, Vert M. Structure-property relationships in the case of the degradation of massive poly(α-hydroxy acids) in aqueous media [J]. Journal of Materials Science: Materials in Medicine, 1990, 1(4): 198-206.

[6] Wang Y, Pan J Z, Han X X, et al. A phenomenological model for the degradation of biodegradable polymers [J]. Biomaterials, 2008, 29(23): 3393-3401.

[7] Han X X, Pan J Z. A model for simultaneous crystallisation and biodegradation of biodegradable polymers [J]. Biomaterials, 2009, 30(3): 423-430.

[8] Siepmann J, Faisant N, Benoit J P. A new mathematical model quantifying drug release from bioerodiblemicroparticles using Monte Carlo simulations [J]. Pharmaceutical Research, 2002, 19(12): 1885-1893.

[9] Siepmann J, Siepmann F, Florence A T. Local controlled drug delivery to the brain: mathematical modeling of the underlying mass transport mechanisms [J]. International Journal of Pharmaceutics, 2006, 314(2): 101-119.

[10] WANG Xiao-peng (王小鵬), CHEN Tian-ning (陳天寧). Degradation and erosion model of biodegradable PLGA (可降解高聚物PLGA 降解溶蝕的仿真模型) [J]. Polymer Materials Science and Engineering (高分子材料科學(xué)與工程), 2012, 28(2): 174-178.

[11] YI Juan (易雋). Computer simulation of a polylactic acid degradation (聚乳酸降解的計(jì)算機(jī)模擬) [D]. Hangzhou (杭州)：Zhejiang University (浙江大學(xué)), 2008.

[12] Zygourakis K, Markenscoff P A. Computer-aided design of bioerodible devices with optimal release characteristics: a cellular automata approach [J]. Biomaterials, 1996, 17(2): 125-135.

[13] Bertrand N, Leclair G, Hildgen P. Modeling drug release from bioerodible microspheres using a cellular automaton [J]. International Journal of Pharmaceutics, 2007, 343(1): 196-207.

[14] Chao G, Xiaobo S, Chenglin C, et al. A cellular automaton simulation of the degradation of porous polylactide scaffold: I. Effect of porosity [J]. Materials Science and Engineering: C, 2009, 29(6): 1950-1958.

[15] Han X X, Pan J Z. Polymer chain scission, oligomer production and diffusion: a two-scale model for degradation of bioresorbablepolyesters [J]. Actabiomaterialia, 2011, 7(2): 538-547.

[16] Chen Y, Zhou S, Li Q. Mathematical modeling of degradation for bulk-erosive polymers: applications in tissue engineering scaffolds and drug delivery systems [J]. Actabiomaterialia, 2011, 7(3): 1140-1149.

Trans-Scale Modeling on Degradation of Biodegradable Polymers

GAO Xiao-hao1,2, ZHANG Tao-hong1,2, YANG Zhi-yong1,2, SUN Le-ran1,2, FAN Chun-mei1,2, ZHANG De-zheng1,2
(1. Department of Computer, School of Computer and Communication Engineering, University of Science and Technology Beijing, Beijing 100083, China；
2. Beijing Key Laboratory of Knowledge Engineering for Materials Science, Beijing 100083, China)

Biodegradable polymers have attracted great interests for their wide applications. Modeling of degradation process is limited within conditions like microscopic, mesoscopic or macroscopic scales, or continuous/discrete and stochastic/certain processes due to the complexity of the degradation process. A new Cellular Monte Carlo Automata method (CMCA) based on Monte Carlo method and Cellular Automata algorithm was studied, which integrated discrete and continuous, stochastic and certain conditions under different scales. The continuous problem was divided by cellular, the stochastic event was implemented by probability, and the evolution was updated by cellular iterations to reveal the degradation process from micro, mesoscopic to macroscopic dimensions. The concrete steps and the executive algorithm are explained in detail and the polyactic acid plate is used as an example. The calculation results prove the applicability of the model.

Monte Carlo method； cellularautomata； biomaterial； biodegradation polymer modeling；cellular Monte carlo Automata method

O63；TP391.9

A

10.3969/j.issn.1003-9015.2016.06.026

1003-9015(2016)06-1419-08

2015-11-04；

：2016-02-23。

中央高?；究蒲袠I(yè)務(wù)費(fèi)專(zhuān)項(xiàng)資金(FRF-BR-15-059A)；留學(xué)回國(guó)人員科研啟動(dòng)基金第48批；國(guó)家科技支撐計(jì)劃(2013BAI13B06)。

高曉浩(1989-)，男，湖北黃岡人，北京科技大學(xué)碩士生。

張?zhí)壹t，E-mail：waterswordzth@163.com