We help IT Professionals succeed at work.

cube root in TurboPascal?

syvapump
syvapump asked
on
This question is simple - how to take a cube root of X (or whatever) in TurboPascal... Also i would appreciate if someone could tell me how the gradation works in TP. Thank you!
Comment
Watch Question

dbruntonQuid, Me Anxius Sum?  Illegitimi non carborundum.
BRONZE EXPERT

Commented:
From

http://www.maths.cam.ac.uk/undergrad/tripos/catam/catsl/node5.html

The following program calculates the cube root of numbers entered from the keyboard, and uses a function routine to assist with the calculation. The logic of this is considerably more complicated than the first program, but should present little difficulty to a mathematics student. In brief, it involves a trial-and-error method of making the function x3 - r as small as possible. Note that the main program comes after the function - in Pascal, dependent routines (procedures and functions) are always defined before the main program.

   
    PROGRAM Cuberoot;
    USES CTbase, CTcntrl, CTgraph, CTmaths, CTmisc;  {  CATSL units }
    VAR
      r, a, b, c, fa, fb, fc: real;
      nc: integer;        
    FUNCTION f (x: real): real;
    BEGIN
       f := x * x * x - r;                     { LHS is the function name   }
    END;
    BEGIN                                      { beginning of main program  }
      nc := 0;                                 { set counter for the main
                                                               loop to zero }
      REPEAT                                   { main (outer) loop          }
        nc := nc + 1;                          { add 1 to counter           }
        Write ('enter number   ');             { screen message             }
        Read (r);                              { read number from keyboard  }
        a := 0;                                { a and b bracket the root   }
        b := (r / abs(r)) * (1 + abs(r));      { abs returns absolute value }
        fa := f(a);                            { find f(a) and store        }
        fb := f(b);                            { find f(b) and store        }
        REPEAT  
          c := (a + b) / 2;                    { midpoint                   }
          fc := f(c);                          { find f(c) and store        }
          IF fc *  fa >= 0 THEN
            a := c                             { no semicolon before else   }
          ELSE
            b := c;
        UNTIL Abs(a - b)  <  1E-8 ;
        Writeln('Cube root is ', c : 12 : 5 );
     UNTIL nc = 10;
    END.

And what do you mean by gradation?

Commented:
These may help:

function Power(Number,Exponent: Double): Double;
begin
  Result:=Exp(Exponent*Ln(Number));
end;

function Root(Number,Exponent: Double): Double;
begin
  Result:=Exp((1/Exponent)*Ln(Number));
end;

What "graduation" do you mean?
My name is MudSystems Engineer

Commented:
>>What "graduation" do you mean?

SIN(x * (PI/180))
COS(x * (PI/180))

Commented:
syvapump, what's going on?

Commented:
>>>>What "graduation" do you mean?

>>SIN(x * (PI/180))
>>COS(x * (PI/180))
Isn't the thing above is called "Radian"?

Commented:
Ok, here is the deal:

for positive powers and real x, we can have:
{change the data type yourself if you like, but keep the output type the same as parameter type}

function power(x,exponent; real): real;
var sign: integer;
begin
  if x <> 0 then
    begin
      sign := ord(x=abs(x)) * 2 - 1;
      power := exp(ln(abs(x)) * exponent) * sign;
    end
  else
    power := 0;
end;

for complex numbers, we can have:
{this is the strict case for cube root, modify the thing for other powers yourself if you like}
{also, assume that the function above is defined in the following disscussion}

type rectcplx = record   {a + bi}
                   Re: real;
                   Im: real;
                end;
     anglcplx = record   {r * exp(i * theta)}
                   r: real;
                   theta: real;
                 end;

procedure r2a(x: rectcplx; var result: anglcplx);
{convert from rectangular to polar}
var L: real;  {absolute of x}
  begin
    result.theta := arctan(x.Im/x.Re);
    result.r := sqrt(sqr(x.Im) + sqr(x.Re));
    {make sure the angle is in the right quadrant}
    if (Im <0) and (Re> 0) then
      result.theta := result.theta + pi

    else if (Im <0) and (Re < 0) then
      result.theta := result.theta + pi

    else if (Im >0) and (Re< 0) then
      result.theta := 2 * pi + result.theta;

  end;

procedure a2r(x: anglcplx; var result: rectcplx);
{convert from polar to rectangular}
  begin
    result.Re := cos(x.theta) * x.r;
    result.Im := sin(x.theta) * x.r;
  end;

procedure fmod(x: real; divisor: real): real
{for calculate the remaining angles (restrict theta from 0 to 2 * pi)}
begin
  fmod = x - (trunc(x/divisor) * divisor);
end;


procedure cplxcuberoot(x: rectcplx; var root1: rectcplx; var root2: rectcplx; var root3: rectcplx);
  var t1,t2: anglcplx; {temporary spaces}
  begin
    t1 = r2a(x);
    t2.r = power(t1.r,(1/3)); {function defined in the first section about real numbers}
    t2.theta = t1.theta * (1/3);
    a2r(t2,root1);
    t2.theta := fmod(t2.theta + 2*pi*(1/3), 2*pi);
    a2r(t2,root2);
    t2.theta := fmod(t2.theta + 2*pi*(1/3), 2*pi);
    a2r(t2,root3);
    {for more roots you need more variable parameter to store them, so it is not demonstrated}
  end;
dbruntonQuid, Me Anxius Sum?  Illegitimi non carborundum.
BRONZE EXPERT

Commented:
k61824

It is normal at EE to post any proposed answers as comments and let the questioner decide whose answer is the best.

Can you change your answer to a comment?

Posting the question as an answer locks it away from other experts to comment.
All,
I am unlocking this question in preparation for cleanup.  I will return in 7 days to finalize this question.  Please leave any recommendations for the final state of this question, I will take all recommendations into consideration.  Failing any feedback, I may decide in 7 days to delete or PAQ this question with no refund.  Thanks.

SpideyMod
Community Support Moderator @Experts Exchange

Commented:
While some fellow experts have posted code, I believe that my solution was the first effective and mathematically correct one (dbrunton's method is a quite slow approach method which seems to come from some math studies of someone).
dbruntonQuid, Me Anxius Sum?  Illegitimi non carborundum.
BRONZE EXPERT

Commented:
Well, no.  Your solution doesn't quite work.

function Root(Number, Exponent: Double): Double;
begin
  Result := Exp((1/Exponent)*Ln(Number));
end;

In your posted function above you have Result which won't return anything.  

In the correction below I have changed Result to Root.  This works correctly.

***   ***   ***   ***

This completed code test the function

program X;
{$N+}
var
 number : double;

function Root(TheNumber, Exponent: Double): Double;
begin
  Root := Exp((1/Exponent)*Ln(TheNumber));
end;

begin
  Writeln('Give us a number');
  Readln(number);
  Writeln('The cube root of ', number:4:2, ' is ', Root(number, 3):4:2);
end.

***   ***   ***

And the points SHOULD go to AvonWyss.  :-)

I haven't tested it for negative numbers by the way.
Commented:
Of course, you are right. I wrote and tested that thing in Delphi, and forgot to fix the "result" variable... thanks for pointing that out.
Answered.  Thanks dbrunton for increasing the value of the PAQ.  If you go here, I've placed a 20 point Q for the extra value added:

http://www.experts-exchange.com/Programming/Programming_Languages/Pascal/Q_20586093.html

SpideyMod
Community Support Moderator @Experts Exchange