
12

Hi. I hope that I'm not embarrassing myself by asking the following.
I believe that
3 2/3 1
y =  (2 x) + 
4 2
is an even function, yet
x = 1:0.04:1;
plot (x, 0.75 * (2 .* x).^(2 / 3) + 0.5)
is not symmetrical about the yaxis. What is wrong with my thinking?
TIA.
John.

Octave is freely available under the terms of the GNU GPL.
Octave's home on the web: http://www.octave.orgHow to fund new projects: http://www.octave.org/funding.htmlSubscription information: http://www.octave.org/archive.html


John,
Could this be it?
x^(2/3) is equal to (x^(1/3))^2 which isn't even.
x^(2/3) is not equal to (x^2)^(1/3) which is even.
Steve
On Mar 19 16:20PM, John B. Thoo wrote:
> Hi. I hope that I'm not embarrassing myself by asking the following.
>
> I believe that
>
> 3 2/3 1
> y =  (2 x) + 
> 4 2
>
> is an even function, yet
>
> x = 1:0.04:1;
> plot (x, 0.75 * (2 .* x).^(2 / 3) + 0.5)
>
> is not symmetrical about the yaxis. What is wrong with my thinking?
>
> TIA.
> John.
>
>
>
> 
> Octave is freely available under the terms of the GNU GPL.
>
> Octave's home on the web: http://www.octave.org> How to fund new projects: http://www.octave.org/funding.html> Subscription information: http://www.octave.org/archive.html> 
>

Octave is freely available under the terms of the GNU GPL.
Octave's home on the web: http://www.octave.orgHow to fund new projects: http://www.octave.org/funding.htmlSubscription information: http://www.octave.org/archive.html


I would have thought that y = x^b is only symmetrical about yaxis if b = 0,
2, 4, etc, i.e. b = even integer.
Henry
on 3/19/05 4:20 PM, John B. Thoo at [hidden email] wrote:
> Hi. I hope that I'm not embarrassing myself by asking the following.
>
> I believe that
>
> 3 2/3 1
> y =  (2 x) + 
> 4 2
>
> is an even function, yet
>
> x = 1:0.04:1;
> plot (x, 0.75 * (2 .* x).^(2 / 3) + 0.5)
>
> is not symmetrical about the yaxis. What is wrong with my thinking?
>
> TIA.
> John.
>
>
>
> 
> Octave is freely available under the terms of the GNU GPL.
>
> Octave's home on the web: http://www.octave.org> How to fund new projects: http://www.octave.org/funding.html> Subscription information: http://www.octave.org/archive.html> 
>

Octave is freely available under the terms of the GNU GPL.
Octave's home on the web: http://www.octave.orgHow to fund new projects: http://www.octave.org/funding.htmlSubscription information: http://www.octave.org/archive.html


BEGIN PGP SIGNED MESSAGE
Hash: SHA1
John B. Thoo wrote:
 Hi. I hope that I'm not embarrassing myself by asking the following.

 I believe that

 3 2/3 1
 y =  (2 x) + 
 4 2

 is an even function, yet

 x = 1:0.04:1;
 plot (x, 0.75 * (2 .* x).^(2 / 3) + 0.5)

 is not symmetrical about the yaxis. What is wrong with my thinking?

 TIA.
 John.
x^(2/3) is not an even function,
octave:1> [ (+8)^(2/3) ; (8)^(2/3) ]
ans =
~ 4.0000 + 0.0000i
~ 2.0000 + 3.4641i
so nor is the function which you are plotting.
 
Geraint Bevan
http://homepage.ntlworld.com/geraint.bevanBEGIN PGP SIGNATURE
Version: GnuPG v1.2.4 (GNU/Linux)
iEYEARECAAYFAkI8z/sACgkQcXV3N50QmNOi1QCdEvZ44CO0bbl0h+TQvgFpD9uy
YQgAn2a4A1x4m7ZdST8vG9IBRgBu1lCW
=7F+W
END PGP SIGNATURE

Octave is freely available under the terms of the GNU GPL.
Octave's home on the web: http://www.octave.orgHow to fund new projects: http://www.octave.org/funding.htmlSubscription information: http://www.octave.org/archive.html


Hmmm...
Well, actually the function really *is* even if you think of it
this way:
y = 3/4*((2*x)^2)^3 + 1/2
Here's the problem: how do you know that the exponent is a
rational fraction, so that you can reinterpret it? Octave
rightfully doesn't. It evidently uses complex analysis to find
nonintegral powers. Try this
octave:16> x = (1)^(1/3)
x = 0.50000 + 0.86603i
octave:17> x^3
ans = 1.0000e+00 + 3.6370e16i
So octave calculates a complex third root of unity, rather than
the real root x = 1 that might spring to mind. Makes sense,
since you can generate the other roots from the complex
primitive root of unity.
Tom Shores
On Sunday 20 March 2005 01:20 am, Geraint Paul Bevan wrote:
> BEGIN PGP SIGNED MESSAGE
> Hash: SHA1
>
> John B. Thoo wrote:
>  Hi. I hope that I'm not embarrassing myself by asking the
>  following.
> 
>  I believe that
> 
>  3 2/3 1
>  y =  (2 x) + 
>  4 2
> 
>  is an even function, yet
> 
>  x = 1:0.04:1;
>  plot (x, 0.75 * (2 .* x).^(2 / 3) + 0.5)
> 
>  is not symmetrical about the yaxis. What is wrong with my
>  thinking?
> 
>  TIA.
>  John.
>
> x^(2/3) is not an even function,
>
> octave:1> [ (+8)^(2/3) ; (8)^(2/3) ]
> ans =
>
> ~ 4.0000 + 0.0000i
> ~ 2.0000 + 3.4641i
>
> so nor is the function which you are plotting.
>
>  
> Geraint Bevan
> http://homepage.ntlworld.com/geraint.bevan>
> BEGIN PGP SIGNATURE
> Version: GnuPG v1.2.4 (GNU/Linux)
>
> iEYEARECAAYFAkI8z/sACgkQcXV3N50QmNOi1QCdEvZ44CO0bbl0h+TQvgFpD
>9uy YQgAn2a4A1x4m7ZdST8vG9IBRgBu1lCW
> =7F+W
> END PGP SIGNATURE
>
>
>
> 
> Octave is freely available under the terms of the GNU GPL.
>
> Octave's home on the web: http://www.octave.org> How to fund new projects: http://www.octave.org/funding.html> Subscription information: http://www.octave.org/archive.html> 

