本發(fā)明提供一種三維溫度場(chǎng)的動(dòng)態(tài)數(shù)據(jù)分析方法，具體涉及一種基于時(shí)空動(dòng)態(tài)耦合的三維溫度傳感數(shù)據(jù)分析方法，屬于工業(yè)工程領(lǐng)域。

背景技術(shù)：

溫度場(chǎng)的動(dòng)態(tài)數(shù)據(jù)分析技術(shù)在工程領(lǐng)域發(fā)揮著重要作用，對(duì)提高產(chǎn)品質(zhì)量、提升系統(tǒng)性能等工程任務(wù)提供重要信息。近年來(lái)，這項(xiàng)技術(shù)引起科研與工程人員的廣泛關(guān)注，已經(jīng)在生態(tài)、氣象、衛(wèi)生保健、糧食倉(cāng)儲(chǔ)等工程領(lǐng)域得到應(yīng)用。溫度場(chǎng)的動(dòng)態(tài)數(shù)據(jù)分析技術(shù)能為復(fù)雜系統(tǒng)的監(jiān)控提供全面、準(zhǔn)確的信息，達(dá)到提高系統(tǒng)性能或服務(wù)質(zhì)量的效果。同時(shí)，該技術(shù)有助于工程結(jié)構(gòu)的優(yōu)化設(shè)計(jì)，實(shí)現(xiàn)節(jié)約資源、降低成本的目標(biāo)。

溫度場(chǎng)的動(dòng)態(tài)數(shù)據(jù)分析技術(shù)旨在實(shí)現(xiàn)對(duì)溫度場(chǎng)的準(zhǔn)確估計(jì)。傳統(tǒng)的溫度場(chǎng)估計(jì)方法是一種基于熱傳遞原理的仿真方法。文獻(xiàn)[1]～[2]中記載，該方法考慮影響溫度場(chǎng)的變化的外部因素，如環(huán)境因素、內(nèi)部傳熱機(jī)制，通過(guò)給定初始溫度和邊界條件，建立三維熱力學(xué)模型，實(shí)現(xiàn)對(duì)三維動(dòng)態(tài)溫度場(chǎng)的估計(jì)。然而，溫度場(chǎng)的變化不僅受到外部因素的影響，還受到內(nèi)部因素及其他多種不確定因素的影響。溫度場(chǎng)可以被視為一個(gè)隨著時(shí)空變化的復(fù)雜的傳熱系統(tǒng)。傳統(tǒng)的溫度場(chǎng)估計(jì)方法僅能刻畫在理想狀態(tài)下溫度場(chǎng)變化的輪廓和趨勢(shì)。對(duì)于溫度場(chǎng)中經(jīng)常出現(xiàn)的由內(nèi)部因素引起的局部溫度變化，該方法不再適用。因此，傳統(tǒng)的溫度場(chǎng)估計(jì)方法誤差較大，無(wú)法提供高精度的溫度場(chǎng)信息。

隨著時(shí)代的發(fā)展，無(wú)線傳感器技術(shù)被廣泛應(yīng)用于工程領(lǐng)域(如文獻(xiàn)[3])。該技術(shù)通過(guò)無(wú)線傳感器采集溫度數(shù)據(jù)，并利用這些數(shù)據(jù)建立溫度場(chǎng)估計(jì)模型。由于傳感器配置成本的限制，目前在溫度場(chǎng)中安裝的傳感器數(shù)目較少，通過(guò)傳感器只能獲取一少部分稀疏的觀測(cè)數(shù)據(jù)。觀測(cè)數(shù)據(jù)的稀疏性會(huì)導(dǎo)致溫度場(chǎng)估計(jì)模型精確度大大降低。同時(shí)，由傳感器采集的溫度數(shù)據(jù)存在測(cè)量誤差，同樣會(huì)降低溫度場(chǎng)估計(jì)模型的精確度。因此，僅依靠無(wú)線傳感器技術(shù)與數(shù)學(xué)模型相結(jié)合的方法無(wú)法實(shí)現(xiàn)對(duì)溫度場(chǎng)的精準(zhǔn)估計(jì)。

此外，作為一個(gè)復(fù)雜的傳熱系統(tǒng)，三維溫度場(chǎng)相鄰位置或相鄰時(shí)刻的溫度數(shù)據(jù)存在時(shí)間和空間的相關(guān)性，且溫度數(shù)據(jù)的時(shí)、空相關(guān)性是相互作用的，即數(shù)據(jù)間既存在時(shí)間相關(guān)性、空間相關(guān)性，又存在時(shí)空相關(guān)性?，F(xiàn)有的研究對(duì)溫度場(chǎng)這一特性的描述多是在時(shí)間和空間相關(guān)性是相互獨(dú)立的假設(shè)條件下進(jìn)行的，沒(méi)有考慮數(shù)據(jù)的時(shí)空相關(guān)性，難以實(shí)現(xiàn)對(duì)三維動(dòng)態(tài)溫度場(chǎng)的時(shí)空相關(guān)性進(jìn)行考量和估計(jì)。

