Code covered by the BSD License  

Highlights from
Variable Precision Integer Arithmetic

  • demo_vpi
  • base2vpi(B,base) bin2vpi: converts an integer in an arbitrary base into vpi (decimal) form
  • bin2vpi(B) bin2vpi: converts a binary representation of an integer into vpi (decimal) form
  • binomfactors(n,k) binomfactors: list all factors of the binomial coefficient nchoosek(n,k)
  • catdigits(N,M) catdigits: concatenates the digits of N and M into an aggregate number
  • createPrimesList createPrimesList - For users of older matlab releases, this function will generate a compatible _primeslist_ file
  • factorialfactors(n) factorialfactors: efficient computation of the prime factors of factorial(n)
  • fibonacci(n) fibonacci: vpi tool to efficiently compute the n'th Fibonacci number and the n'th Lucas number
  • getprimeslist loads the primeslist file, and decompresses it, returning the list of primes up to 2^26
  • ispalindrome(N) ispalindrome: test if the number N (vpi or numeric, or a digit string as a vector) is a palindrome
  • iszero(INT) vpi/iszero: test to see if a numeric object is zero
  • legendresymbol(a,p) legendresymbol: computes the legendre symbol (a/p) for prime p
  • lineardiophantine(A,B,C) lineardiophantine: solve the linear Diophantine equation, A*x + B*y = C
  • mersenne(p) mersenne: identify whether 2^p-1 is a Mersenne prime, using the Lucas-Lehmer test
  • minv(a,p)
  • modfibonacci(n,modulus) fibonacci: compute the n'th Fibonacci number and the n'th Lucas number, all modulo a given value
  • modrank(A,p) modrank: compute the rank of an integer array, modulo p
  • modroot(a,p)
  • modsolve(A,rhs,p)
  • nextprime(N,direction,kprimes) nextprime: finds the next larger prime number directly above (or below) N
  • numberOfPartitions(N) numberOfPartitions: compute the number of partitions of the positive integer n
  • powermod(a,d,n) vpi/powermod: Compute mod(a^d,n)
  • quadraticresidues(N) quadraticresidues: returns a list of the possible quadratic residues of the integer N
  • quotient(numerator,denominato... quotient: divides two integers, computing a quotient and remainder
  • subfactorial(N) subfactorial: The subfactorial of an integer (or integers) N, known as !N
  • totient(N) vpi/totient: the number of positive integers less than N that are coprime to N
  • vpi(N) vpi: Creator function for a variable precision integer
  • View all files
from Variable Precision Integer Arithmetic by John D'Errico
Arithmetic with integers of fully arbitrary size. Arrays and vectors of vpi numbers are supported.

mersenne(p)
function [tf,S] = mersenne(p)
% mersenne: identify whether 2^p-1 is a Mersenne prime, using the Lucas-Lehmer test
% usage: tr = mersenne(p)
%
% This is a simple test to apply, working fairly
% easily for Mersenne primes up to p at least a 
%
% http://mathworld.wolfram.com/Lucas-LehmerTest.html
% http://en.wikipedia.org/wiki/Mersenne_prime
%
% arguments: (input)
% p  - a (relatively small) prime
%
% arguments: (output)
% tf - boolean, true if 2^p-1 is prime
%
% Example:
%
%
%  See also: isprime
%  
% 
%  Author: John D'Errico
%  e-mail: woodchips@rochester.rr.com
%  Release: 1.0
%  Release date: 1/23/09

% ensure that p is a vpi
dp = double(p);

% Save the terms from the Lucas-Lehmer
% recurrence relation, just in case they
% are of interest.
S = cell(1,dp-1);

if ~isprime(dp)
  warning('p must be prime for 2^p-1 to be a mersenne prime')
  tf = false;
  return
end
p = vpi(p);

% compute the candidate Mersenne prime
mp = vpi(2)^p-1;

% S0 = 4
S = repmat(vpi(4),1,dp-1);

% and loop, from n = 1 to (p-2)
h = waitbar(0,'Stop bothering me. I''m thinking, can''t you see?');
for n = 1:(dp-2)
  waitbar(n/dp,h)
  S(n+1) = mod(S(n)*S(n) - 2,mp);
end
delete(h)

% mp is prime IFF this last term is zero
tf = (0 == S(dp-1));



Contact us at files@mathworks.com