[comp.lang.c] # TO THE NTH POWER

CJH101@psuvm.psu.edu (Carl J. Hixon) (11/01/90)

I appologize for bothering you computer wizards with such an elementary
question but, I'm floundering.  Why am I unable to find an opperator which
raises a number to a power. (The fortran equivalent of x**n)  Is there such
an opporator or do I need to write my own function?  It seems like a terrible
oversite to me that such a common operation was overlooked.

Any help would be greatly appreciated!  I'm not sure I could even write the
function, I'm currently writing my first C program ever...

Thankyou,

Carl

edgincd2@mentor.cc.purdue.edu (Chris Edgington *Computer Science Major*) (11/02/90)

In article <90305.005050CJH101@psuvm.psu.edu>, CJH101@psuvm.psu.edu (Carl J. Hixon) writes:
> I appologize for bothering you computer wizards with such an elementary
> question but, I'm floundering.  Why am I unable to find an opperator which
> raises a number to a power. (The fortran equivalent of x**n)  Is there such

I don't know if this is what you are looking for, but this is a neat little
trick to take a number to a power.

     Answer = exp(ln(Root)*Exponent);

Chris Edgington

randy@cs.tamu.edu (Randy Hutson) (11/02/90)

In article <15984@mentor.cc.purdue.edu-> edgincd2@mentor.cc.purdue.edu (Chris Edgington *Computer Science Major*) writes:
->In article <90305.005050CJH101@psuvm.psu.edu->, CJH101@psuvm.psu.edu (Carl J. Hixon) writes:
->-> I appologize for bothering you computer wizards with such an elementary
->-> question but, I'm floundering.  Why am I unable to find an opperator which
->-> raises a number to a power. (The fortran equivalent of x**n)  Is there such
->
->I don't know if this is what you are looking for, but this is a neat little
->trick to take a number to a power.
->
->     Answer = exp(ln(Root)*Exponent);
->
->Chris Edgington

I use an even neater trick!
	
	Answer = pow(Root, Exponent);

mirani-rajiv@cs.yale.edu (Rajiv Mirani) (11/02/90)

In article <15984@mentor.cc.purdue.edu> edgincd2@mentor.cc.purdue.edu (Chris Edgington *Computer Science Major*) writes:
>In article <90305.005050CJH101@psuvm.psu.edu>, CJH101@psuvm.psu.edu (Carl J. Hixon) writes:
>> I appologize for bothering you computer wizards with such an elementary
>> question but, I'm floundering.  Why am I unable to find an opperator which
>> raises a number to a power. (The fortran equivalent of x**n)  Is there such
>
>I don't know if this is what you are looking for, but this is a neat little
>trick to take a number to a power.
>
>     Answer = exp(ln(Root)*Exponent);
>
>Chris Edgington

Assuming, of course, that Root is positive.

--Rajiv

kistler@hp-and.HP.COM (Jim Kistler) (11/02/90)

Try pow

catfood@NCoast.ORG (Mark W. Schumann) (11/02/90)

In article <9750@helios.TAMU.EDU> randy@cs.tamu.edu (Randy Hutson) writes:
>In article <15984@mentor.cc.purdue.edu-> edgincd2@mentor.cc.purdue.edu (Chris Edgington *Computer Science Major*) writes:
>->In article <90305.005050CJH101@psuvm.psu.edu->, CJH101@psuvm.psu.edu (Carl J. Hixon) writes:
>->-> I appologize for bothering you computer wizards with such an elementary
>->-> question but, I'm floundering.  Why am I unable to find an opperator which
>->-> raises a number to a power. (The fortran equivalent of x**n)  Is there such
>->
>->I don't know if this is what you are looking for, but this is a neat little
>->trick to take a number to a power.
>->
>->     Answer = exp(ln(Root)*Exponent);
>->
>->Chris Edgington
>
>I use an even neater trick!
>	
>	Answer = pow(Root, Exponent);
Randy's answer presumes you have the 'pow' function in your library.
It *is* part of the ANSI standard, so if your compiler is ANSI, you're
groovin'.

If you are interested in an
integer version of the same function, try this:

int powi (int root, int exponent)

{

int i;
int result = 1;

  for (i = exponent; i > 1; i--) result *= root;
  return result;

  }

Now you've gotta be careful about overflow here or use long ints.
If you needed an integer version I hope this helps.


-- 
============================================================
Mark W. Schumann  3111 Mapledale Avenue, Cleveland 44109 USA
Domain: catfood@ncoast.org
UUCP:   ...!mailrus!usenet.ins.cwru.edu!ncoast!catfood
============================================================

mcdonald@aries.scs.uiuc.edu (Doug McDonald) (11/02/90)

In article <1990Nov1.232830.17131@NCoast.ORG> catfood@NCoast.ORG (Mark W. Schumann) writes:
>In article <9750@helios.TAMU.EDU> randy@cs.tamu.edu (Randy Hutson) writes:
>>In article <15984@mentor.cc.purdue.edu-> edgincd2@mentor.cc.purdue.edu (Chris Edgington *Computer Science Major*) writes:
>>->In article <90305.005050CJH101@psuvm.psu.edu->, CJH101@psuvm.psu.edu (Carl J. Hixon) writes:
>>->-> I appologize for bothering you computer wizards with such an elementary
>>->-> question but, I'm floundering.  Why am I unable to find an opperator which
>>->-> raises a number to a power. (The fortran equivalent of x**n)  


The original poster said he was looking for an **operator**, not
a function. Now, a function (inlined or not) may or may not make
him happy.

Doug McDonald

rory@maccs.dcss.mcmaster.ca (Rory Jacobs) (11/02/90)

In article <1990Nov1.232830.17131@NCoast.ORG> you write:
>In article <9750@helios.TAMU.EDU> randy@cs.tamu.edu (Randy Hutson) writes:
>>In article <15984@mentor.cc.purdue.edu-> edgincd2@mentor.cc.purdue.edu (Chris Edgington *Computer Science Major*) writes:
>>->In article <90305.005050CJH101@psuvm.psu.edu->, CJH101@psuvm.psu.edu (Carl J. Hixon) writes:

[
  Lots of methods for computing a raised to the power n deleted
]

Oh well, I *mailed* everything said in the above articles to the original
poster.  But the last one (iterative computation of a^n) makes me post the
following version:

int power2 (x,n)
/* note that x, p, sqp if made type double this would allow
 * x^n for real numbers x and integers n.  This just does integers  
 * x raised to the intger powers 
 */

int x;
int n;
{
   int p = 1.0;
   int sqp = x;
   int i=n;

   while (i>0) {
      if ((i % 2) == 1) p *= sqp;
      sqp *= sqp;
      i = i / 2;
   }
   return(p);
}

This calculates x^n and run in O(log2 n) time...

  Rory "Bug Slayer" Jacobs
  rory@maccs.dcss.mcmaster.ca

ok@goanna.cs.rmit.oz.au (Richard A. O'Keefe) (11/02/90)

In article <90305.005050CJH101@psuvm.psu.edu>, CJH101@psuvm.psu.edu (Carl J. Hixon) writes:
> Why am I unable to find an opperator which  raises a number to a power.
Because you're looking for an operator.  Look for a function.
Specifically, look for pow().  In case you're wondering, look in the
manual for your C compiler, in the index under "pow" or "math", or look
in a C textbook.  On a UNIX system, the command "man 3 intro" may help.
On a VMS system, HELP will get you there eventually.

-- 
The problem about real life is that moving one's knight to QB3
may always be replied to with a lob across the net.  --Alasdair Macintyre.

dolf@idca.tds.PHILIPS.nl (Dolf Grunbauer) (11/02/90)

In article <1990Nov1.232830.17131@NCoast.ORG> catfood@NCoast.ORG (Mark W. Schumann) writes:
<int powi (int root, int exponent)
<{
<int i;
<int result = 1;
<  for (i = exponent; i > 1; i--) result *= root;
Oeps, I think this one ---^ should be a 0

<  return result;
<  }

A quicker version (throw in a few registers if still not fast enough):
int powi (int root, int exponent)
{
   int result = 1;

   while (exponent > 0)
   {
      while (!(exponent & 1))
      {
         root     *= root;
         exponent /= 2;
      }
      result *= root;
      exponent--;
   }

   return result;
}
<Now you've gotta be careful about overflow here or use long ints.
<If you needed an integer version I hope this helps.
This is still true for my version as well, and also exponent must be >= 0;
-- 
   _ _ 
  / U |  Dolf Grunbauer  Tel: +31 55 433233 Internet dolf@idca.tds.philips.nl
 /__'<   Philips Information Systems        UUCP     ...!mcsun!philapd!dolf
88  |_\  If you are granted one wish do you know what to wish for right now ?

ok@goanna.cs.rmit.oz.au (Richard A. O'Keefe) (11/02/90)

In article <15984@mentor.cc.purdue.edu>, edgincd2@mentor.cc.purdue.edu (Chris Edgington *Computer Science Major*) writes:
> I don't know if this is what you are looking for, but this is a neat little
> trick to take a number to a power.

>      Answer = exp(ln(Root)*Exponent);

*Please* don't do that.  If you want to raise a number x to an integral
power n, use pow(x, (double)n).  That can handle negative numbers x
(x < 0), while exp(ln(x)*n) can't.  If you want to raise a number x to
a non-integral power y, use pow(x, y).  Why?  Because it can be a _lot_
more accurate.  (If you want to know why, read Cody&Waite.)

-- 
The problem about real life is that moving one's knight to QB3
may always be replied to with a lob across the net.  --Alasdair Macintyre.

gary@hpavla.AVO.HP.COM (Gary Jackoway) (11/02/90)

Mark W. Schumann writes:

> If you are interested in an integer version of the same function, try this:

> int powi (int root, int exponent)
> {
> int i;
> int result = 1;
>  for (i = exponent; i > 1; i--) result *= root;
>  return result;
>  }

> Now you've gotta be careful about overflow here or use long ints.
> If you needed an integer version I hope this helps.

If you're going to post C code, you might TEST it!  Clearly it should
be "i > 0" in the for statement, not "i > 1".
Also, this function gives the wrong result for negative powers.
(Granted, you rarely want to pass this function negative powers,
but it ought to work correctly.)
Further, there is a substantially faster way to do this:

int powi (int root, int exponent)
{
    int i;
    int result = 1;
    int pow = root;

    if (exponent < 0) {
       if (root==1 || root==-1)
          return powi(root,-exponent);
       return 0;
    }
    for (i = exponent; i > 0; i = i>>1) {
	if (i&1) result *= pow;
	pow *= pow;
    }
    return result;
}

The first "if {}" makes sure we always get the right value for negative
exponents.  The "for" loop runs for only as many iterations as the highest
bit in the exponent.  Note that the actual loop is only one line longer
than the slower method.

Gary Jackoway

prk@planet.bt.co.uk (Peter Knight) (11/02/90)

CJH101@psuvm.psu.edu (Carl J. Hixon) writes:

>I appologize for bothering you computer wizards with such an elementary
>question but, I'm floundering.  Why am I unable to find an opperator which
>raises a number to a power. (The fortran equivalent of x**n)  Is there such
>an opporator or do I need to write my own function?  It seems like a terrible
>oversite to me that such a common operation was overlooked.

>Any help would be greatly appreciated!  I'm not sure I could even write the
>function, I'm currently writing my first C program ever...

>Thankyou,

>Carl

The likely offered solutions in C, involving ln()/exp() or the pow() functions
are not going to be a capable as Fortran's x**n.  These functions cannot cope
with Fortran legal operations like (-1e6)**10.  Usually, you will have to
cope with the sign yourself.  The other 'gotcha' you are going to encounter
is that doubles rather than floats have a greater range of values (!), ie


	double d;
	float f;

	d=1e300;
	f=d;		/*	This is a error on 99% of machines	*/

Peter Knight
BT Research

catfood@NCoast.ORG (Mark W. Schumann) (11/03/90)

In article <522@ssp9.idca.tds.philips.nl> dolf@idca.tds.philips.nl (Dolf Grunbauer) writes:

> [I had written]:
>int powi (int root, int exponent)
>{
>int i;
>int result = 1;
>  for (i = exponent; i > 1; i--) result *= root;

Dolf says:

>Oeps, I think this one ---^ should be a 0

Nope.  I intended to leave out the iteration where you'd only be multiplying
by one.

>A quicker version (throw in a few registers if still not fast enough):
>int powi (int root, int exponent)
>{
>   int result = 1;
>
>   while (exponent > 0)
>   {
>      while (!(exponent & 1))
>      {
>         root     *= root;
>         exponent /= 2;
>      }
>      result *= root;
>      exponent--;
>   }
>
>   return result;
>}

Spiffy, but it does depend on (exponent & 1) being the same as saying
"exponent is odd."  Most implementations support this, though.
Again, neither solution supports negative exponents.

-- 
============================================================
Mark W. Schumann  3111 Mapledale Avenue, Cleveland 44109 USA
Domain: catfood@ncoast.org
UUCP:   ...!mailrus!usenet.ins.cwru.edu!ncoast!catfood
============================================================

pfalstad@phoenix.Princeton.EDU (Paul John Falstad) (11/03/90)

In lots of articles people write:
>[
>  Lots of methods for computing a raised to the power n deleted
>]

Just out of curiosity, is there an efficient algorithm for finding
square root of an integer (or the integral part of the square root, if
the number is not a perfect square) using only integer arithmetic?
How about the cube root?

"No" is a valid answer.

--
Paul Falstad, pfalstad@phoenix.princeton.edu PLink:HYPNOS GEnie:P.FALSTAD
I would bring back hanging, and go into rope.  I would cut off the more
disreputable parts of the body and use the space for playing fields.

eager@ringworld.Eng.Sun.COM (Michael J. Eager) (11/03/90)

In article <90305.005050CJH101@psuvm.psu.edu> CJH101@psuvm.psu.edu (Carl J. Hixon) writes:
>Why am I unable to find an opperator which
>raises a number to a power. (The fortran equivalent of x**n)  

There is no exponentiation operator.  (Real programmers use integers :-) )

Look at the pow() function.

-- Mike Eager

enag@ifi.uio.no (Erik Naggum) (11/03/90)

In article <1990Nov2.182217.13958@NCoast.ORG> catfood@NCoast.ORG (Mark W. Schumann) writes:

   Spiffy, but it does depend on (exponent & 1) being the same as saying
   "exponent is odd."  Most implementations support this, though.
   Again, neither solution supports negative exponents.

Which are the implementations that doesn't?

For all positive integers, you can prove that odd integers have the
least significant bit equal to one, and even integers have the least
significant bit equal to zero, when represented in binary.  Hint: odd,
positive integers have a remainder of 1 when divided by 2.

For a ones-complement architectures, this even holds for negative
integers.  Twos-complement architectures will return "odd" for an even
negative number in this simple test.  Was that what you were thinking
of?

--
[Erik Naggum]	Naggum Software
		Holmens gt 3, 2nd floor			<== NEW!
		BOX 1570 VIKA / 0118 OSLO / NORWAY	<== NEW!
--

catfood@NCoast.ORG (Mark W. Schumann) (11/04/90)

In article <ENAG.90Nov3145415@hild.ifi.uio.no> enag@ifi.uio.no (Erik Naggum) writes:
>In article <1990Nov2.182217.13958@NCoast.ORG> catfood@NCoast.ORG (Mark W. Schumann) writes:
>
>   Spiffy, but it does depend on (exponent & 1) being the same as saying
>   "exponent is odd."  Most implementations support this, though.
>   Again, neither solution supports negative exponents.
>
>Which are the implementations that doesn't?

There may be none.  You are simply not allowed to assume *anything*
about the internal representation of integers if your programs are to
be 100% portable ANSI.  I realize this is an extremely conservative
application of the principle, but "all the machines we know about"
is not the same as "all the machines that could ever possibly be."

-- 
============================================================
Mark W. Schumann  3111 Mapledale Avenue, Cleveland 44109 USA
Domain: catfood@ncoast.org
UUCP:   ...!mailrus!usenet.ins.cwru.edu!ncoast!catfood
============================================================

henry@zoo.toronto.edu (Henry Spencer) (11/04/90)

In article <3809@idunno.Princeton.EDU> pfalstad@phoenix.Princeton.EDU (Paul John Falstad) writes:
>Just out of curiosity, is there an efficient algorithm for finding
>square root of an integer (or the integral part of the square root, if
>the number is not a perfect square) using only integer arithmetic?

It's not difficult using Newton's Method, i.e. successive approximations,
if you have good integer division.  Like this:

	sqrtx = ( sqrtx + x/sqrtx ) / 2

although you'll have to watch roundoff effects to avoid the possibility
of an infinite loop when the square root isn't exact.  You also need an
initial approximation, although just x/2 is a good start.

>How about the cube root?

	curtx = ( 2*curtx + x/(curtx*curtx) ) / 3

if I'm not mistaken.  Same caveat about roundoff error and infinite loops.
-- 
"I don't *want* to be normal!"         | Henry Spencer at U of Toronto Zoology
"Not to worry."                        |  henry@zoo.toronto.edu   utzoo!henry

dan@b11.ingr.com (Dan Webb) (11/04/90)

In article <ENAG.90Nov3145415@hild.ifi.uio.no>, enag@ifi.uio.no (Erik Naggum) writes:
> For all positive integers, you can prove that odd integers have the
> least significant bit equal to one, and even integers have the least
> significant bit equal to zero, when represented in binary.  Hint: odd,
> positive integers have a remainder of 1 when divided by 2.
> 
> For a ones-complement architectures, this even holds for negative
> integers.  Twos-complement architectures will return "odd" for an even
> negative number in this simple test.  Was that what you were thinking
> of?

I think you've got that backwards.  In 2's complement, -1 has all bits set,
and it's certainly an odd number.

I think the bottom line is that in 2's complement notation, the LSBit means
an odd number.  Period.

----------------------
Dan Webb
Intergraph Corp.
...!uunet!ingr!b11!dan

donn@brat.UUCP (Donn Pedro) (11/05/90)

In article <90305.005050CJH101@psuvm.psu.edu>, CJH101@psuvm.psu.edu (Carl J. Hixon) writes:
> I appologize for bothering you computer wizards with such an elementary
> question but, I'm floundering.  Why am I unable to find an opperator which
> raises a number to a power. (The fortran equivalent of x**n)  Is there such
> an opporator or do I need to write my own function?  It seems like a terrible
> oversite to me that such a common operation was overlooked.
> 
> Any help would be greatly appreciated!  I'm not sure I could even write the
> function, I'm currently writing my first C program ever...
> 
> Thankyou,
> 
> Carl

Straight from K&R.

Page 25.


/*power: raise base to n-th power; n>= 0 */

int power(int base, int n)
{
	int i, p;

	p = 1;
	for (i = 1; i <= n; ++i)
		p = p * base;
	return p;
}


This is written as a function.  On page 24 you will find the same 
stuff in a main().

      Donn F Pedro ....................a.k.a. uunet!nwnexus!mcgp1!brat!donn
         else:  {the known world}!uunet!nwnexus!mcgp1!brat!donn 

ok@goanna.cs.rmit.oz.au (Richard A. O'Keefe) (11/05/90)

In article <prk.657559483@pegasus>, prk@planet.bt.co.uk (Peter Knight) writes:
> The likely offered solutions in C, involving ln()/exp() or the pow() functions
> are not going to be a capable as Fortran's x**n.  These functions cannot cope
> with Fortran legal operations like (-1e6)**10.

This claim is false.  C's pow() function is precisely as capable as
Fortran's ** operator.  In particular,
	#include <math.h>
	...
	x = pow(-1.0e6, (double)10);
will set x to 1.0e60.  pow() checks for integral exponents.
-- 
The problem about real life is that moving one's knight to QB3
may always be replied to with a lob across the net.  --Alasdair Macintyre.

karl@ima.isc.com (Karl Heuer) (11/05/90)

In article <90305.005050CJH101@psuvm.psu.edu> CJH101@psuvm.psu.edu (Carl J. Hixon) writes:
>Why am I unable to find an opperator which raises a number to a power.  (The
>fortran equivalent of x**n)

If you look up "exponentiation" in the index to K&R, you'll be directed to a
page that mentions the non-existence of such an operator, and an example
function that implements the integer-to-integer version.  (And, if the K&R was
2nd edition, a mention of the existence of the real-to-real version in the
standard library.)

As to *why*, the usual explanation is that C operators tend to mirror common
hardware instructions, with anything more complex (such as string copy) being
implemented as a library function instead.  This explanation isn't entirely
satisfactory, and in fact I believe C should have had such an operator; but
I'll continue that thread in alt.lang.cfutures.

In article <1990Nov1.232830.17131@NCoast.ORG> catfood@NCoast.ORG (Mark W. Schumann) writes:
>	int result = 1;
>	for (i = exponent; i > 1; i--) result *= root;

Fencepost error.  Should be "i > 0".

Karl W. Z. Heuer (karl@ima.isc.com or uunet!ima!karl), The Walking Lint

rh@smds.UUCP (Richard Harter) (11/05/90)

In article <1990Nov3.220554.19404@zoo.toronto.edu>, henry@zoo.toronto.edu (Henry Spencer) writes:
> In article <3809@idunno.Princeton.EDU> pfalstad@phoenix.Princeton.EDU (Paul John Falstad) writes:
> >Just out of curiosity, is there an efficient algorithm for finding
> >square root of an integer (or the integral part of the square root, if
> >the number is not a perfect square) using only integer arithmetic?

> It's not difficult using Newton's Method, i.e. successive approximations,
> if you have good integer division.  Like this:

 	sqrtx = ( sqrtx + x/sqrtx ) / 2

> although you'll have to watch roundoff effects to avoid the possibility
> of an infinite loop when the square root isn't exact.  You also need an
> initial approximation, although just x/2 is a good start.

Actually x/2 is not a good start -- the problem is that the convergence
will be geometric with power .5 (i.e. you only gain one bit in each
iteration) until you get close to the actual root.  Once you are sufficiently
close the convergence rate does become quadratic.  What you want to do is
(in effect) get the nearest power of four and take the square root of that.
A reasonably good way to get a starting value (given the integer arithmetic
restriction) is:

	for (y=x,sqrtx=1;y>0;) {ysave = y;y =/ 16;sqrtx =* 4;}
	if (ysave > 4) sqrtx *= 2;

One can fiddle with this (and the Newton formula) to get maximal efficiency
if there is a real need; however a good starting point should suffice for
all but the most cycle hungry.


-- 
Richard Harter, Software Maintenance and Development Systems, Inc.
Net address: jjmhome!smds!rh Phone: 508-369-7398 
US Mail: SMDS Inc., PO Box 555, Concord MA 01742
This sentence no verb.  This sentence short.  This signature done.

peter@ficc.ferranti.com (Peter da Silva) (11/06/90)

In article <226@smds.UUCP> rh@smds.UUCP (Richard Harter) writes:
> 	for (y=x,sqrtx=1;y>0;) {ysave = y;y =/ 16;sqrtx =* 4;}

Wow.

This man has been programming in C for a *long* time. I'm impressed.

(quiz: why do I come to this conclusion?)
-- 
Peter da Silva.   `-_-'
+1 713 274 5180.   'U`
peter@ferranti.com

catfood@NCoast.ORG (Mark W. Schumann) (11/06/90)

In article <5240@ima.ima.isc.com> karl@ima.isc.com (Karl Heuer) writes:

> catfood@NCoast.ORG (Mark W. Schumann) writes:
>>	int result = 1;
>>	for (i = exponent; i > 1; i--) result *= root;

>Fencepost error.  Should be "i > 0".

Uh, yes, that's right.  Sorry, NetPeople.

-- 
============================================================
Mark W. Schumann  3111 Mapledale Avenue, Cleveland 44109 USA
Domain: catfood@ncoast.org
UUCP:   ...!mailrus!usenet.ins.cwru.edu!ncoast!catfood
============================================================

fredg@marob.masa.com (Fred Goldrich) (11/06/90)

In article <P1X6QN1@xds13.ferranti.com> peter@ficc.ferranti.com (Peter da Silva) writes:
>In article <226@smds.UUCP> rh@smds.UUCP (Richard Harter) writes:
>> 	for (y=x,sqrtx=1;y>0;) {ysave = y;y =/ 16;sqrtx =* 4;}
>
>Wow.
>
>This man has been programming in C for a *long* time. I'm impressed.
>
>(quiz: why do I come to this conclusion?)
>-- 
>Peter da Silva.   `-_-'
>+1 713 274 5180.   'U`
>peter@ferranti.com

I think it's because you are reading "sqrtx =* 4" as ancient-C style
for the assignment "sqrtx *= 4"; my guess is that is was a typo.
-- 
	Fred Goldrich
	{att,philabs,rutgers,cmcl2,hombre}!phri!marob!maestro!fred
	(formerly dasys1!maestro)

pgd@bbt.se (11/07/90)

In article <P1X6QN1@xds13.ferranti.com> peter@ficc.ferranti.com (Peter da Silva) writes:
>In article <226@smds.UUCP> rh@smds.UUCP (Richard Harter) writes:
>> 	for (y=x,sqrtx=1;y>0;) {ysave = y;y =/ 16;sqrtx =* 4;}
                                            ^^          ^^
>
>Wow.
>
>This man has been programming in C for a *long* time. I'm impressed.
>
You mean have not been programming in C for a *long* time, don't you? :-)

rh@smds.UUCP (Richard Harter) (11/07/90)

In article <P1X6QN1@xds13.ferranti.com>, peter@ficc.ferranti.com (Peter da Silva) writes:
> In article <226@smds.UUCP> rh@smds.UUCP (Richard Harter) writes:
> > 	for (y=x,sqrtx=1;y>0;) {ysave = y;y =/ 16;sqrtx =* 4;}

> Wow.

> This man has been programming in C for a *long* time. I'm impressed.

> (quiz: why do I come to this conclusion?)

Er, well, yes.  The moral of the story is that your short term memory
can fail you when you post in the wee hours of the morning.  It makes
a nice quiz question.  Color me embarassed.

Wonder how many people caught it?

-- 
Richard Harter, Software Maintenance and Development Systems, Inc.
Net address: jjmhome!smds!rh Phone: 508-369-7398 
US Mail: SMDS Inc., PO Box 555, Concord MA 01742
This sentence no verb.  This sentence short.  This signature done.

blambert@lotus.com (Brian Lambert) (11/07/90)

In article <P1X6QN1@xds13.ferranti.com> peter@ficc.ferranti.com (Peter da Silva) writes:
>In article <226@smds.UUCP> rh@smds.UUCP (Richard Harter) writes:
>> 	for (y=x,sqrtx=1;y>0;) {ysave = y;y =/ 16;sqrtx =* 4;}
>
>Wow.
>
>This man has been programming in C for a *long* time. I'm impressed.
>
>(quiz: why do I come to this conclusion?)
>-- 
>Peter da Silva.   `-_-'
>+1 713 274 5180.   'U`
>peter@ferranti.com

It's because programs that take up fewer lines of vertical space also
consume less CPU time.

--------------------------------
Brian Lambert
Lotus Development Corporation 
blambert@lotus.com
-- 
--------------------------------
Brian Lambert
Lotus Development Corporation 

gwyn@smoke.brl.mil (Doug Gwyn) (11/07/90)

In article <90305.005050CJH101@psuvm.psu.edu> CJH101@psuvm.psu.edu (Carl J. Hixon) writes:
>Why am I unable to find an opperator which raises a number to a power.

Because the designers of C decided that it was better to leave such a
potentially expensive operation out of the language repertoire and instead
relegate it to a math library function, pow().

engbert@cs.vu.nl (Engbert Gerrit IJff) (11/07/90)

In article <P1X6QN1@xds13.ferranti.com>,
	peter@ficc.ferranti.com (Peter da Silva) writes:
) In article <226@smds.UUCP> rh@smds.UUCP (Richard Harter) writes:
) > 	for (y=x,sqrtx=1;y>0;) {ysave = y;y =/ 16;sqrtx =* 4;}
) 
) Wow.
) 
) This man has been programming in C for a *long* time. I'm impressed.
) 
) (quiz: why do I come to this conclusion?)
) -- 
) Peter da Silva.   `-_-'
) +1 713 274 5180.   'U`
) peter@ferranti.com
That is an easy quiz. Look at the =/ and =* operators.

