Message Boards Message Boards

1
|
14174 Views
|
3 Replies
|
4 Total Likes
View groups...
Share
Share this post:

Combine Plot Graphics

Posted 12 years ago
I have one Graphics3D and one ParametricPlot3D, which I need to combine. Unfortunately, the Show[] Command does not work here, maybe because I use a Manipulate[] Command. The Prolog-> also does not work.

Plot1, the Graphics 3D, is a set of 3 Points (the initial positions of Plot2):
(* Plot 1 *)
P1 = {0, 0, 0}; P2 = {-1, 0, 2}; P3 = {-3, 4, -2};
Plot1 = Graphics3D[Point[{P1, P2, P3}],
PlotRange -> {{-5, 5}, {-5, 5}, {-3, 7}}, AspectRatio -> 1];

Plot 2, the ParametricPlot3D, is a 3D simulation of the three-body-problem. I need to have the Plot1 plotted right into Plot2:

(Edit: I dont know if the Code gets lost due conversion to plain text. If so, see .nb File at http://yukterez.ist.org/3kp/3k3D.nb

  (* Plot 2 *)
  G = 667384/
     10^16;(* Gravitationskonstante *)
  M1 = 100000000; M2 = 70000000; \
  M3 = 50000000;(* Massen 1,2,3 in kg *)
  v1x = 1/70; v1y = -1/130; v1z =
    1/20;(* v0 von M1 in m/sek *)
  v2x = -1/150; v2y = 1/100; v2z =
   1/90;(* v0 von M2 in m/sek *)
 v3x = 1/200; v3y =
  1/110; v3z = -1/170;(* v0 von M3 in m/sek *)
 
 P1 = {0, 0, 0}; (* Position 1 XYZ *)
 
 P2 = {-1, 0, 2}; (* Position 2 XYZ *)
 
 P3 = {-3, 4, -2}; (* Position 3 XYZ *)
 
 (* Container *)
 Funktion[{{x10_, y10_, z10_}, {x20_, y20_, z20_}, {x30_, y30_,
     z30_}}, {{vx10_, vy10_, vz10_}, {vx20_, vy20_, vz20_}, {vx30_,
     vy30_, vz30_}}, {m1_, m2_, m3_}, T_, plotType : ("x" | "v"),
   plotOptions___] := Module[{nds, Tmax, prolog, funcToPlot},
    
    (* Differentialgleichung *)nds = NDSolve[{
       
       x1'[t] == vx1[t], y1'[t] == vy1[t], z1'[t] == vz1[t],
       
       x2'[t] == vx2[t], y2'[t] == vy2[t], z2'[t] == vz2[t],
       
       x3'[t] == vx3[t], y3'[t] == vy3[t], z3'[t] == vz3[t],
       
       
       vx1'[
         t] == -((G  m2 (x1[t] - x2[t]))/
          Sqrt[((x1[t] - x2[t])^2 + (y1[t] - y2[t])^2 + (z1[t] -
              z2[t])^2)^3]) - (G  m3 (x1[t] - x3[t]))/
         Sqrt[((x1[t] - x3[t])^2 + (y1[t] - y3[t])^2 + (z1[t] -
             z3[t])^2)^3],
       
       
        vy1'[
         t] == -((G  m2 (y1[t] - y2[t]))/
          Sqrt[((x1[t] - x2[t])^2 + (y1[t] - y2[t])^2 + (z1[t] -
              z2[t])^2)^3]) - (G m3 (y1[t] - y3[t]))/
         Sqrt[((x1[t] - x3[t])^2 + (y1[t] - y3[t])^2 + (z1[t] -
             z3[t])^2)^3],
       
       
        vz1'[
         t] == -((G m2 (z1[t] - z2[t]))/
          Sqrt[((x1[t] - x2[t])^2 + (y1[t] - y2[t])^2 + (z1[t] -
              z2[t])^2)^3]) - (G m3 (z1[t] - z3[t]))/
         Sqrt[((x1[t] - x3[t])^2 + (y1[t] - y3[t])^2 + (z1[t] -
             z3[t])^2)^3],
       
       
       vx2'[t] == (G m1  (x1[t] - x2[t]))/
         Sqrt[((x1[t] - x2[t])^2 + (y1[t] - y2[t])^2 + (z1[t] -
             z2[t])^2)^3] - (G  m3 (x2[t] - x3[t]))/
         Sqrt[((x2[t] - x3[t])^2 + (y2[t] - y3[t])^2 + (z2[t] -
             z3[t])^2)^3],
       
       
        vy2'[t] == (G m1  (y1[t] - y2[t]))/
         Sqrt[((x1[t] - x2[t])^2 + (y1[t] - y2[t])^2 + (z1[t] -
             z2[t])^2)^3] - (G  m3 (y2[t] - y3[t]))/
         Sqrt[((x2[t] - x3[t])^2 + (y2[t] - y3[t])^2 + (z2[t] -
             z3[t])^2)^3],
       
       
        vz2'[t] == (G m1  (z1[t] - z2[t]))/
         Sqrt[((x1[t] - x2[t])^2 + (y1[t] - y2[t])^2 + (z1[t] -
             z2[t])^2)^3] - (G m3 (z2[t] - z3[t]))/
         Sqrt[((x2[t] - x3[t])^2 + (y2[t] - y3[t])^2 + (z2[t] -
             z3[t])^2)^3],
       
       
        vx3'[t] == (G m1  (x1[t] - x3[t]))/
         Sqrt[((x1[t] - x3[t])^2 + (y1[t] - y3[t])^2 + (z1[t] -
             z3[t])^2)^3] + (G m2  (x2[t] - x3[t]))/
         Sqrt[((x2[t] - x3[t])^2 + (y2[t] - y3[t])^2 + (z2[t] -
             z3[t])^2)^3],
       
       
        vy3'[t] == (G m1  (y1[t] - y3[t]))/
         Sqrt[((x1[t] - x3[t])^2 + (y1[t] - y3[t])^2 + (z1[t] -
             z3[t])^2)^3] + (G m2  (y2[t] - y3[t]))/
         Sqrt[((x2[t] - x3[t])^2 + (y2[t] - y3[t])^2 + (z2[t] -
             z3[t])^2)^3],
       
       
        vz3'[t] == (G m1  (z1[t] - z3[t]))/
         Sqrt[((x1[t] - x3[t])^2 + (y1[t] - y3[t])^2 + (z1[t] -
             z3[t])^2)^3] + (G m2  (z2[t] - z3[t]))/
         Sqrt[((x2[t] - x3[t])^2 + (y2[t] - y3[t])^2 + (z2[t] -
             z3[t])^2)^3],
       
       
      x1[0] == x10, y1[0] == y10, z1[0] == z10,
      
      x2[0] == x20, y2[0] == y20, z2[0] == z20,
      
      x3[0] == x30, y3[0] == y30, z3[0] == z30,
      
      vx1[0] == vx10, vy1[0] == vy10, vz1[0] == vz10,
      
      vx2[0] == vx20, vy2[0] == vy20, vz2[0] == vz20,
      
      vx3[0] == vx30, vy3[0] == vy30, vz3[0] == vz30},
    
     {x1, x2, x3, y1, y2, y3, z1, z2, z3, vx1, vx2, vx3, vy1, vy2,
      vy3, vz1, vz2, vz3}, {t, 0, T}];
   
   If[Head[nds] =!= NDSolve, Tmax = nds[[1, 1, 2, 1, 1, 2]];
   
    funcToPlot =
     If[plotType ===
        "x", {{x1[t], y1[t], z1[t]}, {x2[t], y2[t], z2[t]}, {x3[t],
         y3[t], z3[t]}}, {{vx1[t], vy1[t], vz1[t]}, {vx2[t], vy2[t],
         vz2[t]}, {vx3[t], vy3[t], vz3[t]}}] /. nds[[1]];
   
    (* Plot Specifications *)
   
    prolog = {PointSize[0.3], Point[{P1, P2, P3}]};
   
    ParametricPlot3D[Evaluate[funcToPlot], {t, 0, Tmax},
     PlotStyle -> {Red, Blue, Green},
    
     (* Plot Range und Seitenverhältnis *)
    
     PlotRange -> {{-5, 5}, {-5, 5}, {-3, 7}}, AspectRatio -> 1,
     (* Prolog -> prolog, *) MaxRecursion -> ControlActive[3, 9],
     Axes -> True, plotOptions], Text["Yukterez Mod. of demonstrations.wolfram.com/PlanarThreeBodyProblem"]]] // Quiet