Octave is freely available under the terms of the GNU GPL.
Octave's home on the web: http://www.octave.orgHow to fund new projects: http://www.octave.org/funding.htmlSubscription information: http://www.octave.org/archive.html


Postscript to my earlier comment: try this
octave:44> x=1:.04:1;
octave:45> y = 3/4*(2*x).^(2/3)+1/2;
octave:46> plot(x,y)
octave:47> hold on
octave:48> plot(x,real(y))
octave:49> plot(x,imag(y))
and it will be clear what you're seeing.
Tom Shores
On Sunday 20 March 2005 01:20 am, Geraint Paul Bevan wrote:
> BEGIN PGP SIGNED MESSAGE
> Hash: SHA1
>
> John B. Thoo wrote:
>  Hi. I hope that I'm not embarrassing myself by asking the
>  following.
> 
>  I believe that
> 
>  3 2/3 1
>  y =  (2 x) + 
>  4 2
> 
>  is an even function, yet
> 
>  x = 1:0.04:1;
>  plot (x, 0.75 * (2 .* x).^(2 / 3) + 0.5)
> 
>  is not symmetrical about the yaxis. What is wrong with my
>  thinking?
> 
>  TIA.
>  John.
>
> x^(2/3) is not an even function,
>
> octave:1> [ (+8)^(2/3) ; (8)^(2/3) ]
> ans =
>
> ~ 4.0000 + 0.0000i
> ~ 2.0000 + 3.4641i
>
> so nor is the function which you are plotting.
>
>  
> Geraint Bevan
> http://homepage.ntlworld.com/geraint.bevan>
> BEGIN PGP SIGNATURE
> Version: GnuPG v1.2.4 (GNU/Linux)
>
> iEYEARECAAYFAkI8z/sACgkQcXV3N50QmNOi1QCdEvZ44CO0bbl0h+TQvgFpD
>9uy YQgAn2a4A1x4m7ZdST8vG9IBRgBu1lCW
> =7F+W
> END PGP SIGNATURE
>
>
>
> 
> Octave is freely available under the terms of the GNU GPL.
>
> Octave's home on the web: http://www.octave.org> How to fund new projects: http://www.octave.org/funding.html> Subscription information: http://www.octave.org/archive.html> 

Octave is freely available under the terms of the GNU GPL.
Octave's home on the web: http://www.octave.orgHow to fund new projects: http://www.octave.org/funding.htmlSubscription information: http://www.octave.org/archive.html


Umm... make that
y = 3/4*((2*x)^2)^(1/3) + 1/2
So add to my last plotting suggestions
octave:5> y = 3/4*((2*x).^2).^(1/3)+1/2;
octave:6> plot(x,y)
Tom Shores
On Saturday 19 March 2005 08:36 pm, Thomas Shores wrote:
> Hmmm...
> Well, actually the function really *is* even if you think of
> it this way:
>
> y = 3/4*((2*x)^2)^3 + 1/2
>
> Here's the problem: how do you know that the exponent is a
> rational fraction, so that you can reinterpret it? Octave
> rightfully doesn't. It evidently uses complex analysis to
> find nonintegral powers. Try this
>
> octave:16> x = (1)^(1/3)
> x = 0.50000 + 0.86603i
> octave:17> x^3
> ans = 1.0000e+00 + 3.6370e16i
>
> So octave calculates a complex third root of unity, rather
> than the real root x = 1 that might spring to mind. Makes
> sense, since you can generate the other roots from the
> complex primitive root of unity.
>
> Tom Shores
>
> On Sunday 20 March 2005 01:20 am, Geraint Paul Bevan wrote:
> > BEGIN PGP SIGNED MESSAGE
> > Hash: SHA1
> >
> > John B. Thoo wrote:
> >  Hi. I hope that I'm not embarrassing myself by asking
> >  the following.
> > 
> >  I believe that
> > 
> >  3 2/3 1
> >  y =  (2 x) + 
> >  4 2
> > 
> >  is an even function, yet
> > 
> >  x = 1:0.04:1;
> >  plot (x, 0.75 * (2 .* x).^(2 / 3) + 0.5)
> > 
> >  is not symmetrical about the yaxis. What is wrong with
> >  my thinking?
> > 
> >  TIA.
> >  John.
> >
> > x^(2/3) is not an even function,
> >
> > octave:1> [ (+8)^(2/3) ; (8)^(2/3) ]
> > ans =
> >
> > ~ 4.0000 + 0.0000i
> > ~ 2.0000 + 3.4641i
> >
> > so nor is the function which you are plotting.
> >
> >  
> > Geraint Bevan
> > http://homepage.ntlworld.com/geraint.bevan> >
> > BEGIN PGP SIGNATURE
> > Version: GnuPG v1.2.4 (GNU/Linux)
> >
> > iEYEARECAAYFAkI8z/sACgkQcXV3N50QmNOi1QCdEvZ44CO0bbl0h+TQvgF
> >pD 9uy YQgAn2a4A1x4m7ZdST8vG9IBRgBu1lCW
> > =7F+W
> > END PGP SIGNATURE
> >
> >
> >
> > 
> > Octave is freely available under the terms of the GNU
> > GPL.
> >
> > Octave's home on the web: http://www.octave.org> > How to fund new projects:
> > http://www.octave.org/funding.html Subscription
> > information: http://www.octave.org/archive.html> > 
> >
>
> 
> Octave is freely available under the terms of the GNU GPL.
>
> Octave's home on the web: http://www.octave.org> How to fund new projects: http://www.octave.org/funding.html> Subscription information: http://www.octave.org/archive.html> 

