# How to highlight the intersections between two surfaces in a 3D plot

I would like to see the intersection of a 3D function g with the plane z=0.

En[a_, ru_] := (ru^(3/2) - 2*ru^(1/2) + a)/(ru^(3/4)*
Sqrt[ru^(3/2) - 3*ru^(1/2) + 2*a]);
L[a_, ru_] := (ru^(2) - 2*a*ru^(1/2) + a^2)/(ru^(3/4)*
Sqrt[ru^(3/2) - 3*ru^(1/2) + 2*a]);
\[Rho][a_, r_] := r^2 + a^2 + 2*a^2/r;
\[CapitalDelta][a_, r_] := r^2 - 2*r + a^2;
ra[a_, ru_] := (2*(a*En[a, ru] - L[a, ru])^2)/(ru^2*(1 - En[a, ru]^2));
Rr[a_, ru_, r_] := -(1 - En[a, ru]^2)*r^4 +
2*r^3 - (a^2*(1 - En[a, ru]^2) + L[a, ru]^2)*r^2 +
2*r*(a*En[a, ru] - L[a, ru])^2;
\[Upsilon][a_, ru_, r_] :=
Sqrt[Rr[a, ru, r]/(r*\[CapitalDelta][a, r]) +
L[a, ru]^2/\[Rho][a, r]];
\[Gamma][a_, ru_, r_] := 1/Sqrt[1 - \[Upsilon][a, ru, r]^2];
Nr[r_, a_] := Sqrt[\[CapitalDelta][a, r]/\[Rho][a, r]];
Nf[r_, a_] := -2*a/(r*\[Rho][a, r]);
gtphi[a_, RS_] := -2*a/RS;
gphiphi[a_, RS_] := \[Rho][a, RS];
gtt[a_, RS_] := -(1 - 2/RS);
\[CapitalOmega]max[a_,
RS_] := (-gtphi[a, RS] +
Sqrt[gtphi[a, RS]^2 - gphiphi[a, RS]*gtt[a, RS]])/
gphiphi[a, RS];
bb[a_, RS_, \[CapitalOmega]_] := -(gtphi[a, RS] +
gphiphi[a, RS]*\[CapitalOmega])/(gtt[a, RS] +
gtphi[a, RS]*\[CapitalOmega]);
\[Beta][r_, RS_, a_, \[CapitalOmega]_] :=
ArcCos[bb[a, RS, \[CapitalOmega]]*
Nr[r, a]/(Sqrt[\[Rho][a, r]]*(1 +
bb[a, RS, \[CapitalOmega]]*Nf[r, a]))];
Ar[a_, ru_, r_, RS_, \[CapitalOmega]_] := \[Gamma][a, ru, r] -
Sqrt[Rr[a, ru, r]]*Sin[\[Beta][r, RS, a, \[CapitalOmega]]]/r -
L[a, ru]*Cos[\[Beta][r, RS, a, \[CapitalOmega]]]/Sqrt[\[Rho][a, r]];
g[a_, ru_, r_,
RS_, \[CapitalOmega]_] := (2*a*L[a, ru] -
r*\[Rho][a, r]*En[a, ru])*(1/\[Gamma][a, ru, r] -
Ar[a, ru, r, RS, \[CapitalOmega]]) +
((r - 2)*L[a, ru] +
2*a*En[a, ru])*(Cos[\[Beta][r, RS, a, \[CapitalOmega]]] -
L[a, ru]*
Ar[a, ru, r, RS, \[CapitalOmega]]/Sqrt[\[Rho][a, r]]) + (Sqrt[
Rr[a, ru, r]]*\[CapitalDelta][a, r]/
r)*(Sin[\[Beta][r, RS, a, \[CapitalOmega]]] -
Sqrt[Rr[a, ru, r]/\[CapitalDelta][a, r]]*
Ar[a, ru, r, RS, \[CapitalOmega]]/r);

Z1[a_] := 1 + (1 - a^2)^(1/3)*((1 + a)^(1/3) + (1 - a)^(1/3));
Z2[a_] := Sqrt[3*a^2 + Z1[a]^2];
rISCO[a_] := 3 + Z2[a] - Sqrt[(3 - Z1[a])*(3 + Z1[a] + 2*Z2[a])];
rIBCO[a_] := 2 - a + 2*Sqrt[1 - a];
rH[a_] := 1 + Sqrt[1 - a^2];
barR[a_] := rIBCO[a] + 2/5*(rIBCO[a] + rISCO[a]);
RS = 6;
a = 0.2;
ru = barR[a] - 0.001;
\[CapitalOmega]max[a, RS];


After defined the function g, I do the following plot:

Plot3D[{g[a, ru, r, RS, \[CapitalOmega]], 0}, {r, ru,
ra[a, ru]}, {\[CapitalOmega], 0, \[CapitalOmega]max[a, RS]},
PlotRange -> {{ru, ra[a, ru]}, {0, 0.14}, {-20, 20}},
AxesStyle -> Directive[20, Black],
AxesLabel -> {Style[Rotate["\!$$\* StyleBox[\"r\",\nFontColor->GrayLevel[0]]$$ (M)", 0 Degree]],
Style[Rotate[
"\!$$\*SubscriptBox[\(\[CapitalOmega]$$, $$\[FivePointedStar]$$]\
\) (\!$$\*SuperscriptBox[\(M$$, $$-1$$]\))", 65 Degree]],
Style[Rotate["g", 0 Degree]]}, ImageSize -> 800,
PlotStyle -> Opacity[0.5], MeshFunctions -> {#3 &},
MeshStyle -> {Red}]


The result of the plot is the following

Now how can I highlight the intersection between the function g (yellow surface) and the plane z=0 (blue surface).

• just add the option Mesh -> {{0}}?
– kglr
Oct 21 '19 at 23:33
• Thank you very much! It works!
– VDF
Oct 21 '19 at 23:41

You can add the option Mesh -> {{0}}

Plot3D[{g[a, ru, r, RS, Ω], 0}, {r, ru, ra[a, ru]}, {Ω, 0, Ωmax[a, RS]},
PlotRange -> {{ru, ra[a, ru]}, {0, 0.14}, {-20, 20}}, PlotStyle -> Opacity[0.5],
Mesh -> {{0}}, MeshFunctions -> {#3 &}, MeshStyle -> {Red}]


to get

If you need the additional mesh lines you can style the line corresponding to z=0 differently:

Plot3D[{g[a, ru, r, RS, Ω], 0}, {r, ru, ra[a, ru]}, {Ω, 0, Ωmax[a, RS]},
PlotRange -> {{ru, ra[a, ru]}, {0, 0.14}, {-20, 20}}, PlotStyle -> Opacity[0.5],
Mesh -> {Append[{#, Orange} & /@ FindDivisions[{-20, 20}, 10], {0,
Directive[Red, Thick]}]},
MeshFunctions -> {#3 &}]