數(shù)學(xué)建模全國(guó)一等獎(jiǎng)?wù)撐膬?chǔ)油罐的變位識(shí)別與罐容表標(biāo)定

1、2010高教社杯全國(guó)大學(xué)生數(shù)學(xué)建模競(jìng)賽承 諾 書(shū)我們仔細(xì)閱讀了中國(guó)大學(xué)生數(shù)學(xué)建模競(jìng)賽的競(jìng)賽規(guī)則。我們完全明白，在競(jìng)賽開(kāi)始后參賽隊(duì)員不能以任何方式（包括電話(huà)、電子郵件、網(wǎng)上咨詢(xún)等）與隊(duì)外的任何人（包括指導(dǎo)教師）研究、討論與賽題有關(guān)的問(wèn)題。我們知道，抄襲別人的成果是違反競(jìng)賽規(guī)則的, 如果引用別人的成果或其他公開(kāi)的資料（包括網(wǎng)上查到的資料），必須按照規(guī)定的參考文獻(xiàn)的表述方式在正文引用處和參考文獻(xiàn)中明確列出。我們鄭重承諾，嚴(yán)格遵守競(jìng)賽規(guī)則，以保證競(jìng)賽的公正、公平性。如有違反競(jìng)賽規(guī)則的行為，我們將受到嚴(yán)肅處理。我們參賽選擇的題號(hào)是（從A/B/C/D中選擇一項(xiàng)填寫(xiě)）： A 我們的參賽報(bào)名號(hào)為（如果賽區(qū)設(shè)置

2、報(bào)名號(hào)的話(huà)）： 所屬學(xué)校（請(qǐng)?zhí)顚?xiě)完整的全名）： 重慶大學(xué) 參賽隊(duì)員 (打印并簽名) ：1. 王建 2. 丁超 3. 王昌贏 指導(dǎo)教師或指導(dǎo)教師組負(fù)責(zé)人 (打印并簽名)： 龔劬 日期： 年 月 日賽區(qū)評(píng)閱編號(hào)（由賽區(qū)組委會(huì)評(píng)閱前進(jìn)行編號(hào)）：2010高教社杯全國(guó)大學(xué)生數(shù)學(xué)建模競(jìng)賽編 號(hào) 專(zhuān) 用 頁(yè)賽區(qū)評(píng)閱編號(hào)（由賽區(qū)組委會(huì)評(píng)閱前進(jìn)行編號(hào)）：賽區(qū)評(píng)閱記錄（可供賽區(qū)評(píng)閱時(shí)使用）：評(píng)閱人評(píng)分備注全國(guó)統(tǒng)一編號(hào)（由賽區(qū)組委會(huì)送交全國(guó)前編號(hào)）：全國(guó)評(píng)閱編號(hào)（由全國(guó)組委會(huì)評(píng)閱前進(jìn)行編號(hào)）：儲(chǔ)油罐的變位識(shí)別與罐容表標(biāo)定摘要本文對(duì)儲(chǔ)油罐的變位識(shí)別和罐容表的標(biāo)定問(wèn)題進(jìn)行了深入探討，建立了儲(chǔ)油量和油位高度以及變位參數(shù)之

3、間關(guān)系的數(shù)學(xué)模型，主要應(yīng)用了mtalab進(jìn)行求解。針對(duì)問(wèn)題一，我們利用積分的方法推導(dǎo)出小橢圓儲(chǔ)油罐在無(wú)變位和發(fā)生縱向傾斜變位時(shí)的一般公式。討論了在儲(chǔ)油罐發(fā)生縱向傾斜變位后對(duì)罐容表的影響，定義了平均影響率(變位前后儲(chǔ)油量之差絕對(duì)值的平均值占總罐體容積的比例)作為評(píng)價(jià)罐體變位對(duì)罐容表的影響程度的大小的指標(biāo)，求出。并分別給出了小橢圓儲(chǔ)油罐在無(wú)變位和在縱向傾斜變位角取的罐容表。表1 小橢圓儲(chǔ)油罐罐容表(縱向變位)油位高度0.000.010.020.030.040.05儲(chǔ)油體積00.00170.00350.00630.01000.01480.02071.151.161.171.17131.181.191

4、.203.91033.93383.95603.95883.97663.99554.10174.1101針對(duì)問(wèn)題二，將儲(chǔ)油罐分為5個(gè)區(qū)域分別進(jìn)行討論，考慮到在球冠處的體積表達(dá)式過(guò)于復(fù)雜，我們省略了球冠處的一小部分體積，進(jìn)行了近似求解，得出了罐內(nèi)儲(chǔ)油量與油位高度以及變位參數(shù)之間的一般關(guān)系的數(shù)學(xué)模型。在利用儲(chǔ)油罐的實(shí)際測(cè)量值估計(jì)變位參數(shù)時(shí)，我們建立了最小二乘擬合模型，得到了最佳的變位參數(shù)為：縱向傾斜變位，橫向偏轉(zhuǎn)變位。并據(jù)此對(duì)儲(chǔ)油罐的罐容表進(jìn)行了標(biāo)定(見(jiàn)表3)。在模型驗(yàn)證中，我們又采用蒙特卡洛模擬的方法對(duì)在問(wèn)題二的模型中忽略的部分球冠體積進(jìn)行了模擬計(jì)算。得到用問(wèn)題二模型中求出的總儲(chǔ)油量與模擬得出的總

