MATLAB error using fzero function

My following code generates the plot of V and D values in figure 1. In the graph, the parabolas and straight lines intersect, and I need to find the roots from the plot. So I tried to use fzero function, but the error appeared: Operands to the logical AND (&&) and OR (||) operators must be convertible to logical scalar values. Use the ANY or ALL functions to reduce operands to logical scalar values.

Error in fzero (line 326) elseif ~isfinite(fx) || ~isreal(fx)
Error in HW1 (line 35) x=fzero(fun,1);

My code is:

clear all; close all
W = 10000; %[N]
S = 40; %[m^2]
AR = 7;
cd0 = 0.01;
k = 1 / pi / AR;
clalpha = 2*pi;
Tsl=800; 
figure(1);hold on; xlabel('V');ylabel('D')
for h=0:1:8; 
    i=0;
    for alpha = 1:0.25:12
        i=i+1;
        rho(i)=1.2*exp(-h/10.4);
        cl(i) = clalpha * alpha * pi/180;
        V(i) = sqrt(2*W/rho(i)/S/cl(i));
        L(i) = 0.5 * rho(i) * V(i) * V(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);
        clcd(i) = cl(i)/cd(i);
        p(i) = D(i)*V(i);
        ang(i) = alpha;
        T(i)=Tsl*(rho(i)/1.2).^0.75;
        
        
    end
    figure(1); plot(V,D); hold on 
    plot(V,T);
end

fun = @(V) 0.5*V.*V.*rho.*S.*cd-T;
x=fzero(fun,1);

Operands to the logical AND (&&) and OR (||) operators must be convertible to logical scalar values. Use the ANY or ALL functions to reduce operands to logical scalar values.

Error in fzero (line 326)
elseif ~isfinite(fx) || ~isreal(fx)

Probably, I should not use the fzero function, but the task is to find the roots of V from a plot (figure 1). There are parabolas and straight lines respectively.

Answer 1

[Star Strider]

I have no idea what intersections you want to find, so this is a guess.

It should at least get you started —

% clear all; close all
W = 10000; %[N]
S = 40; %[m^2]
AR = 7;
cd0 = 0.01;
k = 1 / pi / AR;
clalpha = 2*pi;
Tsl=800; 
figure(1);hold on; xlabel('V');ylabel('D')
hv=0:1:8;
alphav = 1:0.25:12;
for k1 = 1:numel(hv)
    h = hv(k1);

    i=0;
    for alpha = 1:0.25:12
        i=i+1;
        rho(i)=1.2*exp(-h/10.4);
        cl(i) = clalpha * alpha * pi/180;
        V(i) = sqrt(2*W/rho(i)/S/cl(i));
        L(i) = 0.5 * rho(i) * V(i) * V(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);
        clcd(i) = cl(i)/cd(i);
        p(i) = D(i)*V(i);
        ang(i) = alpha;
        T(i)=Tsl*(rho(i)/1.2).^0.75;
        
        
    end
    
    Vm(k1,:) = V;                                                           % Save Results
    Dm(k1,:) = D;                                                           % Save Results
    Tm(k1,:) = T;                                                           % Save Results
    
    xix = find(diff(sign(0.5*V.*V.*rho.*S.*cd-T)));                         % Approximate Indices Of Intersections
    if ~isempty(xix)                                                        % Skip If Empty
    for k3 = 1:numel(xix)
        idxrng = max(xix(k3)-1,1) : min(numel(V),xix(k3)+1);
        xv(k3) = interp1(D(idxrng)-T(idxrng), V(idxrng), 0);
        yv(k3) = interp1(V(idxrng),D(idxrng), xv(k3));
        Xm{k1,k3} = xv(k3);                                                 % Intersection X-Matrix
        Ym{k1,k3} = yv(k3);                                                 % Intersection Y-Matrix
    end
    end
    
    figure(1)
    plot(V,D, 'DisplayName',sprintf('D_{%d}',k1))
    hold on 
    plot(V,T, 'DisplayName',sprintf('T_{%d}',k1))
    plot(xv, yv, 's', 'DisplayName',sprintf('Intersections %d',k1))         % Plot Intersections (Optional)
end
legend('Location','eastoutside','NumColumns',2)

Intersections = cell2table(cat(2,Xm,Ym), 'VariableNames',{'X1','X2','Y1','Y2'})

