Help:=proc(): print(`RG3C(n), ZKP3(n,G,C,Opt)`): end:
with(combinat):
RG3C:=proc(n) local ra,i,j,C,ra1,G:
ra:=rand(1..3):
C:=[seq(ra(),i=1..n)]:
ra1:=rand(0..1):
 G:={}:
for i from 1 to n do
 for j from i+1 to n do
   if C[i]<>C[j] then
    if ra1()=1 then
      G:=G union {{i,j}}:
    fi:
   fi:
 od:
od:
n,G,C:
end:


ZKP3:=proc(n,G,C,Opt) local c,L,i,j,sig,BiV,BiC,BjV,BjC,B,S:
#box is a pair [color,vertex]
L:=[seq(seq([c,i],i=1..n),c=1..3)]:
sig:=randperm(3*n):
L:=[seq(L[sig[i]],i=1..3*n)]:

for i from 1 to 3*n do
  for j from i+1 to 3*n do
    BiV:=L[i][2]:
    BiC:=L[i][1]:
    BjV:=L[j][2]:
    BjC:=L[j][1]:
   if C[BiV]=BiC and C[BjV]=BjC and member({BiV,BjV},G) then
     B[i,j]:=1:
    else
     B[i,j]:=0:
   fi:
 od:
od:
#L,op(B):

if Opt=1 then
   RETURN([op(B),[seq(L[i][2],i=1..nops(L))]]):
else
 S:={}:
 for i from 1 to 3*n do
   for j from i+1 to 3*n do
    if L[i][1]=L[j][1] then
       S:=S union {B[i,j]}:
    fi:
   od:
od:
RETURN([seq(L[i][2],i=1..nops(L))],S):
fi:
end: