//implementation of few simple arithmetic procedures under big numbers
//not very efficient but working properly
//author: Michal Pilipczuk
//XIV Secondary High School, Warsaw, Poland

unit bignumber;
interface

type mstring=record
    zn:longint;
    war:string;
end;

type wuh=record
    aa:mstring;
    bb:mstring;
end;

procedure initialize;
//initializing whole unit, should be done before using any functions below

function wch(fu:longint):char;
function wlo(c:char):longint;
//functions for exporting longints into chars and vice versa

function przez2(a:string;spod:longint):string;
//dividing by 2 in any system (based on 'spod')

function zmienw2(b:string;spod:longint):string;
//changing number into binar system from system with base 'spod' (Note: spod should be even)

function mn(x,y:string):boolean;
//gives true if x<=y

function wyczysc(ab:string):string;
//clears zeros at the beginning of the number 

function bezkon(ab:string):string;
//clears zeros at the end

function bez(tu:string):string;
//clears last digit

function od(g,h:string):string;
//difference between g and h (Note: they should be in binar system and g should be greater or equal)

function nwduj(c,d:string):string;
//counting gcd of two numbers (Note: they should be in binar system)

function modulo(d,r:string):string;
//counting d mod r

function dod(a,b:string;pod:longint):string;
//adding a and b in system with base 'spod'

function zmien(co:string;dow:longint):string;
//changing binar number into any other system (based on 'dow')

function rand(max:string):string;
//random, but not very efficient

function mul(w,q:string):string;
//multiplication of w and q in binar system

function poteg(a,b,r:string):string;
//power a^b, everything modulo r (Note: a,b,r should be binar)

function jacobi(g,d:string):longint;
//counts jacobi symbol: g is up d is down, d should be odd

function sprawdz(qur:string):boolean;
//says if a string is a proper number

function longintuj(a:string;base:longint):longint;
//makes a longint from a string in a given positional system(base)

function millero(pod,po,pw:string):boolean;
//function used in counting miller's test

function mill(dosp:string;ile:longint):boolean;
//using miller's test finds out if dosp is prime. The number of tests is ile.

function divid(a,b:string):string;
//computes a/b

function megadod(a,b:mstring):mstring;
//adds two mstrings

function megamul(a,b:mstring):mstring;
//multiplies two mstrings

function rand_prim(bo:string):string;
//gives a random prime number lower or equal to bo

function megamodu(ga,gb:mstring):mstring;
//a modulo operation for mstrings

function odw(wej:wuh;paa:mstring):wuh;
//used in odwrot

function odwrot(t,p:string):string;
//computes t^(-1) mod p

function wyp(ll:mstring):string;
//gives mstring as it should be written

function strassen(p:string;ile:longint):boolean;
//using strassen's test finds out if dosp is prime. The number of tests is ile.

implementation

const MAX=10000;
//maximum number of digits of numbers used

var zera:array[0..MAX] of string;
//array of zeros: zera[i] contains string with i zeros :))))

procedure initialize;
var ze:string;
ku:longint;
begin
    randomize;
    ze:='';
    zera[0]:='';
    for ku:=1 to MAX do
    begin
	ze:=ze+'0';
	zera[ku]:=ze;
    end;
end;

function wch(fu:longint):char;
begin
    wch:=chr(fu+48);
end;

function wlo(c:char):longint;
begin
    wlo:=ord(c)-48;
end;

function przez2(a:string;spod:longint):string;
var k,i,u:longint;
pom:string;
begin
    k:=0;
    pom:='';
    for i:=1 to length(a) do
    begin
	u:=wlo(a[i])+k*spod;
	k:=u mod 2;
	if not((i=1) and ((u shr 1)=0)) then
	pom:=pom+wch(u shr 1);
    end;
    przez2:=pom;
end;

function zmienw2(b:string;spod:longint):string;
begin
    if b='' then zmienw2:='' else
    zmienw2:=zmienw2(przez2(b,spod),spod)+wch(wlo(b[length(b)]) mod 2);
end;

