function miscmixbound(m0,n0)
global n m
randomize
n = n0;
m = m0;

P0 = ones(m^n,1); P0 = P0/sum(P0);
zz = P0*0 + 1e-12; oo = P0*0+1;

global pstr
pstr = 'miscmix'

if 0 % require sym. stationarity
  Aeq = is_stat(m,n);
  Aeq = [Aeq;oo'];
  Beq = zeros(size(Aeq,1),1); Beq(end)=1;
else
  Aeq = oo';
  Beq = 1;
end

global maxr
maxr = 0;
while 1
  %Py = rand(m,1); Py = Py/sum(Py);
  %P = randgivenstat(m,n,Py);
  P = rand(m^n,1); P = P/sum(P);
  
  x = [P(:)];
  
  y = fmincon(@radrat,x,[],[],Aeq,Beq,zz,oo);
  
  fprintf('%s: n=%d;  maxr = %1.9f \n',pstr,n,maxr);
end

return



function y = radrat(P)
global m n maxr
global pstr
P = makepos(P);
U = prodmeas(P,m,n);

[hh,Ht,H] = gethhn(P,m,n);
%hh = hh + eye(n);
%P2n = getK1(P,m,n);
%[hh,Ht,H2n] = gethhn(P2n,m,n-1);
%lhs = TV(P-U) + 2*H2n;
%rhs = 2*H;

%lhs = TV(P-U);
%rhs = max(diag(hh,1));

md = 0;
for k=1:n-1
  md = md + max(diag(hh,k));
end

%lhs = sum(sum(hh));
lhs = md;
rhs = H;

r = lhs/rhs;
y = -r;
if r > maxr
  global optP
  optP = P;
  maxr = r;
  fprintf('%s: n=%d;  maxr = %1.9f \n',pstr,n,maxr);
end


return