关于喷嘴计算的程序问题指教
自己编写的程序，运行后出现下列问题，请高手指教，程序见附件Optimization terminated: first-order optimality is less than options.TolFun.
f =
0 -1.4086
Warning: Default trust-region dogleg method of FSOLVE cannot
handle non-square systems; using Gauss-Newton method instead.
> In fsolve at 248
In swirlnozzle at 35
f =
0 -1.4086
f =
0 -0.3986
f =
0 -0.1218
f =
0 0.6373
f =
0 -0.1440
f =
0 -0.1218
f =
0 -0.1218
f =
0 0.1749
f =
0 0.0159
f =
1.0e-003 *
0 -0.8681
f =
1.0e-003 *
0 -0.8681
f =
1.0e-003 *
0 -0.8680
f =
1.0e-004 *
0 0.5751
f =
1.0e-003 *
0 0.4645
f =
1.0e-009 *
0 -0.1381
f =
1.0e-009 *
0 -0.1381
f =
1.0e-007 *
0 0.5120
Optimization terminated: directional derivative along
search direction less than TolFun and infinity-norm of
gradient less than 10*(TolFun+TolX).
f =
0 -0.1066
Warning: Default trust-region dogleg method of FSOLVE cannot
handle non-square systems; using Gauss-Newton method instead.
> In fsolve at 248
In swirlnozzle at 45
f =
0 -0.1066
f =
0 0.0257
f =
0 0.0993
f =
1.0e-003 *
0 -0.3986
f =
1.0e-003 *
0 -0.3986
f =
1.0e-003 *
0 -0.3986
f =
1.0e-004 *
0 -0.6689
f =
1.0e-004 *
0 0.8004
f =
1.0e-011 *
0 -0.2870
f =
1.0e-011 *
0 -0.2870
f =
1.0e-007 *
0 0.1555
Optimization terminated: directional derivative along
search direction less than TolFun and infinity-norm of
gradient less than 10*(TolFun+TolX).
单位：mm
排放口直径=6.2,修正后的入口直径=6.2旋流半径=4.97,喷射角=48.5 0.5616
m选的合适，不需要再调整，m=
0.8000
参数合适，不需要调整
旋流室直径=16.2>>
[ 本帖最后由 eight 于 2008-3-10 16:39 编辑 ]
下次要贴程序,请直接贴出来,不要用附件的形式
function swirlnozzle% xFsolve.m
% 喷嘴计算，采用《简式旋流喷嘴的设计》文章中的方法，作者“姚增权．马进”
clear all
clc
% -------------------------------------------------------------------------
% 初始参数设定
Kv = 0.8667; % 流量系数Kv
Q = Kv/3600; % 流量，单位m3/s
i = 3; % 旋流入口数目
alpha = 60; % 雾化角，alpha单位为度
m = 0.8; % 旋流腔半径与排放口直径之比,m取0.5～3
a = 2*pi*alpha/360; % a将alpha单位转为弧度
Pl = 1e5; % 压力为100kPa
density = 1000; % 水的密度
miu = 1.0034e-6; % 水的运动粘度
x0 = 0.5;
% -------------------------------------------------------------------------
theta = fsolve(@NonlinEqs,x0,[],a); % 第一步计算填充系数theta
K = sqrt(2)*(1-theta)/(theta*sqrt(theta)); % 计算几何常数K
Cd = theta*sqrt(theta/(2-theta)); % 计算流量系数Cd
% -------------------------------------------------------------------------
Dout = sqrt(4*Q/(pi*Cd*sqrt(2*Pl/density))); % 第二步计算排放口直径Dout
% -------------------------------------------------------------------------
R = m*Dout; % 第三步初步选定旋流半径R
Dp = sqrt(2*R*Dout/(i*K)); % 计算喷嘴入口直径Dp
% -------------------------------------------------------------------------
Re = 4*Q/(pi*miu*sqrt(i)*Dp); % 第四步计算雷诺数Re
% -------------------------------------------------------------------------
s = 10^(25.8/(log10(Re))^2.58-2); % 第五步计算摩擦系数
Ky = R*(Dout/2)/(i*(Dp/2)^2+(s/2)*R*(R-Dout/2)); % 计算修正后的几何常数
% -------------------------------------------------------------------------
% 迭代计算填充系数theta1、流量系数Cd1、排放口直径Dout1、旋流半径R1、喷嘴入口直径Dp1
theta1 = fsolve(@NonlinEqs1,x0,[],Ky);
Cd1 = theta1^1.5/sqrt(2-theta1);
Dout1 = sqrt(4*Q/(pi*Cd1*sqrt(2*Pl/density)));
R1 = m*Dout1;
K1 = sqrt(2)*(1-theta1)/(theta1^1.5);
Dp1 = sqrt(2*R1*Dout1/(i*K1));
% -------------------------------------------------------------------------
% 计算修正后的入口直径Dp2
Dp2 = Dp1/sqrt(0.9); % 修正后的入口直径Dp2
% -------------------------------------------------------------------------
S = fsolve(@NonlinEqs2,x0,[],Cd1,K1); % 计算空气芯半径
spray = 2*atan((2*Cd1*K1)/sqrt((1+S)^2-4*(Cd1*K1)^2));% 计算喷射角
% -------------------------------------------------------------------------
disp(' 单位：mm')
fprintf('排放口直径=%.1f,修正后的入口直径=%.1f',Dout1*1000,Dp2*1000)
fprintf('旋流半径=%.2f,喷射角=%.1f',R1*1000,spray/(2*pi)*360)
disp(Cd1)
% -------------------------------------------------------------------------
% 最后参数校验
% 旋流腔半径与排放口直径之比m的校验
if(R1>(Dp2/2))
disp('m选的合适，不需要再调整，m=');
disp(m);
else disp('m选的不合适，需要再调整')
end
% 旋流室直径Dk的校验
Dk = 2*R1+Dp2;
if(or(Dk<(2.5*Dout1),Dk>(3.3*Dout1)))
disp('参数需要调整')
else disp('参数合适，不需要调整')
fprintf('旋流室直径=%.1f',Dk*1000)
end
% -------------------------------------------------------------------------
% 程序结束
% -------------------------------------------------------------------------
% ------------------------------------------------------------------
function f = NonlinEqs(x,a)
f(1) = tan(a/2)-(2^1.5*(1-x)/(sqrt(x)*(1+sqrt(1-x))));
% -------------------------------------------------------------------
function f = NonlinEqs1(x,Ky)
f(2) = Ky-sqrt(2)*(1-x)/(x^1.5)
% -------------------------------------------------------------------
function f = NonlinEqs2(x,Cd1,K1)
f(2) = Cd1-(sqrt(1-(Cd1*K1)^2)-x*sqrt(x^2-(Cd1*K1)^2)-(Cd1*K1)^2*log((1+sqrt(1-(Cd1*K1)^2))/(x+sqrt(x^2-(Cd1*K1)^2))))
还有请具体说一下是出现什么问题,我看只有警告.
[ 本帖最后由 sigma665 于 2008-3-9 08:56 编辑 ] 这个警告会不会影响结果
回复 3楼 的帖子
Warning: Default trust-region dogleg method of FSOLVE cannothandle non-square systems; using Gauss-Newton method instead.
用了高斯牛顿法,方法不一样,可能精度会不一样
页:
[1]