參考文獻(xiàn)

[1]d.wangandx.zhang,“apredictionmethodforinteriortemperatureofgrainstorageviadynamicsmodel:asimulationstudy”,proceed-ingsofieeeinternationalconferenceofautomationscienceandengineering,pp.1477-1483,2015.

[2]c.jia,d.sunandc.cao,“finiteelementpredictionoftransienttemperaturedistributioninagrainstoragebin”,journalofagriculturalengineeringresearch,vol.76,no.4,pp.323–330,2000.

[3]y.ding,e.a.elsayed,s.kumara,j.lu,f.niuandj.shi,“distributedsensingforqualityandproductivityimprovements”,ieeetransactionsonautomationscienceandengineering,vol.3,no.4,pp.344–359,2006.

技術(shù)實(shí)現(xiàn)要素：

針對(duì)上述現(xiàn)有技術(shù)存在的不足，本發(fā)明提供一種基于時(shí)空動(dòng)態(tài)耦合的三維溫度傳感數(shù)據(jù)分析方法，在混合效應(yīng)模型的框架下，將反映傳熱機(jī)理的溫度場(chǎng)物理模型與基于時(shí)空相關(guān)性的統(tǒng)計(jì)模型相結(jié)合，實(shí)現(xiàn)對(duì)三維動(dòng)態(tài)溫度場(chǎng)的估計(jì)。通過(guò)本發(fā)明能夠解決目前在工程領(lǐng)域中存在的由于傳感器數(shù)據(jù)不足無(wú)法獲得準(zhǔn)確、全面溫度場(chǎng)信息的問(wèn)題，可以得到溫度場(chǎng)的全面而準(zhǔn)確的信息，為工程領(lǐng)域中監(jiān)控、決策等措施提供依據(jù)。

本發(fā)明提供的技術(shù)方案如下：

一種基于時(shí)空動(dòng)態(tài)耦合的三維溫度傳感數(shù)據(jù)分析方法，利用反映傳熱機(jī)理的溫度場(chǎng)物理模型與基于時(shí)空相關(guān)性的統(tǒng)計(jì)模型，建立混合效應(yīng)模型，對(duì)三維溫度傳感數(shù)據(jù)進(jìn)行分析，獲得三維空間位置上的動(dòng)態(tài)溫度值(即任一空間位置或任一時(shí)刻溫度值)；包括如下步驟：

1)利用溫度場(chǎng)物理模型與基于時(shí)空相關(guān)性的統(tǒng)計(jì)模型，建立混合效應(yīng)模型框架，三維溫度場(chǎng)某時(shí)空點(diǎn)的溫度值包括表示溫度場(chǎng)的全局溫度變化的均值函數(shù)項(xiàng)、表示溫度場(chǎng)的局部溫度變化的局部變化項(xiàng)和表示由隨機(jī)或不可控因素引起的溫度變化的隨機(jī)噪聲項(xiàng)；

建立混合效應(yīng)模型框架如下：

假設(shè)y(s,t)表示三維溫度場(chǎng)在時(shí)空點(diǎn)(s,t)的溫度值，s和t分別表示空間和時(shí)間的自變量，混合效應(yīng)模型的框架為：

y(s,t)＝μ(s,t)+w(s,t)+∈(s,t)(式1)

其中，μ(s,t)表示在時(shí)空點(diǎn)(s,t)處的均值函數(shù)項(xiàng)，用來(lái)刻畫溫度場(chǎng)的全局溫度變化情況；w(s,t)表示在時(shí)空點(diǎn)(s,t)處的局部變化項(xiàng)，用來(lái)刻畫溫度場(chǎng)的局部溫度變化情況；∈(s,t)表示在時(shí)空點(diǎn)(s,t)處的隨機(jī)噪聲項(xiàng)，用來(lái)刻畫由隨機(jī)或不可控因素引起的溫度變化情況，通常假設(shè)其為白噪聲。下面，將分別對(duì)均值函數(shù)項(xiàng)μ(s,t)和局部變化項(xiàng)w(s,t)建立模型。

2)均值函數(shù)項(xiàng)μ(s,t)的建模，實(shí)現(xiàn)方法為：

b1.建立均值模型：

