Groups | Search | Server Info | Keyboard shortcuts | Login | Register [http] [https] [nntp] [nntps]
Groups > comp.soft-sys.math.mathematica > #2103 > unrolled thread
| Started by | Daniel Lichtblau <danl@wolfram.com> |
|---|---|
| First post | 2011-05-04 23:47 +0000 |
| Last post | 2011-05-04 23:47 +0000 |
| Articles | 1 — 1 participant |
Back to article view | Back to comp.soft-sys.math.mathematica
Re: envelope of curves Daniel Lichtblau <danl@wolfram.com> - 2011-05-04 23:47 +0000
| From | Daniel Lichtblau <danl@wolfram.com> |
|---|---|
| Date | 2011-05-04 23:47 +0000 |
| Subject | Re: envelope of curves |
| Message-ID | <ipsoj7$8be$1@smc.vnet.net> |
Dan wrote:
> I would like to find the envelope of a set of curves. This evelope may
> be formed by moving the center of one curve along the path of the
> other curve. It doesn't matter which curve is held fixed and which is
> translated as can be seen in the two plots. I am pretty sure there is
> a fairly straightforward algorithm for this, but it is escaping me at
> the moment. One of the curves, drawn here in grey, is always an
> ellipse with the major and minor axes aligned with x and y. The other
> curve, drawn in red, may be like an ellipse but it may be some other
> shape. In the actual problem it is a fairly involved calculation to
> determine it, but it is always smooth and convex. If there is a closed
> form solution for the elliptical case, that would be useful. On the
> other hand, a numerical solution that can be expressed as an
> interpolated curve would probably be the most general solution for me.
>
> With[{a=5,b=1,c=0.25,d=1.5,theta=-Pi/6,n=50},
> Module[{rot,c1,c2},
> rot={{Cos[theta],-Sin[theta]},{Sin[theta],Cos[theta]}};
> c1[psi_]:={a Sin[psi],b Cos[psi]};
> c2[psi_]:=rot.{c Sin[psi],d Cos[psi]};
> Column[{
> Show[{
> ParametricPlot[c1[psi],{psi,-Pi,Pi},Axes->False,
> PlotStyle->Gray,PlotRange->All],
> ParametricPlot[c2[psi],{psi,-Pi,Pi},Axes->False,
> PlotStyle->Red],
> Table[
> ParametricPlot[c1[phi]+c2[psi],{psi,-Pi,Pi},
> Axes->False,
> PlotStyle->Directive[Red,Opacity[0.2]]
> ],{phi,-Pi,Pi,2Pi/n}
> ]
> },ImageSize->400],
> Show[{
> ParametricPlot[c1[psi],{psi,-Pi,Pi},Axes->False,
> PlotStyle->Gray,PlotRange->All],
> ParametricPlot[c2[psi],{psi,-Pi,Pi},Axes->False,
> PlotStyle->Red],
> Table[
> ParametricPlot[c1[psi]+c2[phi],{psi,-Pi,Pi},
> Axes->False,
> PlotStyle->Directive[Gray,Opacity[0.2]]
> ],{phi,-Pi,Pi,2Pi/n}
> ]
> },ImageSize->400]
> }]
> ]
> ]
>
> Thanks in advance.
> Dan
>
Not sure I understand what it is you want, but I'll try to explain what
is computed below. We take your family of parametric curves
c1[psi]+c2[phi], set up polynomials in {x,y} in terms of these, add trig
identities.
We then "extremize" via a Lagrange multiplier. This gives two new
polynomials and one new variable. We now "polynomialize" by making the
trigs into new variables (keep in mind we have identities to relate them
in pairs).
The next thing is to eliminate all variables other than {x,y}. In doing
so we obtain a polynomial. It is the implicit form of the extremal curve.
Here is the code I use for this.
grad[expr_, a_, b_] := {D[expr, a], D[expr, b]}
vals = {a -> 5, b -> 1, c -> 1/4, d -> 3/2, theta -> -Pi/6};
rot = {{Cos[theta], -Sin[theta]}, {Sin[theta], Cos[theta]}};
c1[psi_] := {a Sin[psi], b Cos[psi]};
c2[psi_] := rot.{c Sin[psi], d Cos[psi]};
coordpolys = {x, y} - (c1[phi] + c2[psi]);
reprules = {Cos[psi] -> cps, Sin[psi] -> sps, Cos[phi] -> cph,
Sin[phi] -> sph};
trigpolys = {Cos[psi]^2 + Sin[psi]^2 - 1, Cos[phi]^2 + Sin[phi]^2 - 1};
vars = {x, y};
gradpolys =
grad[coordpolys[[1]], phi, psi] -
lambda*grad[ coordpolys[[2]], phi, psi];
allpolys = Join[coords, trigpolys, gradpolys] /. vals /. reprules
elimvars = Complement[Variables[allpolys], vars];
(* Here is the polynomial system *)
Out[92]= {-((3*cps)/4) - 5*sph - (Sqrt[3]*sps)/8 +
x, -cph - (3*Sqrt[3]*cps)/4 +
sps/8 + y, -1 + cps^2 + sps^2, -1 + cph^2 + sph^2, -5*cph -
lambda*sph,
-((Sqrt[3]*cps)/8) + (3*sps)/4 - lambda*(cps/8 + (3*Sqrt[3]*sps)/4)}
(* Now compute the implicit polynomial. Can be done via GroebnerBasis,
eliminating all variables but {x,y}. *)
In[98]:= imp = First[GroebnerBasis[allpolys, vars, elimvars]]
Out[98]= -632169101025 - 15438252990*x^2 + 2327819519*x^4 +
50943632*x^6 -
3041536*x^8 + 178880770180*Sqrt[3]*x*y -
18989375300*Sqrt[3]*x^3*y +
365862560*Sqrt[3]*x^5*y + 3906560*Sqrt[3]*x^7*y +
6477954371014*y^2 -
580839892690*x^2*y^2 + 16449887888*x^4*y^2 - 158016512*x^6*y^2 +
328983054260*Sqrt[3]*x*y^3 - 17352995520*Sqrt[3]*x^3*y^3 +
196725760*Sqrt[3]*x^5*y^3 - 1341456141001*y^4 +
117801783152*x^2*y^4 -
2198334976*x^4*y^4 - 66752924000*Sqrt[3]*x*y^5 +
2511488000*Sqrt[3]*x^3*y^5 + 26694841200*y^6 - 3731788800*x^2*y^6 +
873600000*Sqrt[3]*x*y^7 - 243360000*y^8
This polynomial is irreducible (does not factor). We can plot it as below.
ct = ContourPlot[imp, {x, -8, 8}, {y, -4, 4}, Contours -> {0},
ContourShading -> None, PlotPoints -> 50]
It looks like an oblong curved rhombus, with a necktie inside (because
it is elegant, one would presume). The necktie is something of an
artifact, where the extremal equations are satisfied but all the same
other parts overlay it.
You can get an idea of how it all fits together as below.
n = 50;
plts = Table[
ParametricPlot[c1[phi] + c2[psi] /. vals, {psi, -Pi, Pi},
Axes -> False,
PlotStyle -> Directive[Red, Opacity[0.2]]], {phi, -Pi, Pi,
2 Pi/n}];
Show[{ct, plts}]
Can do similarly with the other family.
plts2 = Table[
ParametricPlot[c1[psi] + c2[phi] /. vals, {psi, -Pi, Pi},
Axes -> False,
PlotStyle -> Directive[Gray, Opacity[0.2]]], {phi, -Pi, Pi,
2 Pi/n}];
For more information about this method of computing an envelope, see
http://library.wolfram.com/infocenter/Conferences/7516/
Daniel Lichtblau
Wolfram Research
Back to top | Article view | comp.soft-sys.math.mathematica
csiph-web