Intersections = 9×4 table
      X1           X2           Y1           Y2     
    ______    ____________    ______    ____________

    55.445    {0×0 double}       800    {0×0 double}
    55.654    {0×0 double}    744.34    {0×0 double}
    55.885    {0×0 double}    692.55    {0×0 double}
    55.914    {[ 21.1789]}    644.37    {[644.3652]}
    55.866    {[ 23.3517]}    599.53    {[599.5325]}
    55.586    {[ 25.8488]}    557.82    {[557.8191]}
    54.967    {[ 28.7803]}    519.01    {[519.0081]}
    53.781    {[ 32.3659]}     482.9    {[482.8973]}
    51.479    {[ 37.2235]}     449.3    {[449.2990]}

% Vm
% Dm
% Tm

% fun = @(V) 0.5*V.*V.*rho.*S.*cd-T;
% x=fzero(fun,1);

See the documentation on interp1 to understand the intersections calculations.

.

[Alina Abdikadyr]

Thank you!

[Star Strider]

My pleasure!

If my Answer helped you solve your problem, please Accept it!

.

[Alina Abdikadyr]

Hello!

I modified the code and tried to use interp1 command. In the graph the intersections are displayed correctly, but in the table, only two intersections are displayed.

clear all; close all

W = 10000;

S = 40;

AR = 7;

cd0 = 0.01;

k = 1 / pi / AR;

clalpha = 2*pi;

Tsl=800;

figure(1);hold on; xlabel('V');ylabel('D')

for rho = 1.2:-0.15:0.6

k1=numel(rho);

h=rho(k1);

i=0;

for alpha = 1:0.25:15

i=i+1;

cl(i) = clalpha * alpha * pi/180;

V(i) = sqrt(2*W/rho/S/cl(i));

L(i) = 0.5 * rho * V(i) * V(i) * S * cl(i);

cd(i) = cd0 + k * cl(i) * cl(i);

D(i) = 0.5 * rho * V(i) * V(i) * S * cd(i);

clcd(i) = cl(i)/cd(i);

p(i) = D(i)*V(i);

ang(i) = alpha;

T(i)=Tsl*(rho/1.2).^0.75;

end

Vm(k1,:) = V; % Save Results

Dm(k1,:) = D; % Save Results

Tm(k1,:) = T; % Save Results

xix = find(diff(sign(0.5*V.*V.*rho.*S.*cd-T))); % Approximate Indices Of Intersections

if ~isempty(xix) % Skip If Empty

for k3 = 1:numel(xix)

idxrng = max(xix(k3)-1,1) : min(numel(V),xix(k3)+1);

xv(k3) = interp1(D(idxrng)-T(idxrng), V(idxrng), 0);

yv(k3) = interp1(V(idxrng),D(idxrng), xv(k3));

Xm{k1,k3} = xv(k3); % Intersection X-Matrix

Ym{k1,k3} = yv(k3); % Intersection Y-Matrix

end

end

figure(1)

plot(V,D, 'DisplayName',sprintf('D_{%d}',k1))

hold on

plot(V,T, 'DisplayName',sprintf('T_{%d}',k1))

xlim([0 100])

ylim([0 1000])

plot(xv, yv, 's', 'DisplayName',sprintf('Intersections %d',k1)) % Plot Intersections (Optional)

end

%legend('Location','eastoutside','NumColumns',2)

Intersections = cell2table(cat(2,Xm,Ym), 'VariableNames',{'Xmax','Xmin','Y1','Y2'})

[Star Strider]

It was modified a bit too much, and eliminated the ‘k1’ loop count variable. That was the essential problem.

Creating:

rhov = 1.2:-0.15:0.6;

and then using:

for k1 = 1:numel(rhov);

rho = rhov(k1);

<REST OF THE LOOP CODE>

end

restored the original code architecture.

It now works as the orignal version did —