function mn(x,y:string):boolean;
var flag:boolean;
j:longint;
begin
    if length(x)<length(y) then mn:=true else
    begin
	if length(x)>length(y) then mn:=false
	else
	begin
	    flag:=true;
	    for j:=1 to length(x) do
	    begin
		flag:=(ord(x[j])<ord(y[j]));
		if not(x[j]=y[j]) then break;
	    end;
	    mn:=flag;
	end;
    end;
end;

function wyczysc(ab:string):string;
var poma:string;
flag:boolean;
kur:longint;
begin
    flag:=false;
    poma:='';
    for kur:=1 to length(ab) do
    begin
	if not(ab[kur]='0') then flag:=true;
	if flag then poma:=poma+ab[kur];
    end;
    wyczysc:=poma;
end;

function bezkon(ab:string):string;
var poma:string;
flag:boolean;
kur:longint;
begin
    flag:=false;
    poma:='';
    for kur:=length(ab) downto 1 do
    begin
	if not(ab[kur]='0') then flag:=true;
	if flag then poma:=ab[kur]+poma;
    end;
    bezkon:=poma;
end;

function bez(tu:string):string;
var pom:string;
o:longint;
begin
    pom:='';
    for o:=1 to length(tu)-1 do
	pom:=pom+tu[o];
    bez:=pom;
end;

function od(g,h:string):string;
var k,l,wu:longint;
pom,hu:string;
begin
    hu:=zera[length(g)-length(h)]+h;
    k:=0;
    pom:='';
    for l:=length(g) downto 1 do
    begin
	wu:=wlo(g[l])-wlo(hu[l])-k;
	if wu<0 then
	begin
	    k:=1;
	    wu:=wu+2;
	end
	else
	begin
	    k:=0;
	end;
	pom:=wch(wu)+pom;
    end;
    pom:=wyczysc(pom);
    od:=pom;
end;


function nwduj(c,d:string):string;
begin
    if c=d then nwduj:=c
    else
    begin
	if mn(c,d) then nwduj:=nwduj(d,c)
	else
	begin
	    if ((c[length(c)]='0') and (d[length(d)]='0')) then nwduj:=nwduj(bez(c),bez(d))+'0'
	    else
	    begin
		if c[length(c)]='0' then nwduj:=nwduj(bez(c),d)
		else
		begin
		    if d[length(d)]='0' then nwduj:=nwduj(c,bez(d))
		    else
			nwduj:=nwduj(d,od(c,d));
		end;
	    end;
	end;
    end;
end;

function modulo(d,r:string):string;
begin
    if ((d=r) or (d='')) then modulo:=''
    else
    begin
	if mn(d,r) then modulo:=d
	else
	begin
	    if length(r)<(length(d)-1) then modulo:=modulo((od(d,(r+zera[length(d)-length(r)-1]))),r)
	    else modulo:=modulo(od(d,r),r);
	end;
    end;
end;

function dod(a,b:string;pod:longint):string;
var bu,pom:string;
i,k,l:longint;
begin
    if mn(a,b) then dod:=dod(b,a,pod)
    else
    begin
	bu:=zera[length(a)-length(b)]+b;
	k:=0;
	pom:='';
	for i:=length(a) downto 1 do
	begin
	    l:=wlo(a[i])+wlo(bu[i])+k;
	    if l>(pod-1) then
	    begin
		l:=l-pod;
		k:=1;
	    end else k:=0;
	    pom:=wch(l)+pom;
	end;
	if k=1 then pom:='1'+pom;
	dod:=pom;
    end;
end;

function zmien(co:string;dow:longint):string;
var akt,sum:string;
v:longint;
begin
    akt:='1';
    sum:='';
    for v:=length(co) downto 1 do
    begin
	if co[v]='1' then sum:=dod(sum,akt,dow);
	akt:=dod(akt,akt,dow);
    end;
    zmien:=sum;
end;

function rand(max:string):string;
var t:longint;
pom:string;
begin
    pom:='';
    for t:=1 to length(max) do
    begin
	pom:=pom+wch(random(2));
    end;
    pom:=wyczysc(pom);
    if ((mn(max,pom)) or (pom='') or (max=pom)) then rand:=rand(max) else 
	rand:=pom;
end;

