代碼以及過程學習來源于：西北工業大學 肖華用老師課程 代碼在文章最后

附上鏈接：?9.5 插值與擬合模型(二)-----水塔流量問題_嗶哩嗶哩_bilibili

也可以去中國mooc觀看肖老師視頻

問題描述：美國某洲的各用水管理機構要求各社區提供以每小時多少加侖計的用水率以及每天總的用水量,但許多社區并沒有測量水流入或流出水塔水量的設備,他們只能每小時測量水塔中的水位，精度在0.5%以內,更為重要的是,無論什么時候,只要水塔中的水位下降到某一最低水位L時，水泵就啟動向水塔重新充水至某一最高水位H，但也無法得到水泵的供水量的測量數據。水泵每天向水塔充水一次或兩次，每次約兩小時。水塔是一個垂直圓形柱體，高為40英尺，直徑57英尺。?，當水位升至35.50英尺時停止充水。試估計在任何時刻，甚至包括水泵正在工作期間內，水從水塔流出的流量f(t)，并估計一天的總用水量。

表1?小鎮某天的水塔水位(0.01英尺)

?表2?各時刻水體積表

要從水的流速得到水的流量，可以使用直接積分法，對0到24小時的水的流速進行積分計算即可得到一天的總用水量。或者使用分部積分法，對第一次充水前，第一次充水到第二次期間，第二次充水后，這三段進行分部積分，然后將兩種積分方法進行比較。

計算各時刻點的水流量(加侖/小時)：

?

25個時刻處的水流量采用差分的方法得到，共分三段分別處理

差分公式為：

對每段前兩點采用向前三點差分公式：

?

對每段最后兩點采用向后三點差分公式：

?

對每段中間點采用中心差分公式：

?

?表3?各時刻水流量

?使用數值積分的方式，求解一天的總用水量

方法一：直接積分法

?方法二：分段計算

?

clear,clc c=0.3048;%英尺和米的換算 p=1.0/3.785;%升和加侖的換算 d=57*c;%直徑 h=31.75*c; %v=pi*d*d*h/4*1000*p; data=[0,3175;3316,3110;6635,3054;10619,2994;13937,2947;17921,2892;21240,2850;25223,2797;28543,2752;32284,2697;39435,3550;43318,3445;46636,3350;49953,3260;53936,3167;57254,3087;60574,3012;64554,2927;68535,2842;71854,2767;75021,2697;82649,3550;85968,3475;89953,3397;93270,3340];%數據導入 t=data(:,1)/3600;%計算時間 v=pi*d*d*data(:,2)/100*c/4*1000*p;%計算體積 %計算差分 n=length(v); f=zeros(n,1);%儲存差分值 %計算第一段 n1=10; for i=1:n1if i<=2 %前兩點用向前三點差分f(i)=abs(-3*v(i)+4*v(i+1)-v(i+2))/(2*(t(i+1)-t(i)));elseif i<=n1-2%中心差分公式f(i)=abs(-v(i+2)+8*v(i+1)-8*v(i-1)+v(i-2))/(12*(t(i+1)-t(i)));elseif i>=n1-1 %后兩點用向后三點差分f(i)=abs(3*v(i)-4*v(i-1)+v(i-2))/(2*(t(i)-t(i-1)));end end %計算第二段 n2=21; for i=n1+1:n2if i<=n1+2f(i)=abs(-3*v(i)+4*v(i+1)-v(i+2))/(2*(t(i+1)-t(i)));elseif i<=n2-2f(i)=abs(-v(i+2)+8*v(i+1)-8*v(i-1)+v(i-2))/(12*(t(i+1)-t(i)));elseif i>=n2-1f(i)=abs(3*v(i)-4*v(i-1)+v(i-2))/(2*(t(i)-t(i-1)));end end %計算第三段 n3=25; for i=n2+1:n3if i<=n2+2 f(i)=abs(-3*v(i)+4*v(i+1)-v(i+2))/(2*(t(i+1)-t(i)));elseif i<=n3-2f(i)=abs(-v(i+2)+8*v(i+1)-8*v(i-1)+v(i-2))/(12*(t(i+1)-t(i)));elseif i>=n3-1f(i)=abs(3*v(i)-4*v(i-1)+v(i-2))/(2*(t(i)-t(i-1)));end end plot(t,f,'r*');%原始圖 tmin=min(t);tmax=max(t); tt=tmin:0.1:tmax;%獲取離散的時間點，以便做樣條曲線 ff=spline(t,f,tt);%計算三次樣條插值 hold on plot(tt,ff,'b');%畫樣條曲線 xlabel('時間(小時)'); ylabel('水塔流量圖'); hold off dt=0.05; t2=0.5:dt:24.5;%獲取離散時間點，以便積分 nn=length(t2); f2=spline(t,f,t2); %計算24小時用水量，使用復化梯形公式 s=(f2(1)+f2(nn)+2*sum(f2(2:nn-1)))*dt/2; fprintf('(全部積分法)1天總水量s=%8.2f\n',s); %第一次水塔增加的水 v10=v(11)-v(10); dt1=t(11)-t(10); %第一次充本其間流出的水 tp=t(10):dt:t(11); nn=length(tp); yp=spline(t,f,tp); %計算三次樣條插值 v11=(yp(1)+yp(nn)+2*sum(yp(2:nn-1)))*dt/2; v1=v10+v11;%第一次充水總量 p1=v1/dt1; %第一次充水的平均水流量 %第一次充水前總流量 vv1=v(1)-v(10); %兩次充水間總流量 vv2=v(11)-v(21); % t=83649到85968其間流量 vv3=v(22)-v(23); %第一次充水期間流量 ta=t(10):dt:t(11); %獲得離散的時間點，用于積分 nn=length(ta); fa=spline(t,f,ta); s1=(fa(1)+fa(nn)+2*sum(fa(2:nn-1)))*dt/2; %第二次充水期間流量 tb=t(21):dt:t(22); %獲得離散的時間點用于積分 nn=length(tb); fb=spline(t,f,tb); s2=(f(1)+fb(nn)+2*sum(fb(2:nn-1)))*dt/2; %t=85968到86400流量 tc=t(23):dt:24; %獲得離散的時間點，用于積分 nn=length(tc); fc=spline(t,f,tc); s3=(fc(1)+fc(nn)+2*sum(fc(2:nn-1)))*dt/2; ss=vv1 +vv2+vv3+s1+s2+s3; fprintf('部分積分法)1天總水流量ss= %8.2f\n',ss); err=abs((s-ss)/s); fprintf('兩種計算法總量相對誤差%6.2f%%\n',err* 100);

總結

以上是生活随笔為你收集整理的matlab实现水塔水流量问题的全部內容，希望文章能夠幫你解決所遇到的問題。

如果覺得生活随笔網站內容還不錯，歡迎將生活随笔推薦給好友。