5、儲(chǔ)油量一致度達(dá)到了99%，誤差非常小，驗(yàn)證了我們所建立的模型的合理性和準(zhǔn)確性。關(guān)鍵詞 平均影響率 最小二乘參數(shù)估計(jì)法 蒙特卡洛模擬一 問(wèn)題重述通常加油站都有若干個(gè)儲(chǔ)存燃油的地下儲(chǔ)油罐，并且一般都有與之配套的“油位計(jì)量管理系統(tǒng)”，采用流量計(jì)和油位計(jì)來(lái)測(cè)量進(jìn)/出油量與罐內(nèi)油位高度等數(shù)據(jù)，通過(guò)預(yù)先標(biāo)定的罐容表（即罐內(nèi)油位高度與儲(chǔ)油量的對(duì)應(yīng)關(guān)系）進(jìn)行實(shí)時(shí)計(jì)算，以得到罐內(nèi)油位高度和儲(chǔ)油量的變化情況。許多儲(chǔ)油罐在使用一段時(shí)間后，由于地基變形等原因，使罐體的位置會(huì)發(fā)生縱向傾斜和橫向偏轉(zhuǎn)等變化（以下稱(chēng)為變位），從而導(dǎo)致罐容表發(fā)生改變。按照有關(guān)規(guī)定，需要定期對(duì)罐容表進(jìn)行重新標(biāo)定。問(wèn)題一為了掌握罐體變位后對(duì)罐容表

6、的影響，利用小橢圓型儲(chǔ)油罐（兩端平頭的橢圓柱體），分別對(duì)罐體無(wú)變位和傾斜角為a=4.1的縱向變位兩種情況做了實(shí)驗(yàn)。請(qǐng)建立數(shù)學(xué)模型研究罐體變位后對(duì)罐容表的影響，并給出罐體變位后油位高度間隔為1cm的罐容表標(biāo)定值。問(wèn)題二對(duì)于實(shí)際的儲(chǔ)油罐，試建立罐體變位后罐內(nèi)儲(chǔ)油量與油位高度及變位參數(shù)（縱向傾斜角度a和橫向偏轉(zhuǎn)角度b ）之間的一般關(guān)系。利用罐體變位后在進(jìn)/出油過(guò)程中的實(shí)際檢測(cè)數(shù)據(jù)，根據(jù)你們所建立的數(shù)學(xué)模型確定變位參數(shù)，并給出罐體變位后油位高度間隔為10cm的罐容表標(biāo)定值。然后進(jìn)一步用實(shí)際檢測(cè)數(shù)據(jù)來(lái)分析檢驗(yàn)?zāi)Ｐ偷恼_性與方法的可靠性。二 模型假設(shè)1、假設(shè)題中所給數(shù)據(jù)均為儲(chǔ)油罐內(nèi)壁測(cè)量值；2、不考慮由于

7、溫度、壓強(qiáng)變化等原因而引起儲(chǔ)油罐的體積變化；3、油位探針被固定在儲(chǔ)油罐上，其上油浮子能夠準(zhǔn)確測(cè)量油位高度；三 符號(hào)說(shuō)明 橢圓的半長(zhǎng)軸長(zhǎng) 橢圓的半短軸長(zhǎng) 儲(chǔ)油罐的總長(zhǎng) 油位探針到油罐底部左側(cè)的距離 儲(chǔ)油罐的總體積 儲(chǔ)油罐的縱向傾斜角度 儲(chǔ)油罐的橫向偏轉(zhuǎn)角度 油位高度 球冠體的半徑圓柱體的底面半徑 注：未說(shuō)明符號(hào)在文中用到時(shí)注明四 問(wèn)題一的解答小橢圓儲(chǔ)油罐罐體變位前后都可以應(yīng)用積分的方法求出罐體的儲(chǔ)油量和油位高度之間的關(guān)系。對(duì)于縱向傾斜的小橢圓儲(chǔ)油罐，考慮分段求出其儲(chǔ)油量和油位高度之間的關(guān)系，從而得到重新標(biāo)定后的罐容表。4.1 小橢圓儲(chǔ)油罐無(wú)變位時(shí)的模型由于此時(shí)的橢圓無(wú)變位，考慮先對(duì)二維橢圓進(jìn)行積

8、分。為方便表示油位高度，建立如圖所示的坐標(biāo)系，橢圓的半長(zhǎng)軸長(zhǎng)為，半短軸長(zhǎng)為，則橢圓方程為圖1 對(duì)橢圓的積分示意圖在方向上取橢圓面中的一微元積分得到油的側(cè)面積儲(chǔ)油罐內(nèi)油的體積為查積分表得到 (1)利用matlab計(jì)算得到經(jīng)驗(yàn)證兩種方法得到的體積公式完全等價(jià)，(1)式即為小橢圓儲(chǔ)油罐無(wú)變位時(shí)的儲(chǔ)油量和油位高度關(guān)系的模型。根據(jù)此模型，我們可以求出小橢圓儲(chǔ)油罐無(wú)變位時(shí)罐容表標(biāo)定值(油位高度間隔取，結(jié)果見(jiàn)附錄一)。4.2 小橢圓儲(chǔ)油罐縱向傾斜變位時(shí)的模型儲(chǔ)油罐縱向傾斜之后，油位計(jì)在油位過(guò)高或者過(guò)低時(shí)將不起作用(如圖2所示的和區(qū)域)，考慮到傾斜角變化一般不會(huì)很大，所以我們可以將儲(chǔ)油罐按液面高低分成五個(gè)部分