Plot2 = Manipulate[Funktion[
   (* Positionen xy *){P1, P2, P3},
   (* Geschwindigkeiten xy *){{v1x, v1y, v1z}, {v2x, v2y, v2z}, {v3x,
     v3y, v3z}},
   (* Massen *){M1, M2, M3},
   (* Regler *)T, xv,
   ImageSize -> {320, 320}],
  {{xv, "x", "Position/Geschwindigkeit"},
   {"x" -> "Position", "v" -> "Geschwindigkeit"}},
  {{T, 1, "Zeit"}, 1, 200}, ControlPlacement -> {Bottom}]


(* http://yukterez.ist.org/3kp/3k3D.html *)
POSTED BY: Simon Tyran
3 Replies
May be you ran the Newton's eqilibrium equation formulation  including inverse square law already to verify angular momentum conservation, the ellipse orbit etc.. Suggest taking mass of earth = 1, radius =1 , or better in astronomical units and adjust constants like G accordingly to obtain  reasonable values of distances and times.
Posted 12 years ago
Thanks a lot! I updated the .nb File to simulate sun, earth and moon to test the physics of the code. The Three body simulation now runs with the code
  (* Yukterez 2013, Mathematica Syntax *)
  
  Au = 149597870690; (* sun earth distance in m *)
  Mu = 384400000; (* moon earth distance in m *)
  mSol = 1988*10^26; (* sun mass in kg *)
  mE = 5972*10^21; (* earth mass in kg *)
  vE = Sqrt[G*mSol/Au]; (* orbital velocitvy earth to sun in m/sek *)
  mM = 7346*10^18; (* moon mass in kg *)
  vM1 = Sqrt[G*mSol/(Au + Mu)]; (* v orbit moon to sun in m/sek *)
 vM2 = Sqrt[G*mE/Mu]; (* v orbit moon to earth in m/sek *)
 
 G = 667384/10^16;(* gravitational constant *)
 M1 = mSol; M2 = mE; M3 = mM;(* mass 1,2,3 in kg *)
 v1x = 20; v1y = 210; v1z = 90;(* v0 of M1 in m/sek *)
 v2x = 20 + vE; v2y = 210; v2z = 90;(* v0 of M2 in m/sek *)
 v3x = 20 + vM1; v3y = 210; v3z = 90 + vM2;(* v0 of M3 in m/sek *)
 P1 = {0, 0, 0}; (* Position 1 XYZ *)
 P2 = {0, Au, 0}; (* Position 2 XYZ *)
 P3 = {0, Au + Mu, 0}; (* Position 3 XYZ *)
 
 Plot1 = Graphics3D[Point[{P1, P2, P3}],
    PlotRange -> {{-1.5 Au, 1.5 Au}, {-1.5 Au, 1.5 Au}, {-1.5 Au,
       1.5 Au}}, AspectRatio -> 1];
 
 (* Container *)
 Funktion[{{x10_, y10_, z10_}, {x20_, y20_, z20_}, {x30_, y30_,
     z30_}}, {{vx10_, vy10_, vz10_}, {vx20_, vy20_, vz20_}, {vx30_,
     vy30_, vz30_}}, {m1_, m2_, m3_}, T_, plotType : ("x" | "v"),
   plotOptions___] := Module[{nds, Tmax, prolog, funcToPlot},
    
    (* Differential Equation *)
    
    nds = NDSolve[{
       
       x1'[t] == vx1[t], y1'[t] == vy1[t], z1'[t] == vz1[t],
       
       x2'[t] == vx2[t], y2'[t] == vy2[t], z2'[t] == vz2[t],
       
       x3'[t] == vx3[t], y3'[t] == vy3[t], z3'[t] == vz3[t],
       
       
       vx1'[
         t] == -((G  m2 (x1[t] - x2[t]))/
          Sqrt[((x1[t] - x2[t])^2 + (y1[t] - y2[t])^2 + (z1[t] -
              z2[t])^2)^3]) - (G  m3 (x1[t] - x3[t]))/
         Sqrt[((x1[t] - x3[t])^2 + (y1[t] - y3[t])^2 + (z1[t] -
             z3[t])^2)^3],
       
       
        vy1'[
         t] == -((G  m2 (y1[t] - y2[t]))/
          Sqrt[((x1[t] - x2[t])^2 + (y1[t] - y2[t])^2 + (z1[t] -
              z2[t])^2)^3]) - (G m3 (y1[t] - y3[t]))/
         Sqrt[((x1[t] - x3[t])^2 + (y1[t] - y3[t])^2 + (z1[t] -
             z3[t])^2)^3],
       
       
        vz1'[
         t] == -((G m2 (z1[t] - z2[t]))/
          Sqrt[((x1[t] - x2[t])^2 + (y1[t] - y2[t])^2 + (z1[t] -
              z2[t])^2)^3]) - (G m3 (z1[t] - z3[t]))/
         Sqrt[((x1[t] - x3[t])^2 + (y1[t] - y3[t])^2 + (z1[t] -
             z3[t])^2)^3],
       
       
       vx2'[t] == (G m1  (x1[t] - x2[t]))/
         Sqrt[((x1[t] - x2[t])^2 + (y1[t] - y2[t])^2 + (z1[t] -
             z2[t])^2)^3] - (G  m3 (x2[t] - x3[t]))/
         Sqrt[((x2[t] - x3[t])^2 + (y2[t] - y3[t])^2 + (z2[t] -
             z3[t])^2)^3],
       
       
        vy2'[t] == (G m1  (y1[t] - y2[t]))/
         Sqrt[((x1[t] - x2[t])^2 + (y1[t] - y2[t])^2 + (z1[t] -
             z2[t])^2)^3] - (G  m3 (y2[t] - y3[t]))/
         Sqrt[((x2[t] - x3[t])^2 + (y2[t] - y3[t])^2 + (z2[t] -
             z3[t])^2)^3],
       
       
        vz2'[t] == (G m1  (z1[t] - z2[t]))/
         Sqrt[((x1[t] - x2[t])^2 + (y1[t] - y2[t])^2 + (z1[t] -
             z2[t])^2)^3] - (G m3 (z2[t] - z3[t]))/
         Sqrt[((x2[t] - x3[t])^2 + (y2[t] - y3[t])^2 + (z2[t] -
             z3[t])^2)^3],
       
       
        vx3'[t] == (G m1  (x1[t] - x3[t]))/
         Sqrt[((x1[t] - x3[t])^2 + (y1[t] - y3[t])^2 + (z1[t] -
             z3[t])^2)^3] + (G m2  (x2[t] - x3[t]))/
         Sqrt[((x2[t] - x3[t])^2 + (y2[t] - y3[t])^2 + (z2[t] -
             z3[t])^2)^3],
       
       
        vy3'[t] == (G m1  (y1[t] - y3[t]))/
         Sqrt[((x1[t] - x3[t])^2 + (y1[t] - y3[t])^2 + (z1[t] -
             z3[t])^2)^3] + (G m2  (y2[t] - y3[t]))/
         Sqrt[((x2[t] - x3[t])^2 + (y2[t] - y3[t])^2 + (z2[t] -
             z3[t])^2)^3],
       
      
       vz3'[t] == (G m1  (z1[t] - z3[t]))/
        Sqrt[((x1[t] - x3[t])^2 + (y1[t] - y3[t])^2 + (z1[t] -
            z3[t])^2)^3] + (G m2  (z2[t] - z3[t]))/
        Sqrt[((x2[t] - x3[t])^2 + (y2[t] - y3[t])^2 + (z2[t] -
            z3[t])^2)^3],
      
      
      x1[0] == x10, y1[0] == y10, z1[0] == z10,
      
      x2[0] == x20, y2[0] == y20, z2[0] == z20,
      
      x3[0] == x30, y3[0] == y30, z3[0] == z30,
      
      vx1[0] == vx10, vy1[0] == vy10, vz1[0] == vz10,
      
      vx2[0] == vx20, vy2[0] == vy20, vz2[0] == vz20,
      
      vx3[0] == vx30, vy3[0] == vy30, vz3[0] == vz30},
    
     {x1, x2, x3, y1, y2, y3, z1, z2, z3, vx1, vx2, vx3, vy1, vy2,
      vy3, vz1, vz2, vz3}, {t, 0, T}];
   
   If[Head[nds] =!= NDSolve, Tmax = nds[[1, 1, 2, 1, 1, 2]];
   
    funcToPlot =
     If[plotType ===
        "x", {{x1[t], y1[t], z1[t]}, {x2[t], y2[t], z2[t]}, {x3[t],
         y3[t], z3[t]}}, {{vx1[t], vy1[t], vz1[t]}, {vx2[t], vy2[t],
         vz2[t]}, {vx3[t], vy3[t], vz3[t]}}] /. nds[[1]];
   
    (* Plot Specifications *)
    prolog = {PointSize[0.3], Point[{P1, P2, P3}]};
   
    ParametricPlot3D[Evaluate[funcToPlot], {t, 0, Tmax},
     PlotStyle -> {Red, Blue, Green},
    
     (* Plot Range and Aspect Ratio *)
     PlotRange -> {{-1.5 Au, 1.5 Au}, {-1.5 Au, 1.5 Au}, {-1.5 Au,
        1.5 Au}}, AspectRatio -> 1,(* Prolog->prolog *)
     MaxRecursion -> ControlActive[3, 9], Axes -> True, plotOptions],
    Text["Yukterez Mod."]]] // Quiet