W = 10000;
S = 40;
AR = 7;
cd0 = 0.01;
k = 1 / pi / AR;
clalpha = 2*pi;
Tsl=800; 
rhov = 1.2:-0.15:0.6;
figure(1);hold on; xlabel('V');ylabel('D')
for k1 = 1:numel(rhov);
    rho = rhov(k1);
    i=0;
    for alpha = 1:0.25:15
        i=i+1;
        cl(i) = clalpha * alpha * pi/180;
        V(i) = sqrt(2*W/rho/S/cl(i));
        L(i) = 0.5 * rho * V(i) * V(i) * S * cl(i);
        cd(i) = cd0 + k * cl(i) * cl(i);
        D(i) = 0.5 * rho * V(i) * V(i) * S * cd(i);
        clcd(i) = cl(i)/cd(i);
        p(i) = D(i)*V(i);
        ang(i) = alpha;
        T(i)=Tsl*(rho/1.2).^0.75;
    end
    Vm(k1,:) = V;                                                           % Save Results
    Dm(k1,:) = D;                                                           % Save Results
    Tm(k1,:) = T;                                                           % Save Results
    
    xix = find(diff(sign(0.5*V.*V.*rho.*S.*cd-T)));                         % Approximate Indices Of Intersections
    if ~isempty(xix)                                                        % Skip If Empty
    for k3 = 1:numel(xix)
        idxrng = max(xix(k3)-1,1) : min(numel(V),xix(k3)+1);
        xv(k3) = interp1(D(idxrng)-T(idxrng), V(idxrng), 0);
        yv(k3) = interp1(V(idxrng),D(idxrng), xv(k3));
        Xm{k1,k3} = xv(k3);                                                 % Intersection X-Matrix
        Ym{k1,k3} = yv(k3);                                                 % Intersection Y-Matrix
    end
    
    end
   figure(1)
    plot(V,D, 'DisplayName',sprintf('D_{%d}',k1))
    hold on 
    plot(V,T, 'DisplayName',sprintf('T_{%d}',k1))
    xlim([0 100])
    ylim([0 1000])
    plot(xv, yv, 's', 'DisplayName',sprintf('Intersections %d',k1))         % Plot Intersections (Optional)
end

% Q = cat(1,Xm,Ym)
legend('Location','eastoutside','NumColumns',1)

Intersections = cell2table(cat(2,Xm,Ym), 'VariableNames',{'Xmax','Xmin','Y1','Y2'})

Intersections = 5×4 table
     Xmax      Xmin       Y1        Y2  
    ______    ______    ______    ______

    55.445    16.019       800       800
    55.735    18.193    723.76    723.76
    55.914    21.162    644.74    644.74
    55.614    25.551    562.34    562.34
    53.429    33.237    475.68    475.68

The rest of the code is essentialy unchanged.

.

[Alina Abdikadyr]

Thank you very much!

Could you please help me how could I find a minimum point from each parabola?

[Star Strider]

As always, my pleasure!

The easiest approach would simply to use the min function in this instance.

Adding:

[Dmin,idx] = min(D);

minD(k1,:) = Dmin;

VminD(k1,:) = V(idx);

and the appropriate table entries using the addvars function returns their values and ‘V’ coordinate values in the ‘Intersections’ table —

W = 10000;
S = 40;
AR = 7;
cd0 = 0.01;
k = 1 / pi / AR;
clalpha = 2*pi;
Tsl=800; 
rhov = 1.2:-0.15:0.6;
figure(1);hold on; xlabel('V');ylabel('D')
for k1 = 1:numel(rhov);
    rho = rhov(k1);
    i=0;
    for alpha = 1:0.25:15
        i=i+1;
        cl(i) = clalpha * alpha * pi/180;
        V(i) = sqrt(2*W/rho/S/cl(i));
        L(i) = 0.5 * rho * V(i) * V(i) * S * cl(i);
        cd(i) = cd0 + k * cl(i) * cl(i);
        D(i) = 0.5 * rho * V(i) * V(i) * S * cd(i);
        clcd(i) = cl(i)/cd(i);
        p(i) = D(i)*V(i);
        ang(i) = alpha;
        T(i)=Tsl*(rho/1.2).^0.75;
    end
    Vm(k1,:) = V;                                                           % Save Results
    Dm(k1,:) = D;                                                           % Save Results
    Tm(k1,:) = T;                                                           % Save Results
    
    xix = find(diff(sign(0.5*V.*V.*rho.*S.*cd-T)));                         % Approximate Indices Of Intersections
    if ~isempty(xix)                                                        % Skip If Empty
    for k3 = 1:numel(xix)
        idxrng = max(xix(k3)-1,1) : min(numel(V),xix(k3)+1);
        xv(k3) = interp1(D(idxrng)-T(idxrng), V(idxrng), 0);
        yv(k3) = interp1(V(idxrng),D(idxrng), xv(k3));
        Xm{k1,k3} = xv(k3);                                                 % Intersection X-Matrix
        Ym{k1,k3} = yv(k3);                                                 % Intersection Y-Matrix
    end
    
    end
   figure(1)
    plot(V,D, 'DisplayName',sprintf('D_{%d}',k1))
    hold on 
    plot(V,T, 'DisplayName',sprintf('T_{%d}',k1))
    xlim([0 100])
    ylim([0 1000])
    plot(xv, yv, 's', 'DisplayName',sprintf('Intersections %d',k1))         % Plot Intersections (Optional)
    [Dmin,idx] = min(D);
    minD(k1,:) = Dmin;
    VminD(k1,:) = V(idx);
