建模作品(储油罐)

通常许多大型液体储罐（如加油站的地下储油罐、酒厂的酒罐、食用油的储油罐等）一般都有与之配套的“液位计量管理系统”，采用流量计和液位计来测量进/出液体量与罐内液位高度等数据，通过预先标定的罐容表（即罐内液位高度与储存液体量的对应关系）进行实时计算，以得到罐内液位高度和储液量的变化情况。

许多大型液体储罐在使用一段时间后，由于自然原因或人为原因而导致罐体的位置会发生纵向倾斜和横向偏转等变化，从而导致罐容表发生改变，按照有关规定，需要定期对罐容表进行重新标定。

本文我们利用数学模型方法以研究解决加油站储油罐的变位识别与罐容表标定的问题为例，来讨论大型液体储藏罐变位对罐容表的影响，以及新罐容表的制定。

（1）为了掌握罐体变位后对罐容表的影响，我们先考虑小椭圆型储油罐（两端平头的椭圆柱体），分别对罐体无变位和倾斜角为α=4.10的纵向变位两种情况做了实验，建立数学模型研究罐体变位后对罐容表的影响，并给出了罐体变位后油位高度间隔为1cm 的新罐容表标定值。

（2）对于典型的储油罐，其主体为圆柱体，两端为球冠体，我们建立罐体变位后标定罐容表的数学模型，即罐内储油量与油位高度及变位参数（纵向倾斜角度α和横向偏转角度β）之间的一般关系。利用罐体变位后在进/出油过程中的实际检测数据，建立的数学模型确定变位参数，并给出罐体变位后油位高度间隔为10cm的罐容表标定值。最后利用实际检测数据来分析我们所建立模型的正确性与方法的可靠性

一、摘要

本文通过研究问题中所给出的小椭圆型储油罐的结构，分析出了罐体在无变位情况下罐内油位高度与储油量之间的关系，建立了其数学模型.在此基础上，本文通过讨论在纵向倾斜角 4.1,0αβ=?=时，找出了罐内油位高度与等体积油量下的等效高度的关系，从而建立了罐体变位后油位高度与储油量的数学模型.并通过问题中所给出的观测数据，我们分别拟合出了罐体在变位前后其理论值的函数曲线与观测值的函数曲线.通过分析比较，我们利用最小二乘法建立了其误差随油位高度变化的函数，通过此函数以确保我们建立的数学模型更加准确无误，使之更具有说服力.最后利用建立的模型我们分析出了罐体变位后对罐容表的影响结果，并计算出变位后的罐容表标定值.

本文利用实际储油罐的示意图，先建立了罐体在无变位的环境下（即0,0αβ==时），标定罐容表的数学模型.在此模型的基础之上，本文通过讨论在两种变位方式（即纵向倾斜和横向偏转）下，观测到的油位高度与其同等体积下的等效高度之间的关系，建立了等效高度与纵向倾斜角度α和横向偏转角度β之间的一般关系，于是我们就得出了罐内储油量与油位高度及变位参数（即倾斜角度α和横向偏转角度β）之间的一般关系.进而利用问题中给出的实际检测数据，建立出一个理论计算值和观测值的偏差关于参数,αβ的数学模型，显然，偏差最小取得最小值时就是我们估计出的参数值.所以我们将其转化成了最优化求解的问题，针对此问题，我们采用SNTO 算法，即序贯数论算法，确定出了非常接近实际的参数值，利用此参数值，我们计算出了问题中所要求给出的罐容表标定值，最后我们利用实际检测数据检测的结果与估计出的值进行比较进一步证明了模型的正确性与方法的可靠性.

本文还对该油位计量管理系统的偶然误差和系统误差做出了定量的分析，并给出了误差与油位高度之间的变化关系，进一步对建立的模型进行修正.

二 符号说明

1V 是无变位情况下小椭圆型储油罐的储油量. H 是油位高度.

L 是小椭圆型储油罐的罐长.

1R 是小椭圆型储油罐的椭圆截面的短半轴. 2R 是小椭圆型储油罐的椭圆截面的长半轴. 1()S y 是椭圆截面微分成小矩形的面积.

1dV 是整个椭圆柱体微分成的一个小长方体的体积. 2V 是变位情况下小椭圆型储油罐的储油量.

