✅作者简介：热爱科研的Matlab仿真开发者，修心和技术同步精进，matlab项目合作可私信。
🍎个人主页：Matlab科研工作室
🍊个人信条：格物致知。
更多Matlab仿真内容点击👇
⛄ 内容介绍
%justw 水库调洪演算自编程序
%水量平衡方程:(Q1+Q2)*dt/2-(q1+q2)*dt/2=V2-V1
%水库蓄泄方程:q=m1*L*(z-z0)^1.5+m2*A*(z-z0)^0.5
%zm--演算结果:[初始水位m 最高水位m 水位上涨m 初始库容1e4m3 最大库容1e4m3 调洪库容1e4m3 入库洪峰m3/s 最大泄量m3/s 洪峰削减m3/s]
%Qt--洪水过程线:{[时间h 流量m3/s]}
%Vz--水位库容曲线:[水位m 库容1e4m3]
%fq--蓄泄参数:[流量系数 堰顶净宽m/洞口面积m2 堰顶高程m/洞心高程m 溢洪道(1.5)/泄洪洞(0.5)]
%z0--初始参数:[初始水位m 时段长度h 洪水线性插值<1>/样条插值(2)]
⛄ 部分代码
function zm=justw(Qt,Vz,fq,z0)
if nargin<1, Qt=fin('Qt',''); end
if nargin<2, Vz=fin('Vz',''); end
if nargin<3, fq=fin('fq',''); end
if nargin<4, z0=fin('z0',''); end
Qt=chk('Qt',Qt);
Vz=chk('Vz',Vz);
fq=chk('fq',fq);
z0=chk('z0',z0);
if isempty(z0), z0=[min(fq(:,3)) 0.1 1]; end
save justw Qt Vz fq z0;
method={'l';'s'};
zm=zeros(length(Qt),9);
qt=cell(size(Qt));
for i=1:length(Qt)
t=union(Qt{i}(:,1),[0:z0(2):max(Qt{i}(:,1))]');
t(diff(t).*3600<1)=[];
Q=interp1([0;Qt{i}(:,1)],[0;Qt{i}(:,2)],t,method{z0(3)});
z=repmat(z0(1),size(t));
for j=1:length(t)-1, z(j+1)=funz2(Q(j:j+1),t(j:j+1),z(j),Vz,fq); end
q=funq(fq,z);
zm(i,[1 2 4 5 7 8])=[z0(1) max(z) interp1(Vz(:,1),Vz(:,2),[z0(1) max(z)],'s') max(Q) max(q)];
zm(i,[3 6 9])=[diff(zm(i,[1 2])) diff(zm(i,[4 5])) diff(zm(i,[8 7]))];
qt{i}=[t Q q];
end
for i=1:length(qt)
subplot(length(qt),1,i);
plot(qt{i}(:,1),qt{i}(:,2),qt{i}(:,1),qt{i}(:,3));
end
str=sprintf('\n水库调洪演算\n计算程序:justw\n计算时间:%s\n计算结果......\n',date);
str=[str sprintf('\n 初始水位 最高水位 水位上涨 初始库容 最大库容 调洪库容 入库洪峰 最大泄量 洪峰削减')];
str=[str sprintf('\n m m m 万m3 万m3 万m3 m3/s m3/s m3/s')];
str=[str sprintf('\n%12.3f %12.3f %12.3f %12.2f %12.2f %12.2f %12.2f %12.2f %12.2f',zm')];
disp(str);
fid=fopen('justwout.txt','w+t');
fprintf(fid,'%s\n',str);
fclose(fid);
function y=fin(opt,str)
%fin 输入原始数据
switch opt
case 'Qt', str2='Qt--洪水过程:{[时间h 流量m3/s]} <>:';
case 'Vz', str2='Vz--水位库容曲线:[水位m 库容1e4m3] <>:';
case 'fq', str2='fq--蓄泄参数:[流量系数 堰顶净宽m/洞口面积m2 堰顶高程m/洞心高程m 溢洪道(1.5)/泄洪洞(0.5)] <>:';
case 'z0', str2='z0--初始参数:[初始水位m 时段长度h<0.1> 洪水线性插值<1>/样条插值(2)] <>:';
end
y=input([str str2]);
y=chk(opt,y);
function y=chk(opt,x)
%chk 检查原始数据
str='\n数据错误，请检查后重新......\n';
switch opt
case 'Qt'
if iscell(x) & ~isempty(x)
x=x(:);
n=zeros(size(x));
for i=1:length(x)
if size(x{i},2)~=2, x{i}=x{i}'; end
if size(x{i},2)==2
if isequal(sort(x{i}(:,1)),x{i}(:,1)) & all(x{i}(:,2)>=0), n(i)=1; end
end
end
if all(n), str=''; end
end
%funz 迭代
dq=mean(Q)-funq(fq,z1);
dz=[0 0];
if dq>1e-4, dz(2)=0.1; elseif dq<-1e-4, dz(2)=-0.1; end
while abs(diff(dz))>0.0005
dz(1)=dz(2);
dz(2)=(2.*dq.*dz(1).*diff(t).*3600)./(2.*diff(interp1(Vz(:,1),Vz(:,2).*1e4,[z1 z1+dz(1)],'s'))+diff(funq(fq,[z1 z1+dz(1)])).*diff(t).*3600);
end
z2=z1+dz(2);
⛄ 运行结果
⛄ 参考文献
[1]邹进, 张勇传. 一种多目标决策问题的模糊解法及在洪水调度中的应用[J]. 水利学报, 2003.