end

% Q = cat(1,Xm,Ym)
legend('Location','eastoutside','NumColumns',1)

Intersections = cell2table(cat(2,Xm,Ym), 'VariableNames',{'Xmax','Xmin','Y1','Y2'});
Intersections = addvars(Intersections, minD,VminD, 'After','Y2')

Intersections = 5×6 table
     Xmax      Xmin       Y1        Y2      minD     VminD 
    ______    ______    ______    ______    _____    ______

    55.445    16.019       800       800    426.5      29.9
    55.735    18.193    723.76    723.76    426.5    31.964
    55.914    21.162    644.74    644.74    426.5    34.526
    55.614    25.551    562.34    562.34    426.5    37.821
    53.429    33.237    475.68    475.68    426.5    42.285

All the minima appear to have the same magnitude, however their ‘V’ coordinates (‘VminD’) differ.

.

[Alina Abdikadyr]

Thank you very much!

I'm trying to modify the code.

Sorry for asking so many questions, but if in the intersections table, I get {0×0 double} values.

And I try to convert cell arrays using cell2mat, but the error appears.

Vmax Vmin Pmax Pmin

______ ____________ _____ ______________

56.013 {0×0 double} 25000 {0×0 double }

55.822 {0×0 double} 22198 {0×0 double }

55.235 {0×0 double} 19710 {0×0 double }

53.92 {0×0 double} 17501 {0×0 double }

51.119 {[ 23.8756]} 15539 {[1.5539e+04]}

39.638 {[ 39.6377]} 13798 {[1.3798e+04]}

My code is:

clear all; close all

W = 10000;
S = 40;
AR = 7;
cd0 = 0.005;
k = 1 / pi / AR;
clalpha = 2*pi;
Psl=25000; 
figure(1);hold on; xlabel('V');ylabel('P')
hv=0:1.6484:8.242;
for k1 = 1:numel(hv)
    h = hv(k1);
    i=0;
    for alpha = 1:0.25:15
        i=i+1;
        rho(i)=1.225*exp(-h/10.4);
        cl(i) = clalpha * alpha * pi/180;
        V(i) = sqrt(2*W/rho(i)/S/cl(i));
        L(i) = 0.5 * rho(i)* V(i) * V(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);
        clcd(i) = cl(i)/cd(i);
        P(i) = D(i)*V(i);
        ang(i) = alpha;
        Ph(i)=Psl*(rho(i)/1.225).^0.75;
    end
    Vm(k1,:) = V;                                                          
    Pm(k1,:) = P;                                                         
    Pmh(k1,:) = Ph;                                                          
    
    xix = find(diff(sign(D.*V-Ph)));                         
    if ~isempty(xix)                                                        
    for k3 = 1:numel(xix)
        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));
        Xm{k1,k3} = xv(k3);                                                 % Intersection X-Matrix
        Ym{k1,k3} = yv(k3);                                                 % Intersection Y-Matrix
    end
    end
   figure(1)
    plot(V,P)
    hold on 
    plot(V,Ph)
    xlim([0 90])
    ylim([0 45000])
    plot(xv, yv, 's', 'DisplayName',sprintf('Intersections %d',k1))         % Plot Intersections (Optional)
end
%legend('Location','eastoutside','NumColumns',2)
Intersections = cell2table(cat(2,Xm,Ym), 'VariableNames',{'Vmax','Vmin','Pmax','Pmin'})
Xm1=cell2mat(Xm);

[Star Strider]

You omitted the ‘Pmax’ and ‘Pmin’ calculations.

Adding them back —