9、，來(lái)求其儲(chǔ)油量和油位高度之間的關(guān)系。我們討論的是小橢圓儲(chǔ)油罐縱向傾斜變位為逆時(shí)針旋轉(zhuǎn)，如圖2。對(duì)于儲(chǔ)油罐順時(shí)針旋轉(zhuǎn)變位(即為負(fù)值)時(shí)的情況與此非常類(lèi)似，在此不再詳細(xì)討論。圖2 儲(chǔ)油罐分區(qū)示意圖 對(duì)區(qū)域的討論 在區(qū)域，其油位低于油位探針的油浮子，所以油位計(jì)量系統(tǒng)中顯示油位高度為零。當(dāng)油位計(jì)剛開(kāi)始有示數(shù)時(shí)，計(jì)算其儲(chǔ)油體積。將區(qū)域放大得到圖3圖3 區(qū)域的放大圖圖中，從原點(diǎn)紙面向里為軸，利用三重積分可以得到其中為油位探針到儲(chǔ)油罐左側(cè)的距離 積分得到 (2) 對(duì)區(qū)域的討論由區(qū)域很容易得到區(qū)域的儲(chǔ)油量和油位高度的變化關(guān)系，直接給出結(jié)論：所以 (3) 對(duì)區(qū)域的討論圖4 區(qū)域示意圖在小橢圓儲(chǔ)油罐無(wú)變位模型中我

10、們已經(jīng)求出了的計(jì)算公式，同區(qū)域中的積分原理可以計(jì)算出，我們就可以得到此時(shí)的油量體積為 (4)其中 對(duì)區(qū)域的討論由區(qū)域4和區(qū)域2的相似性，將(3)式中的換為，將換為，并用總體積減去即為區(qū)域4的儲(chǔ)油體積和油位高度的變化關(guān)系。其中為小橢圓儲(chǔ)油罐的總體積化簡(jiǎn)并積分可得 (5) 對(duì)區(qū)域的討論在此區(qū)域中油浮子到達(dá)油位探針頂點(diǎn)，無(wú)法進(jìn)一步測(cè)量油位高度。無(wú)法測(cè)量的總體積為： (6) 綜合各區(qū)域的罐容表標(biāo)定的數(shù)學(xué)模型綜上所述，我們得到了儲(chǔ)油量和油位高度、縱向傾斜角之間的分段函數(shù)關(guān)系式：表2 分段函數(shù)關(guān)系區(qū)域油位高度儲(chǔ)油量12345根據(jù)儲(chǔ)油量和油位高度的分段函數(shù)關(guān)系我們可以得到罐體縱向傾斜變位()后油位高度間隔為

11、的罐容表標(biāo)定值(見(jiàn)附錄一)。當(dāng)時(shí)，各區(qū)域油位高度及體積變化范圍為：表3 各區(qū)域油位高度和儲(chǔ)油量變化范圍區(qū)域油位高度儲(chǔ)油量123454.3 罐體變位后對(duì)罐容表的影響為了能更加準(zhǔn)確地刻畫(huà)罐體的縱向傾斜變位對(duì)罐容表的影響，我們分別對(duì)罐體變位前后的理論值和測(cè)量值進(jìn)行多方面的比較。 罐體變位前理論值與測(cè)量值比較根據(jù)附件一中所給數(shù)據(jù)，我們計(jì)算出在附件中所給的油位高度下理論值和測(cè)量值，并畫(huà)出其曲線(xiàn)。圖5 罐體變位前的曲線(xiàn)對(duì)比通過(guò)對(duì)比我們發(fā)現(xiàn)，對(duì)于任意，儲(chǔ)油量的理論值和實(shí)際值始終成如(7)式的比例關(guān)系。 (7) 罐體變位后理論值與測(cè)量值比較同樣根據(jù)附件一中所給數(shù)據(jù)，我們計(jì)算出在附件中所給的油位高度下理論值和測(cè)

12、量值，并畫(huà)出其曲線(xiàn)，如圖6。圖6 罐體變位后的曲線(xiàn)對(duì)比從圖6中可以看出測(cè)量值仍然始終小于理論值，進(jìn)一步求得理論值與測(cè)量值之差的變化范圍為0.0454,0.0910，測(cè)量值的相對(duì)誤差范圍為1.56%，5.18%。 罐體縱向傾斜變位前后理論值比較圖7 罐體變位前后的曲線(xiàn)對(duì)比圖8 同一高度下儲(chǔ)油量的理論值與測(cè)量值之差變化關(guān)系圖9 儲(chǔ)油量的測(cè)量值的相對(duì)誤差隨油位高度的變化關(guān)系由以上各圖可以清晰地看出縱向傾斜變位后，使得在同一個(gè)油位高度下，變位后比變位前的儲(chǔ)油量減小。但是這樣仍不夠直觀，我們需要找到一個(gè)指標(biāo)來(lái)定量刻畫(huà)罐體變位后對(duì)罐容表的影響。從圖9中可以看出，當(dāng)油位高度較小時(shí)(0.1附近)，變位后相對(duì)于