2H 是变位情况下小椭圆型储油罐的等效油位高度. 1L 是油位尺在罐底的垂点距罐侧面的距离. 0D 是下尺点内竖直径.

L '是油罐变位后特定情况下的等效底边长. i v 是每组观测的储油量. i h 是每组观测的油位高度.

三 模型的假设

3.1、假设题目所给的数据真实可靠； 3.2、假设油的温度恒定不变；

3.3、假设油的黏稠度可以忽略不计； 3.4、假设油罐内壁厚度可以忽略不计；

3.5、题目中给定的角度不变即未发生二次位变；

3.6、假设倾斜方向如图中所示且一直按照一个方向倾斜； 3.7、假设发生倾斜的过程当中只发生纵向倾斜； 3.8、假设油罐内部管子体积可以忽略不计； 3.9、假设所有的数据是浮标稳定后测得的；

四 问题的分析与模型的建立

对于小型椭圆型储油罐分析与模型的建立

4.1 分析：

我们把椭圆截面微分成许多个小矩形，这些小矩形的面积用1()S y 表示，再将整个椭圆柱体微分成一个个小长方体，这些小长方体的体积用1dV 表示，最后将其积分.

其推导过程如下：

1

11122

222122221111

111222111111()()211

2[arcsin ]2

1

[(arcsin(1)]2

h

R dV S y dy

S y xL x y R R x R R h V L

L R R R R R h H R R H V L H R R R R R ππ-==+==

==

+=-=

--+?

根据以上分析，我们建立了在罐体无变位的情况下罐内储油量1V 与罐内油位高度H 的数学模型：

2211111211

()[(arcsin(1)]2

R H V H L H R R R R R π=--+ （1）

其中2R 为椭圆的长半轴，1R 为短半轴，L 为椭圆柱体的长度，()erf H 为误差值与罐内油位高度H 的函数.

4.2分析：

如图1所示，直圆筒为椭圆筒，以液面高端截面为参考面，分三个区段计算罐内储油量的体积2V ，罐的倾斜角为 4.1α=?，圆筒长为L ，椭圆筒横截面椭圆的长半轴2R ，短半轴1R ，下尺点F 至液面高端的距离为1L

.

由图2可知：

1

111111

112,,2,,,(2)cos ,,cos 2(2)cos cos R WP H QG H FG R GO L CO L FW R H FP R WP FP FW R H αα

α

α

======-==-=--

则

2112tan cos H H R αα

=- （2）

令

11tan N H L α=+

椭圆方程为2222211X Y R R +=,即2

212

R y R x R =-.液高为22()H CD H =,即112()y R H =--,亦即直线AB 的直线方程.将1y 代入椭圆方程得：

11212(2)R

x H R H R =-，微元面积：

2

211122,()R dA ydx y y y R H R x R -==-=---

液体截面积为：

11210022()x x S ydx H R dx ?==-+??

??

2

1221222()arcsin 2R R x H R x R R ??=-+????

1212

(R H R R ??

=-???

2

212R ???????

21122()H R R H =

-+-??

121

R R +所求体积为：

21121

)V H R LR R =

-+ (3)

式（3）的使用条件是11

tan 2N R L α≤+,即水平液高21H R ≤.

当tan N L α≤，则21

2H N =,底长用cot L N α'=代替L ;

当11tan tan 2L N R L αα<≤+时，则21

tan 2

H N L α=-.

液位继续升高，当12N R ≤(即水平液高21H R >)的情况)时,用2122H R H '=-来代替

（3）式中的2H .则这部分容积的计算式为：

1212121)V R R L R H LR R π=--+??

(4)

其中式（4）中的水平液高21

tan 2

H N L α=-.

当液位升至量油孔处，即12N R >,而10H D ≤(下尺点内竖直径)时，用（112R H -）

代替1H ，用（1L L -）代替1L ，则1112()tan N R H L L α=-+-,此时水平液高21

2

H N =，

底长cot L N α'=，所求体积为：

1221121)cot V R R L H R R R N πα=--+??

(5)

通过以上分析计算，我们建立出罐内储油量2V 与观测高度2H 的数学模型： 1）当21H R ≤时：

221121

)V H R LR R =-+ （6）

在这个范围内，当t an N L α≤，则21

2

H N =

,底长用cot L N α'=代替L ；当11t a n t a n 2L N R L αα<≤+

