Seed7 
 FAQ 
 Manual 
 Screenshots 
 Examples 
 Algorithms 
 Download 
 Links 

 Algorithms 
 Sorting 
 Searching 
 Date & Time 
 String 
 Mathematics 
 File 
 Puzzles 
 Others 
 Algorithms 
Mathematics
 previous   up   next 

Determine wether a given integer number is prime

const func boolean: is_prime (in integer: number) is func
  result
    var boolean: result is FALSE;
  local
    var integer: count is 2;
  begin
    if number = 2 then
      result := TRUE;
    elsif number > 2 then
      while number rem count <> 0 and count * count <= number do
        incr(count);
      end while;
      result := number rem count <> 0;
    end if;
  end func;

Function to calculate the prime factors of a number

const func array integer: factorise (in var integer: number) is func
  result
    var array integer: result is 0 times 0;
  local
    var integer: checker is 2;
  begin
    while checker * checker <= number do
      if number rem checker = 0 then
        result &:= [](checker);
        number := number div checker;
      else
        incr(checker);
      end if;
    end while;
    if number <> 1 then
      result &:= [](number);
    end if;
  end func;

Verify if a given Mersenne number is prime

A Mersenne number is a number that is one less than a power of two (2**p-1). The Lucas-Lehmer test allows to check if a given Mersenne number is prime: For the prime p, the Mersenne number 2**p-1 is prime if and only if S(p-1) rem 2**p-1 = 0 where S(1)=4 and S(n)=S(n-1)**2-2.

const func boolean: mersenne_number_is_prime (in integer: p) is func
  result
    var boolean: result is FALSE;
  local
    var bigInteger: m_p is 0_;
    var bigInteger: s is 4_;
    var integer: i is 0;
  begin
    m_p := 2_ ** p - 1_;
    for i range 2 to pred(p) do
      s := (s ** 2 - 2_) rem m_p;
    end for;
    result := s = 0_;
  end func;

Determine the greatest common divisor of two positive integer numbers

const func integer: gcd (in var integer: a, in var integer: b) is func
  result
    var integer: result is 0;
  local
    var integer: help is 0;
  begin
    while a <> 0 do
      help := b rem a;
      b := a;
      a := help;
    end while;
    result := b;
  end func;

Determine the least common multiple of two positive integer numbers

const func integer: lcm (in integer: a, in integer: b) is
  return a * b div gcd(a, b);

Compute PI with 1000 decimal digits using John Machin's formula

$ include "seed7_05.s7i";
  include "bigint.s7i";
  include "bigrat.s7i";

# John Machin's formula from 1706:
# PI = 16 * arctan(1 / 5) - 4 * arctan(1 / 239)

# Taylor series for arctan:
# arctan(x) = sum_n_from_0_to_inf((-1) ** n / succ(2 * n) * x ** succ(2 * n))

# Taylor series of John Machin's formula:
# PI = sum_n_from_0_to_inf(16 * (-1) ** n / (succ(2 * n) *   5 ** succ(2 * n)) -
#                           4 * (-1) ** n / (succ(2 * n) * 239 ** succ(2 * n)))

const func bigRational: compute_pi_machin is func
  result
    var bigRational: sum is 0_ / 1_;
  local
    var integer: n is 0;
  begin
    for n range 0 to 713 do
      sum +:= 16_ * (-1_) ** n / (bigInteger conv succ(2 * n) *   5_ ** succ(2 * n)) -
               4_ * (-1_) ** n / (bigInteger conv succ(2 * n) * 239_ ** succ(2 * n));
    end for;
  end func;

 const proc: main is func
   begin
     writeln(compute_pi_machin digits 1000);
   end func;

Compute PI with 1000 decimal digits using the Bailey/Borwein/Plouffe formula

$ include "seed7_05.s7i";
  include "bigint.s7i";
  include "bigrat.s7i";

# In 1997, David H. Bailey, Peter Borwein and Simon Plouffe published a
# paper (Bailey, 1997) on a new formula for PI as an infinite series:

# PI = sum_n_from_0_to_inf(16 ** (-n) *
#      (4 / (8 * n + 1) - 2 / (8 * n + 4) - 1 / (8 * n + 5) - 1 / (8 * n + 6)))

const func bigRational: compute_pi_bailey_borwein_plouffe is func
  result
    var bigRational: sum is 0_ / 1_;
  local
    var integer: n is 0;
    var bigInteger: k is 0_;
  begin
    for n range 0 to 825 do
      k := bigInteger conv n;
      sum +:= 1_ / 16_ ** n *
          (4_ / (8_ * k + 1_) - 2_ / (8_ * k + 4_) - 1_ / (8_ * k + 5_) - 1_ / (8_ * k + 6_));
    end for;
  end func;

 const proc: main is func
   begin
     writeln(compute_pi_bailey_borwein_plouffe digits 1000);
   end func;

Write PI with 1000 decimal digits using Newtons formula

$ include "seed7_05.s7i";

# Newtons formula for PI is:
# PI / 2 = sum_n_from_0_to_inf(n! / (2 * n + 1)!!)

# This can be written as:
# PI / 2 = 1 + 1/3 * (1 + 2/5 * (1 + 3/7 * (1 + 4/9 * (1 + ... ))))