function mul(w,q:string):string;
var pom:string;
k:longint;
begin
    if mn(w,q) then mul:=mul(q,w)
    else
    begin
	pom:='';
	for k:=length(q) downto 1 do
	    if q[k]='1' then pom:=dod(pom,(w+zera[length(q)-k]),2);
	mul:=pom;
    end;
end;

function poteg(a,b,r:string):string;
var akt,pocz:string;
k:longint;
begin
    if ((a=r) or (a='')) then poteg:=''
    else
    begin
	if not(mn(a,r)) then poteg:=poteg(modulo(a,r),b,r)
	else
	begin
	    akt:=a;
	    pocz:='1';
	    for k:=length(b) downto 1 do
	    begin
		if b[k]='1' then
		begin
		    pocz:=modulo((mul(pocz,akt)),r);
		    end;
		akt:=modulo((mul(akt,akt)),r);
	    end;
	    poteg:=pocz;
	end;
    end;
end;

function resztymod8(a:string):longint;
var u:longint;
begin
    if a='1' then exit(1);
    if a='11' then exit(-1);
    if a='101' then exit(-1);
    if a='111' then exit(1);
    u:=length(a);
    resztymod8:=resztymod8(wyczysc(a[u-2]+a[u-1]+a[u]));
end;

function resztymod4(a:string):longint;
var u:longint;
begin
    if a='1' then exit(1);
    if a='11' then exit(-1);
    u:=length(a);
    resztymod4:=resztymod4(wyczysc(a[u-1]+a[u]));
end;

function faf(oi,pu:string):longint;
begin
    if (resztymod4(oi)+resztymod4(pu))<0 then faf:=(-1) else faf:=1;
end;

function jacobi(g,d:string):longint;
begin
    if g='' then jacobi:=0
    else
    begin
	if g='1' then jacobi:=1
	else
	begin
	    if ((d=g) or (mn(d,g))) then jacobi:=jacobi(modulo(g,d),d)
	    else
	    begin
		if g[length(g)]='0' then jacobi:=resztymod8(d)*jacobi(bez(g),d)
		else
		begin
		    jacobi:=(faf(g,d)*jacobi(d,g));
		end; 
	    end;
	end;
    end;
end;

function sprawdz(qur:string):boolean;
var f,h:longint;
flag:boolean;
begin
    flag:=true;
    for f:=1 to length(qur) do
    begin
	h:=wlo(qur[f]);
	if ((h<0) or (h>9)) then 
	begin
	    flag:=false;
	    break;
	end;
    end;
    sprawdz:=flag;
end;

function longintuj(a:string;base:longint):longint;
var f,g,h:longint;
begin
    f:=0;
    h:=1;
    for g:=length(a) downto 1 do
    begin
	f:=f+wlo(a[g])*h;
	h:=h*base;
    end;
    longintuj:=f;
end;

function millero(pod,po,pw:string):boolean;
var w1,w2:string;
flag2:boolean;
l:longint;
begin
    w1:='1';
    flag2:=true;
    for l:=1 to length(pw) do
    begin
	w2:=modulo(mul(w1,w1),po);
	if ((w2='1') and (not((w1='1') or (w1=pw)))) then
	begin
	    flag2:=false;
	    break;
	end;
	if pw[l]='1' then w1:=modulo(mul(w2,pod),po) else w1:=w2;
    end;
    millero:=((flag2)and(w1='1'));
end;

function mill(dosp:string;ile:longint):boolean;
var ggg,fff:string;
t:longint;
flag:boolean;
begin
    if dosp[length(dosp)]='0' then mill:=false else
    begin
	ggg:=od(dosp,'1');
	flag:=true;
	for t:=1 to ile do
	begin
	    fff:=rand(ggg);
	    flag:=millero(fff,dosp,ggg);
	    if not(flag) then break;
	end;
	mill:=flag;
    end;
end;

function divid(a,b:string):string;
var wynf,za:string;
begin
    if b='' then divid:='error' else
    begin
	wynf:=zera[length(a)];
	while (mn(b,a) or (a=b)) do
	begin
	    za:=b+zera[length(a)-length(b)];
	    if mn(a,za) then 
	    begin
		za:=b+zera[length(a)-length(b)-1];
		wynf[length(wynf)-length(a)+1+length(b)]:='1';
	    end
	    else
	    begin
		wynf[length(wynf)-length(a)+length(b)]:='1';
	    end;
	    a:=od(a,za);
	end;
	divid:=wyczysc(wynf);
    end;