13、變位前的相對(duì)誤差幾乎達(dá)到了60%以上，但是此時(shí)的儲(chǔ)油量的差別并不大，鑒于此，我們定義平均影響率：來(lái)刻畫(huà)罐體變位后對(duì)罐容表的影響?？梢郧蟪鲈诳v向傾斜變位時(shí)，4.87%。五 問(wèn)題二的解答如圖8實(shí)際的儲(chǔ)油罐示意圖，其主體為圓柱體，兩端為球冠體，在儲(chǔ)油罐無(wú)變位時(shí)我們計(jì)算其各部分體積。圖10 實(shí)際儲(chǔ)油罐示意圖圓柱體積計(jì)算公式帶入數(shù)據(jù)得到：主體圓柱體積為56.5487一端球缺體積計(jì)算公式帶入數(shù)據(jù)得到：兩端球缺總體積為，則儲(chǔ)油罐的總體積為64.66455.1考慮不發(fā)生變位時(shí)儲(chǔ)油量和油位高度的關(guān)系圓柱內(nèi)的油體積隨油位高度變化關(guān)系： (8)一端球缺內(nèi)油體積隨油位高度變化關(guān)系： (9)其中為圓柱體底面半徑為球缺對(duì)

14、應(yīng)的半徑為球內(nèi)小圓半徑總的儲(chǔ)油罐內(nèi)的油量對(duì)油位高度的變化關(guān)系為：5.2只考慮橫向偏轉(zhuǎn)變位時(shí)的儲(chǔ)油量和油位高度的關(guān)系只考慮儲(chǔ)油罐橫向偏轉(zhuǎn)為時(shí)：對(duì)實(shí)際的油位高度沒(méi)有影響，但此時(shí)的油位探針已經(jīng)隨儲(chǔ)油罐發(fā)生偏轉(zhuǎn)(如下圖)圖11 只考慮橫向偏轉(zhuǎn)示意圖由油浮子測(cè)量得到的油位高度與實(shí)際油位高度的關(guān)系為： (10)所以對(duì)罐容表影響轉(zhuǎn)換公式：5.3只考慮縱向傾斜變位時(shí)的儲(chǔ)油量和油位高度的關(guān)系只考慮縱向傾斜變位時(shí)，我們利用問(wèn)題一中的思想，將儲(chǔ)油罐分成5個(gè)區(qū)域，分別計(jì)算儲(chǔ)油量和油位高度的關(guān)系(如圖12所示)。圖12 實(shí)際儲(chǔ)油罐分區(qū)示意圖對(duì)區(qū)域1的討論油位探針測(cè)得的油位高度始終為零，分成三部分來(lái)計(jì)算這部分的體積，如區(qū)

15、域1的放大圖(13)圖13 實(shí)際儲(chǔ)油罐分區(qū)示意圖區(qū)域1的總體積應(yīng)為：對(duì)圓柱體部分進(jìn)行三重積分得到：由球缺部分的體積隨油位高度的計(jì)算公式(9)得：由于部分的體積精確計(jì)算非常復(fù)雜，而且儲(chǔ)油罐的縱向傾斜角度一般不會(huì)高于所我們考慮將這部分體積省略，進(jìn)行近似計(jì)算。由于傾斜角較小，所以區(qū)域3占據(jù)了儲(chǔ)油罐的絕大部分，而在區(qū)域3中這種近似計(jì)算的誤差將由于左右兩個(gè)球冠的省略體積一正一負(fù)而有所減小，所以由此得到區(qū)域1的總體積公式： (11) 對(duì)區(qū)域2的討論在區(qū)域2中，油位高度的變化范圍。各部分儲(chǔ)油體積和油位高度的變化關(guān)系:, 總儲(chǔ)油量和油位高度變化關(guān)系： (12) 對(duì)區(qū)域3的討論在區(qū)域3中，油位高度的變化范圍。圓

16、柱體部分的儲(chǔ)油量：其中球缺部分的近似儲(chǔ)油量：，總儲(chǔ)油量和油位高度變化關(guān)系： (13) 對(duì)區(qū)域4的討論在區(qū)域4中，油位高度的范圍。圓柱體部分的儲(chǔ)油量：球缺部分的總體積：，總儲(chǔ)油量和油位高度變化關(guān)系： (14) 對(duì)區(qū)域5的討論在區(qū)域4中，油位高度始終等于。各部分體積分別為：，總儲(chǔ)油量和油位高度變化關(guān)系：5.4綜合考慮儲(chǔ)油罐縱向傾斜和橫向偏轉(zhuǎn)根據(jù)如上討論，我們可以得出結(jié)論，可以直接把5.3中各區(qū)域的公式用式(10)進(jìn)行變換即可得到綜合考慮了儲(chǔ)油罐縱向傾斜和橫向偏轉(zhuǎn)的一般關(guān)系式。我們依然考慮將儲(chǔ)油罐分成五個(gè)區(qū)域分別求解。 對(duì)區(qū)域1的討論當(dāng)儲(chǔ)油罐未發(fā)生偏轉(zhuǎn)時(shí)，油位高度為，而當(dāng)儲(chǔ)油罐發(fā)生橫向偏轉(zhuǎn)后，就可能

17、使得油位探針測(cè)得的示數(shù)變?yōu)榱?。圖14 對(duì)區(qū)域一的討論示意圖即當(dāng)時(shí)，儲(chǔ)油罐發(fā)生橫向偏轉(zhuǎn)后，油位探針測(cè)得的油位高度恰好為零。本屬于第二區(qū)域的部分橫向偏轉(zhuǎn)變位后轉(zhuǎn)為了變位后的第一區(qū)域。在積分時(shí)只需要將替換在第一區(qū)域的積分即可得到的一般關(guān)系式。在此區(qū)域依然恒等于零，的(因此區(qū)域，所以只是求其最大體積)一般關(guān)系式為 (15) 對(duì)區(qū)域2的討論由式(22)可以得到在第二區(qū)域須滿(mǎn)足條件計(jì)算可得再由(12)式： (16)然后只需要做下式的變換即可得到一般關(guān)系式： 對(duì)區(qū)域3的討論同區(qū)域2中討論，應(yīng)該滿(mǎn)足得到：由(13)式： (16)然后只需要做下式的變換即可得到一般關(guān)系式： 對(duì)區(qū)域4的討論應(yīng)該滿(mǎn)足關(guān)系： 對(duì)區(qū)域5