Octave is freely available under the terms of the GNU GPL.
Octave's home on the web: http://www.octave.orgHow to fund new projects: http://www.octave.org/funding.htmlSubscription information: http://www.octave.org/archive.html


First, thanks to Steve T, Henry M, Geraint B, and Thomas S for their
replies.
Thomas, your plotting example
octave:5> y = 3/4*((2*x).^2).^(1/3)+1/2;
octave:6> plot(x,y)
gives what I expected.
I learnt something tonight. When I teach (high school) algebra, I tell
my students that for rational m/n in lowest terms,
a^(m/n) = (a^m)^(1/n) = ( a^(1/n) )^m
as long as a^(1/n) is real. And, for n odd, there is no reason at
this level to think that a^(1/n) is not real.
In my example, 2/3 is clearly(?) a rational exponent with n=3 odd,
so I never expected that Octave would treat a^(1/3) as complex.
Certainly food for thought for me.
So, thanks again to all for your help.
Oh, in case you're interested in how this came up for me, I was
plotting normals to the curve y = x^2,
> x = 1:0.04:1; t = 1:1:1;
> for i = 1:51
> plot (t, t ./ (2 * x (i)) + (1 + 2 * x (i).^2) ./ 2, "b")
> endfor
and the envelope of the normals (the evolute of the parabola?) is given
by the curve I had asked about,
y = (3/4) (2x)^(2/3) + 1/2.
Anyhow, thanks again.
John.
On Mar 19, 2005, at 12:48 PM, Thomas Shores wrote:
> Umm... make that
> y = 3/4*((2*x)^2)^(1/3) + 1/2
>
> So add to my last plotting suggestions
>
> octave:5> y = 3/4*((2*x).^2).^(1/3)+1/2;
> octave:6> plot(x,y)
>
> Tom Shores
>
> On Saturday 19 March 2005 08:36 pm, Thomas Shores wrote:
>> Hmmm...
>> Well, actually the function really *is* even if you think of
>> it this way:
>>
>> y = 3/4*((2*x)^2)^3 + 1/2
>>
>> Here's the problem: how do you know that the exponent is a
>> rational fraction, so that you can reinterpret it? Octave
>> rightfully doesn't. It evidently uses complex analysis to
>> find nonintegral powers. Try this
>>
>> octave:16> x = (1)^(1/3)
>> x = 0.50000 + 0.86603i
>> octave:17> x^3
>> ans = 1.0000e+00 + 3.6370e16i
>>
>> So octave calculates a complex third root of unity, rather
>> than the real root x = 1 that might spring to mind. Makes
>> sense, since you can generate the other roots from the
>> complex primitive root of unity.
>>
>> Tom Shores
>>
>> On Sunday 20 March 2005 01:20 am, Geraint Paul Bevan wrote:
>>> BEGIN PGP SIGNED MESSAGE
>>> Hash: SHA1
>>>
>>> John B. Thoo wrote:
>>>  Hi. I hope that I'm not embarrassing myself by asking
>>>  the following.
>>> 
>>>  I believe that
>>> 
>>>  3 2/3 1
>>>  y =  (2 x) + 
>>>  4 2
>>> 
>>>  is an even function, yet
>>> 
>>>  x = 1:0.04:1;
>>>  plot (x, 0.75 * (2 .* x).^(2 / 3) + 0.5)
>>> 
>>>  is not symmetrical about the yaxis. What is wrong with
>>>  my thinking?
>>> 
>>>  TIA.
>>>  John.
>>>
>>> x^(2/3) is not an even function,
>>>
>>> octave:1> [ (+8)^(2/3) ; (8)^(2/3) ]
>>> ans =
>>>
>>> ~ 4.0000 + 0.0000i
>>> ~ 2.0000 + 3.4641i
>>>
>>> so nor is the function which you are plotting.
>>>
>>>  
>>> Geraint Bevan
>>> http://homepage.ntlworld.com/geraint.bevan>>>
>>> BEGIN PGP SIGNATURE
>>> Version: GnuPG v1.2.4 (GNU/Linux)
>>>
>>> iEYEARECAAYFAkI8z/sACgkQcXV3N50QmNOi1QCdEvZ44CO0bbl0h+TQvgF
>>> pD 9uy YQgAn2a4A1x4m7ZdST8vG9IBRgBu1lCW
>>> =7F+W
>>> END PGP SIGNATURE
>>>
>>>
>>>
>>> 
>>>  Octave is freely available under the terms of the GNU
>>> GPL.
>>>
>>> Octave's home on the web: http://www.octave.org>>> How to fund new projects:
>>> http://www.octave.org/funding.html Subscription
>>> information: http://www.octave.org/archive.html>>> 
>>> 
>>
>> 
>> Octave is freely available under the terms of the GNU GPL.
>>
>> Octave's home on the web: http://www.octave.org>> How to fund new projects: http://www.octave.org/funding.html>> Subscription information: http://www.octave.org/archive.html>> 
>

Octave is freely available under the terms of the GNU GPL.
Octave's home on the web: http://www.octave.orgHow to fund new projects: http://www.octave.org/funding.htmlSubscription information: http://www.octave.org/archive.html


