How to convert a system of parametric equations to a normal equation?

For example, I have a system of parametric equations (R is a constant number) :

{ x == p +  R Cos[k],
y == Cos[p] +  R Sin[k],
k == ArcTan[ 1/Sin[p] ]   }


Now I am going to find the normal equation of $\; y = f(x)\;$ without other variables. (Theoretically it is possible because there are totally 4 variables and 3 equations.)

How can I do this in Mathematica?

-
Is px a variable, or is it p * x ? What about py ? –  Artes Jul 20 '12 at 13:59
@Artes they are variables, sorry for confusing. I have already eliminated py, and rename px to p. –  Skyler Jul 20 '12 at 14:03
But there are 5 variables in your 3 equations. Otherwise, perhaps you can use Eliminate. –  yulinlinyu Jul 20 '12 at 14:26
@yulinlinyu there are x,y,p,k 4 variables in the 3 equations. R is a constant number. I wonder how to do elimination in Mathematica. –  Skyler Jul 20 '12 at 14:35

First, for the sake of simplicity let's define eqs - the system of our interest :

eqs = { x == p + R Cos[k],
y == Cos[p] + R Sin[k],
k == ArcTan[ 1/Sin[p]]  };


Equations y = y(x)

For this system we can find an explicit equation $\;y = y(x)\;$ only assuming R == 0, otherwise we could find only implicit solutions.

Solve[ eqs /. R -> 0, y, {p}, MaxExtraConditions -> All]

{{y -> ConditionalExpression[Cos[x], k == ArcTan[Csc[x]]]}}


See e.g. Inverse function theorem for the general issue of invertibility, let's restrict to :

eqs[[1]] /. (Rule @@@ eqs)[[3]]

x == p + R/Sqrt[1 + Csc[p]^2]


When R is given, one could invert $\;x = x(p)\;$ in a neighbourhood of $p_{0}$ to get $\;p = p(x)\;$ if $\frac{dx}{dp}|_{p=p_{0}} \neq 0$, i.e. only for special ranges of p where $\frac{dx}{dp}$ doesn't vanish. While it would work well on the symbolic level using InverseFunction, but generically there is no way to get explicit formulae. We can animate dependence of the function $\;x = x(p)\;$ : p + R Cos[ArcTan[1/Sin[p]]] and its inverse function $\;p = p(x)\;$ on the running parameter -1 < R < 1 (then $x(p)$ is globally invertible). We add the symmetry axis p == x.

Animate[ Plot[ { p + R Cos[ArcTan[1/Sin[p]]],
p,
InverseFunction[ # + R Cos[ ArcTan[ 1/Sin[#]]] &]@p }, {p, -Pi, Pi},

PlotPoints -> 80, PlotRange -> {-3.3, 3.3},
PlotStyle -> {Thick, Dashed, Thick}],
{R, -1, 1} ]


Equations y = y(k)

We can get an equation for y in an implicit form using Eliminate and then solve the equation with respect to y :

Eliminate[ eqs, {x, p}, InverseFunctions -> True]
Reduce[%, y]

(-1 + y^2 - 2 R y Sin[k] + R^2 Sin[k]^2) Tan[k]^2 == -1

Tan[k] != 0 &&
(y == Cot[k]^2 (R Sin[k] Tan[k]^2 - Sqrt[-Tan[k]^2 + Tan[k]^4]) ||
y == Cot[k]^2 (R Sin[k] Tan[k]^2 + Sqrt[-Tan[k]^2 + Tan[k]^4])    )


We've got solutions in terms of k only. This is the best one can do for your system,
otherwise to get an equation for y in terms of x we should invert $\;k = k(x)\;$ to get $\;y = y(x)\;$

First @ Solve[ k == ArcTan[ Csc[ p]], p]

{p -> ArcCsc[Tan[k]]}


We could do that analogically as in the case $\;p = p(x)\;$, i.e using InverseFunction.

Equations y = y(p)

You could as well try another way of eliminating variables while solving equations (as pointed out by Searke) i.e. Solve[ eqns, vars, elims], and using MaxExtraConditions -> All to get also non-generic solutions, here we eliminate k to get a solution in terms of p :

Solve[ eqs, y, k, MaxExtraConditions -> All]

{{y -> ConditionalExpression[(R Csc[p] + Cos[p] Sqrt[1 + Csc[p]^2])/ Sqrt[1 + Csc[p]^2],
Sqrt[1 + Csc[p]^2] != 0 &&
p - x + R/Sqrt[1 + Csc[p]^2] == 0]}}


The second argument in ConditionalExpression is the condition which has to be satisfied.

Nevertheless you could do the following to express y in terms of p and x :

Solve[ eqs, {y, k}]

{{y -> -Sqrt[-p^2 + R^2 + 2 p x - x^2] + Cos[p], k -> ArcTan[Csc[p]]},
{y -> Sqrt[-p^2 + R^2 + 2 p x - x^2] + Cos[p], k -> ArcTan[Csc[p]]}}


or

 Reduce[eqs, y, Backsubstitution -> True]

x == (R + p Sqrt[1 + Csc[p]^2])/Sqrt[1 + Csc[p]^2] && k == ArcTan[Csc[p]] &&
y == (R Csc[p] + Cos[p] Sqrt[1 + Csc[p]^2])/Sqrt[1 + Csc[p]^2]


There is another approach to get x in terms of y and p :

Eliminate[eqs, {x, k}, InverseFunctions -> True] //
Simplify[#, x == p + R Cos[k]] &

Reduce[%, x]

R^2 + 2 y Cos[p] == (p - x)^2 + y^2 + Cos[p]^2

x == p - Sqrt[R^2 - y^2 + 2 y Cos[p] - Cos[p]^2] ||
x == p + Sqrt[R^2 - y^2 + 2 y Cos[p] - Cos[p]^2]


Above, in all cases we have to use ReplaceAll (i.e. /.) and InverseFunction wherever appear p or k to get $\;p = p(x)\;$ or $\;k = k(x)\;$.

-

The Mathematica documentation has a tutorial on eliminating variables that is helpful:

http://reference.wolfram.com/mathematica/tutorial/EliminatingVariables.html

Solve[{x == p + R*Cos[k], y == Cos[p] + R*Sin[k],
k == ArcTan[1/Sin[p]]}, y, {p}]

y -> Cot[k]^2 (R Sin[k] Tan[k]^2 - Sqrt[-Tan[k]^2 + Tan[k]^4])


Solving for y in terms of just x turns out to be more difficult because it is difficult to expression k in terms of x:

Solve[{x == p + R*Cos[k], y == Cos[p] + R*Sin[k],
k == ArcTan[1/Sin[p]]}, x, {y, p}]

{x -> ArcCsc[Tan[k]] + R Cos[k]}


The equation above is difficult to invert. There isn't a way to solve for k in terms of x except to invent a function.

Solve[{x == ArcCsc[Tan[k]] + R Cos[k]}, k]
{{}}


Also, there's a ton of inverse trig functions everywhere. Issues with branch cuts are likely.

-