18、的討論在區(qū)域5中油位高度始終保持為。5.5 變位參數(shù)的確定和罐容表的計(jì)算用最小二乘參數(shù)估計(jì)法確定參數(shù)最小二乘參數(shù)估計(jì)法基本思想：根據(jù)的關(guān)系表達(dá)式求得幾組油量高度，計(jì)算出相鄰高度油量的體積之差 通過(guò)與附件的實(shí)際儲(chǔ)油量進(jìn)行比較，通過(guò)對(duì)進(jìn)行等間距的窮舉最終求得理論值與實(shí)際值的差值的平方和，當(dāng)取得最小值，此時(shí)的即為所求的最佳值。即求解如下最小二乘擬合模型算法描述：輸入： 組顯示油高。輸出： 縱橫向偏角的值Step1: 根據(jù)的高度值以及關(guān)系式，求得關(guān)于的表達(dá)式；Step2: 根據(jù)附錄找出實(shí)際對(duì)應(yīng)的出油量；Step3: 對(duì)進(jìn)行等間距窮舉，同時(shí)計(jì)算出，當(dāng)取得最小值時(shí)，求得的值。算法結(jié)束用最小二乘參數(shù)估計(jì)法得

19、到的變位參數(shù)為：，角度都符合實(shí)際情況。5.5.2 實(shí)際儲(chǔ)油罐罐容表的制定估計(jì)出實(shí)際儲(chǔ)油罐的縱向傾斜變位參數(shù)和橫向旋轉(zhuǎn)變位參數(shù)后，我們就可以根據(jù)所建立的儲(chǔ)油量和油位高度以及變位參數(shù)的一般模型計(jì)算得到罐容表。表4 實(shí)際儲(chǔ)油罐罐容表(縱向變位，橫向偏轉(zhuǎn))油位高度00.10.20.22240.30.40.5儲(chǔ)油體積00.04780.35581.05821.27892.20433.67645.40050.60.70.80.911.11.21.37.33459.447111.712814.109216.616519.216221.89124.62431.41.51.61.71.81.922.127.400

20、130.202933.017235.827938.619441.376544.083446.72372.22.32.42.52.62.72.82.949.280751.736454.071956.26658.295260.131461.739263.06652.929363.386264.003464.66455.6利用所給數(shù)據(jù)進(jìn)行進(jìn)一步的檢驗(yàn)我們根據(jù)推導(dǎo)的到的關(guān)系式對(duì)附件2中的第一次補(bǔ)充進(jìn)油后的300組數(shù)據(jù)進(jìn)行了驗(yàn)證，得到在相同油位高度差下，理論出油值值與實(shí)際出油值的關(guān)系圖：圖15 理論儲(chǔ)油量和實(shí)際儲(chǔ)油量的差別對(duì)比圖從圖可以看出所得的理論值與實(shí)際值基本重合，表示所算出的和的值是較為合理的。六

21、 模型的檢驗(yàn)與評(píng)價(jià)6.1 模型的檢驗(yàn)（蒙特卡洛模擬方法）在實(shí)際儲(chǔ)油罐罐容表模型的建立和求解過(guò)程中，我們對(duì)球冠體內(nèi)傾斜的部分燃油的體積進(jìn)行了近似計(jì)算，忽略了一小部分球冠體體積。鑒于此，我們考慮通過(guò)用計(jì)算機(jī)模擬對(duì)該部分的體積進(jìn)行模擬計(jì)算，觀察近似計(jì)算值與精確模擬數(shù)值的吻合情況，同時(shí)也多我們建立的模型進(jìn)行驗(yàn)證。模擬過(guò)程中的主要步驟：Step1: 劃分空間，確定被忽略區(qū)域Q的空間限制范圍，建立空間限制函數(shù)表達(dá)式。并尋求一包含該區(qū)域Q的最小長(zhǎng)方體。建立坐標(biāo)系，確定Q所在的區(qū)域范圍；Step2: 均勻做點(diǎn)，在長(zhǎng)方體內(nèi)分別從三個(gè)坐標(biāo)軸依次等間距的產(chǎn)點(diǎn)，記錄落入該區(qū)域的點(diǎn)以及生成的點(diǎn)的總數(shù)，計(jì)算該長(zhǎng)方體區(qū)域的

22、總體積；Step3：統(tǒng)計(jì)落在該區(qū)域的點(diǎn)的個(gè)數(shù)，求該部分體積，計(jì)算公式為：。在模擬中對(duì)不同區(qū)域分別進(jìn)行求解所忽略部分的體積,再與所得到的理論值相加可得實(shí)際測(cè)量的精確值。根據(jù)蒙特卡洛模擬得到的被省略部分的體積，我們可以畫(huà)出實(shí)際體積和近似體積的差別圖。圖16蒙特卡洛模擬驗(yàn)證圖由圖中可以看到，模型所推導(dǎo)出的結(jié)果與實(shí)際精確的數(shù)值吻合的非常好，模型對(duì)數(shù)據(jù)求解的精度很高，驗(yàn)證了模型的穩(wěn)定性與可行性。6.2 模型的評(píng)價(jià)優(yōu)點(diǎn)：1、模型是由簡(jiǎn)單到復(fù)雜一步步建立的，增強(qiáng)可讀性，邏輯性強(qiáng)；2、在建模過(guò)程中將一些給定的數(shù)據(jù)參數(shù)化，使得模型更具一般性，可應(yīng)用范圍變廣；3、通過(guò)計(jì)算機(jī)模擬的方法較好地對(duì)近似計(jì)算所得結(jié)果做出了