# This algorithm puts 2 * 1000 on the right side and computes everything from inside out.

const integer: SCALE is 10000;
const integer: MAXARR is 3500;
const integer: ARRINIT is 2000;

const proc: main is func
  local
    var integer: i is 0;
    var integer: j is 0;
    var integer: carry is 0;
    var array integer: arr is MAXARR times ARRINIT;
    var integer: sum is 0;
  begin
    for i range MAXARR downto 1 step 14 do
      sum := 0;
      for j range i downto 1 do
        sum := sum*j + SCALE*arr[j];
        arr[j] := sum rem pred(j*2);
        sum := sum div pred(j*2);
      end for;
      write(carry + sum div SCALE lpad0 4);
      carry := sum rem SCALE;
    end for;
    writeln;
  end func;

Determine the truncated square root of a big integer number

The bigInteger type already has an integer square root which can be called with 'sqrt(number)'. A hand written version would be:

const func bigInteger: intSqrt (in var bigInteger: number) is func
  result
    var bigInteger: result is 0_;
  local
    var bigInteger: res2 is 0_;
  begin
    if number > 0_ then
      res2 := number;
      repeat
        result := res2;
        res2 := (result + number div result) div 2_;
      until result <= res2;
    end if;
  end func;

Exponentiation function for integer numbers

The integer type already has an exponentiation operator which can be called with 'base ** exponent'. A hand written version would be:

const func integer: intPow (in var integer: base, in var integer: exponent) is func
  result
    var integer: result is 0;
  begin
    if exponent < 0 then
      raise(NUMERIC_ERROR);
    else
      if odd(exponent) then
        result := base;
      else
        result := 1;
      end if;
      exponent := exponent div 2;
      while exponent <> 0 do
        base *:= base;
        if odd(exponent) then
          result *:= base;
        end if;
        exponent := exponent div 2;
      end while;
    end if;
  end func;

Multiply integers using the peasant multiplication

const func integer: peasantMult (in var integer: a, in var integer: b) is func
  result
    var integer: result is 0;
  begin
    while a <> 0 do
      if odd(a) then
        result +:= b;
      end if;
      a := a div 2;
      b *:= 2;
    end while;
  end func;

Matrix addition

const type: matrix is array array float;

const func matrix: (in matrix: left) + (in matrix: right) is func
  result
    var matrix: result is matrix.value;
  local
    var integer: i is 0;
    var integer: j is 0;
  begin
    if length(left) <> length(right) or
        length(left[1]) <> length(right[1]) then
      raise RANGE_ERROR;
    else
      result := length(left) times length(left[1]) times 0.0;
      for i range 1 to length(left) do
        for j range 1 to length(left[i]) do
          result[i][j] := left[i][j] + right[i][j];
        end for;
      end for;
    end if;
  end func;

Matrix multiplication

const type: matrix is array array float;

const func matrix: (in matrix: left) * (in matrix: right) is func
  result
    var matrix: result is matrix.value;
  local
    var integer: i is 0;
    var integer: j is 0;
    var integer: k is 0;
    var float: accumulator is 0.0;
  begin
    if length(left[1]) <> length(right) then
      raise RANGE_ERROR;
    else
      result := length(left) times length(right[1]) times 0.0;
      for i range 1 to length(left) do
        for j range 1 to length(right) do
          accumulator := 0.0;
          for k range 1 to length(left) do
            accumulator +:= left[i][k] * right[k][j];
          end for;
          result[i][j] := accumulator;
        end for;
      end for;
    end if;
  end func;

Random number generator

This random generator delivers random numbers between 0 and 2**32-1. To get different pseudorandom sequence seed can be initialized with another value.

var bigInteger: seed is 987654321_;

const func bigInteger: rand_32 is func
  result
    var integer: result is 0;
  begin
    seed := (seed * 1073741825_ + 1234567891_) rem 2_ ** 64;
    result := seed div 2_ ** 32;
  end func;

Recursive fibonacci function

const func integer: fib (in integer: number) is func
  result
    var integer: result is 1;
  begin
    if number > 2 then
      result := fib(pred(number)) + fib(number - 2);
    elsif number = 0 then
      result := 0;
    end if;
  end func;

Display 100 numbers of the fibonacci sequence

$ include "seed7_05.s7i";
  include "bigint.s7i";
 
const proc: main is func
  local
    var integer: i is 0;
    var bigInteger: a is 0_;
    var bigInteger: b is 1_;
    var bigInteger: c is 0_;
  begin
    writeln(a);
    for i range 1 to 100 do
      writeln(b);
      c := a;
      a := b;
      b +:= c;
    end for;
  end func;

The tak function

const func integer: tak (in integer: x, in integer: y, in integer: z) is func
  result
    var integer: result is 0;
  begin
    if y >= x then
      result := z;
    else
      result := tak(tak(pred(x), y, z),
                    tak(pred(y), z, x),
                    tak(pred(z), x, y));
    end if;
  end func;

The ackermann function

const func integer: ackermann (in integer: m, in integer: n) is func
  result
    var integer: result is 0;
  begin
    if m = 0 then
      result := succ(n);
    elsif n = 0 then
      result := ackermann(pred(m), 1);
    else
      result := ackermann(pred(m), ackermann(m, pred(n)));
    end if;
  end func;


 previous   up   next