function [p,d]=ps23a(h)
% function [p,d]=ps23a(h)
% 
% given  h (n-by-1, h>0, default h=rand(3,1),
% minimize p=r*prod(h) subject to d>0,
%  r*v(1)^2+v(1)*v(n+1)+d'*(h.*v(2:n+1).*v(2:n+1)-v(2:n+1).*v(1:n)) > 0
if nargin<1, h=rand(2,1); end
h=h(:);
n=length(h);

r=msspoly('r');          % r
v=msspoly('v',n+1);      % v
d=msspoly('d',n);        % d
q=msspoly('Q',nchoosek(n+2,2));        % coefficients of Q

pr = mssprog;            % initialize program
pr.free = r;             % register r
pr.pos = d;              % register d
pr.psd = q;              % register q
Q = mss_v2s(q);          % symmetric Q
pr.eq = r*v(1)^2+v(1)*v(n+1)+ ...
    d'*(h.*v(2:n+1).*v(2:n+1)-v(2:n+1).*v(1:n))-v'*Q*v;                  
pr.sedumi = r;           % minimize r
p=prod(h)*pr({r});       % get p
d=pr({d});               % get d