溫度場(chǎng)的變化通常是由環(huán)境因素引起的。本發(fā)明中，充分考慮環(huán)境因素的影響，在笛卡爾坐標(biāo)系下，對(duì)均值函數(shù)項(xiàng)μ(s,t)建立三維非穩(wěn)態(tài)傅里葉傳熱模型：

式2中，μ(x,y,z,t)表示在笛卡爾坐標(biāo)系下的均值函數(shù)項(xiàng)，其中空間坐標(biāo)s＝{x,y,z}，x、y和z分別表示在x、y和z方向的坐標(biāo)；ρ表示物質(zhì)的密度，c表示物質(zhì)的比熱容，λx，λy和λz分別表示在三維笛卡爾坐標(biāo)系下物質(zhì)沿x，y和z方向的熱導(dǎo)率。

b2.均值模型的求解

采用有限差分法，利用式3對(duì)式2求解：

其中，(i,j,k,m)表示坐標(biāo)(x,y,z,t)對(duì)應(yīng)的網(wǎng)格位置，δx，δy，δz和δt分別表示空間、時(shí)間方向的網(wǎng)格間隔。為保證這種數(shù)值解法的精確度，參數(shù)值的選取應(yīng)保證給定初始溫度及邊界條件，采用有限差分法即可實(shí)現(xiàn)對(duì)均值函數(shù)的求解。

3)局部變化項(xiàng)w(s,t)的建模，實(shí)現(xiàn)方法為：

充分考慮溫度場(chǎng)的時(shí)空相關(guān)性，采用高斯隨機(jī)場(chǎng)與克里金模型相結(jié)合的方法來(lái)刻畫溫度場(chǎng)的局部溫度變化。

c1.采用高斯隨機(jī)場(chǎng)確定模型框架

假設(shè)在時(shí)刻tm(m＝1,2,…,m)，局部溫度變化具有空間平穩(wěn)性，采用高斯隨機(jī)場(chǎng)模型描述si點(diǎn)的溫度與其相鄰位置溫度的關(guān)系，即在給定si點(diǎn)相鄰位置溫度值的條件下，si點(diǎn)的溫度與服從于正態(tài)分布：

其中，sj～si表示si和sj在高斯隨機(jī)場(chǎng)中處于相鄰位置，ωij(tm)表示在tm時(shí)刻的空間權(quán)重參數(shù)，σ²(si,tm)表示條件方差。

c2.采用克里金模型確定權(quán)重參數(shù)

采用克里金模型，通過(guò)式5確定權(quán)重參數(shù)ωij(tm)：

其中，c(si,sj)表示si與sj的協(xié)方差矩陣，c(sj,sj)表示sj之間的協(xié)方差矩陣。協(xié)方差矩陣由協(xié)方差函數(shù)確定，協(xié)方差函數(shù)的表達(dá)式如式6：

