7
$\begingroup$

How to reverse formula

$y(x)=x (\frac{1}{sin \frac{\pi x}{2}})^\alpha$

i.e. express it as

$x = x(y)$

in Mathematica?

I did this way

Clear[α, x, y]
Solve[y == x (1/Sin[Pi x/2])^α, x]

and it answered

During evaluation of In[39]:= Solve::nsmet: This system cannot be solved with the methods available to Solve. >>

Out[40]= Solve[y == x Csc[(π x)/2]^α, x]

How to know, what prevent equation from solvation?

I have plotted formula for some values

Plot[Table[
  x (1/Sin[Pi x/2])^α, {α, {-1, -0.5, 0, 0.5, 1}}], {x, 
  0, 1}]

and found nothing criminal in my expected domain

enter image description here

I have tried to add domain to the formula, but it didn't gave me answer anyway

Solve[y == x (1/Sin[Pi x/2])^α, x, 
 x >= 0 && x <= 1 && α >= -1 && α <= 1]

During evaluation of In[41]:= Solve::bdomv: Warning: x>=0&&x<=1&&α>=-1&&α<=1 is not a valid domain specification. Mathematica is assuming it is a variable to eliminate. >>

During evaluation of In[41]:= Solve::ivar: x>=0&&x<=1&&α>=-1&&α<=1 is not a valid variable. >>

Out[41]= Solve[y == x Csc[(π x)/2]^α, x, 
 x >= 0 && x <= 1 && α >= -1 && α <= 1]

How to set conditions correctly?

How to force to give computable result like as series?

UPDATE

Using InverseFunction I wrote:

MyFun2[x_, α_] := x*(1/Sin[Pi x/2])^α
MyFunInverse2 = InverseFunction[MyFun2, 1, 2]

which gave apparently good result:

Plot[Table[MyFunInverse2[y, α], {α, -1, 1, 0.25}], {y, 
  0, 1}]

enter image description here

but trying to produce computable expression failed:

Series[MyFunInverse2[y, α], {y, 0, 5}]

Out[73]= SeriesData[y, 0, {
InverseFunction[MyFun2, 1, 2][0, α], 
Derivative[1, 0][
InverseFunction[MyFun2, 1, 2]][0, α], 
  Rational[1, 2] Derivative[2, 0][
InverseFunction[MyFun2, 1, 2]][0, α], 
  Rational[1, 6] Derivative[3, 0][
InverseFunction[MyFun2, 1, 2]][0, α], 
  Rational[1, 24] Derivative[4, 0][
InverseFunction[MyFun2, 1, 2]][0, α], 
  Rational[1, 120] Derivative[5, 0][
InverseFunction[MyFun2, 1, 2]][0, α]}, 0, 6, 1]

UPDATE 2

I need to port function to C-like programming language, so final result should consist of "computable" operations. Is it possible?

How does Mathematica itself plots the result? Does just perform optimization algorithm inside? Is it really impossible to do otherwise, for example do with multivariate series?

UPDATE 3

If I take 3 terms of a series

Fwd[x_, \[Alpha]_] := x (1/Sin[Pi x/2])^\[Alpha]
In[115]:= Normal[Series[Fwd[x, \[Alpha]], {x, 0, 3}]]
Fwd1 = Function[{x, \[Alpha]}, %]

I see that result is quite applicable for me

enter image description here

nevertheless I can't reverse this (now much simpler) formula with neither Solve, InverseFunction and InverseSeries.

Isn't it really possible to get some short numeric series with 2-3-4 terms in the case???

$\endgroup$
6
  • 1
    $\begingroup$ Some (perhaps it's better to say most) expressions cannot be inverted analytically no matter how smooth they are. You can define a numerical inverse using InverseFunction if all parameters are known. $\endgroup$ Apr 22, 2015 at 15:34
  • 1
    $\begingroup$ "Invert" is a better term than "reverse" here. $\endgroup$
    – Szabolcs
    Apr 22, 2015 at 15:43
  • $\begingroup$ @SjoerdC.deVries so how to invert with InverseFunction having that $\alpha$ is an argument? $\endgroup$
    – Suzan Cioc
    Apr 22, 2015 at 15:46
  • 1
    $\begingroup$ If all you need is a power series for the inverse function, you can calculate the power series for $x(y)$ and then use InverseSeries to find the power series for $y(x)$. Note that $\alpha$ has to be an explicit rational number for this to work, though. $\endgroup$ Apr 22, 2015 at 15:48
  • $\begingroup$ What range of $\alpha$ are you expecting? $\endgroup$ Apr 22, 2015 at 18:25

2 Answers 2

5
$\begingroup$

If all you're interested in is the inverse power series, then don't calculate the InverseFunction; instead, use the InverseSeries function:

fseries[α_] := Series[x (1/Sin[Pi x/2])^α, {x, 0, 5}]
Table[InverseSeries[fseries[α], y], {α, -1, 1, 1/2}]

enter image description here

Note that $\alpha$ has to be a rational number for this to work. In particular, if you replace the 1/2 with 0.5 in the Table command above, it doesn't work:

fseries[α_] := Series[x (1/Sin[Pi x/2])^α, {x, 0, 5}]
Table[InverseSeries[fseries[α], y], {α, -1, 1, 0.5}]

enter image description here

$\endgroup$
2
  • $\begingroup$ So, how to solve for variable alpha? $\endgroup$
    – Suzan Cioc
    Apr 22, 2015 at 17:06
  • $\begingroup$ InverseSeries[fseries[1], y] gives a correct result under v10.2, but under 11.3 I get strange result: ComposeSeries[SeriesData[ SeriesDump`z$65839, 0, {1, 0, Rational[-7, 480] Pi^2}, 1, 5, 1]... $\endgroup$ Jun 23, 2018 at 16:56
5
$\begingroup$

You can get very close to the solution in three iterations of Newton's method:

f[x_, a_] := x/Sin[Pi/2 x]^a

fInvNewton[y_?NumericQ, a_?NumericQ, iter_: 3] :=
 Nest[# - (f[#, a] - y)/Derivative[1, 0][f][#, a] &, 
  If[y > 2/Pi, y, ((Pi/2)^a y)^(1/(1 - a))], iter]

Example:

Plot[Evaluate@{fInvNewton[y, 0.5, 3], 
   InverseFunction[f, 1, 2][y, 0.5]}, {y, 0, 1}]

enter image description here

LogPlot[Evaluate@
  Abs[fInvNewton[y, 0.5, 3] - InverseFunction[f, 1, 2][y, 0.5]], {y, 
  0, 1}]

enter image description here

Here's some pseudocode for the method:

def fInverseNewton(y, a):
  # choose initial value
  if y > math.pi/2
    x = y
  else:
    x = ((math.pi/2)**a * y)**(1/(1-a))
  # three iterations of netwton's method
  for i in range(3):
    p = math.pi/2*x
    c = math.cos(p)
    s = math.sin(p)
    x = (a*x*p*c - y*s**(a+1)) / (a*p*c - s)
  return x
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

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