Hi, Steve. I'm sorry for dragging this out, but I'm not yet
understanding.
Why is
x^(2/3) equal to (x^(1/3))^2
but
x^(2/3) not equal to (x^2)^(1/3)?
Also, why is (x^(1/3))^2 not even? If f(x) = (x^(1/3))^2, then
isn't
f(x) = ((x)^(1/3))^2 = ( x^(1/3))^2 = (x^(1/3))^2 = f(x)
(for real x)? I must not be understanding something rudimentary.
TIA.
John.
(a numerical computation and Octave novice)
On Mar 19, 2005, at 5:03 PM, Steve C. Thompson wrote:
> John,
>
> Could this be it?
>
> x^(2/3) is equal to (x^(1/3))^2 which isn't even.
>
> x^(2/3) is not equal to (x^2)^(1/3) which is even.
>
>
> Steve
>
>
> On Mar 19 16:20PM, John B. Thoo wrote:
>> Hi. I hope that I'm not embarrassing myself by asking the following.
>>
>> I believe that
>>
>> 3 2/3 1
>> y =  (2 x) + 
>> 4 2
>>
>> is an even function, yet
>>
>> x = 1:0.04:1;
>> plot (x, 0.75 * (2 .* x).^(2 / 3) + 0.5)
>>
>> is not symmetrical about the yaxis. What is wrong with my thinking?
>>
>> TIA.
>> John.
>>
>>
>>
>> 
>> Octave is freely available under the terms of the GNU GPL.
>>
>> Octave's home on the web: http://www.octave.org>> How to fund new projects: http://www.octave.org/funding.html>> Subscription information: http://www.octave.org/archive.html>> 
>>
>

Octave is freely available under the terms of the GNU GPL.
Octave's home on the web: http://www.octave.orgHow to fund new projects: http://www.octave.org/funding.htmlSubscription information: http://www.octave.org/archive.html


BEGIN PGP SIGNED MESSAGE
Hash: SHA1
John B. Thoo wrote:
> Also, why is (x^(1/3))^2 not even? If f(x) = (x^(1/3))^2, then isn't
> f(x) = ((x)^(1/3))^2 = ( x^(1/3))^2 = (x^(1/3))^2 = f(x)
It is not true to say that (x)^(1/3) is the same as (x^(1/3)). That
is one possibility but there are two more. If n is an integer, x^n has
a unique value but x^(1/n) has n possible values i.e. the square root
x^(1/2) has two possible values, the cube root x^(1/3) has three
possible values, etc.
Consider that (x)^(1/n) can be rewritten as ((1)*(+x))^(1/n) which
is ((1)^(1/n))*((+x)^(1/n))
For symmetry about the yaxis, you therefore require that:
(1)^(1/n) = 1.
However, for n=3, y=(1)^(1/n) has three solutions:
y = 1
y = 0.5+i*sqrt(3)/2
y = 0.5i*sqrt(3)/2
You obviously want to use the first of these solutions to get symmetry
about the yaxis but Octave doesn't know that unless you tell it.
 
Geraint Bevan
http://homepage.ntlworld.com/geraint.bevanBEGIN PGP SIGNATURE
Version: GnuPG v1.2.4 (GNU/Linux)
iEYEARECAAYFAkI9rkcACgkQcXV3N50QmNM56QCeJl6XCh0ZnD+WgYz06N3HAhmK
qD8AnirasH4pDgipzksWCGLJlcUAGzoM
=Bt7P
END PGP SIGNATURE

Octave is freely available under the terms of the GNU GPL.
Octave's home on the web: http://www.octave.orgHow to fund new projects: http://www.octave.org/funding.htmlSubscription information: http://www.octave.org/archive.html


On Sun, 20 Mar 2005, Geraint Paul Bevan wrote:
> BEGIN PGP SIGNED MESSAGE
> Hash: SHA1
>
> John B. Thoo wrote:
>
>> Also, why is (x^(1/3))^2 not even? If f(x) = (x^(1/3))^2, then isn't
>> f(x) = ((x)^(1/3))^2 = ( x^(1/3))^2 = (x^(1/3))^2 = f(x)
>
>
> It is not true to say that (x)^(1/3) is the same as (x^(1/3)). That is
> one possibility but there are two more. If n is an integer, x^n has a
> unique value but x^(1/n) has n possible values i.e. the square root
> x^(1/2) has two possible values, the cube root x^(1/3) has three
> possible values, etc.
>
> Consider that (x)^(1/n) can be rewritten as ((1)*(+x))^(1/n) which is
> ((1)^(1/n))*((+x)^(1/n))
>
> For symmetry about the yaxis, you therefore require that:
> (1)^(1/n) = 1.
>
> However, for n=3, y=(1)^(1/n) has three solutions:
> y = 1
> y = 0.5+i*sqrt(3)/2
> y = 0.5i*sqrt(3)/2
>
> You obviously want to use the first of these solutions to get symmetry
> about the yaxis but Octave doesn't know that unless you tell it.
I am not a mathematician, but I believe that you are missing something.
The equation y^3 + 1 = 0 has three solutions, as noted above. The
equation y^3  1 = 0 also has three solutions:
y = 1
y = 0.5i*sqrt(3)/2
y = 0.5+i*sqrt(3)/2
How does Octave "know what we want" (1) in the latter case, but it does
not know what we want (1) in the former case? I think the answer has to
do with numerical analysis. The binary representation of 1/3 is never
exactly 1/3, but it is only at the exact value of 1/3 that (1)^(1/3)
equals 1. This is not a problem for (+1)^(1/3).
For (1)^(333/1000), Octave gives a very similar answer to what it gives
for (1)^(1/3), which seems appropriate. Octave doesn't seem to know or
care that there are 1000 possible solutions in the complex plane!
Mike

Octave is freely available under the terms of the GNU GPL.
Octave's home on the web: http://www.octave.orgHow to fund new projects: http://www.octave.org/funding.htmlSubscription information: http://www.octave.org/archive.html