其中sp和sq分別表示p和q點(diǎn)處的空間坐標(biāo)，在笛卡爾坐標(biāo)系下，sp＝{xp,yp,zp}和sq＝{xq,yq,zq}；為在tm時(shí)刻下協(xié)方差函數(shù)的參數(shù)，將η＝{{η(t1),η(t2),…,{η(tm)}記為協(xié)方差參數(shù)。

c3.協(xié)方差參數(shù)求解

在初始時(shí)刻t1，采用極大似然估計(jì)方法求解協(xié)方差參數(shù)。為了刻畫溫度場(chǎng)的時(shí)間相關(guān)性，采用貝葉斯估計(jì)方法對(duì)t2到tm時(shí)刻的協(xié)方差參數(shù)進(jìn)行更新。假設(shè)已知t1到tm-1時(shí)刻的協(xié)方差參數(shù)及溫度數(shù)據(jù)(m＝2,…,m)，通過(guò)式7得到tm時(shí)刻協(xié)方差參數(shù)的預(yù)測(cè)分布：

p(η(tm)|y(s,t1:m-1))＝∫p(η(tm)|η(tm-1))p(η(tm-1)|y(s,t1:m-1))dη(tm)(式7)

其中，p(η(tm)|y(s,t1:m-1))表示在已知t1到tm-1時(shí)刻溫度數(shù)據(jù)y(s,t1:m-1)的條件下，tm時(shí)刻的協(xié)方差參數(shù)η(tm)的概率密度；p(η(tm)|η(tm-1))表示在已知tm-1時(shí)刻的協(xié)方差參數(shù)η(tm-1)的條件下，tm時(shí)刻的協(xié)方差參數(shù)η(tm)的概率密度；p(η(tm-1)|y(s,t1:m-1))表示在已知t1到tm-1時(shí)刻溫度數(shù)據(jù)y(s,t1:m-1)的條件下，tm-1時(shí)刻的協(xié)方差參數(shù)η(tm-1)的概率密度。

再通過(guò)式8，利用tm時(shí)刻的溫度數(shù)據(jù)實(shí)現(xiàn)協(xié)方差參數(shù)的更新，即求協(xié)方差參數(shù)的后驗(yàn)概率密度分布p(η(tm)|y(s,t1:m))：

其中，

p(y(s,tm)|y(s,t1:m-1))＝∫p(y(s,tm)|η(tm))p(η(tm)|y(s,t1:m-1))dη(tm)(式9)

p(y(s,tm)|η(tm))表示在tm時(shí)刻的協(xié)方差參數(shù)η(tm)的條件下，y(s,tm)的概率密度。由于參數(shù)的分布是非正態(tài)的，因此采用粒子濾波方法實(shí)現(xiàn)貝葉斯估計(jì)的迭代計(jì)算。

4)時(shí)空溫度場(chǎng)估計(jì)，實(shí)現(xiàn)方法為：

參數(shù)估計(jì)之后，給定觀測(cè)數(shù)據(jù)，通過(guò)式10實(shí)現(xiàn)對(duì)溫度場(chǎng)中任何位置的溫度值的估計(jì)：

其中，在tm時(shí)刻點(diǎn)的溫度估計(jì)值。

與現(xiàn)有技術(shù)相比，本發(fā)明的有益效果是：

本發(fā)明提供一種基于時(shí)空動(dòng)態(tài)耦合的三維溫度傳感數(shù)據(jù)分析方法，利用反映傳熱機(jī)理的溫度場(chǎng)物理模型與基于時(shí)空相關(guān)性的統(tǒng)計(jì)模型，建立混合效應(yīng)模型，對(duì)三維溫度傳感數(shù)據(jù)進(jìn)行分析，獲得三維溫度場(chǎng)中任一空間位置的動(dòng)態(tài)溫度值(即任一空間位置或任一時(shí)刻溫度值)。通過(guò)本發(fā)明所提供的技術(shù)方案，將反映傳熱機(jī)理的溫度場(chǎng)物理模型與基于時(shí)空相關(guān)性的統(tǒng)計(jì)模型相結(jié)合，建立混合效應(yīng)模型，實(shí)現(xiàn)對(duì)三維動(dòng)態(tài)溫度場(chǎng)的估計(jì)，為溫度場(chǎng)的實(shí)時(shí)監(jiān)控提供準(zhǔn)確、全面的信息，有助于實(shí)現(xiàn)溫度場(chǎng)中傳感器的最優(yōu)配置，達(dá)到降低成本、節(jié)約能源的效果。

附圖說(shuō)明

圖1是本發(fā)明提供的基于時(shí)空動(dòng)態(tài)耦合的三維溫度傳感數(shù)據(jù)分析方法的流程框圖。

圖2是本發(fā)明實(shí)施例的協(xié)方差參數(shù)估計(jì)結(jié)果；

其中，橫坐標(biāo)為時(shí)間，單位為天；縱坐標(biāo)為協(xié)方差參數(shù)值；(a)為參數(shù)ax的估計(jì)值；(b)為參數(shù)ay的估計(jì)值；(c)為參數(shù)az的估計(jì)值；(d)為參數(shù)的估計(jì)值。

圖3是本發(fā)明實(shí)施例的糧食溫度場(chǎng)的估計(jì)結(jié)果；

其中，圖例中數(shù)字是采樣點(diǎn)的空間坐標(biāo)。

具體實(shí)施方式

下面結(jié)合附圖，通過(guò)實(shí)施例進(jìn)一步描述本發(fā)明，但不以任何方式限制本發(fā)明的范圍。

本發(fā)明提供一種基于時(shí)空動(dòng)態(tài)耦合的三維溫度傳感數(shù)據(jù)分析方法，在混合效應(yīng)模型的框架下，將反映傳熱機(jī)理的溫度場(chǎng)物理模型與基于時(shí)空相關(guān)性的統(tǒng)計(jì)模型相結(jié)合，實(shí)現(xiàn)對(duì)三維動(dòng)態(tài)溫度場(chǎng)的估計(jì)。通過(guò)本發(fā)明能夠解決目前在工程領(lǐng)域中存在的由于傳感器數(shù)據(jù)不足無(wú)法獲得準(zhǔn)確、全面溫度場(chǎng)信息的問(wèn)題，可以得到溫度場(chǎng)的全面而準(zhǔn)確的信息，為工程領(lǐng)域中監(jiān)控、決策等措施提供依據(jù)。

圖1是本發(fā)明提供的基于時(shí)空動(dòng)態(tài)耦合的三維溫度傳感數(shù)據(jù)分析方法的流程框圖。實(shí)施例設(shè)定溫度場(chǎng)范圍為(－50℃，50℃)，利用傳感器采集溫度數(shù)據(jù)，采集時(shí)間間隔為7天，采集得到的觀測(cè)數(shù)據(jù)表示為yd(s,t)。利用本發(fā)明方法對(duì)三維溫度傳感數(shù)據(jù)進(jìn)行分析的具體實(shí)施方式如下。

1)根據(jù)研究對(duì)象的溫度場(chǎng)的特征，建立混合效應(yīng)模型。

