You selected design0.pl

title:-
A=['% design0.pl : experimental design using extended Galois field'
,'% By Kenryo Indo (Kanto Gakuen University)'
,'% file: design0.pl ( for SWI-Prolog 5.0.9)'
,'% created: 20-21 Aug 2003'
,'% modified: 15 Jul 2004 (corrected the positions of dynamic declarations)'
,'% Reference:'
,'% Lecture Note of Simulation (by Prof. Iwaro Takahashi, University of Tsukuba, 1987). '
,'% import: matrix01.pl(20 Aug 2003) matrix operations.'
,'%-- the main has prolonged in the last part of the program ---------'
,'% othogonal_table/2.'
,'% start_design/0.'
]
,forall(member(B,A),(nl,write(B))).

:- [matrix01].



% an expermental design problem
%---------------------------------

:- dynamic  number_of_categories/1.
:- dynamic  number_of_factors/1.
:- dynamic  number_of_trials/1.

number_of_categories(2).
number_of_factors(7).
number_of_trials(8).


% experimental design using GF(S).
%---------------------------------

% careate theta(TH,K), K=1, ..., L, a minimal 
% L-dimensional independent vectors on GF(2) 
% which covers all factors above.

design_pattern(S,S):-member(S,[2,3]).
design_pattern(4,2^2).

all_auxiliary_matrix(using_GF(S),dimension(L),aux(AQ)):-
   design_pattern(S,S1),
   findall(Q,
     auxiliary_vector(using_GF(S1),dimension(L),aux(Q)),
   AQ).

minimal_auxiliary_matrix(design(S,M,N),dimension(L),aux(AQ),min(AP)):-
   number_of_categories(S),  
   minimal_dimension(L,design(S,M,N)),
   !,
   all_auxiliary_matrix(using_GF(S),dimension(L),aux(AQ)),
   findall([(Max,Sum)|P],
     (
      member(P,AQ),
      has_a_non_zero_member(P,Sum),
      max(Max,member(Max,P))
     ),
   AP0),
   sort(AP0,AP1),
   findall_up_to(M,P,member([_|P],AP1),AP).

