Finding intersection of 2 parabolic curve fits with different y-axes

I need assistance trying to find the intersection point of the two curves in lines 124-151 (in bold). This particular figure has two different scaled y-axes and I need a way to find where the two curves meet.

clear all
clc

Dp=.5; %input('Diameter of pipe: ');
Di=.5; %input('impeller diameter');
Wi=3000*2*pi/(60); %input('impeller rotational speed');


Temp= 10; %deg C
gamma_w = 62.4; %specific weight of the water from the aquifer
h=100; %height m
l=5000; %length m
h2=40; %second height m
Vol=5000; %cubic  volume m^3
p=999.7;
mu=(1.31*10^-3);
L=5000;
Zb=140;

Kl90elb=0.2;
Kl45bnd=0.2;
Klscv=2;
KlRee=.8;
Klext=1;
Kl=2*Kl90elb+15*Kl45bnd+4*Klscv+KlRee+Klext;

D_1 = 7.75/12; %impeller diameter in feet
omega_1 = 1750*2*pi/60; %rotational speed of pump one in rads/s
Q_1 = 0.002228.*[0;100;200;300;400;500;600;630]; %volume flow rate for pump one in ft^3/s
H_1 = [55;54;50.5;46;39;30;17;13]; %pump one head in feet
eta_1 = [NaN;NaN;.415;.48;.505;.48;.39;NaN]; %efficiency of pump one as a percentage
Pump_1 = [Q_1 H_1 eta_1]; %table for pump one
Power_1 = gamma_w.*Q_1.*H_1; %Pump power for pump one

D_2 = 7.75/12; %impeller diameter in feet
omega_2 = 2850*2*pi/60; %rotational speed of pump one in rads/s
Q_2 = 0.002228.*[0;100;200;300;400;500;600;680]; %volume flow rate for pump two in ft^3/s
H_2 = [151;149.5;155.5;140;130;119;102;87]; %pump two head in feet
eta_2 = [NaN;NaN;NaN;.42;.485;.525;.525;.48]; %efficiency of pump two as a percentage
Pump_2 = [Q_2 H_2 eta_2]; %table for pump two
Power_2 = gamma_w.*Q_2.*H_2; %Pump power for pump two

D_3 = 7.75/12; %impeller diameter in feet
omega_3 = 3500*2*pi/60; %rotational speed of pump one in rads/s
Q_3 = 0.002228.*[0;100;200;300;400;500;600;700]; %volume flow rate for pump three in ft^3/s
H_3 = [225;222;218;212;205;195;179;142]; %pump three head in feet
eta_3 = [NaN;NaN;NaN;.40;.455;.51;.53;.49]; %efficiency of pump three as a percentage
Pump_3 = [Q_3 H_3 eta_3]; %table for pump three
Power_3 = gamma_w.*Q_3.*H_3; %Pump power for pump three

D_4 = 6/12; %impeller diameter in feet
omega_4 = 3500*2*pi/60; %rotational speed of pump one in rads/s
Q_4 = 0.002228.* [0;100;200;300;400;500;580]; %volume flow rate for pump four in ft^3/s
H_4 = [113;110.5;105;95;81;68;39]; %pump four head in feet
eta_4 = [NaN;NaN;NaN;.47;.49;.40;.34]; %efficiency of pump four as a percentage
Pump_4 = [Q_4 H_4 eta_4]; %table for pump four
Power_4 = gamma_w.*Q_4.*H_4; %Pump power for pump four

D_5 = 7/12; %impeller diameter in feet
omega_5 = 3500*2*pi/60; %rotational speed of pump one in rads/s
Q_5 = 0.002228.*[0;100;200;300;400;500;600;640]; %volume flow rate for pump five in ft^3/s
H_5 = [170;165;160;151.5;141;128;111;102]; %pump five head in feet
eta_5 = [NaN;NaN;NaN;.405;.47;.525;.525;.52]; %efficiency of pump five as a percentage
Pump_5 = [Q_5 H_5 eta_5]; %table for pump five
Power_5 = gamma_w.*Q_5.*H_5; %Pump power for pump five

C_Q1 = Q_1.*(1/(omega_1*D_1^3)); %Flow coefficient for Pump 1
C_Q2 = Q_2.*(1/(omega_2*D_2^3)); %Flow coefficient for Pump 2
C_Q3 = Q_3.*(1/(omega_3*D_3^3)); %Flow coefficient for Pump 3
C_Q4 = Q_4.*(1/(omega_4*D_4^3)); %Flow coefficient for Pump 4
C_Q5 = Q_5.*(1/(omega_5*D_5^3)); %Flow coefficient for Pump 5

C_H1 = (32.2.*H_1)/(omega_1^2*D_1^2); %Head coefficient for Pump 1
C_H2 = (32.2.*H_2)/(omega_2^2*D_2^2); %Head coefficient for Pump 2
C_H3 = (32.2.*H_3)/(omega_3^2*D_3^2); %Head coefficient for Pump 3
C_H4 = (32.2.*H_4)/(omega_4^2*D_4^2); %Head coefficient for Pump 4
C_H5 = (32.2.*H_5)/(omega_5^2*D_5^2); %Head coefficient for Pump 5