时，则21tan 2

H N L α=-. 2）当21H R >，且111

tan 22

R L N R α+<≤时：

21212121)V R R L R H LR R π=--+??

（7）

其中21

tan 2

H N L α=-.

3）当液位升至量油孔处，即12N R >,且10H D ≤(下尺点内竖直径)时：

21221121)arcsin cot V R R L H R R R N πα

=--+??

（8）

其中21

2

H N =，cot L N α'=.

对于典型储油罐分析与模型的建立

通过对储油罐在无变位情况下的分析.首先我们先将储油罐内储油量V 分成三个部分，即中间的圆柱部分储油量1V ，和两边球缺部分的储油量2V （由于左右的对称的，所以我们只需要计算其中一个）.

为简单计算我们先计算图示的灰色部分体积V

.

（1） 计算中间圆柱体部分储油量1V

显然，1V LS =，其中S 为其圆柱体截面.因此我们只需要计算截面面积S .

利用微积分，dS =，计算定积分

22

20

arctan S R h R ??== ? ?-?

?

?

所以，

2

1

()arctan

V h LS LR

==

??

（9）

计算两端球缺部分储油量

2

V

同之前一样的思路，先计算水平截面的面积()

S h，然后在垂直方向上对h积分：

2

()

h

V S h dh

=?

下面求()

S h：

如图求的水平截面是一个圆，半径

h

r=，()

S h就是上图阴影部分的面积.

显然，

()h r

r H r H

S h

--

==

?

计算得

2222

1

()(()()arctan

2

S h H h r h h r

π

??=--+-从而

2

()()

h

V h S h dh

=?

23

2

()

26

hr h

V h

ππ

=-+

32332arctan 3H H r r ??

-+-

323arctan 3h h r ??--

32arctan 3r ??+ 于是，得到了()V h 的计算公式

12()()2()V h V h V h =+

分别代入求得：

3

2

22232

34

()33

)arctan 24(3)arctan arctan 33h V h hr LR r R h hr r ππ=-

+????

++ ??????+-+

下面来导出储油罐中储油量2V 和和油位高度h 的关系式（如下图所示，分2种情况）：

()()2

()()()2

g

g V V h R h R V h V V R h h R ?+-??????>??=??--??????≤?? （10）

其中g V 表示储油罐的总容积.

