This question is related to:

How to obtain the coordinates of the intersection line of two surfaces?

Read the link first.

First, I have tried MaxCellMeasure->0.5, MaxCellMeasure->{"Length"->0.5} and MaxCellMeasure->{"Area"->0.5}, and change the number 0.5. It seems that the numbers of point doesn't change.

Secondly, it's very fast to run your code, but it takes about minutes to run my code with my f[x,y,z], g[x,y,z]. And this is my code:

Clear[f,g];
f[x_,y_,z_,t_]:=x^2+y^2+(z-1)^2-t^2/4;
g[x_,y_,z_]:=Piecewise[{{(x-1/2)^2+y^2+z^2-1/16,(x-1/2)^2+y^2-1/16<=0&&z>=0}},z];
reg=DiscretizeRegion[ImplicitRegion[{f[x,y,z,2]==0, g[x, y, z] == 0}, {x, y, z}]];
MeshCoordinates[reg]

Is it because I used Piecewise in g function that it needs more time to run the code?

Third, what I want to do next is below. In fact, I want to plot the projection of the intersection line on xy plane. And this is point projection, as shown below

enter image description here

Connect (0,0,1) and one point on the intersection line, extend it to intersect with xy plane. The yellow pint is want I am looking for. Repeat this for all points on the intersection line, then I can get a closed curve on xy plane. And I need to plot the curve in a 2-dimensional figure, instead of three-dimensional figure here.(I don't want to use Plot3D+Viewpoint to show a 2d curve in 3d coordinate, because the exported pdf file will be large.) I can calculate the yellow points' coordinate (x',y') based on (x,y,z) of red points on the intersection line. x'=x/(1-z),y'=y(1-z). But how to achieve this?

Or maybe I should use another way: use ParametricPlot

ParametricPlot[{x/(1-z),y/(1-z)},{x,y,z}] and x,y,z satisfy f[x,y,z]==g[x,y,z]
  • 3
    "Read the link first" - Not the best approach to asking well-received questions. Please make your question self-contained, and perhaps shorted as well. Note also that you are linking to a deleted question; not everybody will have access to those. Finally, please format your code (see e.g. here). – MarcoB Apr 30 '16 at 20:30

It actually takes me a fraction of a second to run your code (version 10.4 on OS X 10.11.4).

You can use MaxCellMeasure to get more points:

Length@MeshCoordinates@
     DiscretizeRegion[
      ImplicitRegion[{f[x, y, z, 2] == 0, g[x, y, z] == 0}, {x, y, 
        z}], MaxCellMeasure -> {"Length" -> #}] & /@ {0.1, 0.05, 0.005} // AbsoluteTiming
(* {0.368696, {76, 101, 965}} *)

One way to achieve the projection is to use ScalingTransform:

curve = MeshCoordinates[reg];
center = {0, 0, 1};
proj = ScalingTransform[{1, 1, 1}/(1 - #[[3]]), center][#] & /@ curve;
Graphics3D[{Point@curve, Blue,Point@proj, Red, PointSize[Large], Point[{0, 0, 1}]}]

enter image description here

If you want to plot the points in 2D, you can sort the list by the closest neighbor and then plot

ls = proj;
pt = proj[[1]];

line = Prepend[Table[
    ls = DeleteCases[ls, pt];
    pt = First@Nearest[ls, pt],
    {n, 1, Length@ls - 1}
    ], pt];

Graphics[Line@line[[All, 1 ;; 2]], PlotRange -> All]

enter image description here

  • Good, thank you. My mathematica is Version 10 on Win8. And the data in proj is like (x,y,z) right. What if I only want to plot (x,y)? – Qi Zhong Apr 30 '16 at 20:53
  • Version 10.0 student edition on Win8. I am confused why it takes so long. – Qi Zhong Apr 30 '16 at 21:00
  • Length@MeshCoordinates@DiscretizeRegion[ImplicitRegion[{f[x, y, z, 2] == 0, g[x, y, z] == 0}, {x, y, z}], MaxCellMeasure -> {"Length" -> #}] & /@ {0.1, 0.5} // AbsoluteTiming ({240.931,{157,157}} *) if g=x+y+z, (*{1.4330,{115,115}}) the maxcellmeasure didn't work. – Qi Zhong Apr 30 '16 at 21:21
  • Perfect. I still need to figure out why maxcellmeasure doesn't work in my computer and why it takes so long time to sun the code. Is there any hint to solve this? – Qi Zhong Apr 30 '16 at 21:32
  • I can reproduce the slowness in version 10.0 and 10.1 and I would consider it as a bug. It looks like it has been fixed started from version 10.2. – xslittlegrass Apr 30 '16 at 21:56

Your Answer

By clicking "Post Your Answer", you acknowledge that you have read our updated terms of service, privacy policy and cookie policy, and that your continued use of the website is subject to these policies.

Not the answer you're looking for? Browse other questions tagged or ask your own question.