C_P1 = Power_1/(1.94*omega_1^3*D_1^5); %Power coefficient for Pump 1
C_P2 = Power_2/(1.94*omega_2^3*D_2^5); %Power coefficient for Pump 2
C_P3 = Power_3/(1.94*omega_3^3*D_3^5); %Power coefficient for Pump 3
C_P4 = Power_4/(1.94*omega_4^3*D_4^5); %Power coefficient for Pump 4
C_P5 = Power_5/(1.94*omega_5^3*D_5^5); %Power coefficient for Pump 5

Q = [C_Q1' C_Q2' C_Q3' C_Q4' C_Q5']; %All flow data in a cell array
H = [C_H1' C_H2' C_H3' C_H4' C_H5']; %All head data in a cell array
P = [C_P1' C_P2' C_P3' C_P4' C_P5']; %All power data in a cell array
eta = [eta_1' eta_2' eta_3' eta_4' eta_5']; %All efficiency data in a cell array

%Plotting all datapoints in single graph
x1 = Q;
y1 = H;
y2 = eta;
ind = ~isnan(y2);
y2 = y2(ind);
x2 = x1;
x2 = x2(ind);

%now we do the curve fits:
%QHE
%we found them in matlab using the polyfit and printing the coefficients to
%find the curve fit lines, then plotted them with 100 intervals between 0
%and 0.3
QH = polyfit(x1,H,2);
format compact;
QE = polyfit(x2,y2,4);
format compact;
figure('Name','Fluid Flow of My Tears')
plot(x1,y1,'^');
hold on
yyaxis right
plot(x2,y2,'x');
hold on
QHcap= zeros(1,100); %empty vec
Qspace= 0:0.0003:0.03; %100 points for line
QHline= polyval(QH,Qspace);
yyaxis left;
plot1 = plot(Qspace,QHline,'b');
hold on
QEline= polyval(QE,Qspace);
yyaxis right;
plot2 = plot(Qspace,QEline,'r');

%now we do fig 2 with Cp instead of eta

y3=P;
ind = ~isnan(y3);
y3 = y3(ind);


%now we do the curve fits:
%QHP
%we found them in matlab using the polyfit and polyval the coefficients to
%find the curve fit lines, then plotted them with 100 intervals between 0
%and 0.3
QH2 = polyfit(x1,H,2);
QHline= polyval(QH2,Qspace);
QP = polyfit(x1,y3,2);
format compact;
figure(2)
plot(x1,y1,'^');
hold on
yyaxis right
plot(x1,y3,'o');
Qspace= 0:0.0003:0.03; %100 points for line
yyaxis left;
plot3 = plot(Qspace,QHline,'b');
hold on
QPline= polyval(QP,Qspace);
yyaxis right;
plot4 = plot(Qspace,QPline,'r');

[darova]

Use polyxpoly()

[ShruDuebgen59]

i've tried polyxpoly() and it wouldn't work for my code

[ShruDuebgen59]

if you could figure out where how it could work, that would be very much appreciated. Maybe I might've just typed something in wrong.

Answer 1

[John D'Errico]

Um, trivial, sort of? You want to find the intersection of these two polynomial curves as plotted?

They are quadratics, for gods sake, with coefficients in QH2 and QP.

However, those two curves DON'T intersect where you may think they do!!!!!!!!!!!

The blue curve has its axis on the left. The red curve has its y axis on the right. While the curves do intersect, and roots can find those locations, as:

roots(QP - QH2)
ans =
       -0.0415352767284989
         0.033492450277296

Those roots are NOT the x-locations of where those curves APPEAR to cross. Those roots are the locations of where the curves ACTUALLY DO CROSS.

You can't find the intersection of apples and oranges. It simply does not make sense.

Did you ever read the book, "How to Lie With Statistics"? Perhaps you should do so. You are well on your way to write another chapter, if you are not careful. Unless you can explain why it is that you think the intersection of two curves with totally different axes is a meaningful thing to do, then you are doing something highly non-kosher with the mathematics by trying to find that intersection.

Ok, so why do I suggest that what you are TRYING to do is completely meaningless in terms of mathematics? When you plotted the blue curve and data, you let MATLAB auto-scale the axes. It did the same thing for the red curves and data. MATLAB could as easily have chosen another, subtly different scaling for either of the y-axes. Now had MATLAB done exactly that, that is, chosen a different auto-=scaling for one of the axes, the two curves would have now appeared to have intersected in different points!

Do you understand what I am saying? Merely by changing the axis scaling for one or both of the curves, you could in theory get ANY answer you want to get, and make it seem as if that answer makes sense! And this is why I am calling this a lie with statistics. It is something you should NEVER do. Unless, of course, you want to get a government job, at least with SOME governments.

Edit:

Perhaps you need proof.

x = rand(1,100);
y1 = x - 2 + randn(size(x))/10;
y2 = -x/100 + 1 + randn(size(x))/1000;
p1 = polyfit(x,y1,1)
p1 =
          0.99502346004009         -2.01239125455423
      
p2 = polyfit(x,y2,1)
p2 =
      -0.00977846480086918         0.999853116084099

So y1 and y2 are scaled very differently. In fact, but for the noise in the data, the original functions would have been y1 = x - 2, and y2 = -x/100 + 1. So the intersection between the two curves should have been at y1 == y2. Solving for x, we would get

x - 2 = -x/100 + 1

101/100*x = 3

Which implies an intersection at least near 2.97.

roots(p1 - p2)
ans =
          2.99784892541394

In fact, roots seems to agree reasonably well with that prediction. If we plot the two sets of data together, we would see:

plot(x,y1,'ro',x,y2,'bs')

And plot seems to agree with that conclusion. Any intersection of the two curves is way off to the right. Around x=3 would seem just about correct.

But, can we lie with the plot here?

[Ax,H1,H2] = plotyy(x,y1,x,y2);
H1.LineStyle = 'none';
H1.Marker = 'o';
H2.LineStyle = 'none';
H2.Marker = 's';

And of course, now we could comfortably draw the conclusion that the two curves intersect just around x==0.5. Surely that cannot be correct, can it? After all, I just let MATLAB auto-scale the axes. MATLAB would not lie, would it? (Adding in the linear fits will not change anything, as they will pass right through tha data on the respective axis.)

We can change that simply enough. Where do you want the intersection to be? How about near x==0? Easy, peasy. We can change the axis scaling for the second set of axes simply.

Ax(2).YLim = [0.985 1.05];

Very good. I could now comfortably conclude the two lines will cross quite near zero. Or, do you want them to cross near x==1? Just as easy.

ax(2).YLim = [0.95 1.005];

Again, if this is a course in how to lie with statistics, I think I will pass with flying colors. (In fact, I might even want to apply for that government job now.) The conclusion is simple. It is very easy to trick the eye into believing anything you want it to believe, with a careful choice of how we manipulate a plot. But is that really what you are looking to do?

My guess is you were yourself misled by what you saw on the plot, possibly because something in your data analysis is incorrect, creating badly scaled data. If your professor told you to do what you are asking to do, then your professor is mistken. If you are trying to do what you THINK your professor told you to do, then you are mistaken. I'm sorry, but one of those alternatives would be the truth.

[ShruDuebgen59]

The plot given by MATLAB automatically scaled the x-axis like that. This is for a project I am working on for a class. When my professor when he helped out this afternoon, he stated that the power and head rise coefficients need to be functions of the flow coefficient, then plot said functions on one graph and find the intersection point of both curves.

[John D'Errico]

If you are going to be nasty, as you initially did before, then I will act appropriately. Don't attack the messenger, merely for telling you that you did something mathematically wrong.

If you plot two curves with arbtrary scaling, then you CANNOT find the intersection, as I explain in my answer. That action makes absolutely no mathematical sense. Sorry, but it does not.

I showed you how to find the intersectinos of the curves. But those intersection points are NOT where you think they are.

That means either you misunderstood your professor, or you made a mistake in your computations. Period.

[ShruDuebgen59]

I knew the plot was giving the incorrect intersection based on the axes, but it did find the correct intersection point regardless. I just needed to show the curve fits based on what my professor wanted in the project assignment and his hints to me and my group. There was no miscommunication, misleading or misunderstanding from my professor or the assignment. I wasn't trying to mislead my professor either with the graphs, just had trouble comprehending why the script was showing the right coordinates but the wrong plots. Through process of elimination, discussing with my groupmates, and professor, without the constant editing of your long sermon pertaining to the issue, we arrived that the plots were indeed off. We remedied that and made note to our professor that they were off.

No need to go on and on about me trying to mislead when I am just trying to better myself in MATLAB since I'm still new to the program.

[John D'Errico]

Um, you are the one who decided to call me a rather nasty name, for absolutely no good reason. It seems like you were the one who rather than understanding the message you were being given, decided to insult the messenger, rather than try to understand the message.

And, no. The plot you show in your question does NOT show the correct intersectino point. I proved that as fact.

As I said in my answer, perhaps you just misunderstood what you were doing. Your use pf plotyy might easily mislead you. But my point is still entirely valid. If you go down the road you were taking, you could "prove" anything you decide to "prove".

Answer 2

[darova]

Just change boundaries for "x" axis variable

Qspace= linspace(-0.01,0.04,40); %100 points for line
% ...
[x0, y0] = polyxpoly(Qspace,QHline, Qspace,QPline);

Answer 3

[ShruDuebgen59]

[John D'Errico]

Please stop adding answers every time you comment.

Answer 4

[ShruDuebgen59]

Hey darova, did you have both curves set to one y-axis, or two?

[John D'Errico]

Please don't add answers just to make a comment.