23、驗(yàn)證。缺點(diǎn)：在對(duì)附件2中實(shí)際儲(chǔ)油罐的變位參數(shù)進(jìn)行參數(shù)估計(jì)時(shí)，由于儲(chǔ)油量和油位高度以及變位參數(shù)的函數(shù)關(guān)系非常復(fù)雜，用最小二乘參數(shù)估計(jì)實(shí)現(xiàn)很困難，我們僅采用了附件2中一些數(shù)據(jù)進(jìn)行了參數(shù)估計(jì)，使得得到的變?yōu)閰?shù)與實(shí)際值可能相差較遠(yuǎn)，值得進(jìn)一步改進(jìn)的函數(shù)關(guān)系以使得更好地利用最小二乘參數(shù)估計(jì)法。八 模型的改進(jìn)和推廣7.1 模型的改進(jìn)建立模型積分出的函數(shù)關(guān)系很復(fù)雜，用最小二乘參數(shù)估計(jì)法進(jìn)行參數(shù)估計(jì)時(shí)，程序運(yùn)行時(shí)間較長(zhǎng)。對(duì)此改進(jìn)方案有：1、將模型積分出的函數(shù)關(guān)系式通過(guò)一定的變換，近似化簡(jiǎn)函數(shù)關(guān)系，這樣會(huì)減少一定的程序運(yùn)行時(shí)間；2、利用其它更為簡(jiǎn)便的參數(shù)估計(jì)方法進(jìn)行求解。7.2 模型的推廣對(duì)于不同于題中所給的

24、油罐也可用此方法和理論進(jìn)行分析罐內(nèi)油位高度與儲(chǔ)油量的對(duì)應(yīng)關(guān)系，并且進(jìn)行變位識(shí)別和罐容表的重新標(biāo)定；我們用到的參數(shù)估計(jì)方法也可用到其他的實(shí)際生活中參數(shù)估計(jì)問(wèn)題的求解。參考文獻(xiàn)1 林雪松，周婧，林德新，MATLAB7.0應(yīng)用錦集，北京：機(jī)械工業(yè)出版社，20062 王鄭耀 ，臥式加油灌剩余油料體積的計(jì)算，西安交通大學(xué) 理學(xué)院，20043 王連群，李莉，石油油罐體積計(jì)算方法的探討，吉林化工學(xué)院學(xué)報(bào)，第6期：45-50,1989附錄附錄一：部分結(jié)果附表1 小橢圓儲(chǔ)油罐罐容表(無(wú)變位)油位高度00.010.020.030.040.050.06儲(chǔ)油體積00.00530.01490.02740.0420.05

25、860.07680.070.080.090.10.110.120.130.140.09660.11770.140.16360.18820.21390.24050.26810.150.160.170.180.190.20.210.220.29650.32580.35580.38660.41810.45030.48310.51650.230.240.250.260.270.280.290.30.55060.58520.62040.6560.69220.72890.7660.80350.310.320.330.340.350.360.370.380.84150.87990.91860.95780.9

26、9721.03711.07721.11760.390.40.410.420.430.440.450.461.15831.19931.24051.2821.32371.36571.40781.45010.470.480.490.50.510.520.530.541.49261.53531.57811.6211.66411.70721.75051.79380.550.560.570.580.590.60.610.621.83731.88081.92431.96792.01152.05512.09872.14230.630.640.650.660.670.680.690.72.18582.22942

27、.27292.31632.35962.40292.44612.48910.710.720.730.740.750.760.770.782.53212.57492.61752.662.70232.74452.78642.82810.790.80.810.820.830.840.850.862.86962.91082.95182.99253.0333.07313.11293.15240.870.880.890.90.910.920.930.943.19153.23033.26863.30663.34423.38133.41793.45410.950.960.970.980.9911.011.023

28、.48983.52493.55963.59363.6273.65993.6923.72351.031.041.051.061.071.081.091.13.75433.78443.81363.8423.86963.89623.92193.94661.111.121.131.141.151.161.171.183.97013.99254.01364.03334.05154.06814.08284.09521.191.24.10494.1101附表2 小橢圓儲(chǔ)油罐罐容表(縱向變位)油位高度0.000.010.020.030.040.050.06儲(chǔ)油體積00.00170.00350.00630.01

29、000.01480.02070.02790.070.080.090.100.110.120.130.140.03630.04610.05740.07010.08440.10020.11770.13690.150.150.160.170.180.190.200.210.15130.15780.18020.20400.22890.25480.28180.30970.220.230.240.250.260.270.280.290.33850.36810.39850.42960.46140.49390.52710.56090.300.310.320.330.340.350.360.370.59520.

30、63010.66550.70150.73790.77480.81210.84990.380.390.400.410.420.430.440.450.88810.92670.96561.00491.04451.08451.12471.16530.460.470.480.490.500.510.520.531.20611.24721.28851.33001.37181.41381.45591.49830.540.550.560.570.580.590.600.611.54081.58341.62621.66911.71211.75531.79851.84170.620.630.640.650.66

