基于混合計(jì)算架構(gòu)的用于核素遷移控制方程的求解方法
【專(zhuān)利摘要】本發(fā)明公開(kāi)了一種基于混合計(jì)算架構(gòu)的用于核素遷移控制方程的求解方法，其包括以下步驟：讀取參數(shù)階段，讀入計(jì)算參數(shù)、邊界條件和計(jì)算網(wǎng)格文件；生成矩陣階段：生成剛度稀疏矩陣、時(shí)間演化稀疏矩陣和載荷向量，并生成對(duì)應(yīng)的稀疏矩陣，之后剛度稀疏矩陣、時(shí)間演化稀疏矩陣；循環(huán)計(jì)算求解階段，使用有限差分法循環(huán)計(jì)算，循環(huán)過(guò)程中采用Jacobi迭代法求解AX＝B的矩陣方程，得到的解向量X即為該時(shí)間點(diǎn)的核素濃度場(chǎng)數(shù)據(jù)。本發(fā)明不僅可以存儲(chǔ)更高網(wǎng)格體數(shù)規(guī)模下的數(shù)據(jù)矩陣，求解速度更快。
【專(zhuān)利說(shuō)明】基于混合計(jì)算架構(gòu)的用于核素遷移控制方程的求解方法

【技術(shù)領(lǐng)域】
[0001] 本發(fā)明涉及一種方法，特別是涉及一種基于混合計(jì)算架構(gòu)的用于核素遷移控制方 程的求解方法。

【背景技術(shù)】
[0002] 截止2012年，我國(guó)已經(jīng)有15個(gè)核反應(yīng)堆投入運(yùn)行，另有26個(gè)在建，為解決日益臨 近的核廢物處置問(wèn)題，我國(guó)已經(jīng)明確提出在2020年建成核廢物處置地下實(shí)驗(yàn)室，2050年建 成核廢物地下處置庫(kù)。核廢物處置的發(fā)展將直接影響和制約我國(guó)核事業(yè)和國(guó)防力量的規(guī)劃 和發(fā)展。
[0003] 如何有效安全的處置核廢物，不僅僅影響國(guó)家的發(fā)展和未來(lái)，對(duì)于個(gè)人和社會(huì)也 有著直接的影響，最顯著的例子是前蘇聯(lián)的切爾諾貝利核電站事故和今年的日本福島核電 站事故，這些鮮明的災(zāi)難性事故都說(shuō)明：核廢物一旦無(wú)法安全處置，將給人類(lèi)社會(huì)帶來(lái)災(zāi)難 性的后果。
[0004] 目前，人類(lèi)社會(huì)一般采用深埋（地下500-1000米左右）的方式來(lái)處置核廢物（主 要針對(duì)中高放廢物），其核心思想是通過(guò)各種工程手段、自然屏障來(lái)阻止延遲核素從固化體 向周?chē)h(huán)境（生物圈）遷移。
[0005] 目前核素遷移研究常見(jiàn)的研究手段包括：實(shí)驗(yàn)室模擬研究、天然類(lèi)比研究、示蹤研 究、計(jì)算機(jī)模擬研究。
[0006] 然而，在傳統(tǒng)的求解計(jì)算方法中，使用了全矩陣的方式來(lái)存儲(chǔ)和計(jì)算，這種方式在 一定網(wǎng)格規(guī)模下（例如：小于10萬(wàn)個(gè)網(wǎng)格體數(shù)）還是可行的，一旦網(wǎng)格規(guī)模超過(guò)這個(gè)規(guī)模， 一般的計(jì)算機(jī)已經(jīng)非常困難。
[0007] 現(xiàn)有的技術(shù)是讀取計(jì)算參數(shù)和網(wǎng)格數(shù)據(jù)，采用全矩陣記錄和存儲(chǔ)矩陣數(shù)據(jù)，采用 LU分解求解AX = B的矩陣方程，采用基于CPU的串行方式計(jì)算。這種技術(shù)的缺點(diǎn)是：這種 方法在一定規(guī)模的網(wǎng)格體數(shù)據(jù)之下（10萬(wàn)以下），在計(jì)算機(jī)上求解核素在緩沖/回填材料中 擴(kuò)散數(shù)學(xué)模型（即控制方程）的時(shí)間消耗和內(nèi)存消耗還可以接受，一旦超出這個(gè)范圍，由于 內(nèi)存的物理容量限制，會(huì)導(dǎo)致在整個(gè)求解過(guò)程中反復(fù)在內(nèi)存和硬盤(pán)之間交換數(shù)據(jù)，降低求 解速度，且由于CPU計(jì)算能力和設(shè)計(jì)目的的限制，浮點(diǎn)運(yùn)算耗時(shí)會(huì)讓用戶難以接受。