BEGIN PGP SIGNED MESSAGE
Hash: SHA1
Mike Miller wrote:
> I am not a mathematician, but I believe that you are missing something.
> The equation y^3 + 1 = 0 has three solutions, as noted above. The
> equation y^3  1 = 0 also has three solutions:
>
> y = 1
> y = 0.5i*sqrt(3)/2
> y = 0.5+i*sqrt(3)/2
>
> How does Octave "know what we want" (1) in the latter case, but it does
> not know what we want (1) in the former case? I think the answer has
> to do with numerical analysis. The binary representation of 1/3 is never
> exactly 1/3, but it is only at the exact value of 1/3 that (1)^(1/3)
> equals 1. This is not a problem for (+1)^(1/3).
>
> For (1)^(333/1000), Octave gives a very similar answer to what it gives
> for (1)^(1/3), which seems appropriate. Octave doesn't seem to know or
> care that there are 1000 possible solutions in the complex plane!
>
> Mike
I believe that Octave uses the standard pow function within the system's
math library to perform this operation. The cause of the difference in
behaviour can be found in the Octave sources, xpow.cc:
octave_value
xpow (double a, double b)
{
if (a < 0.0 && static_cast<int> (b) != b)
{
// XXX FIXME XXX  avoid apparent GNU libm bug by converting
// A and B to complex instead of just A.
Complex atmp (a);
Complex btmp (b);
return pow (atmp, btmp);
}
else
return pow (a, b);
}
If the base is not negative, or the exponent is an integer, Octave uses
the standard double version of pow. This routine never finds a complex
solution  it either finds a real one or fails:
a.cc:
#include <cmath>
#include <iostream>
int main(void) {
std::cout << "+1^2.0\t" << pow(+1.0,2.0) << std::endl;
std::cout << "1^2.0\t" << pow(1.0,2.0) << std::endl;
std::cout << "+1^0.5\t" << pow(+1.0,0.5) << std::endl;
std::cout << "1^0.5\t" << pow(1.0,0.5) << std::endl;
return 0;
}
$ g++ a.cc && ./a.out
+1^2.0 1
 1^2.0 1
+1^0.5 1
 1^0.5 nan
As you can see from the results above, the double version of pow cannot
cope with a negative base and noninteger exponent. Therefore in this
case Octave uses the Complex (i.e. complex<double>) version of pow. It
is this routine which is finding complex solutions for negative x
instead of the real root that would be required for symmetry about the
yaxis.
 
Geraint Bevan
http://homepage.ntlworld.com/geraint.bevanBEGIN PGP SIGNATURE
Version: GnuPG v1.2.4 (GNU/Linux)
iEYEARECAAYFAkI9xogACgkQcXV3N50QmNPprACfZA5fdTDy2IHLef4XmUeyd+yW
4fkAniXoa1L3VlyCxC5JyiWNsM1UBjua
=eDGK
END PGP SIGNATURE

Octave is freely available under the terms of the GNU GPL.
Octave's home on the web: http://www.octave.orgHow to fund new projects: http://www.octave.org/funding.htmlSubscription information: http://www.octave.org/archive.html


On Sun, 20 Mar 2005, Geraint Paul Bevan wrote:
> $ g++ a.cc && ./a.out
> +1^2.0 1
>  1^2.0 1
> +1^0.5 1
>  1^0.5 nan
>
> As you can see from the results above, the double version of pow cannot
> cope with a negative base and noninteger exponent. Therefore in this
> case Octave uses the Complex (i.e. complex<double>) version of pow. It
> is this routine which is finding complex solutions for negative x
> instead of the real root that would be required for symmetry about the
> yaxis.
That's a very complete answer! Thanks.
Mike

Octave is freely available under the terms of the GNU GPL.
Octave's home on the web: http://www.octave.orgHow to fund new projects: http://www.octave.org/funding.htmlSubscription information: http://www.octave.org/archive.html


