headline:-
 write('% ------------------------------------------- %'),nl,
 write('%   simulating cooperative games by Prolog.   %'),nl,
 write('% ------------------------------------------- %'),nl,
 h0.
h0:-
 write('%  game(G,value,X,P):-  values (worth) of game '),nl,
 write('%  core(G,C,P):-  coa '),nl,
 write('%  imputation(G,C,P):- imputation '),nl,
 write('%  coalitionally_complain(G,C,Z,B):- complain (excess) '),nl,
 write('%  is_more_acceptable_than(G,A,A1,Z,Z1):-  so '),nl,
 write('%  nucleolus(G,A):- so '),nl,
 write('%  coalition_formation(G,[J|Z],Y/N,[A|B],P):- so '),nl,
 write('%  shapley(G,N,V):- so '),nl,
 write('%  h0:-  this.'),nl.
me:-
   write('% file: coop.pl.'),nl,
   write('% created: 7 Feb 2003.'),nl,
   write('% modified: 9 Feb 2003.'),nl,
   write('% author: Kenryo INDO (Kanto Gakuen University) '),nl,
   write('% url: http://www.us.kanto-gakuen.ac.jp/indo/front.html'),nl.
reference:-
   write('% references: '),nl,
   write('%  Muto, S. (2001). An Introduction to Game Theory.'),nl,
   write('%     Nikkei Bunko. pp.161-194.(Japanese) '),
   write('%  Aumann, R.J.(1989). Lectures on Game Theory. '),nl,
   write('%    Stanford University. [T. Maruyama and H. Tateishi,  '),nl,
   write('%    "Game Ron no Kiso", Keiso Shobo, 1991. (translated in Japanese) ]'),
   nl.

:- headline.

%
% -------------------------------------------------  %
% examples of cooperative games
%   with transferable utilities or side payments
% -------------------------------------------------  %
% these examples below cited from Muto(2001), pp.165-8.
% the game c2 with minor modifications in the values 
% of its characteristic function.
%
% game c1: a majority vote.
game(c1,
   form(characteristic), 
   players([a,b,c]), 
   coalitions([[],[a],[b],[c],[a,b],[b,c],[a,c],[a,b,c]])).
game(c1,value,[],0).
game(c1,value,[a],0).
game(c1,value,[b],0).
game(c1,value,[c],0).
game(c1,value,[a,b],1).
game(c1,value,[b,c],1).
game(c1,value,[a,c],1).
game(c1,value,[a,b,c],1).
%
% game c2: selling the asset of a to whom, b or c?
game(c2,
   form(characteristic), 
   players([a,b,c]), 
   coalitions([[],[a],[b],[c],[a,b],[b,c],[a,c],[a,b,c]])).
game(c2,value,[],0).
game(c2,value,[a],0).
game(c2,value,[b],0).
% coa is empty if modified as below
%game(c2,value,[a],1).
%game(c2,value,[b],1).
game(c2,value,[c],0).
game(c2,value,[a,b],2).
game(c2,value,[b,c],0).
game(c2,value,[a,c],5).
game(c2,value,[a,b,c],5).
%
% game c3: cost-sharing problem among 3 cities.
game(c3,
   form(characteristic), 
   players([a,b,c]), 
   coalitions([[],[a],[b],[c],[a,b],[b,c],[a,c],[a,b,c]])).
game(c3,value,[],0).
game(c3,value,[a],0).
game(c3,value,[b],0).
game(c3,value,[c],0).
game(c3,value,[a,b],6).
game(c3,value,[b,c],8).
game(c3,value,[a,c],0).
game(c3,value,[a,b,c],20).
% game c0: cost-sharing problem among 3 cities.
game(c0,
   form(characteristic), 
   players([a,b]), 
   coalitions([[],[a],[b],[a,b]])).
game(c0,value,[],0).
game(c0,value,[a],0).
game(c0,value,[b],0).
game(c0,value,[a,b],1).
%
% -------------------------------------------------  %
% imputation, and core   
% -------------------------------------------------  %
%
col_rat_outcome(G,players(N),payoff(A)):-
   var(A),
   game(G,form(characteristic),players(N),coalitions(C)),
   member(N,C),
   game(G,value,N,V),
   length(N,LN),
   allocation(LN,V,A).

