IMPLEMENTATION MODULE CombLib;
(*------------------------------------------------------------------------*)
(* Combination of numbers / Kombinatorische Routinen *)
(* *)
(* Letzte Veraenderung *)
(* *)
(* 02.01.18, MRi: Erstellen der ersten Version von Combination *)
(* 03.01.18, MRi: Erstellen der ersten Version von NoverK *)
(* 04.01.18, MRi: Erstellen der ersten Version von NextP *)
(* 10.12.18, MRi: Erstellen der ersten Version von Comb *)
(*------------------------------------------------------------------------*)
(* Implementation : Michael Riedl *)
(* Licence : GNU Lesser General Public License (LGPL) *)
(*------------------------------------------------------------------------*)
(* $Id: CombLib.mod,v 1.1 2018/12/10 09:14:39 mriedl Exp mriedl $ *)
FROM LongMath IMPORT exp;
IMPORT LMathLib;
FROM SpezFunkt2 IMPORT dLnGamma;
PROCEDURE NoverK(N,K : INTEGER) : LONGREAL;
(*----------------------------------------------------------------*)
(* Computes the combinatorial coefficient c(n,k). *)
(* c(n,k) is the number of distinct combinations of k objects *)
(* chosen from a set of n distinct objects. Within a combination *)
(* the order does not matter. *)
(* *)
(* c(n,k) := n! / ( (n-k)! * k! ); *)
(* *)
(* This procedure is base on subroutine comnin of john burkardt *)
(*----------------------------------------------------------------*)
VAR Cnk,fack,facn,facnmk : LONGREAL;
BEGIN
IF (N < 0) THEN
Cnk:=0.0;
ELSIF (K = 0) THEN
Cnk:=1.0;
ELSIF (K = 1) THEN
Cnk:=LFLOAT(N);
ELSIF (1 < K) AND (K < N-1) THEN
facn := dLnGamma(LFLOAT(N+1));
fack := dLnGamma(LFLOAT(K+1));
facnmk := dLnGamma(LFLOAT(N-K+1));
Cnk:=LFLOAT(VAL(INTEGER,(exp(facn-fack-facnmk)+0.49)));
ELSIF (K = N-1) THEN
Cnk:=LFLOAT(N);
ELSIF (K = N) THEN
Cnk:=1.0;
ELSE
Cnk:=0.0;
END;
RETURN Cnk;
END NoverK;
PROCEDURE Combination(VAR Comb : ARRAY OF ARRAY OF SHORTCARD;
VAR NComb : INTEGER;
NSet : INTEGER;
IMax : INTEGER;
VAR iFehl : INTEGER);
VAR more : BOOLEAN;
keep,last,m,n : INTEGER;
ncomb : CARDINAL;
BEGIN
ncomb := VAL(INTEGER,NoverK(IMax,NSet));
IF (ncomb >= HIGH(Comb)) OR (VAL(CARDINAL,NSet) >= HIGH(Comb[0])) THEN
iFehl:=1;
RETURN;
END;
iFehl:=0;
NComb := 1;
FOR m:=1 TO NSet DO
Comb[NComb-1,m-1]:=m; (* the first combination: 1,2,3,...,many *)
END;
(* set up an infinite loop in which to pick up the next combination *)
last:=0;
more:=TRUE;
WHILE more DO
more:=FALSE; (* assume that no more combination is to be found *)
m:=NSet;
LOOP
IF (m < 1) THEN EXIT; END;
IF (Comb[NComb-1,m-1] < IMax-NSet+m) THEN
more:=TRUE;
last:=m;
EXIT
END; (* IF *)
DEC(m);
END;
IF more THEN
INC(NComb);
FOR n:=1 TO last-1 DO
Comb[NComb-1,n-1]:=Comb[NComb-2,n-1];
END;
keep:=Comb[NComb-2,last-1];
FOR n:=last TO NSet DO
INC(keep);
Comb[NComb-1,n-1]:=keep;
END;
END; (* IF more *)
END; (* WHILE *)
END Combination;
PROCEDURE NextP( N : INTEGER;
VAR A : ARRAY OF INTEGER) : BOOLEAN;
VAR i,j,k,t : INTEGER;
nextp : BOOLEAN;
BEGIN
IF (N <= 1) THEN RETURN FALSE; END;
i := N - 1;
WHILE ((i # 0) AND( A[i-1] >= A[i])) DO DEC(i); END;
j := i + 1;
k := N;
LOOP
t := A[j-1]; A[j-1] := A[k-1]; A[k-1] := t;
j:=j + 1;
k:=k - 1;
IF (j >= k) THEN
IF (i = 0) THEN nextp:=FALSE; EXIT; END;
j := i+1;
WHILE (A[j-1] < A[i-1]) DO INC(j); END;
t := A[i-1]; A[i-1] := A[j-1]; A[j-1] := t;
nextp:=TRUE;
EXIT;
END; (* IF *)
END; (* LOOP *)
RETURN nextp;
END NextP;
PROCEDURE CombNK(N,K : CARDINAL) : CARDINAL;
VAR i,m,maxKm : CARDINAL;
numer,denom : LONGREAL;
BEGIN
IF (K > N) THEN
RETURN 0;
ELSE
m := N - K;
numer := 1.0;
maxKm := LMathLib.MaxCard(K,m);
FOR i:=N TO maxKm+1 BY - 1 DO
numer := numer*LFLOAT(i);
END;
denom := LMathLib.fact(LMathLib.MinCard(K,m));
RETURN VAL(CARDINAL,INT(numer/denom));
END; (* IF *)
END CombNK;
END CombLib.