研究對(duì)象可以是工程領(lǐng)域中任何與溫度場(chǎng)有關(guān)的物質(zhì)溫度場(chǎng)，例如倉(cāng)儲(chǔ)中的糧食溫度場(chǎng)，某一地區(qū)的環(huán)境溫度場(chǎng)，人體體溫場(chǎng)。混合效應(yīng)模型的框架如下：假設(shè)y(s,t)表示三維溫度場(chǎng)在時(shí)空點(diǎn)(s,t)的溫度值，s和t分別表示空間和時(shí)間的自變量，混合效應(yīng)模型的框架為：

y(s,t)＝μ(s,t)+w(s,t)+∈(s,t)(式1)

其中，μ(s,t)表示在時(shí)空點(diǎn)(s,t)處的均值函數(shù)項(xiàng)，用來(lái)刻畫溫度場(chǎng)的全局溫度變化情況；w(s,t)表示在時(shí)空點(diǎn)(s,t)處的局部變化項(xiàng)，用來(lái)刻畫溫度場(chǎng)的局部溫度變化情況；∈(s,t)表示在時(shí)空點(diǎn)(s,t)處的隨機(jī)噪聲項(xiàng)，用來(lái)刻畫由隨機(jī)或不可控因素引起的溫度變化情況，通常假設(shè)其為白噪聲。下面，將分別對(duì)均值函數(shù)項(xiàng)μ(s,t)和局部變化項(xiàng)w(s,t)建立模型。

2)均值函數(shù)項(xiàng)μ(s,t)的建模

充分考慮研究對(duì)象的物理屬性，建立均值模型并求解。

b1.均值模型的建立

溫度場(chǎng)的變化通常是由環(huán)境因素引起的。本發(fā)明中，充分考慮環(huán)境因素的影響，在笛卡爾坐標(biāo)系下，對(duì)均值函數(shù)項(xiàng)μ(s,t)建立三維非穩(wěn)態(tài)傅里葉傳熱模型：

式2中，μ(x,y,z,t)表示在笛卡爾坐標(biāo)系下的均值函數(shù)項(xiàng)，其中空間坐標(biāo)s＝{x,y,z}，x、y和z分別表示在x、y和z方向的坐標(biāo)；ρ表示物質(zhì)的密度，c表示物質(zhì)的比熱容，λx，λy和λz分別表示在三維笛卡爾坐標(biāo)系下物質(zhì)沿x，y和z方向的熱導(dǎo)率。

b2.均值模型的求解

采用有限差分法對(duì)式2求解，實(shí)現(xiàn)方法為：

其中(i,j,k,m)表示坐標(biāo)(x,y,z,t)對(duì)應(yīng)的網(wǎng)格位置，δx，δy，δz和δt分別表示空間、時(shí)間方向的網(wǎng)格間隔。為保證這種數(shù)值解法的精確度，參數(shù)的選取應(yīng)保證給定初始溫度值及邊界溫度值，采用有限差分法即可實(shí)現(xiàn)對(duì)均值函數(shù)的求解。

均值模型是基于研究對(duì)象的物理意義建立的，不涉及溫度場(chǎng)的觀測(cè)數(shù)據(jù)，用于描述理想狀態(tài)下溫度場(chǎng)的全局變化。

3)局部變化項(xiàng)w(s,t)的建模

充分考慮溫度場(chǎng)的時(shí)空相關(guān)性，結(jié)合高斯隨機(jī)場(chǎng)與克里金模型建立局部溫度變化模型。

c1.高斯隨機(jī)場(chǎng)確定模型框架

假設(shè)在時(shí)刻tm(m＝1,2,…,m)，局部溫度變化具有空間平穩(wěn)性，采用高斯隨機(jī)場(chǎng)模型描述si點(diǎn)的溫度與其相鄰位置溫度的關(guān)系，即在給定si點(diǎn)相鄰位置溫度值的條件下，si點(diǎn)的溫度與服從于正態(tài)分布：

其中，sj～si表示si和sj在高斯隨機(jī)場(chǎng)中處于相鄰位置，ωij(tm)表示在tm時(shí)刻的空間權(quán)重參數(shù)，σ²(si,tm)表示條件方差。