col_rat_outcome(G,players(N),payoff(A)):-
   \+ var(A),
   game(G,form(characteristic),players(N),coalitions(C)),
   length(N,LN),
   length(A,LN),
   member(N,C),
   game(G,value,N,V),
   \+ (member(X,A), X < 0),
   sum(A,V).

individually_complain(G,J/N,RJ-AJ=Z/A,X):-
   col_rat_outcome(G,players(N),payoff(A)),
   nth1(K,N,J),
   nth1(K,A,AJ),
   game(G,value,[J],RJ),
   Z is RJ - AJ,
   (AJ < RJ -> X = yes; X = no).

imputation(game(G),players(N),payoff(A)):-
   % collectively (i.e.,group) rational outcome.
   col_rat_outcome(G,players(N),payoff(A)),
   % individual rationality.
   \+ individually_complain(G,_J/N,_RJ-_AJ=_Z/A,yes).

excess_of_coalition(G,Y/N,RY-AY=Z/A,X):-
   coalitionally_complain(G,Y/N,RY-AY=Z/A,X).

coalitionally_complain(G,Y/N,RY-AY=Z/A,X):-
   imputation(game(G),players(N),payoff(A)),
   game(G,value,Y,RY),
   Y \= N,
   selected_sum(Y/N,_B/A,AY),
   Z is RY - AY,
   (AY < RY -> X = yes; X = no).

core(game(G),players(N),payoff(A)):-
   imputation(game(G),players(N),payoff(A)),
   % coaltional rationality.
   \+ coalitionally_complain(G,_Y/N,_RY-_AY=_Z/A,yes).

%
% sample execution
%------------------------------------------------------
/*
?- core(game(c2),B,C).

B = players([a, b, c])
C = payoff([5, 0, 0]) ;

B = players([a, b, c])
C = payoff([4, 0, 1]) ;

B = players([a, b, c])
C = payoff([3, 0, 2]) ;

B = players([a, b, c])
C = payoff([2, 0, 3]) ;

Yes
?- coalitionally_complain(c3,B/[a,b,c],_=Z/[6,0,14],no).

B = []
Z = 0 ;

B = [a]
Z = -6 ;

B = [b]
Z = 0 ;

B = [c]
Z = -14 ;

B = [a, b]
Z = 0 ;

B = [b, c]
Z = -6 ;

B = [a, c]
Z = -20 ;

No
?- 
*/
%
% -------------------------------------------------  %
% Schmeidler(1969)'s nucleolus 
% -------------------------------------------------  %
%
% nucleolus: 
% lexicographically minimizing the sorted complaining vector.

complain_vector(G,A,Zs):-
   imputation(game(G),players(N),payoff(A)),
   findall(Z,coalitionally_complain(G,_B/N,_=Z/A,_),Zs).
complain_vector_indexed(G,A,Bs):-
   imputation(game(G),players(N),payoff(A)),
   findall((B,Z),coalitionally_complain(G,B/N,_=Z/A,_),Bs).

sorted_complain_vector(G,A,Z):-
   complain_vector(G,A,S0),
   asort(S0,S),
   reverse(S,Z).

is_more_acceptable_than(G,A,A1):-
   is_more_acceptable_than(G,A,A1,_,_).

is_more_acceptable_than(G,A,A1,Z,Z1):-
   sorted_complain_vector(G,A,Z),
   sorted_complain_vector(G,A1,Z1),
   Z @< Z1.

nucleolus(G,A):-
   imputation(game(G),players(_N),payoff(A)),
   \+ is_more_acceptable_than(G,_A1,A).
%
% sample execution
%------------------------------------------------------
/*
?- is_more_acceptable_than(c3,[12,4,4],[6,0,14],B,C).

B = [0, 0, -4, -4, -10, -12, -16]
C = [0, 0, 0, -6, -6, -14, -20] 

Yes
?- nucleolus(A,B).

A = c1
B = [1, 0, 0] ;

A = c1
B = [0, 1, 0] ;

A = c1
B = [0, 0, 1] ;

A = c2
B = [3, 0, 2] ;

A = c3
B = [6, 7, 7] ;

No
?- 
*/

%
% -------------------------------------------------  %
% Shapley(1953)'s value   
% -------------------------------------------------  %
%
contribution(G,J,X,Y,A):-
   game(G,form(characteristic),players(_N),coalitions(C)),
   member(Y,C),
   game(G,value,Y,VY),
   member(J,Y),
   subtract(Y,[J],X),
   game(G,value,X,VX),
   A is VY - VX.