31、0.670.680.691.88511.92841.97192.01532.05882.10222.14562.18910.700.710.720.730.740.750.760.772.23242.27582.31902.36222.40532.44832.49122.53400.780.790.800.810.820.830.840.852.57662.61902.66142.70352.74542.78722.82872.87000.860.870.880.890.900.910.920.932.91102.95182.99233.03253.07243.11193.15123.1900

32、0.940.950.960.970.980.991.001.013.22853.26673.30443.34163.37853.41483.45073.48601.021.031.041.051.061.071.081.093.52083.55513.58873.62183.65413.68593.71693.74711.101.111.121.131.141.151.161.173.77663.80523.83303.85983.88563.91033.93383.95601.17131.181.191.203.95883.97663.99554.10174.1101附錄二：文中所涉及的程序

33、2.1求解小橢圓儲(chǔ)油罐變位前的公式及罐容表One1.msyms y a b h Lv0=L*int(2*a*sqrt(1-(y-b)2/b2),y,0,h);v=simple(v0)pretty(v)one2.mh=0:0.01:1.2;v=5/9*(-13083/10000*(18/25-3/5*h).(1/2).*h+4361/2000*(18/25-3/5*h).(1/2).*h.2+39249/125000*15(1/2)*asin(1/6*30(1/2)*h.(1/2).*h.(1/2)./h.(1/2)*15(1/2)plot(h,v)grid onaxis(0 1.2 0 4.5

34、)xlabel('油位高度/m')ylabel('儲(chǔ)油量/L')one3.m 測(cè)量值和理論值數(shù)據(jù)對(duì)比v1=50:50:2050 2053.83 2103.83 2105.06 2155.06 2205.06 2255.06 2305.06 2355.06 2404.98 2406.83:50:2906.83 2906.91:50:3706.91;v1=(262+v1)/1000;h=159.02 176.14 192.59 208.50 223.93 238.97 253.66 268.04 282.16 296.03309.69 323.15 336.44 3

35、49.57 362.56 375.42 388.16 400.79 413.32 425.76438.12 450.40 462.62 474.78 486.89 498.95 510.97 522.95 534.90 546.82558.72 570.61 582.48 594.35 606.22 618.09 629.96 641.85 653.75 665.67677.63 678.54 690.53 690.82 702.85 714.91 727.03 739.19 751.42 763.70764.16 776.53 788.99 801.54 814.19 826.95 839.

36、83 852.84 866.00 879.32892.82 892.84 906.53 920.45 934.61 949.05 963.80 978.91 994.43 1010.431026.99 1044.25 1062.37 1081.59 1102.33 1125.32 1152.36 1193.49'h=h/1000;v2=5/9*(-13083/10000*(18/25-3/5*h).(1/2).*h+4361/2000*(18/25-3/5*h).(1/2).*h.2+39249/125000*15(1/2)*asin(1/6*30(1/2)*h.(1/2).*h.(1

37、/2)./h.(1/2)*15(1/2);plot(h,v1,'r',h,v2,'k','markersize',3)xlabel('油位高度h/m')ylabel('儲(chǔ)油量V/m3')v1./v2%變位后理論值和實(shí)際值相比較v1=1.01 1.0583 1.118 1.1675 1.2221 1.2791 1.3274 1.3825 1.4335 1.4849 1.536 1.5918 1.6451 1.6973 1.7498 1.7966 1.8487 1.9 1.9527 1.9991 2.0526 2.10

38、32 2.1484 2.2022 2.2527 2.3037 2.3521 2.4027 2.4522 2.4978 2.5486 2.5995 2.6446 2.6962 2.7488 2.7988 2.8472 2.8958 2.942 2.9925 3.0441 3.0896 3.1408 3.1897 3.234 3.2841 3.3345 3.379 3.4235 3.4731 3.5176 3.5694 3.5732;v2=0.9629 1.0129 1.0629 1.1129 1.1629 1.2129 1.2629 1.3128 1.3628 1.4127 1.4627 1.5

39、127 1.5627 1.6127 1.6627 1.7127 1.7627 1.8127 1.8627 1.9127 1.9627 2.0127 2.0627 2.1127 2.1627 2.2127 2.2627 2.3127 2.3627 2.4127 2.4627 2.5127 2.5627 2.6127 2.6627 2.7127 2.7627 2.8127 2.8627 2.9127 2.9627 3.0127 3.0627 3.1127 3.1627 3.2127 3.2627 3.3127 3.3627 3.4127 3.4627 3.5127 3.5147 ;s1=v1-v2

40、;s=(v1-v2)./v1;h=411.29 423.45 438.33 450.54 463.90 477.74 489.37 502.56 514.69 526.84 538.88 551.96 564.40 576.56 588.74 599.56 611.62 623.44 635.58 646.28 658.59 670.22 680.63 693.03 704.67 716.45 727.66 739.39 750.90 761.55 773.43 785.39 796.04 808.27 820.80 832.80 844.47 856.29 867.60 880.06 892

41、.92 904.34 917.34 929.90 941.42 954.60 968.09 980.14 992.41 1006.34 1019.07 1034.24 1035.36;h=h/1000;plot(h,s)grid onxlabel('油位高度/m')ylabel('儲(chǔ)油量V/立方米')figure(2)plot(h,s1)grid onxlabel('油位高度/m')ylabel('儲(chǔ)油量/m³')figure(3)plot(h,v1,'k*',h,v2,'r')grid

42、onxlabel('油位高度/m')ylabel('儲(chǔ)油量/m³')clear allsyms y hva=2*int(0.89*sqrt(1-y2/0.62)*(h+0.4*0.0717-0.6+y)/0.0717),y,(0.6-(h/0.0717+0.4)*0.0717),0.6);for h=0:0.0005:2.05*tan(4.1*pi/180) r=eval(va); hold on plot(h,r) grid onendaxis(0 1.4 0 4.5)syms hv1=2*int(h+0.4*0.0717)-y)/0.0717*(0

43、.89*(1-(y-0.6)/0.6)2)0.5),y,h,(h+0.4*0.0717);v=0.89*0.6*2.45*(h-0.6)/0.36*(h*(1.2-h)0.5+asin(h/0.6-1)+pi/2);v2=2*int(2.45-(h+0.4*0.0717)-y)/0.0717)*(0.89*(1-(y-0.6)/0.6)2)0.5),y,(h-2.05*0.0717),h);v=v1+v-v2;for h=2.05*tan(4.1*pi/180):0.005:1.2-0.4*tan(4.1*pi/180) r=eval(v); hold on plot(h,r,'r.&

44、#39;) grid onendsyms y hva=4.1101-2*int(0.89*sqrt(1-y2/0.62)*(1.2-h+2.05*0.0717-0.6+y)/0.0717),y,(0.6-(1.2-h)/0.0717+2.05)*0.0717),0.6);for h=1.2-0.4*tan(4.1*pi/180):0.0001:1.2 r=eval(va); hold on plot(h,r,'g.') grid onendh=0:0.01:1.2;v=5/9*(-13083/10000*(18/25-3/5*h).(1/2).*h+4361/2000*(18/

45、25-3/5*h).(1/2).*h.2+39249/125000*15(1/2)*asin(1/6*30(1/2)*h.(1/2).*h.(1/2)./h.(1/2)*15(1/2);plot(h,v)xlabel('油位高度/m')ylabel('儲(chǔ)油量/m³')gtext('變位前')gtext('變化后')%變位前后理論值相比較clear allsyms y hva=2*int(0.89*sqrt(1-y2/0.62)*(h+0.4*0.0717-0.6+y)/0.0717),y,(0.6-(h/0.0717+0

46、.4)*0.0717),0.6);for h=0:0.0005:2.05*tan(4.1*pi/180) r=eval(va); hold on plot(h,r) grid onendaxis(0 1.4 0 4.5)syms hv1=2*int(h+0.4*0.0717)-y)/0.0717*(0.89*(1-(y-0.6)/0.6)2)0.5),y,h,(h+0.4*0.0717);v=0.89*0.6*2.45*(h-0.6)/0.36*(h*(1.2-h)0.5+asin(h/0.6-1)+pi/2);v2=2*int(2.45-(h+0.4*0.0717)-y)/0.0717)*(

47、0.89*(1-(y-0.6)/0.6)2)0.5),y,(h-2.05*0.0717),h);v=v1+v-v2;for h=2.05*tan(4.1*pi/180):0.005:1.2-0.4*tan(4.1*pi/180) r=eval(v); hold on plot(h,r,'r.') grid onendsyms y hva=4.1101-2*int(0.89*sqrt(1-y2/0.62)*(1.2-h+2.05*0.0717-0.6+y)/0.0717),y,(0.6-(1.2-h)/0.0717+2.05)*0.0717),0.6);for h=1.2-0.4