c2.克里金模型確定權(quán)重參數(shù)

采用克里金模型確定權(quán)重參數(shù)ωij(tm)，實(shí)現(xiàn)方法為：

其中，ωij(tm)表示si點(diǎn)相對(duì)其所有相鄰位置sj的權(quán)重參數(shù)，c(si,sj)表示si與sj的協(xié)方差矩陣，c(sj,sj)表示sj之間的協(xié)方差矩陣。協(xié)方差矩陣是由協(xié)方差函數(shù)計(jì)算得到，協(xié)方差函數(shù)的表達(dá)式如下：

其中sp和sq分別表示p和q點(diǎn)處的空間坐標(biāo)，在笛卡爾坐標(biāo)系下，sp＝{xp,yp,zp}和sq＝{xq,yq,zq}；為在tm時(shí)刻下協(xié)方差函數(shù)的參數(shù)，將η＝{{η(t1),η(t2),…,{η(tm)}記為協(xié)方差參數(shù)。

c3.協(xié)方差參數(shù)求解

在初始時(shí)刻t1，使用t1時(shí)刻的觀測(cè)數(shù)據(jù)yd(s,t1)，采用極大似然估計(jì)方法求解協(xié)方差參數(shù)。為了刻畫溫度場(chǎng)的時(shí)間相關(guān)性，采用貝葉斯估計(jì)方法對(duì)t2到tm時(shí)刻的協(xié)方差參數(shù)進(jìn)行更新。假設(shè)已知t1到tm-1時(shí)刻的協(xié)方差參數(shù)及觀測(cè)數(shù)據(jù)yd(s,t1:m-1)(m＝2,…,m)，通過(guò)式7得到tm時(shí)刻協(xié)方差參數(shù)的預(yù)測(cè)分布：

p(η(tm)|yd(s,t1:m-1))＝∫p(η(tm)|η(tm-1))p(η(tm-1)|yd(s,t1:m-1))dη(tm)(式7)

其中，p(η(tm)|yd(s,t1:m-1))表示在已知t1到tm-1時(shí)刻觀測(cè)數(shù)據(jù)yd(s,t1:m-1)的條件下，tm時(shí)刻的協(xié)方差參數(shù)η(tm)的概率密度；p(η(tm)|η(tm-1))表示在已知tm-1時(shí)刻的協(xié)方差參數(shù)η(tm-1)的條件下，tm時(shí)刻的協(xié)方差參數(shù)η(tm)的概率密度；p(η(tm-1)|yd(s,t1:m-1))表示在已知t1到tm-1時(shí)刻觀測(cè)數(shù)據(jù)yd(s,t1:m-1)的條件下，tm-1時(shí)刻的協(xié)方差參數(shù)η(tm-1)的概率密度。

再通過(guò)式8，利用tm時(shí)刻的觀測(cè)數(shù)據(jù)yd(s,tm)實(shí)現(xiàn)協(xié)方差參數(shù)的更新，即求協(xié)方差參數(shù)的后驗(yàn)概率密度分布p(η(tm)|yd(s,t1:m))：

其中：

p(yd(s,tm)|yd(s,t1:m-1))＝∫p(yd(s,tm)|η(tm))p(η(tm)|yd(s,t1:m-1))dη(tm)(式9)

p(yd(s,tm)|η(tm))表示在tm時(shí)刻的協(xié)方差參數(shù)η(tm)的條件下，yd(s,tm)的概率密度。由于參數(shù)的分布是非正態(tài)的，因此采用粒子濾波方法實(shí)現(xiàn)貝葉斯估計(jì)的迭代計(jì)算。

4)估計(jì)得到時(shí)空溫度場(chǎng)中任何位置的溫度

協(xié)方差參數(shù)估計(jì)之后，給定觀測(cè)數(shù)據(jù)yd(s,t)，便可實(shí)現(xiàn)對(duì)溫度場(chǎng)中任何位置的溫度的估計(jì)。

其中，在tm時(shí)刻點(diǎn)的溫度估計(jì)值。

下面通過(guò)實(shí)例對(duì)本發(fā)明做進(jìn)一步說(shuō)明。

實(shí)施例：

以中國(guó)華中地區(qū)某一國(guó)家儲(chǔ)備糧庫(kù)為例，建立倉(cāng)儲(chǔ)過(guò)程中三維糧食溫度場(chǎng)的估計(jì)模型。糧溫?cái)?shù)據(jù)是通過(guò)溫度傳感器采集得到，記為觀測(cè)數(shù)據(jù)yd(s,t),s∈s,t∈t，其中s表示觀測(cè)數(shù)據(jù)的空間自變量s的取值范圍，t表示觀測(cè)數(shù)據(jù)的時(shí)間自變量t的取值范圍。在本例中，s的定義如下：糧食溫度場(chǎng)的規(guī)格為26m*46m*6m，傳感器的布局為在長(zhǎng)、寬方向相鄰傳感器的間隔為5m，在高度方向相鄰傳感器的間隔為1.8m。t的定義如下：糧溫?cái)?shù)據(jù)的采樣時(shí)間為2012年1月31日至2013年3月4日，每隔7天采集一組傳感器數(shù)據(jù)，共73組數(shù)據(jù)。