On Sat, Mar 19, 2005 at 11:19:21PM 0800, John B. Thoo wrote:
> First, thanks to Steve T, Henry M, Geraint B, and Thomas S for their
> replies.
>
> Thomas, your plotting example
>
> octave:5> y = 3/4*((2*x).^2).^(1/3)+1/2;
> octave:6> plot(x,y)
>
> gives what I expected.
>
> I learnt something tonight. When I teach (high school) algebra, I tell
> my students that for rational m/n in lowest terms,
>
> a^(m/n) = (a^m)^(1/n) = ( a^(1/n) )^m
>
> as long as a^(1/n) is real. And, for n odd, there is no reason at
> this level to think that a^(1/n) is not real.
>
> In my example, 2/3 is clearly(?) a rational exponent with n=3 odd,
> so I never expected that Octave would treat a^(1/3) as complex.
> Certainly food for thought for me.
>
> So, thanks again to all for your help.
>
> Oh, in case you're interested in how this came up for me, I was
> plotting normals to the curve y = x^2,
>
> > x = 1:0.04:1; t = 1:1:1;
> > for i = 1:51
> > plot (t, t ./ (2 * x (i)) + (1 + 2 * x (i).^2) ./ 2, "b")
> > endfor
>
> and the envelope of the normals (the evolute of the parabola?) is given
> by the curve I had asked about,
>
> y = (3/4) (2x)^(2/3) + 1/2.
>
> Anyhow, thanks again.
>
> John.
>
>
> On Mar 19, 2005, at 12:48 PM, Thomas Shores wrote:
>
> >Umm... make that
> >y = 3/4*((2*x)^2)^(1/3) + 1/2
> >
> >So add to my last plotting suggestions
> >
> >octave:5> y = 3/4*((2*x).^2).^(1/3)+1/2;
> >octave:6> plot(x,y)
> >
> >Tom Shores
> >
> >On Saturday 19 March 2005 08:36 pm, Thomas Shores wrote:
> >>Hmmm...
> >>Well, actually the function really *is* even if you think of
> >>it this way:
> >>
> >>y = 3/4*((2*x)^2)^3 + 1/2
> >>
> >>Here's the problem: how do you know that the exponent is a
> >>rational fraction, so that you can reinterpret it? Octave
> >>rightfully doesn't. It evidently uses complex analysis to
> >>find nonintegral powers. Try this
> >>
> >>octave:16> x = (1)^(1/3)
> >>x = 0.50000 + 0.86603i
> >>octave:17> x^3
> >>ans = 1.0000e+00 + 3.6370e16i
> >>
> >>So octave calculates a complex third root of unity, rather
> >>than the real root x = 1 that might spring to mind. Makes
> >>sense, since you can generate the other roots from the
> >>complex primitive root of unity.
> >>
> >>Tom Shores
> >>
> >>On Sunday 20 March 2005 01:20 am, Geraint Paul Bevan wrote:
> >>>BEGIN PGP SIGNED MESSAGE
> >>>Hash: SHA1
> >>>
> >>>John B. Thoo wrote:
> >>> Hi. I hope that I'm not embarrassing myself by asking
> >>> the following.
> >>>
> >>> I believe that
> >>>
> >>> 3 2/3 1
> >>> y =  (2 x) + 
> >>> 4 2
> >>>
> >>> is an even function, yet
> >>>
> >>> x = 1:0.04:1;
> >>> plot (x, 0.75 * (2 .* x).^(2 / 3) + 0.5)
> >>>
> >>> is not symmetrical about the yaxis. What is wrong with
> >>> my thinking?
> >>>
> >>> TIA.
> >>> John.
> >>>
> >>>x^(2/3) is not an even function,
> >>>
> >>>octave:1> [ (+8)^(2/3) ; (8)^(2/3) ]
> >>>ans =
> >>>
> >>>~ 4.0000 + 0.0000i
> >>>~ 2.0000 + 3.4641i
> >>>
> >>>so nor is the function which you are plotting.
> >>>
> >>> 
> >>>Geraint Bevan
> >>> http://homepage.ntlworld.com/geraint.bevan> >>>
> >>>BEGIN PGP SIGNATURE
> >>>Version: GnuPG v1.2.4 (GNU/Linux)
> >>>
> >>>iEYEARECAAYFAkI8z/sACgkQcXV3N50QmNOi1QCdEvZ44CO0bbl0h+TQvgF
> >>>pD 9uy YQgAn2a4A1x4m7ZdST8vG9IBRgBu1lCW
> >>>=7F+W
> >>>END PGP SIGNATURE
> >>>
> >>>
> >>>
> >>>
> >>> Octave is freely available under the terms of the GNU
> >>>GPL.
> >>>
> >>>Octave's home on the web: http://www.octave.org> >>>How to fund new projects:
> >>> http://www.octave.org/funding.html Subscription
> >>>information: http://www.octave.org/archive.html> >>>
> >>>
> >>
> >>
> >>Octave is freely available under the terms of the GNU GPL.
> >>
> >>Octave's home on the web: http://www.octave.org> >>How to fund new projects: http://www.octave.org/funding.html> >>Subscription information: http://www.octave.org/archive.html> >>
> >
>
>
>
> 
> Octave is freely available under the terms of the GNU GPL.
>
> Octave's home on the web: http://www.octave.org> How to fund new projects: http://www.octave.org/funding.html> Subscription information: http://www.octave.org/archive.html> 
> Hi,
OK, if I understood correctly the very interesting discussion, in
principle (a^b)^c==(a^c)^b (a is a real scalar) if, and only if
b*c is a rational number.
As octave treats evry inputed value as a double defined to some
finite precision, and as within this uncertainty there are both
numbers commesurate with each other and noncommensurate ones,
it has to make some arbitrary decision.
Now consider the folowing simple example:
588~> x=[10:10];
588~> y=x.^(1/3);
589~> plot(x,y) % the plot is asymetric with respect to the
origin
590~> floor(1)
ans = 1
591~> floor(3)
ans = 3
592~> y=x.^(floor(1)/floor(3));
589~> plot(x,y) % I get the same plot, although floor(1) and
floor(3) are integers.
It's a pity, because that mighty have been a way to obtain from
the software the expected behaviour.
Cheers and thanks to all hte participants for the pleasant hour I
spent recalling long forgotten basic math, Avraham

Octave is freely available under the terms of the GNU GPL.
Octave's home on the web: http://www.octave.orgHow to fund new projects: http://www.octave.org/funding.htmlSubscription information: http://www.octave.org/archive.html


First, thanks to Geraint B, Mike M, Avraham for their replies.
On Mar 20, 2005, at 9:09 AM, Geraint Paul Bevan wrote:
> For symmetry about the yaxis, you therefore require that:
> (1)^(1/n) = 1.
>
> However, for n=3, y=(1)^(1/n) has three solutions:
> y = 1
> y = 0.5+i*sqrt(3)/2
> y = 0.5i*sqrt(3)/2
>
> You obviously want to use the first of these solutions to get symmetry
> about the yaxis but Octave doesn't know that unless you tell it.
So, am I now right then in thinking that if z = r * (cos t + i * sin
t), then Octave will return z^(1/n) = r^(1/n) * (cos (t/n) + i * sin
(t/n))?
And on Mar 20, 2005, at 10:52 AM, Geraint Paul Bevan added:
> $ g++ a.cc && ./a.out
> +1^2.0 1
>  1^2.0 1
> +1^0.5 1
>  1^0.5 nan
>
> As you can see from the results above, the double version of pow cannot
> cope with a negative base and noninteger exponent. Therefore in this
> case Octave uses the Complex (i.e. complex<double>) version of pow. It
> is this routine which is finding complex solutions for negative x
> instead of the real root that would be required for symmetry about the
> yaxis.
I think I have a better understanding now. Of course, I know that I
want the real root :O but it never occurred to me that Octave
wouldn't use the real root if there is one. Duh!
So, if you would allow me at least one more daft question, how do I
tell Octave to use the real root so that the plot of y = x^(1/3) is
symmetrical about the origin?
TIA.
John.

