Represent the intersection curve between a sphere and a paraboloid

Questions
, ,
1

I have a simple problem and it seems that I can’t reach a proper solution for it. I need to find and also graphically represent the intersection curve between two surfaces: a sphere of radius R and a paraboloid.
These two figures are represented with the help of the vectors x1,x2,x3 which are written in terms of the spherical coordinates (r,\theta,\phi). Their equations are given by eqs1 and eqs2, respectively, in the code.
The figures are drawn with the help of the function SphericalPlot3D

Now, my idea was to find the solution for the equation

eqs1[r,\theta,\phi]=eqs2[r,\theta,\phi]

expressed in terms of (\theta,\phi), then just add that solution to the Show[] function. However, doing a SphericalPlot for that line term which I’ve obtained, doesn’t solve the problem. That intersection should be en ellipse only.

How can I solve the problem and get an ellipse right where the intersection of the two surfaces should be?
(I found some examples, but those were for cartesian coordinates, and I have to work with the components x1,x2,x3 expressed in spherical ones.)

Thank you in advance!

Here is the code that I’ve written so far:

x1[r_, theta_, fi_] := r*Cos[theta];
x2[r_, theta_, fi_] := r*Sin[theta]*Cos[fi];
x3[r_, theta_, fi_] := r*Sin[theta]*Sin[fi];

(*Constants*)
A1 = 1;
A2 = 3;
A3 = 6;
R = 4;

(*Surface 1 (paraboloid)*)

eqs1[r_, theta_, fi_] := 
  A1*x1[r, theta, fi]^2 + A2*x2[r, theta, fi]^2 + A3*x3[r, theta, fi];
s1 = Solve[eqs1[r, \[Theta], \[Phi]] == R^2, r][[2]];
surface1 = r /. s1;
plot1 = SphericalPlot3D[{surface1}, {\[Theta], 0, \[Pi]}, {\[Phi], 0, 
    2 \[Pi]}, PlotStyle -> {Blue, Opacity[0.3]}];

(*Surface 2 (sphere)*)

eqs2[r_, theta_, fi_] := 
  x1[r, theta, fi]^2 + x2[r, theta, fi]^2 + x3[r, theta, fi]^2;
surface2 = r /. s2;
plot2 = SphericalPlot3D[{surface2}, {\[Theta], 0, \[Pi]}, {\[Phi], 0, 
    2 \[Pi]}, PlotStyle -> {Red}];
Show[plot2, plot1]

(*my attempt to find the curve*)
pts = Solve[eqs1[r, \[Theta], \[Phi]] == eqs2[r, \[Theta], \[Phi]], 
    r][[2]];
line = r /. pts;

SphericalPlot3D[line, {\[Theta], 0, \[Pi]}, {\[Phi], 0, 2 \[Pi]}, 
 PlotStyle -> {Red}]
2

You can draw the intersection with MeshFunctions:

eqs1[r_, theta_, fi_] := 
  A1*x1[r, theta, fi]^2 + A2*x2[r, theta, fi]^2 + 
   A3*x3[r, theta, fi];
s1 = Solve[eqs1[r, \[Theta], \[Phi]] == R^2, r][[2]];
surface1 = r /. s1;
eqs2[r_, theta_, fi_] := 
  x1[r, theta, fi]^2 + x2[r, theta, fi]^2 + x3[r, theta, fi]^2;
s2 = Solve[eqs2[r, \[Theta], \[Phi]] == R^2, r][[2]] // Simplify;
surface2 = r /. s2;
SphericalPlot3D[{surface1, surface2}, {\[Theta], 0, \[Pi]}, {\[Phi], 
  0, 2 \[Pi]}, PlotStyle -> {Directive[Blue, Opacity[0.3]], Red}, 
 MeshFunctions -> {Function[{x, y, z, \[Theta], \[Phi], r}, 
    Evaluate[surface1]]}, Mesh -> {{4}}, BoundaryStyle -> None, 
 MeshStyle -> Thick, PlotPoints -> 50]
3

@Gianluca thank you so much for the answer!!!
Your solution works as intended for me, however, now I took the problem to another level and I am facing an issue that I can’t solve.
In my initial problem, my second surface was a sphere, and the first one a basic paraboloid. The equation of the paraboloid was straightforward.
Now, I still have the sphere, but the paraboloid is parameterized by three numbers (in the code: e, u and v).

