(計(jì)算機(jī)在材料科學(xué)與工程中的應(yīng)用)材料科學(xué)與工過程中典型物理場(chǎng)的數(shù)值模擬

計(jì)算機(jī)在材料科學(xué)與工程中的應(yīng)用第四章材料科學(xué)與工過程中典型

物理場(chǎng)的數(shù)值模擬引言引言引言材料制備和加工主要內(nèi)容：4.1溫度場(chǎng)的數(shù)學(xué)模型及求解

4.2濃度場(chǎng)的數(shù)學(xué)模型及求解4.3應(yīng)力場(chǎng)的數(shù)學(xué)模型簡介微分方程描述物理、化學(xué)、力學(xué)過程求解（定解條件）4.1溫度場(chǎng)的數(shù)學(xué)模型及求解4.1溫度場(chǎng)的數(shù)學(xué)模型及求解4.1溫度場(chǎng)的數(shù)學(xué)模型及求解加熱設(shè)備冷卻過程加熱方式電阻爐、燃料爐、浴爐、流化床爐、真空爐等感應(yīng)加熱、電子束加熱、激光外表處理、離子轟擊加熱等各種冷卻介質(zhì)的冷卻性能和各種冷卻方式高效率、節(jié)能的新加熱方法。各種冷卻介質(zhì)冷卻，各種冷卻方式。材料學(xué)中傳熱學(xué)重要課題24.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.1溫度場(chǎng)的根本知識(shí)傳熱方式——導(dǎo)熱〔熱傳導(dǎo)〕1〕定義指溫度不同的物體各局部或溫度不同的兩物體間直接接觸時(shí)，依靠分子、原子及自由電子等微觀粒子熱運(yùn)動(dòng)而進(jìn)行的熱量傳遞現(xiàn)象。導(dǎo)熱可以在固體、液體、氣體中發(fā)生。2〕導(dǎo)熱的特點(diǎn)：必須有溫差；物體直接接觸；依靠分子、原子及自由電子等微觀粒子熱運(yùn)動(dòng)而傳遞熱量；不發(fā)生宏觀的相對(duì)位移。3〕導(dǎo)熱機(jī)理氣體：氣體分子不規(guī)那么熱運(yùn)動(dòng)時(shí)相互碰撞的結(jié)果。導(dǎo)電固體：自由電子運(yùn)動(dòng)(類似氣體分子)。非導(dǎo)電固體：晶格結(jié)構(gòu)的振動(dòng)(彈性波)。液體：很復(fù)雜。類似于氣體；類似于非導(dǎo)電固體4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.1溫度場(chǎng)的根本知識(shí)導(dǎo)熱根本定律——Fourier(傅里葉)定律單位時(shí)間內(nèi)通過給定截面所傳遞的熱量，正比例于垂直于該截面法線方向上的溫度變化率，而熱量傳遞的方向與溫度升高的方向相反。qx熱流密度；λ－材料熱傳導(dǎo)系數(shù)〔熱導(dǎo)率〕W/(m·K)；負(fù)號(hào)表示傳熱方向與溫度梯度方向相反〔P84〕。熱導(dǎo)率是指具有單位溫度差〔1K〕的單位厚度的物體(1m)，在它的單位面積上(1m2)、每單位時(shí)間(1s)的導(dǎo)熱量(J)。熱導(dǎo)率表示材料導(dǎo)熱能力大??；物性參數(shù)。不同材料的導(dǎo)熱系數(shù)值不同，即使同一種材料導(dǎo)熱系數(shù)值與溫度等因素有關(guān)。金屬材料最高，良導(dǎo)電體，也是良導(dǎo)熱體，液體次之，氣體最小。對(duì)于鐵、碳鋼和低合金鋼，其值隨溫度的增加而下降；對(duì)于高合金鋼〔不銹鋼、耐熱鋼〕，其值隨溫度的增加而增加。

x方向上的溫度梯度4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.1溫度場(chǎng)的根本知識(shí)定義：流體中〔氣體或液體〕溫度不同的各局部之間，由于發(fā)生相對(duì)的宏觀運(yùn)動(dòng)而把熱量由一處傳遞到另一處的現(xiàn)象。伴隨著熱傳導(dǎo)：熱對(duì)流必然同時(shí)伴隨著熱傳導(dǎo)，自然界不存在單一的熱對(duì)流。對(duì)流換熱：流體與溫度不同的固體壁間接觸時(shí)的熱量交換過程。分類：有相變/無相變對(duì)流換熱；自然對(duì)流/強(qiáng)制對(duì)流/沸騰換熱及凝結(jié)換熱傳熱方式——熱對(duì)流1〕定義與特征4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.1溫度場(chǎng)的根本知識(shí)2〕對(duì)流換熱根本定律——牛頓冷卻公式h—對(duì)流換熱系數(shù)加熱時(shí)：冷卻時(shí)：h的物理意義：對(duì)流換熱系數(shù)，比例系數(shù)(外表傳熱系數(shù))。外表傳熱系數(shù)的大小與傳熱過程中的許多因素有關(guān)。它不僅取決于物體的物性、換熱外表的形狀、大小相對(duì)位置，而且與流體的流速有關(guān)。一般地，就介質(zhì)而言：水的對(duì)流換熱比空氣強(qiáng)烈；就換熱方式而言：有相變的強(qiáng)于無相變的；強(qiáng)制對(duì)流強(qiáng)于自然對(duì)流。對(duì)流換熱研究的根本任務(wù)：用理論分析或?qū)嶒?yàn)的方法推出各種場(chǎng)合下外表換熱導(dǎo)數(shù)的關(guān)系式。4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.1溫度場(chǎng)的根本知識(shí)1〕定義

