PROGRAM d10r4(input,output);
(* driver for routine DBRENT *)
CONST
   tol=1.0e-6;
   eql=1.0e-4;
VAR
   ax,bx,cx,fa,fb,fc,xmin,dbr : real;
   i,iflag,j,nmin : integer;
   amin : ARRAY [1..20] OF real;

(*$I MODFILE.PAS *)
(*$I BESSJ1.PAS *)

FUNCTION dfunc(x: real): real;
BEGIN
   dfunc := -bessj1(x)
END;

(*$I BESSJ0.PAS *)

FUNCTION func(x: real): real;
BEGIN
   func := bessj0(x)
END;

(*$I MNBRAK.PAS *)

(*$I DBRENT.PAS *)

BEGIN
   nmin := 0;
   writeln;
   writeln('minima of the function bessj0');
   writeln('min. #':10,'x':8,
      'bessj0(x)':16,'bessj1(x)':12,'DBRENT':11);
   FOR i := 1 to 100 DO BEGIN
      ax := i;
      bx := i+1.0;
      mnbrak(ax,bx,cx,fa,fb,fc);
      dbr := dbrent(ax,bx,cx,tol,xmin);
      IF (nmin = 0) THEN BEGIN
         amin[1] := xmin;
         nmin := 1;
         writeln(nmin:7,xmin:15:6,func(xmin):12:6,
            dfunc(xmin):12:6,dbr:12:6);
      END ELSE BEGIN
         iflag := 0;
         FOR j := 1 to nmin DO BEGIN
            IF  (abs(xmin-amin[j]) <= eql*xmin)
               THEN iflag := 1
         END;
         IF (iflag = 0) THEN BEGIN
            nmin := nmin+1;
            amin[nmin] := xmin;
            writeln(nmin:7,xmin:15:6,func(xmin):12:6,
               dfunc(xmin):12:6,dbr:12:6)
         END
      END
   END
END.
PROGRAM d5r2(input,output);
(* driver for routine DDPOLY *)
(* polynomial (x-1)**5 *)
CONST
   nc=6;
   nd=5;   (* nd=nc-1 *)
   np=20;
TYPE
   words = PACKED ARRAY [1..15] OF char;
   glcarray = ARRAY [1..nc] OF real;
   glpdarray = ARRAY [1..nd] OF real;
VAR
   i,j : integer;
   x : real;
   c : glcarray;
   pd : glpdarray;
   d : ARRAY [1..nd,1..np] OF real;
   a : ARRAY [1..nd] OF words;
   glntop : integer;
   gla : ARRAY [1..33] OF real;

FUNCTION power(x: real; n: integer): real;
BEGIN
   IF (n=0) THEN power := 1.0 ELSE power := x*power(x,n-1)
END;

(*$I MODFILE.PAS *)
(*$I GAMMLN.PAS *)

(*$I FACTRL.PAS *)

(*$I DDPOLY.PAS *)

BEGIN
   glntop := 0;
   gla[1] := 1.0;
   a[1] := 'polynomial:    '; a[2] := 'first deriv:   ';
   a[3] := 'second deriv:  '; a[4] := 'third deriv:   ';
   a[5] := 'fourth deriv:  ';
   c[1] := -1.0; c[2] := 5.0; c[3] := -10.0;
   c[4] := 10.0; c[5] := -5.0; c[6] := 1.0;
   FOR i := 1 to np DO BEGIN
      x := 0.1*i;
      ddpoly(c,nc,x,pd,nc-1);
      FOR j := 1 to nc-1 DO BEGIN
         d[j,i] := pd[j]
      END
   END;
   FOR i := 1 to nc-1 DO BEGIN
      writeln(' ':7,a[i]);
      writeln('x':12,'DDPOLY':17,'actual':15);
      FOR j := 1 to np DO BEGIN
         x := 0.1*j;
         writeln(x:15:6,d[i,j]:15:6,
           (factrl(nc-1)/factrl(nc-i))*power(x-1.0,nc-i):15:6)
      END;
      writeln('press ENTER to continue...');
      readln
   END
END.
PROGRAM d7r13 (input,output,dfile);
(* driver for routine DES *)
LABEL 5,10,99;
TYPE
   gl32array = ARRAY [1..32] OF integer;
   gl48array = ARRAY [1..48] OF integer;
   gl56array = ARRAY [1..56] OF integer;
   gl64array = ARRAY [1..64] OF integer;
VAR
   glicd,ipc1 : gl56array;
   i,idirec,idum,j,jdum,k,l : integer;
   m,mm,nciphr,newkey,twopwr : integer;
   icmp,iin,iout,key,ip,ipm : gl64array;
   hin,hkey,hout,hcmp : string[17];
   verdct : string[8];
   txt : string[60];
   txt2 : string[6];
   ipc2,iet : gl48array;
   ipp : gl32array;
   is : ARRAY [1..16,1..4,1..8] OF integer;
   ibin : ARRAY [1..4,1..16] OF integer;
   ksflg,cyflg,desflg : boolean;
   dfile : text;