48、*tan(4.1*pi/180):0.0001:1.2 r=eval(va); hold on plot(h,r,'g.') grid onendh=0:0.01:1.2;v=5/9*(-13083/10000*(18/25-3/5*h).(1/2).*h+4361/2000*(18/25-3/5*h).(1/2).*h.2+39249/125000*15(1/2)*asin(1/6*30(1/2)*h.(1/2).*h.(1/2)./h.(1/2)*15(1/2)plot(h,v)2.2v1=5.4400 7.3156 9.3777 11.6020 13.9601 16.43

49、59 19.0078 21.6554 24.3670 27.1208 29.9065 32.7028 35.4998 38.2793 41.0227 43.7196 46.3481 48.8963 51.3434 53.6675 55.8495 57.8667 59.6864 61.2643 ;v01=0.0202 0.0270 0.0399 0.0405 0.0466 0.0520 0.0566 0.0603 0.0630 0.0646 0.0651 0.0645 0.0629 0.0602 0.0565 0.0519 0.0466 0.0406 0.0342 0 .0276 0.0210

50、0.0148 0.00910.0045 0.00019 ;v02=flipud(v01); v0=v01-v02;v0(1)=; v2=v1+v0;hold onplot(v1,'r*')plot(v2)grid on%估計(jì)參數(shù)和detav=243.85/1000; minv=100;h1=2485.73/1000;h2=2474.4/1000;L=8;R=1.5;R0=1.625;H=1;l=2;syms yfor a=linspace(0.001,4.5*pi/180,20); for b=linspace(0.001,4.5*pi/180,20); hb=h1; v1=-6*cos(b)*(9-9*cos(b)2+.+2*cos(b)*hb)(1/2); hb=h2; v2=-6*cos(b)*(9-9*cos(b)2+.+2*cos(b)*hb)(1/2); mv=v1-v2; if minv>abs(mv-detav) minv=ab