【發(fā)明內(nèi)容】

[0008] 本發(fā)明所要解決的技術(shù)問(wèn)題是提供一種基于混合計(jì)算架構(gòu)的用于核素遷移控制 方程的求解方法，其不僅可以存儲(chǔ)更高網(wǎng)格體數(shù)規(guī)模下的數(shù)據(jù)矩陣，節(jié)省了內(nèi)存開(kāi)銷(xiāo)，求解 速度更快。
[0009] 本發(fā)明是通過(guò)下述技術(shù)方案來(lái)解決上述技術(shù)問(wèn)題的：一種基于混合計(jì)算架構(gòu)的用 于核素遷移控制方程的求解方法，其特征在于，其包括以下步驟：
[0010] S1:讀取參數(shù)階段：讀入計(jì)算參數(shù)、邊界條件和計(jì)算網(wǎng)格文件；
[0011] S2 :生成矩陣階段：生成剛度稀疏矩陣、時(shí)間演化稀疏矩陣和載荷向量，并生成對(duì) 應(yīng)的稀疏矩陣，之后剛度稀疏矩陣、時(shí)間演化稀疏矩陣；
[0012] S3 :循環(huán)計(jì)算求解階段：使用有限差分法循環(huán)計(jì)算，循環(huán)過(guò)程中采用Jacobi迭代 法求解AX = B的矩陣方程，得到的解向量X即為該時(shí)間點(diǎn)的核素濃度場(chǎng)數(shù)據(jù)。
[0013] 優(yōu)選地，所述核素在緩沖/回填材料中遷移控制方程為：
[0014]