Octave is freely available under the terms of the GNU GPL.
Octave's home on the web: http://www.octave.orgHow to fund new projects: http://www.octave.org/funding.htmlSubscription information: http://www.octave.org/archive.html


On Mar 20, 2005, at 1:25 PM, John B. Thoo wrote:
> So, if you would allow me at least one more daft question, how do I
> tell Octave to use the real root so that the plot of y = x^(1/3) is
> symmetrical about the origin?
OK, this is how I've done it.
> x1 = 1:0.1:0; x2 = 0:0.1:1;
> plot (x1, (x1).^(1/3), x2, x2.^(1/3)) % odd symmetry
> plot (x1, (x1).^(2/3), x2, x2.^(2/3)) % even symmetry
But is there a more elegant way?
TIA again.
John.

Octave is freely available under the terms of the GNU GPL.
Octave's home on the web: http://www.octave.orgHow to fund new projects: http://www.octave.org/funding.htmlSubscription information: http://www.octave.org/archive.html


BEGIN PGP SIGNED MESSAGE
Hash: SHA1
[hidden email] wrote:
> OK, if I understood correctly the very interesting discussion, in
> principle (a^b)^c==(a^c)^b (a is a real scalar) if, and only if
> b*c is a rational number.
This is not correct. Let b=pi (which is not rational) and the identity
will still hold true.
>From the laws of exponentiation, (a^b)^c is the same as a^(b*c).
Similarly (a^c)^b is the same as a^(c*b). As long as you are dealing
with real numbers rather than, say, matrices, multiplication is
commutative, i.e. b*c=c*b. Hence the identity holds for all real a, b
and c, rational or not.
> 592~> y=x.^(floor(1)/floor(3));
> 589~> plot(x,y) % I get the same plot, although floor(1) and
> floor(3) are integers.
I could be wrong, but I don't think that the floor() function is
necessary to get Octave to treat integers as such:
octave> format long
octave> 1.1
ans = 1.10000000000000
octave> 1.0
ans = 1
octave> floor(1)
ans = 1
Nevertheless, this could not work anyway. Although 1 and 3 would be
integers, the division would have to be performed before the
exponentiation and there is no way that 1/3 can be represented as an
integer, at least not without catastrophic loss of precision!
 
Geraint Bevan
http://homepage.ntlworld.com/geraint.bevanBEGIN PGP SIGNATURE
Version: GnuPG v1.2.4 (GNU/Linux)
iEYEARECAAYFAkI97wYACgkQcXV3N50QmNOSGgCaA/hM+xTy3CMpQai90CmRz5mR
dvUAoIDQtDJ7/3FVCktVFqvIFbvp+3tt
=oB8Y
END PGP SIGNATURE

Octave is freely available under the terms of the GNU GPL.
Octave's home on the web: http://www.octave.orgHow to fund new projects: http://www.octave.org/funding.htmlSubscription information: http://www.octave.org/archive.html


On Sun, 20 Mar 2005 [hidden email] wrote:
> OK, if I understood correctly the very interesting discussion, in
> principle (a^b)^c==(a^c)^b (a is a real scalar) if, and only if
> b*c is a rational number.
I think, mathematically speaking, (a^b)^c==(a^c)^b when b and c are real,
not just for rational b*c. Computationally, however, the two may differ
greatly depending on choice of algorithms. As noted earlier, 1/3 cannot
be represented as a finite binary expansion without some imprecision.
This means that the order matters: ((1)^2) = 1 and ((1)^2)^(1/3) = 1,
approximately, but (1)^(1/3) is not 1, and it is complex, so
((1)^(1/3))^2 is also complex, and not equal to 1.
> As octave treats evry inputed value as a double defined to some finite
> precision, and as within this uncertainty there are both numbers
> commesurate with each other and noncommensurate ones, it has to make
> some arbitrary decision.
I don't know that the decision should be called "arbitrary." It seems
like a good set of choices to me. It probably uses De Moivre's theorem:
http://scholar.hw.ac.uk/site/maths/topic17.aspThis also provides consistency with things like...
log(1)
exp(pi*i)
> Now consider the folowing simple example:
> 588~> x=[10:10];
> 588~> y=x.^(1/3);
> 589~> plot(x,y) % the plot is asymetric with respect to the
> origin
> 590~> floor(1)
> ans = 1
> 591~> floor(3)
> ans = 3
> 592~> y=x.^(floor(1)/floor(3));
> 589~> plot(x,y) % I get the same plot, although floor(1) and
> floor(3) are integers.
>
> It's a pity, because that mighty have been a way to obtain from the
> software the expected behaviour.
When you're working with x^(a/b), and a,b are scalar integers, I think you
are expecting something like this:
when a,b both odd:
sign(x).*(abs(x).^(a/b))
when a is even:
abs(x).^(a/b)
otherwise (a odd, b even):
x.^(a/b) # the usual behavior
altogether using indicators...
abs(rem(a,2))*abs(rem(b,2))*sign(x).*(abs(x).^(a/b)) + \
(1abs(rem(a,2)))*abs(x).^(a/b) + \
abs(rem(a,2))*(1abs(rem(b,2)))*x.^(a/b)
...on a single line:
abs(rem(a,2))*abs(rem(b,2))*sign(x).*(abs(x).^(a/b)) + (1abs(rem(a,2)))*abs(x).^(a/b) + abs(rem(a,2))*(1abs(rem(b,2)))*x.^(a/b)
That might not be what everyone expects for things like (8)^(2/3), and it
isn't what Octave currently does, nor is it what MATLAB does, so I guess
you should just code it like that by hand when needed. I think it works
for integers a,b up to about 16 digits long. Does that seem like a good
strategy?
> Cheers and thanks to all hte participants for the pleasant hour I spent
> recalling long forgotten basic math, Avraham
Same here, believe me. I don't use much of this every day.
Mike