(

2

2

2(23

g V R L r r ππ=+

下面我们分析在两种变位方式只为横向偏转（即纵向无偏转）情况下，观测到的油位高度H '与其同等体积下的无偏转等效高度H 之间的关系：

横向倾斜时，如图所示

通过分析可知

cos (1cos )H H R ββ'=+-

将其代入()V h 得出罐容表在横向偏转时储油量与油位高度的数学模型：

cos ()()2

()cos ()()2

g

g V V H R H R V H V V R H H R ββ?''+-??????>??'=??''--??????≤?? （11）

当罐体在罐体发生变化时，罐内油位高度H ''与等效高度H 的关系： 1）当1H R '≤时：

在这个范围内，当tan N L α≤，则1

2

H N '=,底长用co t L N α'=代替L ；当

1t a n t a n 2L N R L αα<≤+时，则1

tan 2

H N L α'=-.

2）当1H R >，且1

tan 22

R L N R α+<≤时：

用1

tan 2

H N L α'=-代替式（11）中的H '.

3）当液位升至量油孔处，即2N R >,且10H D ≤(下尺点内竖直径)时：

用1

2

H N '=代替式（11）中的H '，用cot L N α'=代替式（11）中的L ,其中

1tan

N H L α''=+ 分别将上述H ''与无偏转等效高度H 的表达式代入模型（6）即可得到

(,,)V V h αβ=

因为问题中观测到的N 组实际监测数据,()i i v h ，都满足含有两个参数,αβ的函数表达式(,,)V h αβ，任意给定一组有意义的参数,αβ计算值和测量值的总偏差是：

1

(,)((,,))N

i i i F V h v αβαβ==-∑

显然，使得计算值和测量值的总偏差(,)F αβ最小的参数,αβ就是实际的参数，于是变成下面的最优化问题：

**0,(,)(,)F Min F αβπ

αβαβ<<=

解决此最优化求解问题，我们采用了SNTO 算法，即序贯数论算法，在参数空间中“均匀”的撒入很多点，然后在计算函数在各点的值.并通过这些点值逐步缩小参数的空间范围，上述过程继续进行，直到找到足够精确的参数为止.

接下来我们用实际监测数据检测了模型的正确性，我们分别闭合了其函数图像，通过图像进一步证明了方法的可靠性.

最后我们用此对参数值代入模型，计算出油位高度为10cm 的罐容表标定值.

五 模型的求解

对于小型椭圆型储油罐模型的求解 5.1由题目中的信息我们知道：

20.89R m =，10.6R m =， 2.45L m =

代入（1）式可得：

1() 3.6342[(0.36arcsin(1.66671)0.5655]V H H H =--+

通过这个模型可以由观测到的罐内油位高度计算出对应的理论上的储油量（见附表一）.通过图像

我们发现存在一定误差，这些误差可能是由于温度，黏稠度，罐壁厚度、罐内管子体积等一系列因素造成的.因此我们假定一个误差值V ?与实际观测的油位高度H '之间的函数关系：

()V erf H '?=

根据得到的数据，我们拟合出了其函数图像（参考文献[1]）：

通过图像我们假定其函数为三次函数，利用最小二乘法原理求得其误差函数模型：

32()0.08410.15070.05820.0017V H H H H '''?=-++-

最后我们得出罐体无变位罐容表中的油位高度H 与储油量0V 之间的模型：

0() 3.6342[(0.36arcsin(1.66671)0.5655]V H H H V =--+-?

利用这个模型我们计算出了罐体无变位的储油表（见附表二）.

5.2由题目中的信息我们知道：

11214.1, 2.45,0.4,0.6,0.89,0.0287L m L m R m R m N H α=?=====+ 将数据代入式（6）得出模型：

(220.6) 1.3083arcsin V H =-+（12） 令tan N L α=，可以计算出10.1469H m =. 1）当100.14H m ≤≤时：

将观测油位高度代入式（12）得出观测油量（见附表三）.

通过带入数据求得2110.50.0144,13.95070.4004H H L H '=+=+.将其代入上述模型得到理论储油量（见附表三）.

令

11

tan 2

N R L α=+

带入数据，解得：

10.6591H m =

2）当10.140.65m H m <≤时：

将观测油位高度代入式（12）得出观测油量（见附表三）.将等效高度210.0591H H =-代入式（12）得到理论储油量（见附表三）.

3）

将题目中数据代入式（7）得：

(224.1080.6)1.3083arcsin V H ?=--+

?

???

（13）

令

12N R =

解得：

1 1.1713H =

当10.65 1.17H <≤时：

将观测油位高度代入式（13）得出观测油量（见附表三）. 210.0591H H =-

将等效高度210.0591H H =-代入5.2.2得出的模型中得到理论储油量（见附表三）. 4）

将数据代入式（8）得：

(224.10827.90140.6)0.534arcsin V H H ?=--+

??

??

（14）

根据图（2）可求得

0 1.2D m =

当11.17 1.2H <≤时，210.50.0144H H =+，并用113.95070.4004L H '=+代替式（8）中的L ,此时，

将观测油位高度代入模型得出观测油量（见附表三）.

将等效高度210.50.0144H H =+代入式（14）得到理论储油量（见附表三）. 利用求得的数学模型，我们计算出了罐体变位后观测到的整体的数据（见附表三）. 最后，通过得到的这些数据，首先我们分析了变位后的罐体对之前的罐容表的影响.所谓影响既是同等油位高度下储油量之差的变化规律，于是我们拟合出：

通过图像我们可以看出其中的变化规律：

H 在0~0.4之间随油位高度的增加影响越来越大，

0.6~0.7之间突然变小之后又突然增大，0.7~1.2之间逐渐变小.

其次通过这些观测数据和理论数据上分析，我们发现之间存在一定的误差，于是我们拟合出了观测数据与理论数据的差值与油位高度变化的图像，如图：

通过观察，我们发现与二次图像接近，于是我们用最小二乘法估计出了其误差函数模型：

2

220.38360.57630.1267V H H ?=-+- 于是罐体变位后罐容表的油位高度H 与储油量V 之间的数学模型就是： 2()()()V H V H V H =-?

我们利用这个模型计算得出罐体变位后油位高度间隔为1cm 的罐容表标定值（见附表四）.

这里特别指出的是当油位高度为0刻度时，罐内储油量的体积在0~124.1207L 之间. 对于典型的储油罐模型的求解

题目中给出的信息： 1.5, 1.625,8,R m r m L m ===

32.3322cos ( 1.5)( 1.5)

()32.3322cos (1.5)( 1.5)

V h h V h V h h ββ?+-??????>?=?--??????≤??

其中

33

()8.2958 1.047280.833314.8617arctan 0.6667 2.6406)arctan 5.7214arctan V h h h h h =-+??+??+- ?

??+ ?

（15）

1）当 1.5h m ≤时：

当0.3143N ≤，则1

2

h N =,底长用25.4517L N '=代替式（6）中的L

当0.34138.1707N <≤时，此时1

tan 2

h N L α=-.

2）当 1.5H m >，且1.67073m N m <≤时，此时0.1707h N =-.

3）当液位升至量油孔处，即3N m >,且 3.0024h m ≤(下尺点内竖直径)时，用

25.4517L N '=代替L ，此时1

2

h N =

将以上三种不同情况下的h 分别代入式（15）中，我们得出了储油量V 与油位高

度H ，竖向倾斜角α横向偏转角β之间的关系：

32.3322(,)cos ( 1.5)( 1.5)(,,)32.3322(,)cos (1.5)

( 1.5)V H H H V H V H H H αβαβαβ?+->?

=?--≤?? （16）

首先我们将,αβ取遍区间[0,]π内以0.01为间隔的所有的点集合，将这些点的取值代入式（16），计算出每组参数对应的(,)F αβ，找出使之取值最小那对参数组合.

经过反复的计算（详细数据见附录），我们找出了使得(,)F αβ在此区间内取到的最小值的参数组合：

2.1, 4.2αβ=?=?

将此值代入到模型（16）中，得出：

32.33220.9986(,2.25)( 1.5)

( 1.5)(,,)32.33220.9986(,2.25)(1.5)

( 1.5)V H H H V H V H H H αβ?+?->?=?-?-≤??

其中：

33

()8.2958 1.047280.833314.8617arctan 0.6667 2.6406)arctan 5.7214arctan V H H H H H =-+??+??+- ???+ ?

（1）当 1.5H m ≤时：

1)当0.3143N ≤时，底长用25.4517L N '=代替L .此时12

H N =

2）当0.34138.1707N <≤时， 1

tan 2

H N L α=-.

（2）当 1.5H m >，且1.67073m N m <≤时，0.1707H N =-.

（3）当液位升至量油孔处，即3N m >,且 3.0024H m ≤(下尺点内竖直径)时，用

25.4517L N '=代替式（9）中的L ，此时1

2

H N =.

利用此模型我们计算出罐体变位后油位高度间隔为10cm 的罐容表标定值（见附表五）.

模型的分析与检测：

利用问题附件中实际检测数据我们比较了在同等进/出油量的情况理论上的数值与实际数值的油位高度变化情况，我们拟合出了其函数图像：

通过对图像的观察，我们发现误差在2%之内，因此进一步验证了我们所建立的模型的正确性与方法的可靠性.

六模型的评价与改进

本文通过设计所需要的模型有效的模拟出储油量与油位高度以及变位参数之间关系的数学模型，具有较强的现实意义，可以在现实生活中广泛推广应用.

对罐内储油量与油位高度以及变位参数之间关系进行预测，拟合效果良好.而且，设定参数后，可以对模型进行实时控制和实时预测，有利于模型的建立与求解.

但由于数据的可得性以及准确性，在选取的数据上存在一定的误差，但都尽量控制在合理的范围内.对于容量测量误差，液位测量误差，温度与油品密度的因素等次要因素均加以忽略，与实际情况存在偏差.

另外，模型模拟实验不可能与实际情况完全一致，而且在模型建立时，做出了些模型假设，也忽略了多种并未量化的影响因子，数据收集误差等因素.因此，不能完全准确的进行分析，模型仍需要修正和完善.

七参考文献

[1] 郑阿奇，MATLAB实用教程，出版地：电子工业出版社，2004

[2] J.F.博南，最优化问题的扰动分析，出版地：科学出版社，2008

[3] 郑阿奇，MATLAB实用教程，出版地：电子工业出版社，2004