figure(1);hold on; xlabel('V');ylabel('P');
rho(i)=1.225*exp(-h/10.4);
cl(i) = clalpha * alpha * pi/180;
V(i) = sqrt(2*W/rho(i)/S/cl(i));
cd(i) = cd0 + k * cl(i) * cl(i);
D(i) = 0.5 * rho(i) * V(i) * V(i) * S * cd(i);
Ph(i)=Psl*(rho(i)/1.225).^0.75;
xix = find(diff(sign(D.*V-Ph)));
idxrng = max(xix(k3)-1,1) : min(numel(V),xix(k3)+1);
xv(k3) = interp1(P(idxrng)-Ph(idxrng), V(idxrng), 0);
yv(k3) = interp1(V(idxrng),P(idxrng), xv(k3));
plot(xv, yv, 's', 'DisplayName',sprintf('Intersections %d',k1))
plot([0 V(ix)], [0 P(ix)], '-', 'LineWidth',1, 'Color',c{k1})
plot(V(ix), P(ix), 'd', 'Color',c{k1})
Table = cell2table(cat(2,Xm,Ym), 'VariableNames',{'V_p_max','V_p_min','Pmax','Pmin'});
Table = addvars(Table, Vstall, 'After','V_p_min');
Table = addvars(Table, Pstall, minP,Vminp, 'After','Pmin')
Table =
V_p_max V_p_min Vstall Pmax Pmin Pstall minP Vminp
_______ ____________ ______ _____ ______________ ______ ______ ______
56.013 {0×0 double} 15.752 25000 {0×0 double } 12261 9283.1 26.626
55.822 {0×0 double} 17.051 22198 {0×0 double } 13273 10049 28.822
55.235 {0×0 double} 18.458 19710 {0×0 double } 14367 10878 31.199
53.92 {0×0 double} 19.98 17501 {0×0 double } 15552 11775 33.772
51.119 {[ 23.8756]} 21.628 15539 {[1.5539e+04]} 16835 12746 36.558
39.638 {[ 39.6377]} 23.412 13798 {[1.3798e+04]} 18223 13797 39.573
Vpmin=[Vstall(1:4);cell2mat(Xm(5,2));cell2mat(Xm(6,2))];
lgdstr = [{'H=SI'} compose('H=%.4f',hv(2:end))];
legend([hp{:}],lgdstr, 'Location','best')