(*$I MODFILE.PAS *)

(*$I DESKS.PAS *)

(*$I DES.PAS *)

FUNCTION hex2int(ch:char) : integer;
(* Coverts character representing hexadecimal number to its integer
value in a machine-independent way. *)
BEGIN
   IF (ch >= '0') AND (ch <= '9') THEN
      hex2int := ord(ch)-ord('0')
   ELSE
      hex2int := ord(ch)-ord('A')+10
END;
FUNCTION int2hex(i:integer) : char;
(* Inverse of hex2int *)
BEGIN
   IF i <= 9 THEN
      int2hex := chr(i+ord('0'))
   ELSE
      int2hex := chr(i-10+ord('A'))
END;

BEGIN
   ksflg := true;
   cyflg := true;
   desflg := true;
   glopen(dfile,'destst.dat');
   readln (dfile,txt);
   writeln (txt);
5:   readln (dfile,txt);
   IF eof(dfile) THEN GOTO 99;
   writeln (txt);
   readln (dfile,nciphr);
   readln (dfile,txt2);
   IF (txt2 = 'encode') THEN idirec := 0;
   IF (txt2 = 'decode') THEN idirec := 1;
10:   writeln ('key':10,'plaintext':20,
      'expected cipher':21,'actual cipher':15);
   mm := 16;
   IF (nciphr < 16) THEN mm := nciphr;
   nciphr := nciphr-16;
   FOR m := 1 to mm DO BEGIN
      readln(dfile,hkey,hin,hcmp);
      FOR i := 1 to 16 DO BEGIN
         j := i+1;
         idum := hex2int(hin[j]);
         jdum := hex2int(hkey[j]);
         FOR k := 1 to 4 DO BEGIN
            l := 4*i+1-k;
            iin[l] := idum MOD 2;
            idum := idum DIV 2;
            key[l] := jdum MOD 2;
            jdum := jdum DIV 2
         END
      END;
      newkey := 1;
      des(iin,key,newkey,idirec,iout);
      hout := '                 ';
      FOR i := 1 to 16 DO BEGIN
         jdum := 0;
         twopwr := 16;
         FOR j := 1 to 4 DO BEGIN
            twopwr := twopwr DIV 2;
            jdum := jdum+twopwr*iout[4*(i-1)+j]
         END;
         hout[i+1] := int2hex(jdum)
      END;
      verdct := 'wrong   ';
      IF (hcmp = hout) THEN verdct := 'o.k.    ';
      writeln(hkey,hin,hcmp,hout,verdct:10)
   END;
   writeln ('press RETURN to continue ...');
   readln;
   IF (nciphr <= 0) THEN GOTO 5;
   GOTO 10;
99:
END.
PROGRAM d7r14(input,output);
(* driver for routines KS and CYFUN in file DESKS.PAS *)
TYPE
   gl32array = ARRAY [1..32] OF integer;
   gl48array = ARRAY [1..48] OF integer;
   gl56array = ARRAY [1..56] OF integer;
   gl64array = ARRAY [1..64] OF integer;
VAR
   glicd,ipc1 : gl56array;
   i,k : integer;
   iout,ir,ipp : gl32array;
   kn,ipc2,iet : gl48array;
   key,ip,ipm : gl64array;
   is : ARRAY [1..16,1..4,1..8] OF integer;
   ibin : ARRAY [1..4,1..16] OF integer;
   ksflg,desflg,cyflg : boolean;
   txt : PACKED ARRAY [1..64] OF char;

(*$I MODFILE.PAS *)
(*$I DESKS.PAS *)

BEGIN
(* test routine KS *)
   desflg := true;
   ksflg := true;
   cyflg := true;
   FOR i := 1 to 64 DO BEGIN
      key[i] := i MOD 2;
      IF (key[i] = 0) THEN txt[i] := '-';
      IF (key[i] = 1) THEN txt[i] := '*'
   END;
   writeln('legend:');
   writeln('-=0   *=1 ':25);
   writeln('master key:');
   write(' ':12);
   FOR i := 1 to 56 DO write(txt[i]);
   writeln;
   writeln('sub-master keys:');
   FOR i := 1 to 16 DO BEGIN
      ks(key,i,kn);
      FOR k := 1 to 48 DO BEGIN
         IF (kn[k] = 0) THEN txt[k] := '-';
         IF (kn[k] = 1) THEN txt[k] := '*'
      END;
      write(i:6,' ':6);
      FOR k := 1 to 48 DO write(txt[k]);
      writeln
   END;
   writeln('press RETURN to continue...');
   readln;
