A basic ODE where Runge Kutta doesn't work?

Click For Summary
The discussion centers on a separable differential equation with a cubic solution, where the user encounters unexpected complex results when applying various Runge Kutta methods. The original equation, dy/dt=-6(6/13)^(1/3)(y-343/468)^(2/3), leads to complex solutions for y(1) when using RK4, despite the expected real solution. It is revealed that the programs default to calculating the principal cube root of negative numbers, which causes these complex outputs. A suggested modification to the equation, by squaring the term before applying the cube root, resolves the issue by ensuring non-negative arguments. The discussion highlights the importance of understanding how numerical software handles fractional powers of negative numbers.
h1a8
Messages
86
Reaction score
4
From a cubic function where y(0)=1, y(1)=0, and where there is a local max at y(5/13) I created a basic separable differential equation problem. I wanted to analyze how well different ordered Runge Kutta methods works in an interval [0,1]. Here it is:

dy/dt=-6(6/13)1/3(y-343/468)2/3 , y(0)=1

This ODE yields the cubic solution of

y=1/468(-12t+5)3+343/468

Now it is clear that y(1)=0

But using the several Runge Kutta programs with various computer software (mathematica, ti-nspire cas, mathstudio, etc.) yields a complex solutions for y(1). For example, using the classical RK4 with h=.1 yields
y(1)=0.718779+.005811i.

I don't see how the programs get a complex solution when all the functions have no even roots. Does anyone what is going on with these programs?
 
Physics news on Phys.org
You have negative numbers in that root for y<343/468. All of the mentioned programs are won't to go to complex numbers ASAP. Those programs will take the principal root. For example, the principal root of (-1)2/3 is (-1+i√3)/2, not 1.
 
Last edited:
Yep.

Try dy/dt=-6 (6/13)1/3 ((y-343/468)(y-343/468))1/3 instead.

By multiplying (y-343/468) with itself first, you make sure that the argument to the power (1/3) is non-negative.
This way it will do what you intended.
 
I like Serena said:
Yep.

Try dy/dt=-6 (6/13)1/3 ((y-343/468)(y-343/468))1/3 instead.

By multiplying (y-343/468) with itself first, you make sure that the argument to the power (1/3) is non-negative.
This way it will do what you intended.

Thanks that works. I tried ((y-343/468)2)1/3 and that works too.
Didn't know the programs calculated the principle cube root of negatives by default even if they were complex.
 
Well, basically they don't have a choice.

What do you think (-1)^0.667 is?
 
I like Serena said:
Try dy/dt=-6 (6/13)1/3 ((y-343/468)(y-343/468))1/3 instead.
Nice trick.

h1a8, this works because of the identity (ab)c=abc.

This identity is valid if a is positive and b and c are real. It is not valid in general. (The same goes with a lot of other exponentiation identities.)
 
Last edited:

Similar threads

  • · Replies 65 ·
3
Replies
65
Views
7K
  • · Replies 7 ·
Replies
7
Views
3K
  • · Replies 4 ·
Replies
4
Views
2K
Replies
2
Views
2K
  • · Replies 3 ·
Replies
3
Views
4K
  • · Replies 31 ·
2
Replies
31
Views
8K
  • · Replies 1 ·
Replies
1
Views
6K
  • · Replies 15 ·
Replies
15
Views
3K
  • · Replies 1 ·
Replies
1
Views
2K
  • · Replies 14 ·
Replies
14
Views
3K