Bert.

gwyn@smoke.brl.mil (Doug Gwyn) (11/07/90)

In article <ENAG.90Nov3145415@hild.ifi.uio.no> enag@ifi.uio.no (Erik Naggum) writes:
-In article <1990Nov2.182217.13958@NCoast.ORG> catfood@NCoast.ORG (Mark W. Schumann) writes:
-   Spiffy, but it does depend on (exponent & 1) being the same as saying
-   "exponent is odd."  Most implementations support this, though.

That's good, because the C standard requires this, for positive integers.

-For a ones-complement architectures, this even holds for negative
-integers.  Twos-complement architectures will return "odd" for an even
-negative number in this simple test.  Was that what you were thinking of?

Other way around.

gwyn@smoke.brl.mil (Doug Gwyn) (11/07/90)

In article <1990Nov3.204327.19057@NCoast.ORG> catfood@NCoast.ORG (Mark W. Schumann) writes:
>You are simply not allowed to assume *anything* about the internal
>representation of integers if your programs are to be 100% portable ANSI.

Hey, that's not true.  Look it up.

twpierce@amherst.bitnet (Tim Pierce) (11/08/90)

In article <P1X6QN1@xds13.ferranti.com>, peter@ficc.ferranti.com (Peter da Silva) writes:
> In article <226@smds.UUCP> rh@smds.UUCP (Richard Harter) writes:
>> 	for (y=x,sqrtx=1;y>0;) {ysave = y;y =/ 16;sqrtx =* 4;}
> 
> Wow.
> 
> This man has been programming in C for a *long* time. I'm impressed.
> 
> (quiz: why do I come to this conclusion?)

