Implement x^(e) in Extended Pascal, without using exponentiation operators

Question

Part of a Pascal ISO 10206 program I am building, requires me to implement a function to exponentiate a real number (x) to Eulers number (e), without using any exponentiation functions already included in Pascal(**,pow,exp...).

I have been trying different algorithms for hours but I cant figure out how to do it without using the already included exponentiation functions.

Any help would be appreciated. Any mathematical algorithm of some sort etc... Thanks in advance.

This seems like rather a waste of time. Makes no sense not to use exp. In any case, there are algorithms published that you could find with some research. What research have you done? — David Heffernan, Jun 01 '17 at 10:59
Take a look at [Taylor series](https://en.wikipedia.org/wiki/Taylor_series#Exponential_function). They can be used to calculate the `Exp()` of a value x (and a lot more). Just use a loop to increment, add or multiply the values in the series. — Rudy Velthuis, Jun 01 '17 at 11:53
FWIW, you have a tag for Free Pascal, but you say you are using Extended Pascal. These are not the same. And are you looking for x^e or for e^x? The latter is the `Exp()` function, but the former can be calculated with it as well. — Rudy Velthuis, Jun 01 '17 at 11:58
@RudyVelthuis Not really. The convergence of the infinite sum is too slow to be useful. — David Heffernan, Jun 01 '17 at 12:35
@David: no one said it should be fast. I would use `Exp()` anyway, but if that is not allowed, Taylor series are one way to do this. — Rudy Velthuis, Jun 01 '17 at 12:37
@RudyVelthuis Only a fool with no mathematical knowledge would attempt to implement `Exp` on top of Taylor series. — David Heffernan, Jun 01 '17 at 12:48
@David: Like dentists. Anyway, found this: https://stackoverflow.com/questions/3518973/floating-point-exponentiation-without-power-function. Looks like a nice trick, although I am not sure how good it is. — Rudy Velthuis, Jun 01 '17 at 13:22
@Rudy There's no reason why a dentist shouldn't know maths. But this mathematician is telling you why Taylor series expansions are not practical for evaluating `Exp` generally. — David Heffernan, Jun 01 '17 at 13:26
@DavidHeffernan yes, I have searched for a lot of algorithms over the internet, but most of them end up having to use some kind of exponent. Ofcourse I know it can be easily done with exp(), but thats the point of the exercise, kind of rebuilding the exponential operators of pascal, which I cant seem to get done since most algorithms I found require me to use some kind of exponentiation. I will try what you guys proposed and post back. — Antón Valladares Poncela, Jun 01 '17 at 14:08
What we proposed? I proposed nothing. I just pointed that Rudy's idea is not useful. — David Heffernan, Jun 01 '17 at 14:19
@David: what if you can divide them by a factor so that they are close to the center point (which is 0, I assume)? (and then adjust later on, of course). I don't have any programming language on this computer, but ISTM that that should work. FWIW, I know some maths, but such topics were not part of the curriculum for dentists. — Rudy Velthuis, Jun 01 '17 at 14:43
@RudyVelthuis Try and work out the details, then get back to me — David Heffernan, Jun 01 '17 at 14:50
Are you sure you mean `x^e` and not `e^x`? The latter is a _much_ more commonly used function. — , Jun 01 '17 at 23:04
@DavidHeffernan: finally got home and worked out the details. — Rudy Velthuis, Jun 01 '17 at 23:12

Rudy Velthuis · Answer 1 · 2017-06-02T10:29:22.353

As others have said, it doesn't make sense not to use Exp() or a function based upon it. But if you really must use something else, and don't want to get too technical/mathematical, then the following should work (the real algorithm is far more complicated and requires much more knowledge of math).

The program uses a combination of the first N terms of Taylor series for the fraction of X and binary exponentiation for the integer part of X. It is probably not very fast, but pretty accurate, even for larger exponents. For comparison, I also display Exp(X). If your Pascal has a Double or Extended type, use those instead of Real.

program SimplePower;

{ Required for Delphi, you can omit it in other Pascals: }
{$APPTYPE CONSOLE}

{ Returns approximate value of e^X using sum of first N terms of Taylor series.
  Works fine with X values between 0 and 1.0 and N ~ 30. }
function Exponential(N: Integer; X: Real): Real;
var
  I: Integer;
begin
  Result := 1.0;
  for I := N - 1 downto 1 do
    Result := 1.0 + X * Result / I;
end;

{ Binary exponentiation of Base by integer Exponent. }
function IntegerExp(Base: Real; Exponent: Integer): Real;
begin
  Result := 1.0;
  while Exponent > 0 do
  begin
    if Odd(Exponent) then
      Result := Result * Base;
    Base := Base * Base;
    Exponent := Exponent shr 1;
  end;
end;

{ Combines IntegerExp function for integral part with 
  Exponential function for fractional part. }
function MyExp(N: Integer; X: Real): Real;
const
  E = 2.7182818284590452353602874713527; { from Google: "e euler" }
var
  Factor: Real;
  Fraction: Real;
begin
  Fraction := Exponential(N, Frac(X));
  Factor := IntegerExp(E, Trunc(X));
  Result := Factor * Fraction;
end;

{ Simple demo: }
const
  N = 30;
  X = 73.4567890242421234;

begin
  Writeln('MyExp(', N, ', ', X:22:18, ') = ', MyExp(N, X):22:18);
  Writeln('Exp(', X:22:18, ')       = ', Exp(X):22:18);
end.

Ref:

I did not do anything for negative exponents, but just know that Exp(-x) = 1/Exp(x). You could amend MyExpwith that knowledge.

I see two downvotes. A comment about the "why" would be nice. — Rudy Velthuis, Jun 02 '17 at 10:21

score 0 · Accepted Answer · answered Jun 03 '17 at 07:30

I used the solution pointed out by @RudyVelthuis in some other post, but modified it a bit. It is based upon that x^0.5 = sqrt(x), which we can use to our advantage. Ill leave the Pascal ISO 10206 code I used attached. Thank you all for your help, specially Rudy.

function MyPow(base,power,precision:real):real;
begin
    if (power<0) then MyPow:=1/MyPow(base,-power,precision)
    else   if (power>=10) then MyPow:=sqr(MyPow(base,power/2,precision/2))
    else   if (power>=1) then MyPow:=base*MyPow(base,power-1,precision)
    else   if (precision>=1) then MyPow:=sqrt(base)
    else MyPow:=sqrt(MyPow(base,power*2,precision*2));
end;

That's cool. But note that the solution I gave is slightly more accurate. In the original, I used Delphi's Extended (80-bit) type, which gave correct results all the time. — Rudy Velthuis, Jun 04 '17 at 11:56

Implement x^(e) in Extended Pascal, without using exponentiation operators

2 Answers2