物體通過電磁波來傳遞能量的方式稱為輻射；因熱的原因發(fā)出熱輻射能的現(xiàn)象稱為熱輻射，高溫物體失去熱量而低溫物體得到熱量，這種傳熱方式叫輻射傳熱。輻射熱流量用斯蒂芬-波爾茲曼定律表示。2〕條件：溫差、發(fā)射電磁波。3〕取決于：兩物體空間位置〔角度系數(shù)〕、外表特性〔黑度〕。物體的溫度越高、輻射能力越強(qiáng)；假設(shè)物體的種類不同、外表狀況不同，其輻射能力不同傳熱方式——熱輻射4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.1溫度場(chǎng)的基本知識(shí)4〕特點(diǎn)輻射換熱是一個(gè)動(dòng)態(tài)過程，處于熱平衡時(shí)，輻射換熱量為0.不需要冷熱物體直接接觸：即：不需要中間介質(zhì)的存在，在真空中就可以傳遞能量，在真空中輻射能的傳遞最有效。在輻射換熱過程中伴隨著能量形式的轉(zhuǎn)換：物體熱力學(xué)能電磁波能物體熱力學(xué)能。是一種雙向熱流同時(shí)存在的換熱過程。無論溫度上下，物體都在不停地相互發(fā)射電磁波能、相互輻射能量；高溫物體輻射給低溫物體的能量大于低溫物體輻射給高溫物體的能量;總的結(jié)果是熱由高溫傳到低溫。物體的輻射能力與其溫度性質(zhì)有關(guān)。區(qū)別于熱傳導(dǎo)、熱對(duì)流的根本特點(diǎn)。4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.1溫度場(chǎng)的根本知識(shí)斯蒂芬-玻爾茲曼定律〔Stefan-Boltzmannlaw〕黑體：能全部吸收投射到其外表輻射能的物體?；蚍Q絕對(duì)黑體?！睟lackbody〕——吸收率等于1的物體黑體的輻射能力與吸收能力在同溫度的物體中最強(qiáng)。1〕黑體向外發(fā)射的輻射能：—黑體表面的絕對(duì)溫度（熱力學(xué)溫度）—絕對(duì)黑體輻射力—斯蒂芬-玻爾茲曼常數(shù)，2〕實(shí)際物體輻射能力〔低于同溫度黑體〕：實(shí)際物體外表的發(fā)射率〔黑度〕，0-1；與物體的種類、外表狀況和溫度有關(guān)。4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.1溫度場(chǎng)的基本知識(shí)溫度場(chǎng)4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.2溫度場(chǎng)的數(shù)學(xué)模型導(dǎo)熱微分方程的建立建立導(dǎo)熱微分方程，可以揭示連續(xù)溫度場(chǎng)隨空間坐標(biāo)和時(shí)間變化的內(nèi)在聯(lián)系。理論根底：傅里葉定律+能量守恒方程溫度場(chǎng)：Temperaturefield。物質(zhì)系統(tǒng)內(nèi)各個(gè)點(diǎn)上溫度的集合稱為溫度場(chǎng)，它是時(shí)間和空間坐標(biāo)的函數(shù)?！矊?dǎo)入微元體的總熱流量〕+〔微元體內(nèi)熱源生成熱〕=〔導(dǎo)出微元體的總熱流量〕+〔微元體內(nèi)能增加〕微元體內(nèi)熱源生成熱:

其中：ρ、c、Q、τ——微元體的密度、比熱容、單位時(shí)間內(nèi)單位體積內(nèi)熱源的生成熱及時(shí)間。導(dǎo)入微元體的總熱流量Q1=Qx+Qy+Qz；導(dǎo)出微元體的總熱流量Q2=Qx++dx+Qy+dy+Qz+dz笛卡爾坐標(biāo)系中三維非穩(wěn)態(tài)導(dǎo)熱微分方程的一般表達(dá)式4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.2溫度場(chǎng)的數(shù)學(xué)模型4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.1溫度場(chǎng)的根本知識(shí)笛卡爾坐標(biāo)系中三維非穩(wěn)態(tài)導(dǎo)熱微分方程物理意義：反映了物體的溫度隨時(shí)間和空間的變化關(guān)系。能量方程是目前溫度場(chǎng)數(shù)值模擬中普遍使用的描述方程，它不僅適用于固體，也適用于流體。其中，ρ為材料的密度〔kg/m3〕；cp為材料的比熱容〔J/(kgK)〕；τ為時(shí)間〔s〕；λ分別為材料沿x，y，z方向的熱導(dǎo)率〔W/(m*K)〕；q為材料內(nèi)部的熱源密度〔W/kg〕。上式中，第一項(xiàng)為體元升溫需要的熱量；右側(cè)第一、二和三項(xiàng)是由x，y和z方向流入體元的熱量；最后一項(xiàng)體元內(nèi)熱源產(chǎn)生的熱量。微分方程的物理意義：體元升溫所需的熱量應(yīng)該等于流入體元的熱量與體元內(nèi)產(chǎn)生的熱量的總和。4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.2溫度場(chǎng)的數(shù)學(xué)模型假設(shè)邊界條件和內(nèi)部熱源密度Q不隨時(shí)間變化那么經(jīng)過一定時(shí)間后物體內(nèi)部各點(diǎn)的溫度將到達(dá)平衡，即有穩(wěn)態(tài)熱傳導(dǎo)方程：根據(jù)系統(tǒng)有無內(nèi)熱源、是否導(dǎo)熱過程為穩(wěn)態(tài)導(dǎo)熱，以及一維、二維和三維的情況，可進(jìn)行相應(yīng)的簡化。三維穩(wěn)態(tài)熱傳導(dǎo)方程為：4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.2溫度場(chǎng)的數(shù)學(xué)模型二維穩(wěn)態(tài)熱傳導(dǎo)方程為：一維非穩(wěn)態(tài)熱傳導(dǎo)方程為：一維穩(wěn)態(tài)熱傳導(dǎo)方程為：二維非穩(wěn)態(tài)熱傳導(dǎo)方程為：4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.2溫度場(chǎng)的數(shù)學(xué)模型討論：1〕直角坐標(biāo)下有內(nèi)熱源的非穩(wěn)態(tài)導(dǎo)熱微分方程為常數(shù)時(shí):其中—稱擴(kuò)散系數(shù)〔熱擴(kuò)散率〕。2〕直角坐標(biāo)下無內(nèi)熱源的穩(wěn)態(tài)導(dǎo)熱微分方程，=0，且時(shí)3〕常物性、穩(wěn)態(tài)、無內(nèi)熱源，假設(shè)λ為常數(shù)時(shí)，且屬穩(wěn)態(tài)，即：4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.2溫度場(chǎng)的數(shù)學(xué)模型4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.2溫度場(chǎng)的數(shù)學(xué)模型圓柱坐標(biāo)下有內(nèi)熱源的非穩(wěn)態(tài)導(dǎo)熱微分方程：圓柱坐標(biāo)：把直角坐標(biāo)的xy平面變換為極坐標(biāo)，而z軸不變。圓柱坐標(biāo)與直角坐標(biāo)之間的變換關(guān)系：4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.2溫度場(chǎng)的數(shù)學(xué)模型球坐標(biāo)下有內(nèi)熱源的非穩(wěn)態(tài)導(dǎo)熱微分方程：球面坐標(biāo)與直角坐標(biāo)之間的變換關(guān)系：4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.2溫度場(chǎng)的數(shù)學(xué)模型初始條件和邊界條件1）初始條件初始條件是指求解問題的初始溫度場(chǎng)，也就是在零時(shí)刻溫度場(chǎng)的分布。它可以是均勻的，此時(shí)有也可以是不均勻的，各點(diǎn)的溫度值已知或者遵從某一函數(shù)關(guān)系。4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.2溫度場(chǎng)的數(shù)學(xué)模型2〕邊界條件邊界條件是指物體外表或者邊界與周圍環(huán)境的熱交換情況，通常有三類重要的邊界條件。〔1〕第一類邊界條件第一類邊界條件是指物體邊界上的溫度分布函數(shù)，表示為或〔2〕第二類邊界條件第二類邊界條件是指邊界上的熱流密度，表示為：或n為物體邊界的外法線方向，并規(guī)定熱流密度的方向與邊界的外法線方向相同。〔3〕第三類邊界條件又稱為對(duì)流邊界條件，是指物體與其周圍環(huán)境介質(zhì)間的對(duì)流傳熱系數(shù)k和介質(zhì)的溫度，表示為：4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.2溫度場(chǎng)的數(shù)學(xué)模型邊界上的溫度分布函數(shù)邊界上的熱流密度對(duì)流傳熱系數(shù)和介質(zhì)溫度4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.2溫度場(chǎng)的數(shù)學(xué)模型4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.3溫度場(chǎng)的有限差分求解導(dǎo)熱問題的求解是對(duì)導(dǎo)熱微分方程在邊界條件和初始條件下積分求解，即解析解。但目前由于數(shù)學(xué)上的困難，在工程實(shí)際中許多問題還不能采用分析解法進(jìn)行求解，如物體的幾何形狀比較復(fù)雜及邊界形狀不規(guī)那么，材料的物性常數(shù)隨溫度變化等。近年來，隨著計(jì)算機(jī)技術(shù)和計(jì)算技術(shù)的迅速開展，數(shù)值方法已經(jīng)得到廣泛應(yīng)用并成為有力的輔助求解工具，已開展了許多用于工程問題求解的數(shù)值計(jì)算方法。數(shù)值積分法、有限差分法和有限元法離散數(shù)學(xué)為根底計(jì)算機(jī)作為工具4.1溫度場(chǎng)的數(shù)學(xué)模型及求解——4.1.3溫度場(chǎng)的有限差分求解有限差分法1〕把原來物體內(nèi)隨時(shí)間和空間連續(xù)分布的溫度問題轉(zhuǎn)化為在時(shí)間和空間領(lǐng)域內(nèi)有限個(gè)離散點(diǎn)的溫度值問題，再用這些離散點(diǎn)上的溫度值去逼近連續(xù)的溫度分布。2〕有限差分是以差分代替微分，差商代替微商，建立以節(jié)點(diǎn)溫度為未知量的線性代數(shù)方程組，然后求解得到各節(jié)點(diǎn)溫度的近似值。4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.3溫度場(chǎng)的有限差分求解網(wǎng)格劃分網(wǎng)格劃分二維穩(wěn)態(tài)無內(nèi)熱源導(dǎo)熱教材P9214.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.3溫度場(chǎng)的有限差分求解4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.3溫度場(chǎng)的有限差分求解4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.3溫度場(chǎng)的有限差分求解4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.3溫度場(chǎng)的有限差分求解4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.3溫度場(chǎng)的有限差分求解二維穩(wěn)態(tài)導(dǎo)熱1差分方程建立構(gòu)成矩陣形式4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.3溫度場(chǎng)的有限差分求解二維穩(wěn)態(tài)問題求解1Eg1.邊界點(diǎn)溫度，且區(qū)域內(nèi)無內(nèi)熱源。利用有限差分方法來計(jì)算節(jié)點(diǎn)1、2、3的溫度。Eg.邊界點(diǎn)溫度，且區(qū)域內(nèi)無內(nèi)熱源。利用有限差分方法來計(jì)算節(jié)點(diǎn)1、2、3的溫度。解方程得4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.3溫度場(chǎng)的有限差分求解4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.3溫度場(chǎng)的有限差分求解實(shí)際工作中遇到的導(dǎo)熱問題通常為非穩(wěn)態(tài)導(dǎo)熱，其特點(diǎn)是溫度不僅隨空間坐標(biāo)的變化而變化，而且還隨時(shí)間的變化而變化。因此．溫度場(chǎng)的分布與時(shí)間和位置兩個(gè)因素有關(guān)。非穩(wěn)態(tài)問題的求解原理、離散化方法和主要求解步驟與穩(wěn)態(tài)問題的求解類似，但由于非穩(wěn)態(tài)導(dǎo)熱中增加了時(shí)問變量，因此，在差分格式、解的特性以及求解方法上都要復(fù)雜一些。如在區(qū)域離散化中，不僅包括空間區(qū)域的離散化，還有時(shí)間區(qū)域的離散化。一塊無限大平板〔如圖4-11所示〕，其一半厚度為L=0.1m，初始溫度t0=1000℃，突然將其插入溫度t∞=20℃的流體介質(zhì)中。平板的導(dǎo)熱系數(shù)λ=34.89W/m℃，密度ρ=7800kg/m3，比熱c=0.712J/kg℃，平板與介質(zhì)的對(duì)流換熱系數(shù)為h=233W/m2?℃，求平板內(nèi)各點(diǎn)的溫度分布。無限大平板非穩(wěn)態(tài)導(dǎo)熱非穩(wěn)態(tài)問題的有限差分格式24.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.3溫度場(chǎng)的有限差分求解1〕數(shù)學(xué)描述，由于平板換熱關(guān)于中心線是對(duì)稱的，僅對(duì)平板一半?yún)^(qū)域進(jìn)行計(jì)算即可。坐標(biāo)x的原點(diǎn)選在平板中心線上，因而一半?yún)^(qū)域的非穩(wěn)態(tài)導(dǎo)熱的數(shù)學(xué)描述為：該數(shù)學(xué)模型的解析解為:其中,為方程的根，。(4.34)(4.35)(4.36)(4.37)(4.38)4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.3溫度場(chǎng)的有限差分求解表4-1給出了在平板外表(x=L)處由上式計(jì)算得到的不同時(shí)刻的溫度值。表4-1平板外表各不同時(shí)刻溫度值時(shí)間（S）12345678910溫度(℃)981.84974.47968.88964.20960.11956.14953.08949.97947.07944.342〕微分方程的離散一維非穩(wěn)態(tài)導(dǎo)熱指的是空間坐標(biāo)是一維的。假設(shè)考慮時(shí)間坐標(biāo)，那么所謂的一維非穩(wěn)態(tài)導(dǎo)熱實(shí)際上是二維問題(見以下圖)，即：有時(shí)間坐標(biāo)τ和空間坐標(biāo)x兩個(gè)變量。但要注意，時(shí)間坐標(biāo)是單向的，就是說，前一時(shí)刻的狀態(tài)會(huì)對(duì)后一時(shí)刻的狀態(tài)有影響，但后一時(shí)刻的狀態(tài)卻影響不到前一時(shí)刻，如圖示出了以x和τ為坐標(biāo)的計(jì)算區(qū)域的離散，時(shí)間從τ=0開始，經(jīng)過一個(gè)個(gè)時(shí)層增加到K時(shí)層和K+1時(shí)層。平面區(qū)域的時(shí)間和空間離散4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.3溫度場(chǎng)的有限差分求解平面區(qū)域的時(shí)間和空間離散對(duì)于i節(jié)點(diǎn)，在K和K+1時(shí)刻可將微分方程寫成下面式子：將上式的左端溫度對(duì)時(shí)間的偏導(dǎo)數(shù)進(jìn)行差別離散為：這兩個(gè)式子的右端差分式完全相同，但在兩個(gè)式子中卻有不同含義。對(duì)1式，右端項(xiàng)相對(duì)i點(diǎn)在K時(shí)刻的導(dǎo)數(shù)是向前差分。而在2式中，右端項(xiàng)的向后差分。是i點(diǎn)在K+1時(shí)刻的導(dǎo)數(shù)(4.39)(4.40)(4.41)(4.42)4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.3溫度場(chǎng)的有限差分求解將式(4.41)和(4.42)分別代入式(4.39)和(4.40)，并將式(4.39)和(4.40)右端關(guān)于x的二階導(dǎo)數(shù)用相應(yīng)的差分〔二階中心差分〕代替，那么可得到以下顯式和隱式兩種不同的差分格式：顯式:(K=0,1,2,………,i=2,3,…,N-1)全隱式：(K=0,1,2,………,i=2,3,…,N-1)以上兩式中的