Because of the archaic operators! (I knew I'd get something right someday!)
Someone else suggested that this might have been a typo -- however, since he
uses postfix twice, (postfix is not the right term, I know, but for lack of a
better one...) in "y =/ 16" and "sqrtx =* 4", this is less likely.

----------------------------------------------------------------------------
Tim Pierce		| "They call television a medium. That is because
twpierce@amherst.bitnet	|  it is neither rare nor well done." - Ernie Kovacs
----------------------------------------------------------------------------

bhoughto@cmdnfs.intel.com (Blair P. Houghton) (11/08/90)

In article <8149@star.cs.vu.nl> engbert@cs.vu.nl (Engbert Gerrit IJff) writes:
>In article <P1X6QN1@xds13.ferranti.com>,
>	peter@ficc.ferranti.com (Peter da Silva) writes:
>) In article <226@smds.UUCP> rh@smds.UUCP (Richard Harter) writes:
>) > 	for (y=x,sqrtx=1;y>0;) {ysave = y;y =/ 16;sqrtx =* 4;}
>) 
>) This man has been programming in C for a *long* time. I'm impressed.
>) (quiz: why do I come to this conclusion?)
>
>That is an easy quiz. Look at the =/ and =* operators.

It's not an easy quiz (though your answer is correct);
it's an inside joke.