Michael B. Miller, Ph.D.
Assistant Professor
Division of Epidemiology and Community Health
and Institute of Human Genetics
University of Minnesota
http://taxa.epi.umn.edu/~mbmiller/
Octave is freely available under the terms of the GNU GPL.
Octave's home on the web: http://www.octave.orgHow to fund new projects: http://www.octave.org/funding.htmlSubscription information: http://www.octave.org/archive.html


BEGIN PGP SIGNED MESSAGE
Hash: SHA1
John B. Thoo wrote:
> On Mar 20, 2005, at 1:25 PM, John B. Thoo wrote:
>
>> So, if you would allow me at least one more daft question, how do I
>> tell Octave to use the real root so that the plot of y = x^(1/3) is
>> symmetrical about the origin?
>
>
> OK, this is how I've done it.
>
>> x1 = 1:0.1:0; x2 = 0:0.1:1;
>> plot (x1, (x1).^(1/3), x2, x2.^(1/3)) % odd symmetry
>> plot (x1, (x1).^(2/3), x2, x2.^(2/3)) % even symmetry
>
> But is there a more elegant way?
>
> TIA again.
> John.
>
There is a C++ function "cbrt" which calculates the real cube root of
its (real) argument. The attached cbrt.cc uses this function to call it
and return the result to Octave.
Compile it with:
$ mkoctfile cbrt.cc
and then use it:
octave> x = [10:0.1:10];
octave> plot (x, cbrt(x))
 
Geraint Bevan
http://homepage.ntlworld.com/geraint.bevanBEGIN PGP SIGNATURE
Version: GnuPG v1.2.4 (GNU/Linux)
iEYEARECAAYFAkI+AzEACgkQcXV3N50QmNMFHgCfQO5sMyDYS6pBa5scQzVR6TWm
G1QAoI2S67PRh51vd1+AHLfvow8a5njj
=M7IP
END PGP SIGNATURE
// cbrt calculates the real cube root of its arguments
// Copyright (C) 2005 Geraint Paul Bevan
// This program is free software; you can redistribute it and/or modify
// it under the terms of the GNU General Public License as published by
// the Free Software Foundation; either version 2 of the License, or
// (at your option) any later version.
// This program is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
// GNU General Public License for more details.
// You should have received a copy of the GNU General Public License
// along with this program; if not, write to the Free Software
// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 021111307 USA
#include <octave/oct.h>
DEFUN_DLD(cbrt, args, ,
"* texinfo *\n\
@deftypefn {Loadable Function} {@var{y} =} cbrt (@var{x})\n\
\n\
Finds the (real) cube root of @var{x}.\
If @var{x} is a vector or matrix, the cube root of each\n\
element is calculated individually.\n\
@end deftypefn")
{
octave_value_list retval;
for (int i = 0; i < args.length(); i++) {
int nrow = args(i).rows();
int ncol = args(i).columns();
if (args(i).is_complex_type()) {
error ("cbrt cannot handle complex arguments");
}
Matrix tmp(nrow,ncol);
for (int row = 0; row < nrow; row++) {
for (int col = 0; col < ncol; col++) {
tmp(row,col) = cbrt(args(i).matrix_value()(row,col));
}
}
retval(i) = tmp;
}
return retval;
}


x1 = 0.50000
octave:23> x1^(1/3)
ans = 0.79370
octave:24> x2=0.5
x2 = 0.50000
octave:25> x2^(1/3)
ans = 0.39685 + 0.68736i
octave:26> abs(ans)
ans = 0.79370
Will the above work?
Henry
on 3/20/05 1:25 PM, John B. Thoo at [hidden email] wrote:
> First, thanks to Geraint B, Mike M, Avraham for their replies.
>
> On Mar 20, 2005, at 9:09 AM, Geraint Paul Bevan wrote:
>
>> For symmetry about the yaxis, you therefore require that:
>> (1)^(1/n) = 1.
>>
>> However, for n=3, y=(1)^(1/n) has three solutions:
>> y = 1
>> y = 0.5+i*sqrt(3)/2
>> y = 0.5i*sqrt(3)/2
>>
>> You obviously want to use the first of these solutions to get symmetry
>> about the yaxis but Octave doesn't know that unless you tell it.
>
> So, am I now right then in thinking that if z = r * (cos t + i * sin
> t), then Octave will return z^(1/n) = r^(1/n) * (cos (t/n) + i * sin
> (t/n))?
>
>
> And on Mar 20, 2005, at 10:52 AM, Geraint Paul Bevan added:
>
>> $ g++ a.cc && ./a.out
>> +1^2.0 1
>>  1^2.0 1
>> +1^0.5 1
>>  1^0.5 nan
>>
>> As you can see from the results above, the double version of pow cannot
>> cope with a negative base and noninteger exponent. Therefore in this
>> case Octave uses the Complex (i.e. complex<double>) version of pow. It
>> is this routine which is finding complex solutions for negative x
>> instead of the real root that would be required for symmetry about the
>> yaxis.
>
> I think I have a better understanding now. Of course, I know that I
> want the real root :O but it never occurred to me that Octave
> wouldn't use the real root if there is one. Duh!
>
> So, if you would allow me at least one more daft question, how do I
> tell Octave to use the real root so that the plot of y = x^(1/3) is
> symmetrical about the origin?
>
> TIA.
> John.
>
>
>
> 
> Octave is freely available under the terms of the GNU GPL.
>
> Octave's home on the web: http://www.octave.org> How to fund new projects: http://www.octave.org/funding.html> Subscription information: http://www.octave.org/archive.html> 
>

Octave is freely available under the terms of the GNU GPL.
Octave's home on the web: http://www.octave.orgHow to fund new projects: http://www.octave.org/funding.htmlSubscription information: http://www.octave.org/archive.html

12