下面建立三維動(dòng)態(tài)糧食溫度場(chǎng)估計(jì)模型：

(1)建立混合效應(yīng)模型的框架

假設(shè)y(s,t)表示三維溫度場(chǎng)在時(shí)空點(diǎn)(s,t)的溫度值，s和t分別表示空間和時(shí)間的自變量，混合效應(yīng)模型的框架為

y(s,t)＝μ(s,t)+w(s,t)+∈(s,t)(式1)

其中，μ(s,t)表示在時(shí)空點(diǎn)(s,t)處的均值函數(shù)項(xiàng)，用來(lái)刻畫溫度場(chǎng)的全局溫度變化情況；w(s,t)表示在時(shí)空點(diǎn)(s,t)處的局部變化項(xiàng)，用來(lái)刻畫溫度場(chǎng)的局部溫度變化情況；∈(s,t)表示在時(shí)空點(diǎn)(s,t)處的隨機(jī)噪聲項(xiàng)，用來(lái)刻畫由隨機(jī)或不可控因素引起的溫度變化情況，通常假設(shè)其為白噪聲。下面，將分別對(duì)均值函數(shù)項(xiàng)μ(s,t)和局部變化項(xiàng)w(s,t)建立模型。

(2)均值函數(shù)μ(s,t)模型的建立與求解

糧食溫度場(chǎng)的變化與外界環(huán)境溫度及糧食自身傳熱有關(guān)。在此，外界環(huán)境溫度對(duì)糧食溫度的影響，在笛卡爾坐標(biāo)系下，對(duì)均值函數(shù)項(xiàng)μ(s,t)建立三維非穩(wěn)態(tài)傅里葉傳熱模型：

式2中，μ(x,y,z,t)表示在笛卡爾坐標(biāo)系下的均值函數(shù)項(xiàng)，其中空間坐標(biāo)s＝{x,y,z}，x、y和z分別表示在x、y和z方向的坐標(biāo)；ρ表示物質(zhì)的密度，c表示物質(zhì)的比熱容，λx，λy和λz分別表示在三維笛卡爾坐標(biāo)系下物質(zhì)沿x，y和z方向的熱導(dǎo)率。本例中的糧食品種為小麥，其密度ρ為750kg/m²，比熱容c為0.15w/(m·k)，熱導(dǎo)率λx，λy，λz為2000j/(kg·k)。

采用有限差分法對(duì)上述傅立葉傳熱模型求解，實(shí)現(xiàn)方法為：

其中，(i,j,k,m)表示時(shí)空坐標(biāo)(x,y,z,t)對(duì)應(yīng)的網(wǎng)格位置，δx，δy，δz和δt分別表示時(shí)空方向的網(wǎng)格間隔，即δx＝0.5m，δy＝0.5m，δz＝0.3m，δt＝24h。給定初始糧食溫度及邊界條件的溫度值(糧倉(cāng)邊界溫度的變化是由外界環(huán)境溫度引起的)，采用有限差分法即可實(shí)現(xiàn)對(duì)均值函數(shù)的求解。

(3)局部變化項(xiàng)w(s,t)的建立與協(xié)方差參數(shù)求解

高斯隨機(jī)場(chǎng)確定模型框架為：假設(shè)在時(shí)刻tm(m＝1,2,…,m)，局部溫度變化具有空間平穩(wěn)性，采用高斯隨機(jī)場(chǎng)模型描述si點(diǎn)的溫度與其相鄰位置溫度的關(guān)系，即在給定si點(diǎn)相鄰位置溫度值的條件下，si點(diǎn)的溫度與服從于正態(tài)分布：

其中，sj～si表示si和sj在糧倉(cāng)中的距離最近的傳感器的位置，ωij(tm)表示在tm時(shí)刻的空間權(quán)重參數(shù)，σ²(si,tm)表示條件方差。

采用克里金模型確定權(quán)重參數(shù)ωij(tm)，實(shí)現(xiàn)方法為：