W = 10000;
S = 40;
AR = 7;
cd0 = 0.005;
k = 1 / pi / AR;
clalpha = 2*pi;
Psl=25000; 
figure(1);hold on; xlabel('V');ylabel('P')
hv=0:1.6484:8.242;
for k1 = 1:numel(hv)
    h = hv(k1);
    i=0;
    for alpha = 1:0.25:15
        i=i+1;
        rho(i)=1.225*exp(-h/10.4);
        cl(i) = clalpha * alpha * pi/180;
        V(i) = sqrt(2*W/rho(i)/S/cl(i));
        L(i) = 0.5 * rho(i)* V(i) * V(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);
        clcd(i) = cl(i)/cd(i);
        P(i) = D(i)*V(i);
        ang(i) = alpha;
        Ph(i)=Psl*(rho(i)/1.225).^0.75;
    end
    Vm(k1,:) = V;                                                          
    Pm(k1,:) = P;                                                         
    Pmh(k1,:) = Ph;                                                          
    
    xix = find(diff(sign(D.*V-Ph)));                         
    if ~isempty(xix)                                                        
    for k3 = 1:numel(xix)
        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));
        Xm{k1,k3} = xv(k3);                                                 % Intersection X-Matrix
        Ym{k1,k3} = yv(k3);                                                 % Intersection Y-Matrix
    end
    end
   figure(1)
    plot(V,P)
    hold on 
    plot(V,Ph)
    xlim([0 90])
    ylim([0 45000])
    plot(xv, yv, 's', 'DisplayName',sprintf('Intersections %d',k1))         % Plot Intersections (Optional)
    [maxP,idx] = max(P);
    [minP,idx] = min(P);
    Pmax(k1,:) = maxP;
    Pmin(k1,:) = minP;
    [maxV,idx] = max(V);
    [minV,idx] = min(V);
    Vmax(k1,:) = maxV;
    Vmin(k1,:) = minV;
end

%legend('Location','eastoutside','NumColumns',2)
Intersections = cell2table(cat(2,Xm,Ym), 'VariableNames',{'Xmax','Xmin','Y1','Y2'});
Intersections = addvars(Intersections, Pmax,Pmin,Vmax,Vmin, 'After','Y2')

Intersections = 6×8 table
     Xmax         Xmin         Y1            Y2          Pmax      Pmin      Vmax      Vmin 
    ______    ____________    _____    ______________    _____    ______    ______    ______

    56.013    {0×0 double}    25000    {0×0 double  }    30859    9283.1    61.008    15.752
    55.822    {0×0 double}    22198    {0×0 double  }    33404     10049     66.04    17.051
    55.235    {0×0 double}    19710    {0×0 double  }    36159     10878    71.486    18.458
     53.92    {0×0 double}    17501    {0×0 double  }    39141     11775    77.382     19.98
    51.119    {[ 23.8756]}    15539    {[1.5539e+04]}    42369     12746    83.764    21.628
    39.638    {[ 39.6377]}    13798    {[1.3798e+04]}    45863     13797    90.673    23.412

% Xm1=cell2mat(Xm);

Please always include the intersection arrays ‘Xm’ and ‘Ym’ in the ‘Intersections’ table. (The ‘Xm’ and ‘Ym’ arrays are the intersection coordinates.) Some ‘P’ and ‘Ph’ vectors have one intersection and some have two intersections. Those that have only one intersection have empty cells for the second intersection. (This is to be expected.)

I assume ‘Pmax’, ‘Pmin’, ‘Vmax’ and ‘Vmin’ are calculated as you want them to be. I corrected the ‘Intersections’ table to include all of them. It is easiest to calculate them in the plotting loop, and then add them with the addvars call.

It would be straightforward to include the maxima and minima of the ‘Ph’ vectors as well, if you want to. Add their calculations in the existing loop, and then add their names to the addvars argument list.

.

[Alina Abdikadyr]

Sorry, but in the table, there are Vmax and Vmin, how could I find a corresponding values of P for calculated Vmax and Vmin

The table gives Pmax and Pmin from a plot, but how to get P values for Vmax and Vmin?

     Xmax         Xmin         Y1            Y2          Pmax      Pmin      Vmax      Vmin 
    ______    ____________    _____    ______________    _____    ______    ______    ______

    56.013    {0×0 double}    25000    {0×0 double  }    30859    9283.1    61.008    15.752
    55.822    {0×0 double}    22198    {0×0 double  }    33404     10049     66.04    17.051
    55.235    {0×0 double}    19710    {0×0 double  }    36159     10878    71.486    18.458
     53.92    {0×0 double}    17501    {0×0 double  }    39141     11775    77.382     19.98
    51.119    {[ 23.8756]}    15539    {[1.5539e+04]}    42369     12746    83.764    21.628
    39.638    {[ 39.6377]}    13798    {[1.3798e+04]}    45863     13797    90.673    23.412

[Star Strider]