findall_up_to(N,Target,Goal,Bag):-
   findall(Target,Goal,Bag0),
   findall(Target,(kth(K,Bag0,Target),K== Dim.

covering_dimension_1(using_GF(S),factors(M),trials(N),dimension(L)):-
   length(F,M),
   nth1(L,F,L),
   findall(P,
     auxiliary_vector(using_GF(S),dimension(L),aux(P)),
   Ps1),
   length(Ps1,D),
   M =< D,
   N >= D.

auxiliary_vector(using_GF(S),dimension(L),aux(P)):-
   member(S,[2,3]),
   number(L),
   anum_seq(A,S),
   bag0(P,A,L).

auxiliary_vector(using_GF(2^2),dimension(L),aux(P)):-
   number(L),
   bag0(P,[0,1,b,1+b],L).

aux(S,L,P):-
   auxiliary_vector(using_GF(S),dimension(L),aux(P)).


%--------------------------------------
% create othogonal table on GF(S)
%--------------------------------------

:- dynamic table1 / 2.

othogonal_table(design(S,M,N),AQ*AP=W):-
   write('abolish tables ? >'),
   (read(y)->abolish(table1/2);true),
   minimal_auxiliary_matrix(design(S,M,N),dimension(_L),aux(AQ),min(AP)),
   matrix_transpose(AP -> APT),
   design_pattern(S,S1),
   matrix_multiply_GF(S1,AQ * APT = Z),
   matrix_mod(W is (Z mod S1)),
   assert(table1(design(S,M,N),AQ*AP=W)),
   write('display tables ? >'),
   (read(y)->pp_table(S1);true).

pp_table(2^2):-
   table1(design(S,M,N),AQ*AP=W),
   nl,write(table1(design(S,M,N))),
   matrix_transpose(AQ -> AQT),
   matrix_transpose(AP -> APT),
   MM is M * 3,
   nl,write_n_times('-',MM+3),
   write_table_GF(AQT),
   nl,write_n_times('-',MM),
   write_table_GF(APT),
   nl,write_n_times('-',MM),
  %forall(member(X,W),(nl,write(X))),
   write_table_GF(W),
   matrix_transpose(W -> WT),
   sum_table_GF(WT,Sum),
   nl,write_n_times('-',MM),
   nl,write(Sum).

pp_table(S):-
   member(S,[2,3]),
   table1(design(S,M,N),AQ*AP=W),
   nl,write(table1(design(S,M,N))),
   matrix_transpose(AQ -> AQT),
   matrix_transpose(AP -> APT),
   MM is M * 3,
   nl,write_n_times('-',MM+3),
   write_table(AQT),
   nl,write_n_times('-',MM),
   write_table(APT),
   nl,write_n_times('-',MM),
   write_table(W),
   matrix_transpose(W -> WT),
   sum_table(WT,Sum),
   nl,write_n_times('-',MM),
   nl,write(Sum).

sum_table(W,Sum):-
   findall(S,(member(Q,W), sum(Q,S)),Sum).

sum_table_GF(W,Sum):-
   findall(S,(member(Q,W), sum_GF(Q,S)),Sum).

write_table(W):-
   forall(member(Q,W), (sum(Q,S),nl,write(Q),tab(2),write((S)))).

write_table_GF(W):-
   forall(member(Q,W), (sum_GF(Q,S),nl,write(Q),tab(2),write((S)))).

write_n_times(Q,N):-
   N1 is N,
   len(A,N1),
   forall(member(Q,A), (write(Q))).




%-----------------------
%   utilities
%-----------------------
has_a_non_zero_member(P,Sum):-
   \+ ( member(AA,P), \+ number(AA) ),
   sum(P,Sum),
   Sum > 0,
   !.

has_a_non_zero_member(P,Sum):-
   ( member(AA,P), \+ number(AA) ),
   sum_GF(P,Sum).

sum_GF([],0).
sum_GF([A|P],Sum):-
   sum_GF(P,S0),
   is_GF(Sum, A+S0).

min(X,Goal):-
   max(-X,Goal).

max(X,Goal):-
  % X: the objective variable,
  % Goal: the objective function and constraints,
   setof((X,Goal),Goal,Z),
   forall(member(A,Z),(nl,write(A))),
   max_0(X,Z).

max_0(X,Z):-
   \+ ( member((AA,_),Z), \+ number(AA) ),
   member((X,_Goal),Z),
   %nl,write((case1,X,_Goal)),
   \+ (
     member((Y,_),Z),
     Y > X
   ).

% modified: 22 Aug 2003  a cut !.
max_0(X,Z):-
   member((AA,Goal),Z),
   \+ number(AA),
   !,
   %nl,write((case2,AA,Goal)),
   sort(Z,Z1),
   rev(Z1,[(X,_)|_]).



%---------------------------------
%   operations on extended GF
%---------------------------------

operation_GF(E):-
   member(X,[0,1,b,1+b]),
   member(Y,[0,1,b,1+b]),
   member(F,[+,-,*]),
   E=..[F,X,Y].



:- dynamic is_GF/2.

operations(on_GF(2^2),RULE):-
  % primitive reduced polynominal: 1+x+x^2=0
  RULE=[
  % (is_GF(A,B):- A is B),
    is_GF(0, 0+0),
    is_GF(0, 0-0),
    is_GF(0, 0*0),
    is_GF(1, 0+1),
    is_GF(1, 0-1),
    is_GF(0, 0*1),
    is_GF(b, 0+b),
    is_GF(b, 0-b),
    is_GF(0, 0*b),
    is_GF(1+b, 0+ (1+b)),
    is_GF(1+b, 0- (1+b)),
    is_GF(0, 0* (1+b)),
    is_GF(1, 1+0),
    is_GF(1, 1-0),
    is_GF(0, 1*0),
    is_GF(0, 1+1),
    is_GF(0, 1-1),
    is_GF(1, 1*1),
    is_GF(1+b, 1+b),
    is_GF(1+b, 1-b),
    is_GF(b, 1*b),
    is_GF(b, 1+ (1+b)),
    is_GF(b, 1- (1+b)),
    is_GF(1+b, 1* (1+b)),
    is_GF(b, b+0),
    is_GF(b, b-0),
    is_GF(0, b*0),
    is_GF(1+b, b+1),
    is_GF(1+b, b-1), % X+1=b-1+1->X=(1+b)
    is_GF(b, b*1),
    is_GF(0, b+b),
    is_GF(0, b-b),
    is_GF(1+b, b*b), % since p.r.p. b^2 = -(1+b)=1+b
    is_GF(1, b+ (1+b)),
    is_GF(1, b- (1+b)),
    is_GF(1, b* (1+b)), % b*(1+b)=b+b^2=b+(1+ b)=1+(b+b)=1
    is_GF(1+b, 1+b+0),
    is_GF(1+b, 1+b-0),
    is_GF(0, (1+b)*0),
    is_GF(b, 1+b+1),
    is_GF(b, 1+b-1),
    is_GF(1+b, (1+b)*1),
    is_GF(1, 1+b+b),
    is_GF(1, 1+b-b),
    is_GF(1, (1+b)*b), % (1+b)*b=b+b^2=b+(1+b)=1
    is_GF(0, 1+b+ (1+b)),
    is_GF(0, 1+b- (1+b)),
    is_GF(b, (1+b)* (1+b)) %(1+b)* (1+b)=1+2b+b^2=1+(1+b)=b
  ].

assert_GF(2^2):-
   operations(on_GF(2^2),RULE),
   forall(member(X,RULE),assertz((X:-!))).

:- assert_GF(2^2).


% matrix multiplication for GF(2^2)
% -----------------------------------------------------------  %

matrix_multiply_GF(2^2,X * Y = Z):-
   !,
   \+ number(X),
   (matrix(X,A);X=A),
   (matrix(Y,B);Y=B),
   matrix_of_size(A,N,L0),
   matrix_of_size(B,L0,M),
   matrix_of_size(Z,N,M),
   matrix_transpose(Y -> C),!,
   findall(ZK,
     (
      kth(K,Z,ZK),
      findall(ZKL,
        (
         kth(L,ZK,ZKL),
         matrix_rows(A,[K],[AK]),
         matrix_rows(C,[L],[CL]),
         eq_sum_product_GF(AK,CL,_,ZKL,_)
        ),
      ZK)
     ),
   Z).

matrix_multiply_GF(_S,X * Y = Z):-
   matrix_multiply(X * Y = Z).

eq_sum_product_GF([],[],[],0,0).
eq_sum_product_GF([P|Q],[A|B],[E|F],V,Vq):-
    len(Q,N),
    len(B,N),
    eq_sum_product_GF(Q,B,F,V1,Vq1),
    is_GF(E, P * A),
    is_GF(V, V1 + E),
    (Vq1=0 -> Vq = EP*EA; Vq = Vq1 + (EP * EA)).




% problems
%--------------------------------------
problem(a,
  number_of_categories(2),
  number_of_factors(7),
  number_of_trials(8)
).

problem(b,
  number_of_categories(3),
  number_of_factors(12),
  number_of_trials(27)
).

problem(c,
  number_of_categories(2),
  number_of_factors(10),
  number_of_trials(16)
).

problem(d,
  number_of_categories(4),
  number_of_factors(5),
  number_of_trials(16)
).



% setting the problem
%---------------------------------

:- dynamic  current_problem/1.


select_problem(A):-
  problem(A,B,C,D),
  initialize_current_problem,
  update_current_problem(A,[B,C,D]).

initialize_current_problem:-
  abolish(current_problem/1),
  abolish(number_of_categories/1),
  abolish(number_of_factors/1),
  abolish(number_of_trials/1).

update_current_problem(A,Target):-
  assert(current_problem(A)),
  forall(member(Y,Target),assert(Y)).



start_design:-
  CREATE_TABLE= othogonal_table(design(S,M,N),_),
  nl,
  write_n_times('*',15),
  write(' welcome '),
  write_n_times('*',15),
  nl,
  write(' this is prolog program of the experimental design problems. '),
  nl,
  write(' available samples as follows:'),
  P=problem(A,
    number_of_categories(S),
    number_of_factors(M),
    number_of_trials(N)
  ),
  forall(P,(nl,tab(2),write((A,categories:S,factors:M,trials:N)))),
  nl,
  nl,
  write(' please select a problem above (default is problem a):'),
  read(U),
  (select_problem(U)->true;
   (nl,write('no such problem.'),start_design0)
  ),
  nl,
  write(' ok. '),
  nl,
  write(' start to create table? :'),
  read(y),
  CREATE_TABLE,
  nl,
  write(' the table has created (refer table1/2).'),
  nl,
  write(' continue design ? :'),
  (read(y)->start_design;!,fail).

start_design.

% auto start.
:- start_design.



%------------------------------------------------



return to front page.