The reason Peter comes to any (discernible) conclusion is
so that he can see his (discernible) name in (phosphorescent) lights.

:-)

Ob. C:

Now that the cat's out of the bag, I'd like to point out to
the newer C programmers that the assignment operators we
know as '+=' et al did used to be usable as '=+' et al.
However, it was very long ago[*] that this was changed, and
now only the '+=' forms are considered correct.  You may
(unsettingly often) find a compiler that will accept the
obsolete form (properly flagging it with a warning that
you're being a crufty old coot :), but ANSI forbids it, and
most compilers will call it an error.  The reason given
most often is that it can confuse such expressions as
`x=-y'; is it `x =- y' or `x = -y'?  Under maximal-munch[**]
it's the former.  Under the rules of precedence it's the
latter.

[*] "...three days after Marconi invented the f*cking thing."
-the great Warren Oates in 'Blue Thunder'
[**] Maximal munch:  if the next n characters are a valid token,
and the next n+m characters are also a valid token, consider
the longer one to be the token.  (Hey, I _did_ say this was
for the new guys.)

				--Blair
				  "Prior art includes Picasso,
				   Michelangelo, _and_ Renoir.
				   Why choose?"

kentd@FtCollins.NCR.com (Kent.Dalton) (11/08/90)

Peter da Silva - peter@ferranti.com sez:

>> 	for (y=x,sqrtx=1;y>0;) {ysave = y;y =/ 16;sqrtx =* 4;}
>
> Wow.
>
> This man has been programming in C for a *long* time. I'm impressed.
>
>(quiz: why do I come to this conclusion?)
Is it because he's using the old, crufty, archaic form of =/ and =* instead of
/= and *= ?

--
/**************************************************************************/
/* Kent Dalton                     * EMail: Kent.Dalton@FtCollins.NCR.COM */
/* NCR Microelectronics            *   CIS: 72320,3306                    */
/* 2001 Danfield Ct. MS470A        *                                      */
/* Fort Collins, Colorado 80525    * "This mind intentionally left blank" */
/* (303)223-5100 X-319             *    All standard disclaimers apply    */
/**************************************************************************/

cik@l.cc.purdue.edu (Herman Rubin) (11/08/90)

In article <1990Nov3.220554.19404@zoo.toronto.edu>, henry@zoo.toronto.edu (Henry Spencer) writes:
> In article <3809@idunno.Princeton.EDU> pfalstad@phoenix.Princeton.EDU (Paul John Falstad) writes:
> >Just out of curiosity, is there an efficient algorithm for finding
> >square root of an integer (or the integral part of the square root, if
> >the number is not a perfect square) using only integer arithmetic?
 
> It's not difficult using Newton's Method, i.e. successive approximations,
> if you have good integer division.  Like this:
 
> 	sqrtx = ( sqrtx + x/sqrtx ) / 2
 
> although you'll have to watch roundoff effects to avoid the possibility
> of an infinite loop when the square root isn't exact.  You also need an
> initial approximation, although just x/2 is a good start.

Assuming this is done in integer arithmetic, there is only one case of an
infinite loop.  Using real arithmetic (of course impossible on a computer),
the iteration is always too large.  Thus it should be rounded down.  Now
it is possible that this makes it so small that the next iteration is
larger.  If one starts with something too large, this can only happen
if x = n(n+2), and this can be detected by the iterate being larger than
the previous value, so that a stopping rule when the iterate does not
decrease, and taking the previous value, is always right, if the initial
value is too large.
-- 
Herman Rubin, Dept. of Statistics, Purdue Univ., West Lafayette IN47907
Phone: (317)494-6054
hrubin@l.cc.purdue.edu (Internet, bitnet)   {purdue,pur-ee}!l.cc!hrubin(UUCP)

adkins@carp.cis.ohio-state.edu (Brian Adkins) (11/08/90)

In article <1990Nov7.020445.12839@lotus.com> blambert@lotus.UUCP (Brian Lambert) writes:
>>> 	for (y=x,sqrtx=1;y>0;) {ysave = y;y =/ 16;sqrtx =* 4;}
>It's because programs that take up fewer lines of vertical space also
>consume less CPU time.
>Brian Lambert

Mr. Lambert, do you really believe that the following two routines 
consume different amounts of CPU time?

for (y=x,sqrtx=1;y>0;) {
  ysave = y;
  y /= 16;
  sqrtx *= 4;
}

for (y=x,sqrtx=1;y>0;) {ysave = y;y /= 16;sqrtx *= 4;}

Please correct me if I am mistaken, but I believe any compiler will generate
the same code for both of these routines.

Brian Adkins

unhd (Jonathan W Miner) (11/10/90)

In article <90305.005050CJH101@psuvm.psu.edu> CJH101@psuvm.psu.edu (Carl J. Hixon) writes:
>I appologize for bothering you computer wizards with such an elementary
>question but, I'm floundering.  Why am I unable to find an opperator which
>raises a number to a power. (The fortran equivalent of x**n)  Is there such
>an opporator or do I need to write my own function?  It seems like a terrible
>oversite to me that such a common operation was overlooked.

There is a standard C function pow(x,y) which returns x raised to the y.
:x

-- 
-----------------------------------------------------------------
Jonathan Miner        | I don't speak for UNH, and UNH does not 
jwm775@unhd.unh.edu   | speak for me! 
(603)868-3416         | Rather be downhill skiing... or hacking... 

rnewman@lotus.lotus.com (Ron Newman ) (11/10/90)

In article <4200@goanna.cs.rmit.oz.au>, ok@goanna.cs.rmit.oz.au (Richard
A. O'Keefe) writes:
|>.....  C's pow() function is precisely as capable as
|> Fortran's ** operator.  In particular,
|> 	#include <math.h>
|> 	...
|> 	x = pow(-1.0e6, (double)10);
|> will set x to 1.0e60.  pow() checks for integral exponents.

Unfortunately, many implementations of pow(), such as in Unix System V
for PC's, don't produce acceptable results for integral exponents.
You get aberrant results like pow(10.0,3.0) != 1000.0, or pow(3.0,2.0)
!= 9.0.

There's really no excuse for this; if both the arguments and the result
of a floating-point math function can be represented exactly, then the
function should give the correct, exact result.   When it doesn't, the
function has been implemented wrong, in one or both of two ways:

   a) The function was implemented in terms of e^x and ln(x), instead
      of using the floating-point hardware's native operations
(typically
      2^x and y*log2(x))

   b) The function was implemented in terms of intermediate functions
      that return double-precision (64-bit) values.  This led to loss
      of precision because the underlying arithmetic was done with
      extended-precision (80-bit) values.