Use the returned indices to refer to them.

I calculated them as:

PmaxV(k1,:) = P(Vxidx);

PminV(k1,:) = P(Vnidx);

where ‘Vxidx’ and ‘Vnidx’ are returned from:

[maxV,Vxidx] = max(V);

[minV,Vnidx] = min(V);

And added them to the ‘Intersections’ table as well —

W = 10000;
S = 40;
AR = 7;
cd0 = 0.005;
k = 1 / pi / AR;
clalpha = 2*pi;
Psl=25000; 
figure(1);hold on; xlabel('V');ylabel('P')
hv=0:1.6484:8.242;
for k1 = 1:numel(hv)
    h = hv(k1);
    i=0;
    for alpha = 1:0.25:15
        i=i+1;
        rho(i)=1.225*exp(-h/10.4);
        cl(i) = clalpha * alpha * pi/180;
        V(i) = sqrt(2*W/rho(i)/S/cl(i));
        L(i) = 0.5 * rho(i)* V(i) * V(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);
        clcd(i) = cl(i)/cd(i);
        P(i) = D(i)*V(i);
        ang(i) = alpha;
        Ph(i)=Psl*(rho(i)/1.225).^0.75;
    end
    Vm(k1,:) = V;                                                          
    Pm(k1,:) = P;                                                         
    Pmh(k1,:) = Ph;                                                          
    
    xix = find(diff(sign(D.*V-Ph)));                         
    if ~isempty(xix)                                                        
    for k3 = 1:numel(xix)
        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));
        Xm{k1,k3} = xv(k3);                                                 % Intersection X-Matrix
        Ym{k1,k3} = yv(k3);                                                 % Intersection Y-Matrix
    end
    end
   figure(1)
    plot(V,P)
    hold on 
    plot(V,Ph)
    xlim([0 90])
    ylim([0 45000])
    plot(xv, yv, 's', 'DisplayName',sprintf('Intersections %d',k1))         % Plot Intersections (Optional)
    [maxP,idx] = max(P);
    [minP,idx] = min(P);
    Pmax(k1,:) = maxP;
    Pmin(k1,:) = minP;
    [maxV,Vxidx] = max(V);
    [minV,Vnidx] = min(V);
    Vmax(k1,:) = maxV;
    Vmin(k1,:) = minV;
    PmaxV(k1,:) = P(Vxidx);
    PminV(k1,:) = P(Vnidx);
end

%legend('Location','eastoutside','NumColumns',2)
Intersections = cell2table(cat(2,Xm,Ym), 'VariableNames',{'Xmax','Xmin','Y1','Y2'});
Intersections = addvars(Intersections, Pmax,Pmin,Vmax,Vmin,PmaxV,PminV, 'After','Y2')

Intersections = 6×10 table
     Xmax         Xmin         Y1            Y2          Pmax      Pmin      Vmax      Vmin     PmaxV    PminV
    ______    ____________    _____    ______________    _____    ______    ______    ______    _____    _____

    56.013    {0×0 double}    25000    {0×0 double  }    30859    9283.1    61.008    15.752    30859    12261
    55.822    {0×0 double}    22198    {0×0 double  }    33404     10049     66.04    17.051    33404    13273
    55.235    {0×0 double}    19710    {0×0 double  }    36159     10878    71.486    18.458    36159    14367
     53.92    {0×0 double}    17501    {0×0 double  }    39141     11775    77.382     19.98    39141    15552
    51.119    {[ 23.8756]}    15539    {[1.5539e+04]}    42369     12746    83.764    21.628    42369    16835
    39.638    {[ 39.6377]}    13798    {[1.3798e+04]}    45863     13797    90.673    23.412    45863    18223

% Xm1=cell2mat(Xm);

.

Answer 2

[John D'Errico]

You cannot use fzero on a vector. A vector of elements is NOT a function. It is just a list of numbers. And while you may think of it as a function, it is not. Essentially, a vector of numbers is just a picture of a function.

Suppose you were not feeling well today, so you decided to visit your doctor. But rather than appearing for an examination, you decided to just send your doctor a picture of yourself. Would that be sufficient to solve your problem? Rather than bring your car in for an oil change, would you just send a picture of your car to the JiffyLube place? Fzero needs a function, not just a list of numbers.

Instead, you will need to use tools to identify where the curves cross, perhaps using methods of interpolation. Or, better yet, you could use the actual functions themselves, now as true functions that fzero can evaluate at any point.