其中，ωij(tm)表示si點(diǎn)相對(duì)其在糧倉(cāng)中的距離最近的傳感器的位置sj的權(quán)重參數(shù)，c(si,sj)表示si與sj的協(xié)方差矩陣，c(sj,sj)表示sj之間的協(xié)方差矩陣。協(xié)方差矩陣是由協(xié)方差函數(shù)計(jì)算得到，協(xié)方差函數(shù)的表達(dá)式如下：

其中sp和sq分別表示p和q點(diǎn)處的空間坐標(biāo)，在笛卡爾坐標(biāo)系下，sp＝{xp,yp,zp}和sq＝{xq,yq,zq}；為在tm時(shí)刻下協(xié)方差函數(shù)的參數(shù)，將η＝{{η(t1),η(t2),…,{η(tm)}記為協(xié)方差參數(shù)。

在初始時(shí)刻t1，使用t1時(shí)刻的觀測(cè)數(shù)據(jù)yd(s,t1)，采用極大似然估計(jì)方法求解協(xié)方差參數(shù)。為了刻畫溫度場(chǎng)的時(shí)間相關(guān)性，采用貝葉斯估計(jì)方法對(duì)t2到tm時(shí)刻的協(xié)方差參數(shù)進(jìn)行更新。假設(shè)已知t1到tm-1時(shí)刻的協(xié)方差參數(shù)及觀測(cè)數(shù)據(jù)yd(s,t1:m-1)(m＝2,…,m)，通過(guò)式7得到tm時(shí)刻協(xié)方差參數(shù)的預(yù)測(cè)分布：

p(η(tm)|yd(s,t1:m-1))＝∫p(η(tm)|η(tm-1))p(η(tm-1)|yd(s,t1:m-1))dη(tm)(式7)

其中，p(η(tm)|yd(s,t1:m-1))表示在已知t1到tm-1時(shí)刻觀測(cè)數(shù)據(jù)yd(s,t1:m-1)的條件下，tm時(shí)刻的協(xié)方差參數(shù)η(tm)的概率密度；p(η(tm)|η(tm-1))表示在已知tm-1時(shí)刻的協(xié)方差參數(shù)η(tm-1)的條件下，tm時(shí)刻的協(xié)方差參數(shù)η(tm)的概率密度；p(η(tm-1)|yd(s,t1:m-1))表示在已知t1到tm-1時(shí)刻觀測(cè)數(shù)據(jù)yd(s,t1:m-1)的條件下，tm-1時(shí)刻的協(xié)方差參數(shù)η(tm-1)的概率密度。

再通過(guò)式8，利用tm時(shí)刻的觀測(cè)數(shù)據(jù)yd(s,tm)實(shí)現(xiàn)協(xié)方差參數(shù)的更新，即求協(xié)方差參數(shù)的后驗(yàn)概率密度分布p(η(tm)|yd(s,t1:m))：

其中：

p(yd(s,tm)|yd(s,t1:m-1))＝∫p(yd(s,tm)|η(tm))p(η(tm)|yd(s,t1:m-1))dη(tm)(式9)

p(yd(s,tm)|η(tm))表示在tm時(shí)刻的協(xié)方差參數(shù)η(tm)的條件下，yd(s,tm)的概率密度。由于參數(shù)的分布是非正態(tài)的，因此采用粒子濾波方法實(shí)現(xiàn)貝葉斯估計(jì)的迭代計(jì)算。協(xié)方差參數(shù)估計(jì)的結(jié)果如圖2所示。

(4)時(shí)空溫度場(chǎng)估計(jì)

協(xié)方差參數(shù)估計(jì)之后，給定觀測(cè)數(shù)據(jù)yd(s,t)，便可實(shí)現(xiàn)對(duì)溫度場(chǎng)中任何位置的溫度的估計(jì)。

其中，在tm時(shí)刻點(diǎn)的溫度估計(jì)值。圖3為糧食溫度變化曲線，列舉了一部分糧食溫度場(chǎng)的估計(jì)結(jié)果，其中，圖例中數(shù)字是采樣點(diǎn)的空間坐標(biāo)。圖3表明，任一空間位置的溫度是隨著時(shí)間發(fā)生動(dòng)態(tài)變化的，且相鄰空間位置的溫度變化規(guī)律相近。

需要注意的是，公布實(shí)施例的目的在于幫助進(jìn)一步理解本發(fā)明，但是本領(lǐng)域的技術(shù)人員可以理解：在不脫離本發(fā)明及所附權(quán)利要求的精神和范圍內(nèi)，各種替換和修改都是可能的。因此，本發(fā)明不應(yīng)局限于實(shí)施例所公開的內(nèi)容，本發(fā)明要求保護(hù)的范圍以權(quán)利要求書界定的范圍為準(zhǔn)。