Plot2 = Manipulate[Show[Funktion[
    (* positions xyz *){P1, P2, P3},
    (* velocities xyz *){{v1x, v1y, v1z}, {v2x, v2y, v2z}, {v3x,
       v3y, v3z}},
    (* mass 1,2,3 *){M1, M2, M3},
    T, xv,
    ImageSize -> {320, 320}], Plot1],
  {{xv, "x", "Position/Velocity"},
   {"x" -> "Position", "v" -> "Velocity"}},

(* Interval *)
  {{T, 1000, "Time"}, 1000, 10^9}, ControlPlacement -> {Bottom}]

(* by Simon Tyran, Yukterez 2013. Inspired by Michael Trott *)
POSTED BY: Simon Tyran
I am working directly off your notebook (which is different from the pasted code):
Import["http://yukterez.ist.org/3kp/3k3D.nb","NotebookObject"]

I think you need to do something like:
Manipulate[
Show[{
  Plot1,
  Funktion[..]
}],
...]
And also remove the PlotRange from Plot1 since it won't work well with the large values of your points (for example Au = 149597870690, which won't fit in the given plot range).

 
POSTED BY: Arnoud Buzing
Reply to this discussion
Community posts can be styled and formatted using the Markdown syntax.
Reply Preview
Attachments
Remove
or Discard

Group Abstract Group Abstract