coalition_formation(G,[],[]/N,[],0):-
   game(G,form(characteristic),players(N),coalitions(_C)).
coalition_formation(G,[J|Z],Y/N,[A|B],P):-
   coalition_formation(G,Z,X/N,B,_Q),
   (X=N -> (!,fail);true),
   contribution(G,J,X,Y,A),
   game(G,value,Y,P).

contribution_to_coalition_formation(G,J,X,K,VJ/V):-
   coalition_formation(G,X,N/N,VX,V),
   nth1(K,X,J),
   nth1(K,VX,VJ).

shapley(G,J/N,Ps,SV):-
   game(G,form(characteristic),players(N),coalitions(_C)),
   member(J,N),
   bagof(VJ,
      X^K^contribution_to_coalition_formation(G,J,X,K,VJ/_V),
   Ps),
   length(Ps,L),
   sum(Ps,B),
   SV is B / L.

shapley(G,N,V):-
   bagof(SV,
      J^Ps^shapley(G,J/N,Ps,SV),
   V),
   (
    imputation(game(G),players(N),payoff(V))
    ->true
    ;
      write(not_an_imputation(V))
   ).
%
% sample execution
%------------------------------------------------------
/*
?- shapley(A,B,C),col_rat_outcome(A,_,payoff(C)).

A = c1
B = [a, b, c]
C = [0.333333, 0.333333, 0.333333] ;

A = c2
B = [a, b, c]
C = [2.83333, 0.333333, 1.83333] ;

A = c3
B = [a, b, c]
C = [5, 9, 6] ;
*/
%
% -----------------------------------------------------------  %
% Arithmetic and so on including probabilistic operators
% -----------------------------------------------------------  %
%
% sum
% -----------------------------------------------------------  %
sum([],0).
sum([X|Members],Sum):-
   sum(Members,Sum1),
  %number(X),
   Sum is Sum1 + X.
%
% product
% -----------------------------------------------------------  %
product([],1).
product([X|Members],Z):-
   product(Members,Z1),
  %number(X),
   Z is Z1 * X.
%
% weighted sum
% -----------------------------------------------------------  %
product_sum([],[],[],0).
product_sum([P|Q],[A|B],[E|F],V):-
    length(Q,N),
    length(B,N),
    product_sum(Q,B,F,V1),
    E is P * A,
    V is V1 + E.
%
% selected sum
% -----------------------------------------------------------  %
selected_sum(Y/N,B/A,RX):-
   findall(AJ,
     (
      member(J,Y),
      nth1(K,N,J),
      nth1(K,A,AJ)
     ),
   B),
   sum(B,RX).
%
% projected sum
% -----------------------------------------------------------  %
projected_sum(M,A,Cols):-
   index_of_tuple(M,B,Cols),
   sum(B,A).
%
% allocation
% -----------------------------------------------------------  %
allocation(N,A,[X|Y]):-
    allocation(N,A,A,[X|Y]).
allocation(0,_,0,[]).
allocation(N,A,B,[X|Y]):-
    integer(A),
    length([X|Y],N),
    allocation(_N1,A,B1,Y),
    % N1 is N - 1,
    % sum(Y,B1),
    K is A - B1 + 1,
    length(L,K),
    nth0(X,L,X),
    B is B1 + X.
% -------------------------------------------------  %
%  some (local) utilities for probabilistic operations
% -------------------------------------------------  %
precision(100).
make_a_prob(N0,P):-
   number(P),
   precision(N0),
   P =< 1,
   P >= 0.
make_a_prob(N0,P):-
   var(P),
   precision(N0),
   N1 is N0 + 1, 
   length(L,N1),
   nth0(K,L,K),
   P is K / N0.
quotient_prob(user,R, P):-
   (var(R)->read(R1);true),
   (
    R1 = Q1/Q0
    ->
    R = Q1/Q0
    ;
    quotient_prob(user,R, P)
   ),
   P is R.
%
conditional_event_probability(E,H,P):-
    event(E),
    event(H), 
    H \= [],
    intersection(E,H,F),
    probability_of_event(_,H,P0),
    (P0 = 0 -> (nl,write('-- measure 0 --'),nl,fail);true),
    probability_of_event(bp1,F,P1),
    P is P1 / P0.