end;

function megadod(a,b:mstring):mstring;
var pomis:mstring;
begin
    if a.zn<b.zn then megadod:=megadod(b,a) else
    begin
	if ((a.zn=1) and (b.zn=1)) then
	begin
	    pomis.zn:=1;
	    pomis.war:=dod(a.war,b.war,2);
	    megadod:=pomis;
	end
	else
	begin
	    if ((a.zn=(-1)) and (b.zn=(-1))) then
	    begin
		pomis.zn:=(-1);
		pomis.war:=dod(a.war,b.war,2);
		megadod:=pomis;
	    end
	    else
	    begin
		if mn(a.war,b.war) then
		begin
		    pomis.zn:=(-1);
		    pomis.war:=od(b.war,a.war);
		    megadod:=pomis;
		end
		else
		begin
		    pomis.zn:=1;
		    pomis.war:=od(a.war,b.war);
		    megadod:=pomis;
		end;
	    end;
	end;
    end;
end;

function megamul(a,b:mstring):mstring;
var pomis:mstring;
begin
        pomis.zn:=a.zn*b.zn;
	pomis.war:=mul(a.war,b.war);
	megamul:=pomis;
end;

function megamodu(ga,gb:mstring):mstring;
var fw:mstring;
begin
    fw.zn:=ga.zn;
    fw.war:=modulo(ga.war,gb.war);
    megamodu:=fw;
end;

function rand_prim(bo:string):string;
var wek:string;
begin
    while (1=1) do
    begin
	wek:=rand(bo);
	if mill(wek,10) then break;
    end;
    rand_prim:=wek;
end;

function wyp(ll:mstring):string;
begin
    if ll.zn<0 then wyp:=('-'+zmien(ll.war,10)) else wyp:=zmien(ll.war,10);
end;

function odw(wej:wuh;paa:mstring):wuh;
var ffk,ffw:wuh;
pa,pb,pwa,pwb,apb:mstring;
begin
    pa:=wej.aa;
    pb:=wej.bb;
    if pb.war='' then
    begin
	pwa.war:='1';
	pwa.zn:=1;
	pwb.war:='';
	pwb.zn:=1;
	ffw.aa:=pwa;
	ffw.bb:=pwb;
	odw:=ffw;
    end
    else
    begin
	apb.war:=divid(pa.war,pb.war);
	apb.zn:=pa.zn*pb.zn*(-1);
	ffw.aa:=pb;
	ffw.bb:=megamodu(pa,pb);
	ffk:=odw(ffw,paa);
	ffw.aa:=ffk.bb;
	ffw.bb:=megadod(ffk.aa,megamul(apb,ffk.bb));
	ffw.aa:=megamodu(ffw.aa,paa);
	ffw.bb:=megamodu(ffw.bb,paa);
	odw:=ffw;
    end;
end;

function odwrot(t,p:string):string;
var up,cona,conb,pi:mstring;
cons:wuh;
begin
    cona.zn:=1;
    cona.war:=p;
    conb.zn:=1;
    conb.war:=t;
    cons.aa:=cona;
    cons.bb:=conb;
    pi.zn:=1;
    pi.war:=p;
    cons:=odw(cons,pi);
    up:=cons.bb;
    if up.zn<0 then up.war:=od(p,up.war);
    odwrot:=up.war;
end;

function strassen(p:string;ile:longint):boolean;
var flag3:boolean;
k,iiii:longint;
ax,boo,uu,ux:string;
begin
    if p[length(p)]='0' then strassen:=false
    else
    begin
	flag3:=true;
	ux:=od(p,'1');
	uu:=bez(ux);
	for iiii:=1 to ile do
	begin	    
	    ax:=rand(p);
	    k:=jacobi(ax,p);
	    boo:=poteg(ax,uu,p);
	    if (not(((k=1)and(boo='1')) or ((k=(-1))and(boo=ux)))) then
	    begin
		flag3:=false;
		break;
	    end;
	end;
	strassen:=flag3;
    end;
end;

end.