In our System V products for PC's, we had to bypass the operating-system
vendors'  math libraries and reimplement pow() (and several other
functions) in 
raw 80x87 assembly code.

/Ron Newman

jim.nutt@stjhmc.fidonet.org (jim nutt) (11/10/90)

In a message of <Nov 05 21:40> Peter da Silva (peter@ficc.ferranti.com ) 
writes:

 PdS> In article <226@smds.UUCP> rh@smds.UUCP (Richard Harter) writes:
>       for (y=x,sqrtx=1;y>0;) {ysave = y;y =/ 16;sqrtx =* 4;}

 PdS> Wow.

 PdS> This man has been programming in C for a *long* time. I'm impressed.

 PdS> (quiz: why do I come to this conclusion?)

answer: because he is using the obsolescent forms of the assignment 
operators.. (/=, *=, etc.)

jim nutt
'the computer handyman' 


 

--  
Uucp: ...{gatech,ames,rutgers}!ncar!asuvax!stjhmc!jim.nutt
Internet: jim.nutt@stjhmc.fidonet.org

ok@goanna.cs.rmit.oz.au (Richard A. O'Keefe) (11/10/90)

In article <1990Nov9.152921@lotus.lotus.com>, rnewman@lotus.lotus.com (Ron Newman ) writes:
> In article <4200@goanna.cs.rmit.oz.au>, ok@goanna.cs.rmit.oz.au (Richard
> A. O'Keefe) writes:
> |> pow() checks for integral exponents.
> 
> Unfortunately, many implementations of pow(), such as in Unix System V
> for PC's, don't produce acceptable results for integral exponents.
> You get aberrant results like pow(10.0,3.0) != 1000.0, or pow(3.0,2.0)
> != 9.0.

SVID release 2, vol 1, p 170:
        double pow(x, y)
        double x, y;
[so much for putting the type on the previous line being "standard"!]
                                         y
        The functions [sic] pow returns x.  If x is zero, y must be
	positive.  If x is negative, y must be an integer.
	The function pow returns 0 and sets errno to EDOM when x is 0
	and y is non-positive, or when x is negative and y is not an
	integer.  In these cases a message indicating DOMAIN error is
	printed on the standard error output.  When the correct value
	for pow would overfluw or underflow, pow returns +/- HUGE or 0
	respectively and sets errno to ERANGE.
    FUTURE DIRECTIONS
	The return value of pow will be negative HUGE_VAL when an
	illegal combination of input arguments is passed to pow.
	[HUGE_VAL is +infinity if you have it, +MAX_DOUBLE otherwise.]

This definition has now been superseded by the ANSI C standard, but I
think it was the definition which was relevant to UNIX PCs.  If the
answer was significantly in error (compared with the code in Cody &
Waite, say) it would have been appropriate to complain to the vendor.
The company I used to work for once got a bug report from a vendor that
we were porting our product to; it was my pleasure to track that bug
down to a mistake in that vendor's version of the UNIX V.2 math library.

So yes, you may have to work around mistakes in pow() just as you may have
to work around mistakes in printf() -- I've had to do that too -- but that
doesn't mean that pow() or printf() was not the appropriate function to try.

-- 
The problem about real life is that moving one's knight to QB3
may always be replied to with a lob across the net.  --Alasdair Macintyre.

karl@ima.isc.com (Karl Heuer) (11/11/90)

In article <4233@goanna.cs.rmit.oz.au> ok@goanna.cs.rmit.oz.au (Richard A. O'Keefe) writes:
>In article <1990Nov9.152921@lotus.lotus.com>, rnewman@lotus.lotus.com (Ron Newman ) writes:
>> Unfortunately, many implementations of pow(), such as in Unix System V
>> for PC's, don't produce acceptable results for integral exponents.
>
>SVID release 2, vol 1, p 170: [Specs for pow().]  If the answer was
>significantly in error (compared with the code in Cody & Waite, say) it would
>have been appropriate to complain to the vendor.

How much is "significantly"?  An error of a single ulp is enough to cause
(int)pow(2.0, 3.0) to return 7 instead of 8 (which it does, on several current
implementations).  So unless the standards bodies (*not* the individual
vendors) are willing to guarantee exactness for this function, using pow() is
an unportable solution to the problem of integer-to-integer exponentiation.

(It's still possible to use
	#define ipow(i,j) ((int)(pow((double)i, (double)j)+0.5))
for positive numbers, and a more complicated adjustment to handle negatives as
well, but now it's more work than writing an ipow() function from scratch.)

Karl W. Z. Heuer (karl@ima.isc.com or uunet!ima!karl), The Walking Lint

ok@goanna.cs.rmit.oz.au (Richard A. O'Keefe) (11/11/90)

In article <1990Nov10.224939.23622@dirtydog.ima.isc.com>, karl@ima.isc.com (Karl Heuer) writes:
> In article <4233@goanna.cs.rmit.oz.au> ok@goanna.cs.rmit.oz.au (Richard A. O'Keefe) writes:
> >SVID release 2, vol 1, p 170: [Specs for pow().]  If the answer was
> >significantly in error (compared with the code in Cody & Waite, say) it would
> >have been appropriate to complain to the vendor.

> How much is "significantly"?

I said "significantly" in order to avoid giving any definite answer
to that question!  It all depends.  (There, I managed to use the phrase
that proves I'm an expert...)

> An error of a single ulp is enough to cause
> (int)pow(2.0, 3.0) to return 7 instead of 8 (which it does, on several current
> implementations).

Well, yes, but so what?  That's an incredibly stupid way to calculate
integer powers in the first place.  What you want is ROUND TO NEAREST,
not TRUNCATE TOWARDS ZERO.  If you have the IEEE "appendix function" drem(),

	int iround(double x) { return x - drem(x, 1.0); }
...
	eight = iround(pow(2.0, 3.0));

> So unless the standards bodies (*not* the individual
> vendors) are willing to guarantee exactness for this function, using pow() is
> an unportable solution to the problem of integer-to-integer exponentiation.

Using pow() is just fine if you ROUND instead of TRUNCATING.
(Oh for the good old Algol 60 days...)

As a matter of genuine curiosity, just how useful is a general integer
power function anyway?  On a 32-bit machine, you can compute 2**n for
0 <= n <= 30 only (might as well use a table), and for larger numbers
the useful range of n is smaller.  (For 10**n, the range is 0 <= n <= 9.
Again, you would be better off with a table.)

I fully appreciate that integer powers are useful in languages like Lisp
without artificial restrictions on integer size, but in __C__?

Does anyone know why C hasn't got a round() function and why ANSI
failed to add one?  Using ANSI functions, we can do

	#include <math.h>
	
	double round(double x)
	    {
		double f = floor(x);
		double c = ceil(x);
		return x-f > c-x ? c : f;
	    }

	#define iround(x) (int)round(x)

Even if the "round to even when there is a tie" rule were insisted on
(as I think it should) it's easy to implement.  (I've posted one before.)
-- 
The problem about real life is that moving one's knight to QB3
may always be replied to with a lob across the net.  --Alasdair Macintyre.

gwyn@smoke.brl.mil (Doug Gwyn) (11/11/90)

In article <4239@goanna.cs.rmit.oz.au> ok@goanna.cs.rmit.oz.au (Richard A. O'Keefe) writes:
>Does anyone know why C hasn't got a round() function and why ANSI
>failed to add one?

Because it is trivial to roll your own:
#define	Round( d )	floor( (d) + 0.5 )	/* requires <math.h> */
and because there was no established existing practice (e.g. the library
Base Document did not specify such a function).

ok@goanna.cs.rmit.oz.au (Richard A. O'Keefe) (11/12/90)

In article <14427@smoke.brl.mil>, gwyn@smoke.brl.mil (Doug Gwyn) writes:
: In article <4239@goanna.cs.rmit.oz.au> ok@goanna.cs.rmit.oz.au (Richard A. O'Keefe) writes:
: >Does anyone know why C hasn't got a round() function and why ANSI
: >failed to add one?
: 
: Because it is trivial to roll your own:
: #define	Round( d )	floor( (d) + 0.5 )	/* requires <math.h> */
: and because there was no established existing practice (e.g. the library
: Base Document did not specify such a function).

I'm afraid it's not that trivial.  That doesn't provide the preferred
round-to-even-when-there-is-a-tie.  It is a _pain_ that the function
nint() isn't portable.  The Base Document didn't specify fsetpos() either,
but that didn't stop ANSI filling in a nasty gap.

What's the NCEWG up to?  Are they likely to address this?
-- 
The problem about real life is that moving one's knight to QB3
may always be replied to with a lob across the net.  --Alasdair Macintyre.

gwyn@smoke.brl.mil (Doug Gwyn) (11/12/90)

In article <4247@goanna.cs.rmit.oz.au> ok@goanna.cs.rmit.oz.au (Richard A. O'Keefe) writes:
>: #define	Round( d )	floor( (d) + 0.5 )	/* requires <math.h> */
>I'm afraid it's not that trivial.  That doesn't provide the preferred
>round-to-even-when-there-is-a-tie.

It won't matter in any reasonable application.

>The Base Document didn't specify fsetpos() either,
>but that didn't stop ANSI filling in a nasty gap.

It's not clear that it's a particularly pressing need.
I would MUCH rather have had a mathematically correct mod function;
there are something like 4 similar C functions and operations, none
of them correctly providing a direct modulo capability.

>What's the NCEWG up to?  Are they likely to address this?

NCEG has several interests.  I haven't heard any major concern over
this particular topic.

rnewman@lotus.lotus.com (Ron Newman ) (11/13/90)

In article <4239@goanna.cs.rmit.oz.au>, ok@goanna.cs.rmit.oz.au (Richard
A. O'Keefe) writes:
|> In article <1990Nov10.224939.23622@dirtydog.ima.isc.com>,
karl@ima.isc.com (Karl Heuer) writes:

|> > An error of a single ulp is enough to cause
|> > (int)pow(2.0, 3.0) to return 7 instead of 8 (which it does, on
several current
|> > implementations).
|> 
|> Well, yes, but so what?  That's an incredibly stupid way to
calculate
|> integer powers in the first place.  What you want is ROUND TO
NEAREST,
|> not TRUNCATE TOWARDS ZERO. ...
|> 
|> Using pow() is just fine if you ROUND instead of TRUNCATING.
|> (Oh for the good old Algol 60 days...)
|> 
|> As a matter of genuine curiosity, just how useful is a general
integer
|> power function anyway?  On a 32-bit machine, you can compute 2**n
for
|> 0 <= n <= 30 only (might as well use a table), and for larger
numbers
|> the useful range of n is smaller.  (For 10**n, the range is 0 <= n <=
9.
|> Again, you would be better off with a table.)
|> 
|> I fully appreciate that integer powers are useful in languages like
Lisp
|> without artificial restrictions on integer size, but in __C__?

The problem isn't just that some implementations of pow(x,y) are 
inaccurate when x and y are BOTH integers.  They are also inaccurate 
when x isn't an integer but y is.  It's just as bad when

     pow(0.5, 3.0) != 0.125, or
     pow(1.25, 3.0) != 1.953125

It's acceptable for pow() to give inexact results when the exact result
is not representable as a double.  It's not acceptable to get an
inexact
result when an exact one is possible.   If the ANSI standard for pow() 
doesn't say this, then it should.

/Ron Newman

robbie.matthews@f1701.n7802.fido.oz.au (Robbie Matthews) (11/16/90)

Original to: blambert@lotus.com
OK. Here is a way of getting an integer root that I haven't seen yet:
 
 for (sqrtx=1, j=1, k=1; j<x; sqrtx+=1, k+=2, j+=k );
 
Try it and see. Only works for positive integers, but as there is no
multiplication or division, works reasonably fast. Especially on a 8 bit
machine.
 Not claiming it's a _better_ way... it's just a way.

--- TMail v1.19
 * Origin: Prophet BBS (8:7802/1701)