使用 scilab 查找交点

finding intersection point using scilab

如何使用 fsolve 函数 (from scilab) 在下图中找到交点?

这是我到目前为止尝试过的方法:

function y=f(x)
    y = 30 + 0 * x;
endfunction


function y= g(x)
    y=zeros(x)
    k1 = find(x >= 5 & x <= 11); 
    if  k1<>[]  then
        y(k1)= -59.535905 +24.763399*x(k1) -3.135727*x(k1)^2+0.1288967*x(k1)^3;
    end;
    k2=find(x >= 11 & x <= 12); 
    if  k2 <> []    then 
        y(k2)=1023.4465 - 270.59543 * x(k2) + 23.715076 * x(k2)^2 - 0.684764 * x(k2)^3; 
    end;
    k3 = find(x >= 12 & x <= 17);    
    if  k3 <> [] then
        y(k3) =-307.31448 + 62.094807 *x(k3) - 4.0091108 * x(k3)^2 + 0.0853523 * x(k3)^3;
    end;
    k4 = find(x >= 17 & x <= 50); 
    if k4 <> [] then 
        y(k4) = 161.42601 - 20.624104 *x(k4) + 0.8567075 * x(k4)^2 - 0.0100559 * x(k4)^3;
    end;
endfunction

t=[5:50];
plot(t, g(t));
plot2d(t, f(t));
deff('res = fct', ['res(1) = f(x)'; 'res(2) = g(x)']);
k1=[5, 45];
xsol1 = fsolve(k1, f, g)

您的原文 post 完全不可读且混乱。我花了一段时间来编辑它并理解你想要达到的目标。不过我会尽力帮助你。让我们一步步来:

  1. 我不确定您为什么要这样使用 find 函数。可能您正在尝试矢量化 g 函数?请考虑 Scilab 默认情况下不会通过数组广播函数。您需要对它们进行矢量化或使用 feval 来执行此操作。请阅读 I have written before. find is a vectorized operation applying on an array, a Boolean operation and a scalar, finding the elements of the array which satisfy the operation. For example from the find page:
beers = ["Desperados", "Leffe", "Kronenbourg", "Heineken"];
find(beers == "Leffe")

returns 2

A = rand(1, 20);
w = find(A < 0.4)

returns 数组 A 中小于 0.4.

的那些元素
  1. 请了解条件,特别是 if, then, elsif, else, end 语句。如果您了解了这一点,您将不会以那种方式使用 find 函数。有时你连续有很多 if,然后尝试使用 select, case, else, end。你的第二个函数可以写成:
function y = g(x)
  if x < 5 | 50 < x then
    error("Out of range");
  elseif x <= 11 then
    y = -59.535905 + 24.763399 * x - 3.135727 * x^2 + 0.1288967 * x^3;
    return;
  elseif x <= 12 then
    y = 1023.4465 - 270.59543 * x + 23.715076 * x^2 - 0.684764 * x^3;
    return;
  elseif x <= 17 then
    y = -307.31448 + 62.094807 * x - 4.0091108 * x^2 + 0.0853523 * x^3;
    return;
  else
    y = 161.42601 - 20.624104 * x + 0.8567075 * x^2 - 0.0100559 * x^3;
  end
endfunction
  1. 现在显然您想找到此曲线上值为 30 的点。虽然有一些方法可以自动找到这些点,但绘图对于找到合适的范围非常有帮助:
t = [5:50];
plot(t, feval(t, g) - 30)

表明这两个解在20 < x1 < 3040 < x < 50的范围内。

  1. 现在,如果我们使用 fsolve 和适当的初始值,它会给我们带来好的结果:
--> deff('[y] = g2(x)', 'y = g(x) - 30');

--> fsolve([25; 45], g2)
 ans  =

   26.67373
   48.396547
  1. fsolve函数的第三个参数是g(x)函数的雅各宾/导数。您要么应该手动计算上述多项式的导数(或使用适当的符号软件,如 Maxima),要么使用 poly 函数将它们定义为多项式。例如,参见 this tutorial。然后区分它们,定义一个新函数,如dgdx