從式(4.43)可見，其右端只涉及K時(shí)刻的溫度，當(dāng)從K=0〔即τ=0時(shí)刻〕開始計(jì)算時(shí)，在K=0時(shí)等號(hào)右端都是值，因而直接可計(jì)算出K=1時(shí)刻各點(diǎn)的溫度〔顯示〕。由K=1時(shí)刻的各點(diǎn)的溫度值，又可以直接利用式(4.43)計(jì)算K=2時(shí)刻的各點(diǎn)的溫度，這樣一個(gè)時(shí)層一個(gè)時(shí)層的往下推，各時(shí)層的溫度都能用式(4.43)直接計(jì)算出來，不要求解代數(shù)方程組。而對(duì)于式(4.44)等號(hào)右端包含了與等號(hào)左端同一時(shí)刻但不同節(jié)點(diǎn)的溫度，因而必須通過求解代數(shù)方程組才能求得這些節(jié)點(diǎn)的溫度值〔隱式〕。(4.43)(4.44)4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.3溫度場(chǎng)的有限差分求解顯式：對(duì)兩時(shí)間層格式，知道n時(shí)刻各空間層上函數(shù)值而推的下一時(shí)刻的值。

而隱式格式中，知道n時(shí)刻某空間層上函數(shù)值，而不能得出待求的下一時(shí)刻各空間層上之值，必須聯(lián)立求解一個(gè)與網(wǎng)格點(diǎn)數(shù)目相同的方程組。

顯式格式可以直接求解，缺點(diǎn)是計(jì)算穩(wěn)定性差。隱式格式的優(yōu)點(diǎn)是穩(wěn)定性好，只是求解麻煩。

4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.3溫度場(chǎng)的有限差分求解3)邊界條件的離散對(duì)于式(4.36)和(4.37)所給出的邊界條件，可以直接用差分代替微分，也可以用元體平衡法給出相應(yīng)的邊界條件，亦有顯式和隱式之分。通常，當(dāng)內(nèi)部節(jié)點(diǎn)采用顯式時(shí)，邊界節(jié)點(diǎn)也用顯式離散；內(nèi)部節(jié)點(diǎn)用隱式時(shí)，邊界節(jié)點(diǎn)亦用隱式。邊界節(jié)點(diǎn)的差分格式是顯示還是隱式，取決于如何與內(nèi)部節(jié)點(diǎn)的差分方程組合。用K+1時(shí)刻相應(yīng)節(jié)點(diǎn)的差分，代替式(4.36)和(4.37)中的微分，可得到邊界節(jié)點(diǎn)的差分方程：(4.45)最終的離散格式顯式：初始值：i=1,2,3,…,N)(i=2,3,…,N-1)(4.46)(4.47)4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.3溫度場(chǎng)的有限差分求解其中K=0,1,2,……(4.48)隱式：初始值：(i=2,3,…,N-1)其中K=0,1,2,……(4.49)(4.50)(4.51)(4.52)在用隱式差分計(jì)算時(shí)，每個(gè)時(shí)層都需要迭代求解代數(shù)方程組〔4.49〕-式〔4.52〕。在每個(gè)時(shí)層計(jì)算時(shí)，都要先假定一個(gè)溫度場(chǎng)〔一般取上一時(shí)層的溫度場(chǎng)為本時(shí)層的初始場(chǎng)〕，然后迭代計(jì)算直至收斂。顯式差分格式：每個(gè)節(jié)點(diǎn)方程可以獨(dú)立求解，但需要考慮穩(wěn)定性。得到：(4.53)4.1溫度場(chǎng)的數(shù)學(xué)模型及求解—4.1.3溫度場(chǎng)的有限差分求解二維不穩(wěn)態(tài)焊接熱傳導(dǎo)問題求解2

