fprintf('A demonstration of the cocktail algorithm for D-optimal designs (approximate theory)\n\n');
fprintf('model: y ~ 1 + exp(-t) +  t * exp(-t) + exp(-s) + s * exp(-2s) + exp(-t-s)\n'); 
fprintf('design variables: 0<t, s<1 discretized\n\n');

n=50;
m=6;
t=(1:n)/n;
s=(1:n)/n;
X=ones(n*n, m);
for(i=1:n)
  for(j=1:n)
    X((i-1)*n+j,1:m)=[1, exp(-t(i)), t(i)*exp(-t(i)), exp(-s(j)), s(j)*exp(-2*s(j)), exp(-t(i)-s(j))];
  end
end
%[s, index]=sort(normrnd(0, 1, n*n, 1));
%X=X(index, 1:m);

w=cocktail(X);

index=1:(n*n);
index=index(w>0);
t=floor((index-1)/n)+1;
s=index-(t-1)*n;
t=t/n;
s=s/n;
fprintf('output:\noptimal support points (t(i), s(i))\n');
[t; s]'
fprintf('optimal weights w(i)\n');
w(w>0)