【權(quán)利要求】
1. 一種基于混合計(jì)算架構(gòu)的用于核素遷移控制方程的求解方法，其特征在于，其包括 以下步驟： 51 :讀取參數(shù)階段：讀入計(jì)算參數(shù)、邊界條件和計(jì)算網(wǎng)格文件； 52 :生成矩陣階段：生成剛度稀疏矩陣、時(shí)間演化稀疏矩陣和載荷向量，并生成對(duì)應(yīng)的 稀疏矩陣，之后剛度稀疏矩陣、時(shí)間演化稀疏矩陣； 53 :循環(huán)計(jì)算求解階段：使用有限差分法循環(huán)計(jì)算，循環(huán)過(guò)程中采用Jacobi迭代法求 解AX = B的矩陣方程，得到的解向量X即為該時(shí)間點(diǎn)的核素濃度場(chǎng)數(shù)據(jù)。
2. 如權(quán)利要求1所述的基于混合計(jì)算架構(gòu)的用于核素遷移控制方程的求解方法，其特 征在于，所述核素在緩沖/回填材料中遷移控制方程為：
其中C為核素濃度，Du為核素在多孔介質(zhì)中的動(dòng)力擴(kuò)散系數(shù)：Du = Dij'+Di/，
α τ, α ^分別為橫向和縱向彌散系數(shù)，j/為滲流場(chǎng)速度大小，' 為j方向的滲流速度分 量，D，= DdT δ m Dd為分子擴(kuò)散系數(shù)，T為多孔介質(zhì)的彎曲率； 這里Dd = D^xp ( σ Δ Τ), Ο。為參考溫度下的分子擴(kuò)散系數(shù)，σ為溫度修正系數(shù)，Δ T為 相對(duì)于參考溫度的溫度差，F(xiàn)為等效衰減系數(shù)，Q為源匯項(xiàng)。
3. 如權(quán)利要求1所述的基于混合計(jì)算架構(gòu)的用于核素遷移控制方程的求解方法，其特 征在于，所述S2包括以下步驟： 521 :估算稀疏矩陣大小，若是從1到計(jì)算網(wǎng)格中四面體數(shù)量的循環(huán)，則進(jìn)入S22,若不 是則進(jìn)入S210 ; 522 :獲取一個(gè)四面體的四個(gè)頂點(diǎn)坐標(biāo)，得到單元節(jié)點(diǎn)坐標(biāo)矩陣； 523 :計(jì)算雅克比矩陣； 524 :計(jì)算雅克比矩陣行列式值，計(jì)算雅克比矩陣逆矩陣； 525 :計(jì)算控制方程中衰變系數(shù)De對(duì)于剛度稀疏矩陣的影響，并寫(xiě)入到剛度稀疏矩陣 中； 526 :計(jì)算控制方程中滲流系數(shù)\對(duì)剛度稀疏矩陣的影響，并寫(xiě)入到剛度稀疏矩陣中； 527 :計(jì)算控制方程中擴(kuò)散系數(shù)對(duì)剛度稀疏矩陣的影響，并寫(xiě)入到剛度稀疏矩陣中； 528 :計(jì)算控制方程中源匯項(xiàng)Q對(duì)剛度稀疏矩陣的影響，并寫(xiě)入到剛度稀疏矩陣中； 529 :計(jì)算時(shí)間演化稀疏矩陣，若是從1到計(jì)算網(wǎng)格中四面體數(shù)量的循環(huán)，則進(jìn)入S22， 若不是則進(jìn)入S210 ; 5210 :若是從1到計(jì)算網(wǎng)格中三角形數(shù)量的循環(huán)則進(jìn)入S211，若不是則進(jìn)入S214 ; 5211 :獲取三角形點(diǎn)坐標(biāo)，得到三角形節(jié)點(diǎn)坐標(biāo)矩陣； 5212 :計(jì)算三角形面積，以及通量值； 5213 :計(jì)算通量對(duì)載荷向量的影響，若是從1到計(jì)算網(wǎng)格中三角形數(shù)量的循環(huán)則進(jìn)入 S211，若不是則進(jìn)入S214 ; 5214 :裝載矩陣結(jié)束。
4. 如權(quán)利要求1所述的基于混合計(jì)算架構(gòu)的用于核素遷移控制方程的求解方法，其特 征在于，所述生成矩陣階段步驟S21為估算稀疏矩陣大小，該步驟的實(shí)際含義為：在一定的 網(wǎng)格體數(shù)環(huán)境下，估算稀疏矩陣的元素個(gè)數(shù)，方便申請(qǐng)相應(yīng)的稀疏矩陣內(nèi)存，在網(wǎng)格體數(shù)為 1800左右的時(shí)候，即需要生成的全矩陣元素為1800*1800時(shí)，其非零元素的個(gè)數(shù)占全部元 素比例不會(huì)超過(guò)1 %，而在日常計(jì)算中，網(wǎng)格數(shù)一般在2000個(gè)以上，所以，此處可以大致取 非零元素比例為1%，即當(dāng)網(wǎng)格體數(shù)為N的時(shí)候，單一稀疏矩陣元素的總個(gè)數(shù)為N__ ZOT。= Ν*Ν*0. 01，在此過(guò)程中，稀疏矩陣的存儲(chǔ)采用行、列、值的格式存儲(chǔ)。
5. 如權(quán)利要求1所述的基于混合計(jì)算架構(gòu)的用于核素遷移控制方程的求解方法，其特 征在于，所述S3包括以下步驟： 531 :根據(jù)衰變系數(shù)，計(jì)算第一類(lèi)邊界條件點(diǎn)的濃度值； 532 :去除剛度稀疏矩陣、時(shí)間演化稀疏矩陣和載荷向量中的第一類(lèi)邊界條件點(diǎn)，得到 新的K n"稀疏矩陣、CnOT稀疏矩陣、FnOT向量； 533 :計(jì)算得到 tempi 和 temp2，templ = Cnew+Knew*inv*0. 5，temp2 = Knew*inv*0. 5_Cnew ; 334:計(jì)算8向量，8 = ?11"封鮮-仏111?2祁'，8'為上一次的8向量，如果為第一次，則8 為〇 ; 535 ：A = tempi ； 536 :使用Jacobi迭代法求解AX = B ; 537 :得到解向量X，將解向量X從GPU遷移至CPU，與第一類(lèi)邊界條件點(diǎn)濃度合并，形成 完整核素濃度場(chǎng)數(shù)據(jù)，寫(xiě)入數(shù)據(jù)文件； 538 :計(jì)算矩陣結(jié)束。
6. 如權(quán)利要求5所述的基于混合計(jì)算架構(gòu)的用于核素遷移控制方程的求解方法，其特 征在于，所述步驟S32至步驟S36均在GPU中處理。
【文檔編號(hào)】G06F19/00GK104156606SQ201410406546
【公開(kāi)日】2014年11月19日 申請(qǐng)日期:2014年8月18日 優(yōu)先權(quán)日:2014年8月18日
【發(fā)明者】劉濤, 趙冬梅, 李孟偉, 吳亞?wèn)|, 韓永國(guó) 申請(qǐng)人:西南科技大學(xué)