Eg2.薄板焊接中移動(dòng)熱源二維焊接離散化模型(Qm-最大熱源密度；h-板厚；r-離開熱源中心距離)簡單傳熱學(xué)計(jì)算機(jī)分析二維不穩(wěn)態(tài)焊接熱傳導(dǎo)問題求解2

Eg2.薄板焊接中移動(dòng)熱源二維焊接離散化模型邊界條件簡單傳熱學(xué)計(jì)算機(jī)分析二維不穩(wěn)態(tài)焊接熱傳導(dǎo)問題求解2根據(jù)二維不穩(wěn)定導(dǎo)熱方程(1)以及初始條件、邊界條件,可建立差分方程(2)根據(jù)差分方程以及各邊界節(jié)點(diǎn)差分方程,在滿足穩(wěn)定條件下就可求出不同時(shí)刻的溫度分布。簡單傳熱學(xué)計(jì)算機(jī)分析二維不穩(wěn)態(tài)焊接熱傳導(dǎo)問題求解2時(shí)間步長DS=0.01s；熱源作用時(shí)間S1=0sSe=5sβ=0.0008cal/cm2s℃；Qm=2500cal/cm；板厚H=1cm；焊速v=0.05cm/s節(jié)點(diǎn)數(shù)N1=N2=10,T0=Ta=20℃,比熱C=0.16cal/g℃;密度ρ=7.82g/cm2;導(dǎo)熱系數(shù)k=0.1cal/cms℃,節(jié)點(diǎn)間距DX=0.2cm輸入?yún)?shù)簡單傳熱學(xué)計(jì)算機(jī)分析二維不穩(wěn)態(tài)焊接熱傳導(dǎo)問題求解2計(jì)算框圖簡單傳熱學(xué)計(jì)算機(jī)分析二維不穩(wěn)態(tài)焊接熱傳導(dǎo)問題求解20.5s溫度場(chǎng)簡單傳熱學(xué)計(jì)算機(jī)分析二維不穩(wěn)態(tài)焊接熱傳導(dǎo)問題求解23.0s溫度場(chǎng)簡單傳熱學(xué)計(jì)算機(jī)分析MatlabPDE二維不穩(wěn)態(tài)焊接熱傳導(dǎo)求解3區(qū)域設(shè)置單擊工具，在窗口拉出一個(gè)矩形，雙擊矩形區(qū)域，在ObjectDialog對(duì)話框輸入Left為0，Bottom為0，Width為2,Height為2。與默認(rèn)的坐標(biāo)相比，圖形小的看不見，所以要調(diào)整坐標(biāo)顯示比例。方法：選擇Options->AxesLimits,把X，Y軸的自動(dòng)選項(xiàng)打開。設(shè)置Options->Application為HeatTransfer

（設(shè)置程序應(yīng)用熱傳輸模型）簡單傳熱學(xué)計(jì)算機(jī)分析MatlabPDE二維不穩(wěn)態(tài)焊接熱傳導(dǎo)求解3邊界條件設(shè)置單擊使邊界變紅色，然后分別雙擊每段邊界，打開BoundaryConditions對(duì)話框，設(shè)置邊界條件（根據(jù)邊界條件）。所有的邊界都為Neumann條件。簡單傳熱學(xué)計(jì)算機(jī)分析MatlabPDE二維不穩(wěn)態(tài)焊接熱傳導(dǎo)求解3方程類型設(shè)置單擊，在PDESpecification對(duì)話框中設(shè)置 方程類型為Parabolic（拋物型），