%
% probability by using allocation
% -----------------------------------------------------------  %
probabilities(0,[]).
probabilities(N,[X|Y]):-
    integer(N),
    length([X|Y],N),
    allocation(N,100,[X|Y]).
% 
% any ratio (weight) can be interpreted into a prob.
scale(W,1/Z,P):-
    findall(Y,(nth1(_K,W,X),Y is X/Z),P).
probabilities(W,N,P):-
    length(W,N),
    sum(W,Z),
    scale(W,1/Z,P).
%
make_a_prob(P,base(M),steps(L)):-
    var(P),
    length(P,M),
    allocation(M,L,W),
    probabilities(W,M,P).
make_a_prob(P,base(M),_):-
    \+ var(P),
    length(P,M),
    \+ (
     member(P1,P),
     (
      var(P1);
      P1 > 1;
      P1 < 0
     )
    ),
    sum(P,1).
%
% expected value
% -----------------------------------------------------------  %
expected_value(W,A,E/100):-
    length(A,N),
    probability(W,N,P),
    product_sum(P,A,_,E).
%
% -----------------------------------------------------------  %
%   Utilities for list operations
% -----------------------------------------------------------  %
%
% index for a tuple.
% -----------------------------------------------------------  %
% 1) only mention for a direct product of sets.
index_of_tuple(B,A,Index):-
   \+ var(B),
   \+ var(A),
   length(B,LN),  % base sets
   length(A,LN),  
   length(Index,LN),
   findall(L,
     (
      nth1(K,B,BJ), %write(a(K,B,BJ)),
      nth1(L,BJ,AJ),%write(b(L,BJ,AJ)),
      nth1(K,A,AJ)  %,write(c(K,A,AJ)),nl
     ),
   Index).
index_of_tuple(B,A,Index):-
   \+ var(B),
   \+ var(Index),
   var(A),
   length(B,LN),  % base sets
   length(Index,LN),
   length(A,LN),  
   findall(AJ,
     (
      nth1(K,B,BJ),
      nth1(K,Index,L),
      nth1(L,BJ,AJ)
     ),
   A).
%
% descending/ascending natural number sequence less than N.
% -----------------------------------------------------------  %
dnum_seq([],N):-N<0,!.
dnum_seq([0],1).
dnum_seq([A|Q],N):-
   A is N - 1,
   length(Q,A),
   dnum_seq(Q,A).
anum_seq(Aseq,N):-dnum_seq(Dseq,N),sort(Dseq,Aseq).
%
% inquire the goal multiplicity
% -----------------------------------------------------------  %
sea_multiple(Goal,Cond,N,M):-
  Clause=..Goal,
  findall(Cond,Clause,Z),length(Z,N),sort(Z,Q),length(Q,M).
%
bag0([],_A,0).
bag0([C|B],A,N):-
   length([C|B],N),
   bag0(B,A,_N1),
   member(C,A).
zeros(Zero,N):-bag0(Zero,[0],N).
ones(One,N):-bag0(One,[1],N).
%
% subset_of/3 : subset-enumeration 
% -----------------------------------------------------------  %
subset_of(A,N,As):-
   length(As,L),
   length(D,L),
   list_projection(D,As,B),
   length(B,N),
   sort(B,A).
%
% a sequence of binary choice for a list:
%--------------------------------------------------
list_projection([],[],[]).
list_projection([X|Y],[_A|B],C):-
   X = 0,
   list_projection(Y,B,C).
list_projection([X|Y],[A|B],[A|C]):-
   X = 1,
   list_projection(Y,B,C).
%
% sort without removal of duplicates
%--------------------------------------------------
asort(A,B):-
   sort(A,C),
   bagof(CK,
     J^K^(
      nth1(J,C,CK),
      nth1(K,A,CK)
     ),
   B).
%
% -----------------------------------------------------------  %
%   Utilities for outputs
% -----------------------------------------------------------  %
%
% write and new line.
% -----------------------------------------------------------  %
wn(X):-write(X),nl.
%
% output to file.
% -----------------------------------------------------------  %
tell_test(Goal):-
  open('tell.txt',write,S),
  tell('tell.txt'),
  Goal,
  current_stream('tell.txt',write,S),
  tell(user),wn(end),
  close(S).
%


%end