Hsim[r_, theta_, fi_, u_, v_] := 
  x2[r, theta, fi]^2 + u*x3[r, theta, fi]^2 + 2 v*x1[r, theta, fi];
eqs1[e_, u_, v_] := NSolve[Hsim[r, \[Theta], \[Phi], u, v] == e, r];
surface1[e_, u_, v_] := r /. eqs1[e, u, v][[2]];

I did the same thing and defined surface1 as being the solution to an equation, but this time the equation is Hsim[r,\theta,\phi,e,u,v]==e, which again I have to solve for r.

The problem is in the MeshFunction. If I try in the same manner MeshFunction -> Function[{x, y, z, r, \[Theta], \[Phi]}, Evaluate[surface1[e, up, v0p]]], I don’t get the intersection curve (ellipse) between the sphere and the paraboloid.
e is just a number that I can modify in a Manipulate, to see how the intersection curve modifies with the increase if this value.
The parameters u,v (denoted in the document with up, and v0p are know beforehand).
What is the issue here?

PS: I even tried putting the analytic form of the surface in the MeshFunction like so:

MeshFunctions -> {Function[{x, y, z, 
    r, \[Theta], \[Phi]}, (0.3431771859697404` (7.3454135`*^7 Cos[\
\[Theta]] + 
         1.96286546`*^8 Sqrt[
          0.14003978386304503` Cos[\[Theta]]^2 + 
           4.` e Cos[\[Phi]]^2 Sin[\[Theta]]^2 + 
           0.3582649439870417` e Sin[\[Theta]]^2 \
Sin[\[Phi]]^2]))/(1.34722129`*^8 Cos[\[Phi]]^2 Sin[\[Theta]]^2 + 
       1.2066554`*^7 Sin[\[Theta]]^2 Sin[\[Phi]]^2) /. {e -> e}]}

But the lines that I get are not ok.

I attached here a simple document with the example. My first manipulate is without the MeshFunction, and that works fine. The second Manipulate is where I try to implement the intersection curves, unsuccessfully.

simplifiedHamiltonian_Pr135.nb (25.6 KB)

4

This seems to work:

eqs1[e_, u_, v_] = Last@Solve[Hsim[r, \[Theta], \[Phi], u, v] == e, r];
surface1[e_, u_, v_] := r /. eqs1[e, u, v];
Manipulate[
 Show[Graphics3D[{Red, Sphere[{0, 0, 0}, 9.5]}], 
  SphericalPlot3D[9.5, {\[Theta], 0, \[Pi]}, {\[Phi], 0, 2 \[Pi]}, 
   PlotStyle -> None, 
   MeshFunctions -> {Function[{x, y, z, \[Theta], \[Phi], r}, 
      Evaluate[Hsim[r, \[Theta], \[Phi], up, v0p]]]}, Mesh -> {{e}}, 
   MeshStyle -> Directive[Yellow, Thick], PlotPoints -> 50, 
   BoundaryStyle -> None],
  SphericalPlot3D[
   surface1[e, up, v0p], {\[Theta], 0, \[Pi]}, {\[Phi], 0, 2 \[Pi]}, 
   PlotStyle -> Directive[Blue, Opacity[0.5]], Mesh -> None, 
   PlotPoints -> 40],
  PlotRange -> 13],
 {e, 0, 50, 2}]
5

Thank you again! This solution works like a charm!

Best Regards

6

I think

s2 = Solve[eqs2[r, \[Theta], \[Phi]] == R^2, r][[2]] // Simplify;

is not appropriate. We are looking for an intersection of the paraboloid and the sphere. All these points are elements of the spherical surface and therefore have a r-value of R, so we must not look for r.

The intersection being a line means that it is a one-dimensional entity, depending only on ONE pararmeter. So I think fixing r to the radius of the sphere we should ask for theta as function of phi or phi as function of theta and insert this in the representation of a point:

For example