rho(密度）為7.82，C（比熱）為0.16， k（導(dǎo)熱系數(shù)）為0.1， Q（熱源）為4000*exp(-3*(x.^2+(y-0.4*t).^2)/0.49)簡單傳熱學(xué)計(jì)算機(jī)分析MatlabPDE二維不穩(wěn)態(tài)焊接熱傳導(dǎo)求解3網(wǎng)格設(shè)置初值和誤差設(shè)置單擊Solve菜單中Parameters…選項(xiàng)，打開SolveParameters對(duì)話框，輸入Time為0:0.5:5，u(t0)為20,其他不變。解方程單擊，或者加密網(wǎng)格，單擊。單擊。簡單傳熱學(xué)計(jì)算機(jī)分析7)整理數(shù)據(jù)單擊Mesh->ExportMesh…輸出pet的數(shù)值，單擊Solve->ExportSolution…輸出u.回到Matlab主窗口執(zhí)行下面兩條命令：u1=[p',u(:,7)]%將節(jié)點(diǎn)坐標(biāo)和其在3s時(shí)的溫度組成新矩陣u2=sortrows(u1,3)%將u1按溫度值大小升序排列。u1=[p',u(:,4)]%將節(jié)點(diǎn)坐標(biāo)和其在1.5s時(shí)的溫度組成新矩陣u2=sortrows(u1,3)%將u1按溫度值大小升序排列。MatlabPDE二維不穩(wěn)態(tài)焊接熱傳導(dǎo)求解31.5s時(shí)溫度場(chǎng)簡單傳熱學(xué)計(jì)算機(jī)分析MatlabPDE二維不穩(wěn)態(tài)焊接熱傳導(dǎo)求解33s時(shí)溫度場(chǎng)簡單傳熱學(xué)計(jì)算機(jī)分析4.2濃度場(chǎng)的數(shù)學(xué)模型及求解4.2濃度場(chǎng)的數(shù)學(xué)模型及求解任何不均質(zhì)的材料，在—定的熱力學(xué)條件下，都將趨向于均勻化。譬如，通過擴(kuò)散退火可以改善因凝固帶來的成分不均勻性．這是在合金中分布不均勻的溶質(zhì)原子從高濃度區(qū)域向低濃度區(qū)域運(yùn)動(dòng)(擴(kuò)散)的結(jié)果。所以固態(tài)中的擴(kuò)散本質(zhì)是在擴(kuò)散力(濃度、電場(chǎng)、應(yīng)力場(chǎng)等的梯度)作用下，原于定向、宏觀的遷移。這種遷移運(yùn)動(dòng)的結(jié)果是使系統(tǒng)的化學(xué)自由能下降。材料的擴(kuò)散現(xiàn)象在工程中廣泛存在，如壓力加工時(shí)的動(dòng)態(tài)恢復(fù)再結(jié)晶，雙全屬板的生產(chǎn)、焊接過程，熱處理中的相變，化學(xué)熱處理以及粉末冶金的燒結(jié)等。擴(kuò)散理論的研究主要方面之—是宏觀規(guī)律的研究，它重點(diǎn)討論擴(kuò)散物質(zhì)的濃度分布與時(shí)間的關(guān)系，根據(jù)不同條件建立一系列的擴(kuò)散方程，并技其邊界條件求解。用計(jì)算機(jī)數(shù)值計(jì)算方法代替?zhèn)鹘y(tǒng)的、復(fù)雜的數(shù)學(xué)物理方程對(duì)濃度場(chǎng)問題進(jìn)行研究已成為開展的趨勢(shì)。4.2濃度場(chǎng)的數(shù)學(xué)模型及求解上坡擴(kuò)散上坡擴(kuò)散是沿著濃度升高的方向進(jìn)行擴(kuò)散，即由低濃度向高濃度方向擴(kuò)散，使?jié)舛劝l(fā)生兩級(jí)分化。例如，奧氏體向珠光體轉(zhuǎn)變的過程中，碳原子由濃度較低的奧氏體向濃度較高的滲碳體擴(kuò)散，就是上坡擴(kuò)散。例如，將碳含量相近的碳鋼〔Wc=0.441%〕與硅鋼〔Wc=0.478%〕對(duì)焊在一起，在1050℃進(jìn)行長時(shí)間(13d)加熱后，焊接面處兩側(cè)的碳濃度分布如圖。硅的存在推高了碳的活度和化學(xué)位，從而驅(qū)使碳從含硅的一側(cè)向另一側(cè)上坡擴(kuò)散。發(fā)生濃度低的向濃度高的方向擴(kuò)散，產(chǎn)生成分的偏聚而不是成分的均勻化。硅鋼側(cè)碳鋼側(cè)4.2濃度場(chǎng)的數(shù)學(xué)模型及求解菲克第一定律穩(wěn)定濃度場(chǎng)模型如圖4-21和4-22所示，x軸上兩單位面積1〔A〕和2〔B〕，間距dx，面上原子濃度為CA、CB，那么平面1到平面2上原子數(shù)n1=C1dx，平面2到平面1上原子數(shù)n2=C2dx，假設(shè)原子平均跳動(dòng)頻率f,dτ時(shí)間內(nèi)跳離平面1的原子數(shù)為n1f·dτ，跳離平面2的原子數(shù)為n2fdτ。4.2濃度場(chǎng)的數(shù)學(xué)模型及求解菲克第一定律(一定時(shí)間內(nèi)，濃度不隨時(shí)間變化dC/dτ=0)，單位時(shí)間內(nèi)通過垂直于擴(kuò)散方向的單位截面積的擴(kuò)散物質(zhì)流量〔擴(kuò)散通量〕與該面積處的濃度梯度成正比。定義：組分i每單位時(shí)間通過單位面積的質(zhì)量傳輸量正比于濃度梯度。定義式：