(* test routine CYFUN *)
   FOR i := 1 to 32 DO BEGIN
      ir[i] := (i MOD 3) MOD 2;
      IF (ir[i] = 0) THEN txt[i] := '-';
      IF (ir[i] = 1) THEN txt[i] := '*'
   END;
   writeln('legend:');
   writeln('-=0   *=1 ':25);
   writeln('input to cipher function:');
   write(' ':12);
   FOR i := 1 to 32 DO write(txt[i]);
   writeln;
   writeln('ciphered output:');
   FOR i := 1 to 16 DO BEGIN
      ks(key,i,kn);
      cyfun(ir,kn,iout);
      FOR k := 1 to 32 DO BEGIN
         IF (iout[k] = 0) THEN txt[k] := '-';
         IF (iout[k] = 1) THEN txt[k] := '*'
      END;
      write(i:6,' ':6);
      FOR k := 1 to 32 DO write(txt[k]);
      writeln
   END
END.
PROGRAM d10r10(input,output);
(* driver for routine DF1DIM *)
CONST
   ndim=3;
TYPE
   glnarray = ARRAY [1..ndim] OF real;
VAR
   i,j,ncom : integer;
   p,xi : glnarray;
   pcom,xicom : glnarray;

PROCEDURE dfunc(x : glnarray; VAR df : glnarray);
VAR
   i : integer;
BEGIN
   FOR i := 1 to 3 DO BEGIN
      df[i] := sqr(x[i]-1.0)
   END
END;

(*$I MODFILE.PAS *)
(*$I DF1DIM.PAS *)

FUNCTION fx(x: real): real;
BEGIN
   fx := df1dim(x)
END;

(*$I SCRSHO.PAS *)

BEGIN
   p[1] := 0.0; p[2] := 0.0; p[3] := 0.0;
   ncom := ndim;
   writeln;
   writeln('Enter vector direction along which to');
   writeln('plot the function. Minimum is in the');
   writeln('direction 1.0 1.0 1.0 - enter x y z:');
   read(xi[1],xi[2],xi[3]);
   writeln;
   FOR j := 1 to ndim DO BEGIN
      pcom[j] := p[j];
      xicom[j] := xi[j]
   END;
   scrsho
END.
PROGRAM d10r11(input,output);
(* driver for routine DFPMIN *)
CONST
   ndim=3;
   ftol=1.0e-6;
   pio2=1.5707963;
TYPE
   glnarray = ARRAY [1..ndim] OF real;
   glndim = glnarray;
   glnbyn = ARRAY [1..ndim,1..ndim] OF real;
VAR
   ncom : integer;
   pcom,xicom : glnarray; 
   angl,fret : real;
   iter,k : integer;
   p : glnarray;

(*$I MODFILE.PAS *)
(*$I BESSJ0.PAS *)

(*$I BESSJ1.PAS *)

FUNCTION fnc(x: glnarray): real;
BEGIN
   fnc := 1.0-bessj0(x[1]-0.5)*bessj0(x[2]-0.5)*bessj0(x[3]-0.5)
END;

PROCEDURE dfnc(x: glnarray; VAR df: glnarray);
BEGIN
   df[1] := bessj1(x[1]-0.5)*bessj0(x[2]-0.5)*bessj0(x[3]-0.5);
   df[2] := bessj0(x[1]-0.5)*bessj1(x[2]-0.5)*bessj0(x[3]-0.5);
   df[3] := bessj0(x[1]-0.5)*bessj0(x[2]-0.5)*bessj1(x[3]-0.5)
END;

(*$I F1DIM.PAS *)

FUNCTION func(x: real): real;
BEGIN
   func := f1dim(x)
END;

(*$I MNBRAK.PAS *)

(*$I BRENT.PAS *)

(*$I LINMIN.PAS *)

(*$I DFPMIN.PAS *)

BEGIN
   writeln('Program finds the minimum of a function');
   writeln('with different trial starting vectors.');
   writeln('True minimum is (0.5,0.5,0.5)');
   FOR k := 0 to 4 DO BEGIN
      angl := pio2*k/4.0;
      p[1] := 2.0*cos(angl);
      p[2] := 2.0*sin(angl);
      p[3] := 0.0;
      writeln;
      writeln('Starting vector: (',
         p[1]:6:4,',',p[2]:6:4,',',p[3]:6:4,')');
      dfpmin(p,ndim,ftol,iter,fret);
      writeln('Iterations:',iter:3);
      writeln('Solution vector: (',
         p[1]:6:4,',',p[2]:6:4,',',p[3]:6:4,')');
      writeln('Func. value at solution',fret:14)
   END
END.
