Solving for level curves of an elliptic paraboloid given by quadric surface equation

Hi there, I have an equation which describes an elliptic paraboloid: Ax^2+Bxy+Cy^2+Dx+Ey+F = Z

(Note, the coefficients A,B,C,D,E and F all satisfy the necessary conditions to make an elliptic paraboloid).

In general, B is not zero, so the cross-section is a rotated ellipse (not centered at zero). I would like to solve for the ellipse cross-section (level curve) at a given height z, and to get the vertices of this ellipse. It would be nice to plot the ellipse, too. I have to do this over and over again, so the fastest way would be appreciated! How can I do this? I'm stumped. I've tried various versions of solve to no avail. Thanks in advance.

Answer 1

[Stephan]

Hi,

here is an approach - free for play with it and improvement... Just to let you start:

% Define a function to play with
a = -2;
b = -2;
c = -2;
d = -2;
e = -2; 
f = 5;
[y, x] = meshgrid(-10:10,-10:10);
fun = (a.*x.^2+b.*x.*y+c.*y.^2+d.*x+e.*y+f);

% We want to get the ellipse at z=-100 
z_cross = - 100

% Values for y to calculate x from with defined z=-100
y_vals = -15:0.25:15;

% Resulting x-values - only the non-complex solution are interesting
x_cross1 = real((d.*(-1.0./2.0)+sqrt(a.*f.*-4.0+a.*z_cross.*4.0+d.^2+b.^2.*y_vals.^2-a.*e.*y_vals.*4.0+b.*d.*y_vals.*2.0-a.*c.*y_vals.^2.*4.0).*(1.0./2.0)-b.*y_vals.*(1.0./2.0))./a);
x_cross2 = real((d.*(-1.0./2.0)-sqrt(a.*f.*-4.0+a.*z_cross.*4.0+d.^2+b.^2.*y_vals.^2-a.*e.*y_vals.*4.0+b.*d.*y_vals.*2.0-a.*c.*y_vals.^2.*4.0).*(1.0./2.0)-b.*y_vals.*(1.0./2.0))./a);
x_cross = [x_cross1 x_cross2];

% Calculate the corresponding y-values for the non-comlpex x-values
y_cross1 = real((e.*(-1.0./2.0)+sqrt(c.*f.*-4.0+c.*z_cross.*4.0+e.^2+b.^2.*x_cross1.^2+b.*e.*x_cross1.*2.0-c.*d.*x_cross1.*4.0-a.*c.*x_cross1.^2.*4.0).*(1.0./2.0)-b.*x_cross1.*(1.0./2.0))./c);
y_cross2 = real((e.*(-1.0./2.0)-sqrt(c.*f.*-4.0+c.*z_cross.*4.0+e.^2+b.^2.*x_cross2.^2+b.*e.*x_cross2.*2.0-c.*d.*x_cross2.*4.0-a.*c.*x_cross2.^2.*4.0).*(1.0./2.0)-b.*x_cross2.*(1.0./2.0))./c);
y_cross = [y_cross1 y_cross2];

% This is not needed - just for controll --> z = -100
z_cross = (a.*x_cross.^2+b.*x_cross.*y_cross+c.*y_cross.^2+d.*x_cross+e.*y_cross+f);

% plot result
mesh(x,y,fun)
hold on
scatter3(x_cross,y_cross,z_cross,'r','LineWidth',2)
hold off

This gives you:

.

What is still to be done?

• calculate more points in the open areas - my solution doesnt reach all points
• the representation as a line plot leads to undesirable connecting lines, which can be seen from the bottom. That would be in need of improvement.
• check if a conversion into parameter representation / other coordinates makes sense
• since the calculation is vectorized i think performance is not so bad - but perhaps it can be done (much) better...

The points of the ellipse are obtained from the calculations of the vectors x_cross and y_cross.

Since you have these points, you can calculate the ellipse equation and thus find an analytical description of your target - then you could represent this as a continuous line ...

I hope for an interesting discussion to this topic ;-)

Best regards

Stephan

[supernoob]

This is great! I ended up doing something slightly different. Your solution is a little more complete as it gives the whole level curve. I settled for the vertices and the center, and used contour() to plot.

For the vertices, I got the lengths of the semimajor axes and the center of the ellipse. For this I used equations (19)-(22) of http://mathworld.wolfram.com/Ellipse.html For plotting, used contour():

%generate plot points on a grid based on least squares polynomial
Z = (a.*X.^2) + (b.*X.*Y) + (c.*Y.^2) + (d.*X) + (e.*Y) + f;

figure(ind);
v = [val, val]; %plot at desired height z=val
contour(X,Y,Z,v, 'Color', 'k', 'Linewidth', 2);