妳的公式中 cp的取值範圍是cv到cs之間
用fzero在這個區間之內找結果就好了
u=2.65e+10;c1b=5.057e+10;
cv=2596.15;
cs=3128.89;
us01=10;
us12=20;
h=1e-9;
b=(0:20)';
A=cell(21,2);
f=@(b,cp)?tan(b*sqrt((cp/cv)^2-1))-(c1b*sqrt((cp/cv)^2-1)*(u*sqrt(1-(cp/cs)^2)+(b/h)*us12-(b/h)*us01))/((u*(b/h)*sqrt(1-(cp/cs)^2)+((b/h)^2)*us01)*us01+c1b^2*((cp/cv)^2-1));
for?n=1:length(b)
A{n,1}=b(n);
for?p0=cv+0.01:cs-0.01
[s?v?flag]=fzero(@(cp)?f(b(n),cp),p0);
if?flag==1
if?isempty(A{n})
A{n,2}=s;
elseif?all(abs(A{n,2}-s)>1e-3)
A{n,2}=[A{n,2},s];
end
end?
end
end
由於每個b的取值得到cp的解的個數不同，所以這裏用cell矩陣來儲存結果
最後結果再A中
A的第壹列是b的取值
第二列是cp的可能取值，如果結果是[]，也就是無解
結果如下，第壹列是b的取值，第二列是對應的cp的可能取值