s3 = Solve[eqs1[R, \[Theta], \[Phi]] == R^2, \[Phi]]
Length[s3]
plot3 = ParametricPlot3D[{x1[R, \[Theta], \[Phi]], 
    x2[R, \[Theta], \[Phi]], x3[R, \[Theta], \[Phi]]} /. 
   s3[[5, 1]], {\[Theta], 0, Pi},
  PlotStyle -> {Thick, Blue}]
Show[plot2, plot3]

or

s4 = Solve[eqs1[R, \[Theta], \[Phi]] == R^2, \[Theta]]
Length[s4]
plot4 = ParametricPlot3D[{x1[R, \[Theta], \[Phi]], 
    x2[R, \[Theta], \[Phi]], x3[R, \[Theta], \[Phi]]} /. 
   s4[[6, 1]], {\[Phi], 0, 2 Pi},
  PlotStyle -> {Thick, Blue}]
Show[plot2, plot4]

As s3 and s4 yield several solutions it seems that the intersection is given piecewise. But put together it seems that both solutions are equivalent, as it should be:

plot5= ParametricPlot3D[
{x1[R, \[Theta], \[Phi]], x2[R, \[Theta], \[Phi]], 
x3[R, \[Theta], \[Phi]]} /. # & /@ Flatten[s3],
{\[Theta], 0, Pi},
PlotStyle -> {Thick, Blue}]

plot6= ParametricPlot3D[{x1[R, \[Theta], \[Phi]], x2[R, \[Theta], \[Phi]], 
x3[R, \[Theta], \[Phi]]} /. # & /@ Flatten[Drop[s4, 3]],
{\[Phi], 0, 2 Pi}, PlotStyle -> {Thick, Blue}]

Show[plot2, plot5, plot6]
7

But perhaps this is it what is asked for

s5 = Flatten[FullSimplify[Solve[surface1 == R^2, \[Theta]]]];
plot7 = ParametricPlot3D[
  {x1[R, \[Theta], \[Phi]], x2[R, \[Theta], \[Phi]], 
      x3[R, \[Theta], \[Phi]]} /. # & /@ s5,
  {\[Phi], 0, Pi},
  PlotStyle -> {Thick, Blue}]

or

s6 = Flatten[FullSimplify[Solve[surface1 == R^2, \[Phi]]]];
plot6 = ParametricPlot3D[{x1[R, \[Theta], \[Phi]], 
      x2[R, \[Theta], \[Phi]], x3[R, \[Theta], \[Phi]]} /. # & /@ s6,
  {\[Theta], 0, Pi}, PlotStyle -> {Thick, Blue}]
8

but those were for cartesian coordinates, and I have to work with the
components x1,x2,x3 expressed in spherical ones

Whatever that means. I understand “spherical coordinates” as a set { r, theta, phi }.

The conversion to cartesian coordinates are the well known x1, x2, x3 given above.

Last not least I think the problem should be tackled in cartesian coordinates.

Consider a paraboloid “pointing” in the negative y-direction:

para = {a u Cos[v], -u^2, b u Sin[v]};

This surface is an image of a (planar) region (subset) of R2 with 0<= u < Infinity and
0<= v <= 2Pi.
All points with para . para == r^2 lie on or in the surface of a sphere with radius r.

Now we can determine those u (as a function of v) which fullfill this condition, giving somewhat complicated expressions ( here given in para1) , but in principle it is quite easy to extract the
spherical coordinates, at least numerically.

The follwing code gives the line of intersection (on the sphere) for different paraboloids

R = 3;
Manipulate[
 para = {a u Cos[v], -u^2, b u Sin[v]}; 
 sol1 = Flatten[Solve[para.para == R^2, u]];
 para1 = para /. # & /@ sol1;
 pl1 = ParametricPlot3D[
   para1,
   {v, 0, 2 Pi},
   PlotStyle -> {Thick, Blue}];
 pl2 = Graphics3D[Sphere[{0, 0, 0}, R]];
 Show[pl2, pl1],
 {a, -5, 5}, {b, -5, 5}]

and here we have an element of para1 and its spherical coordinates:

vec = Re[para1 /. v -> .345][[4]]
th0 = ArcCos[vec[[3]]/R]
phi0 = ArcTan[vec[[1]], vec[[2]]]
R {Cos[phi0] Sin[th0], Sin[phi0] Sin[th0], Cos[th0]}