D—擴(kuò)散系數(shù)。負(fù)號(hào)表示質(zhì)量傳輸?shù)姆较蚺c濃度梯度的方向相反。J—擴(kuò)散通量，g/cm2*s式中負(fù)號(hào)說明擴(kuò)散通量的方向與濃度梯度方向相反。可見，只要存在濃度梯度，就會(huì)引起原子的擴(kuò)散，物體的擴(kuò)散系數(shù)單位：m2/s，物理意義：單位傳質(zhì)量相當(dāng)于單位濃度梯度下的擴(kuò)散傳質(zhì)通量。影響因素：物體的種類，物體的結(jié)構(gòu)，溫度、壓力等。4.2濃度場(chǎng)的數(shù)學(xué)模型及求解菲克第二定律解決溶質(zhì)濃度隨時(shí)間變化的情況，即dc/dt≠0。如圖：單元模型那么單元體積中溶質(zhì)積累速率為：兩個(gè)相距dx垂直x軸的平面組成的微體積，J1、J2為進(jìn)入、流出兩平面間的擴(kuò)散通量，擴(kuò)散中濃度變化為：由Fick第一定律：4.2濃度場(chǎng)的數(shù)學(xué)模型及求解即第二個(gè)面的擴(kuò)散通量為第一個(gè)面注入的溶質(zhì)與在這一段距離內(nèi)溶質(zhì)濃度變化引起的擴(kuò)散通量之和。假設(shè)D不隨濃度變化，那么菲克第二定律在三維直角坐標(biāo)系下的形式：4.2濃度場(chǎng)的數(shù)學(xué)模型及求解滲碳濃度場(chǎng)的有限差求解滲碳是最常用的一種化學(xué)熱處理工藝，工件滲碳后外表具有高硬度和耐磨性，心部具有良好的強(qiáng)韌性。工件滲碳后的碳濃度分布決定了工件的性能，求解菲克第二定律偏微分方程，能預(yù)測(cè)工件滲碳后的碳濃度分布。實(shí)際工件的形狀是多種多樣的，其外表形狀可歸為平面、凸柱面、凹柱面和球面四類。許多工作說明，工件的外表形狀和曲率半徑影響滲碳后的碳濃度分布，對(duì)不同形狀的工件或外表，應(yīng)采用不同的坐標(biāo)求解。數(shù)學(xué)描述：工件內(nèi)任一時(shí)刻任何位置的碳濃度變化規(guī)律由菲克第二定律偏微分方程描述：4.2濃度場(chǎng)的數(shù)學(xué)模型及求解當(dāng)工件外表為平面，可將碳在工件內(nèi)的擴(kuò)散看成只在垂直于工件外表的一維方向上進(jìn)行，忽略擴(kuò)散系數(shù)D隨碳濃度的變化，考慮到氣體滲碳過程主要由碳向工件外表的傳遞和碳在工件內(nèi)部的擴(kuò)散兩局部組成，那么工件內(nèi)的碳濃度分布可歸結(jié)為求以下定解問題，直角坐標(biāo)系下的方程：4.2濃度場(chǎng)的數(shù)學(xué)模型及求解半徑為R，長徑比較大的圓柱形工件，可認(rèn)為碳在工件內(nèi)只沿半徑方向進(jìn)行擴(kuò)散，利用坐標(biāo)變換，歸結(jié)為求以下定解問題，極坐標(biāo)系下的方程：4.2濃度場(chǎng)的數(shù)學(xué)模型及求解半徑為R的球形工件內(nèi)的碳濃度分布可歸結(jié)為求以下定解問題，球坐標(biāo)系下的方程：4.2濃度場(chǎng)的數(shù)學(xué)模型及求解在一維擴(kuò)散條件下，邊界條件為：利用有限差分法求解時(shí)，一般分兩步進(jìn)行：(1)將連續(xù)函數(shù)C=f(x，τ)離散化,將x-τ平面劃分為如圖4-23所示的網(wǎng)格，圖中△x、△τ分別代表距離步長和時(shí)間步長。兩條平行線的交點(diǎn)稱為節(jié)點(diǎn)，并以有限個(gè)節(jié)點(diǎn)上的函數(shù)值C(xi，tn)代替連續(xù)函數(shù)C＝f(x，τ)。為簡便起見，將C(xi，τn)寫為Ci,n,即表示在τn時(shí)刻xi處的濃度值；同理可用Ci,n+1表示在τn+1(τn+△τ)時(shí)刻xi+1處的濃度值。(2)用差分代替微分，對(duì)每個(gè)節(jié)點(diǎn)用差分代替微分，此時(shí)在xi處可作以下代換：〔4.78〕4.2濃度場(chǎng)的數(shù)學(xué)模型及求解其濃度對(duì)空間的二階偏導(dǎo)數(shù)同樣也用二階中心差分替代：xyni+1(i,n)τn+1n-1ii-1有限差分方法的節(jié)點(diǎn)網(wǎng)格材料學(xué)中的濃度場(chǎng)簡介將上兩式代入菲克第二定律，得：〔4.80〕上式稱為crank-niclson格式，它對(duì)任何時(shí)間步長都是穩(wěn)定的。上式的截?cái)嗾`差為[(△x)十(△τ)]，小于其他差分格式。按Crank-Niclson格式描述滲層濃度場(chǎng)的示意圖見以下圖。224.2濃度場(chǎng)的數(shù)學(xué)模型及求解將式4.80代入式4.78，整理得到式中：4.2濃度場(chǎng)的數(shù)學(xué)模型及求解那么有：4.2濃度場(chǎng)的數(shù)學(xué)模型及求解4.2濃度場(chǎng)的數(shù)學(xué)模型及求解4.2濃度場(chǎng)的數(shù)學(xué)模型及求解(3)差分方程組求解方程的矩陣僅在主對(duì)角線及相鄰的兩條對(duì)角線上有非零元素，屬于m+1階三對(duì)角矩陣，采用追趕法用計(jì)算機(jī)可快速求解。將矩陣及擴(kuò)散系數(shù)D，傳遞系數(shù)β的計(jì)算式，以及相應(yīng)的滲碳工藝參數(shù)輸人計(jì)算機(jī)。如某一時(shí)刻n的碳濃度分布(Co，C1，…，Cm-1,Cm)n，就可以計(jì)算出△τ時(shí)間后(n+1)時(shí)刻的碳濃度分布(Co，C1，…，Cm-1,Cm)n+1。同時(shí)，將直讀光譜實(shí)測(cè)的碳濃度分布數(shù)據(jù)輸人計(jì)算機(jī)，與計(jì)算曲線進(jìn)行比較。以下圖為差分法求解擴(kuò)散方程計(jì)算的流程圖。4.2濃度場(chǎng)的數(shù)學(xué)模型及求解4.2濃度場(chǎng)的數(shù)學(xué)模型及求解輸人參數(shù)中，T—溫度(℃)；Co—鋼的原始含碳量(%)；△τ—時(shí)間步長(s)；△x—滲層深度步長(mm)；d--需到達(dá)的滲層深度；τf—滲碳時(shí)間。初始條件是C0i＝Co。為保證計(jì)算的穩(wěn)定性，在計(jì)算中合理選樣步長比例是非常重要的。(4)計(jì)算結(jié)果分析圖4-26為20#碳鋼在RJJ—35—9式滲碳爐中進(jìn)行氣體滲碳實(shí)測(cè)和用crank-niclson差分格式汁算的結(jié)果。滲碳工藝為：①單段滲碳：滲碳溫度T＝920℃，碳勢(shì)Cg＝1.2%，時(shí)間tf＝6h；②二段滲碳：滲碳溫度T＝920℃．碳勢(shì)Cgl＝1.2%，時(shí)間tf1＝2h，碳勢(shì)Cg2＝0.76%，時(shí)間tf2＝2h。4.2濃度場(chǎng)的數(shù)學(xué)模型及求解20#鋼滲碳實(shí)測(cè)點(diǎn)和crank-niclson差分格式計(jì)算曲線計(jì)算結(jié)果說明，在合理選擇步長比例的條件下，采用crank-niclson差分格式求解擴(kuò)散方程，適用于單段滲碳和二段滲碳模擬計(jì)算，計(jì)算結(jié)果與實(shí)驗(yàn)測(cè)得結(jié)果十分吻合。4.2濃度場(chǎng)的數(shù)學(xué)模型及求解在碳濃度的模擬計(jì)算中，擴(kuò)散系數(shù)D、碳的傳遞系數(shù)β和氣氛碳勢(shì)Cg是幾個(gè)重要參數(shù)。碳在奧氏體中的擴(kuò)散系數(shù)D與含碳量和溫度有關(guān)，當(dāng)鋼中含碳量小于1%時(shí)，擴(kuò)散系數(shù)D隨碳濃度的變化不大，只考慮擴(kuò)散系數(shù)D隨溫度的變化。在800～1000℃溫度范圍內(nèi)，碳在奧氏體中的擴(kuò)散系數(shù)D與溫度的近似關(guān)系是：式中R=8.314(J/K.mol)T—絕對(duì)溫度(K)碳的傳遞系數(shù)β是描述滲碳界面反響快慢的反響速度常數(shù)，隨滲碳?xì)夥盏某煞侄兓?。確定β值的方法有箔片法等。氣氛碳勢(shì)Cg表征氣氛滲碳能力的大小。碳勢(shì)指滲碳反響到達(dá)平衡時(shí)碳鋼的含碳量，合金鋼滲碳時(shí)，應(yīng)考慮鋼的化學(xué)成分對(duì)滲碳的影響，合金元素改變奧氏體中碳的活度，在相同的氣氛碳勢(shì)下將得到不同的平衡碳含量。從工程應(yīng)用的角度，可認(rèn)為合金元素改變了氣氛的有效碳勢(shì)，Cr、Mn、Mo等傾向形成Fe3C更穩(wěn)定碳化物的元素增高氣氛的有效碳勢(shì)，Si和Ni等傾向形成比Fe3C不穩(wěn)定碳化物的元素降低氣氛的有效碳勢(shì)。合金鋼滲碳后外表含碳量按Gunnarson公式計(jì)算：式中Ca指合金鋼外表到達(dá)的實(shí)際含碳量，Cc指在相同的碳勢(shì)時(shí)碳鋼外表的含碳量。4.2濃度場(chǎng)的數(shù)學(xué)模型及求解計(jì)算碳濃度的通用程序由兩個(gè)功能模塊構(gòu)成。功能模塊1是根據(jù)選擇的鋼種、工件的形狀和尺寸、滲碳的有關(guān)工藝參數(shù)，計(jì)算工件滲碳后的碳濃度分布。計(jì)算結(jié)果能以圖形或表格的形式顯示或打印?？梢灾苯舆x定平面、凸柱面、凹柱面、球體四種形狀中的一種計(jì)算滲碳后的碳濃度分布，并以圖形的形式顯示。也可以選定同時(shí)計(jì)算四種形狀的工件滲碳后的碳濃度分布，在一張圖上顯示相同滲碳條件時(shí)，四種形狀的工件滲碳后的碳濃度分布曲線。在相同滲碳工藝條件下，計(jì)算不同曲率半徑、不同形狀工件滲碳后工件內(nèi)的碳濃度分布，能了解工件曲率半徑和工件形狀對(duì)滲層深度和滲層碳濃度分布梯度的影響。該模塊預(yù)測(cè)工件滲碳后的滲層深度和滲層碳濃度分布梯度，為工藝人員進(jìn)行滲碳工藝設(shè)計(jì)提供“離線〞幫助。

4.2濃度場(chǎng)的數(shù)學(xué)模型及求解用通用程序計(jì)算了不同鋼種、不同曲率半徑、不同形狀的工件滲碳后的滲碳層深度和碳濃度分布。結(jié)果說明當(dāng)工件曲率半徑較大時(shí)，形狀對(duì)工件滲碳后的滲碳層深度和碳濃度分布影響小，碳濃度分布曲線重合。當(dāng)工件曲率半徑較小時(shí)，形狀對(duì)工件滲碳后的滲碳層深度和碳濃度分布影響較大，滲碳層深度和碳濃度分布梯度差異較大。如圖是20CrNiMo鋼的一個(gè)算例，計(jì)算條件如下：滲碳溫度930℃，強(qiáng)滲碳勢(shì)Cg=1.2，強(qiáng)滲時(shí)間6h，擴(kuò)散碳勢(shì)Cg=0.8，擴(kuò)散時(shí)間4h，β取1.58×10-4mm/s。以0.4%含碳量作為滲碳層深度。從圖看出，相同滲碳工藝條件下，球體滲碳后滲碳層深度最大，碳濃度分布梯度明顯高于其他形狀工件的碳濃度分布。碳濃度分布曲線(R=3mm)4.2濃度場(chǎng)的數(shù)學(xué)模型及求解例：對(duì)一足夠長的碳含量ωc=0.1%的低碳鋼棒材滲碳，滲碳溫度為930℃，設(shè)滲碳開始時(shí)棒材外表碳含量即達(dá)ωc=1%且始終保持這一水平，試求滲碳進(jìn)行4h后外表4×10-4m處的碳濃度C?！蔡荚讦?Fe中的擴(kuò)散系數(shù)1.61×10-12m2/s〕通過誤差函數(shù)解得C=0.157%。轉(zhuǎn)化為求解如下方程組：4.2濃度場(chǎng)的數(shù)學(xué)模型及求解1〕區(qū)域設(shè)置單擊工具，在窗口拉出一個(gè)矩形，雙擊矩形區(qū)域，在ObjectDialog對(duì)話框輸入Left為0，Bottom為0，Width為1e-4,Height為4e-3。與默認(rèn)的坐標(biāo)相比，圖形小的看不見，所以要調(diào)整坐標(biāo)顯示比例。方法：選擇Options->AxesLimits,把X，Y軸的自動(dòng)選項(xiàng)翻開。用PDETool解題步驟如下：2〕設(shè)置邊界條件單擊，使邊界變紅色，然后分別雙擊每段邊界，翻開BoundaryConditions對(duì)話框，設(shè)置邊界條件。在左邊界和右邊界，選擇Neumann，輸入g為0，q為0?！脖硎咀笥疫吔缗c外界絕緣〕下邊界選擇Dirichlet條件，輸入h為1，r為1e-2〔表示下邊界恒為0.01〕。上邊界選擇Dirichlet條件，輸入h為1，r為1e-3〔表示上邊界恒為0.001〕。4.2濃度場(chǎng)的數(shù)學(xué)模型及求解3〕設(shè)置方程類型單擊，翻開PDESpecification對(duì)話框，設(shè)置方程類型為Parabolic〔拋物型〕，c=1.61e-12，a=0,f=0,d=1。4〕網(wǎng)格劃分單擊，或者加密網(wǎng)格，單擊。5〕初值和誤差的設(shè)置