% runcode for p_criterion

randn('seed',5);
rand('seed',5);

% simulated data

path('./p_criterion_with_K',path);

fid = fopen('p_sim.txt','w');
fclose(fid);

repeats = 30;
maxiter = 10000;

k = 3;

feas = [];

for design = 1:4
  
  if design ~= 3
    
    for m = [10000 50000 100000]
      
      for p = [-0.25 -0.75 -1.1 -1.2]
      
        Table = zeros(repeats,4);
      
        for repeat = 1:repeats
        
          [M, K, w0] = real_data(design,m);
          n = sqrt(size(M,1));
      
          K = randn(n,k);
      
          tic
          [w fval2] = multiplicativep(M,p,K,w0,10000,2e-4,maxiter);
          time2 = toc;
          
          feas = [feas; sum(w) min(w) fval2];
      
          tic
          [w fval3] = IPM_pK(M,p,K,1e-10,inf);
          time3 = toc;
          
          feas = [feas; sum(w) min(w) fval3];
      
          Table(repeat,:) = [time2 time3 fval2 fval3];
        
        end
      
        Table = sum(Table,1)/repeats;
      
        fid = fopen('p_sim.txt','a');
        fprintf(fid, '%3d & %3d & %4.2f &  %5.2f  &  %5.2f  &  %13.6g  &  %13.6g\n', design, m, p, Table);
        fclose(fid);
        
      end
      
    end
    
  else 
    
    for m = [100 200 300]
      
      for p = [-0.25 -0.75 -1.1 -1.2]
      
        Table = zeros(repeats,4);
      
        for repeat = 1:repeats
        
          [M, K, w0] = real_data(design,m);
          n = sqrt(size(M,1));
      
          K = randn(n,k);
      
          tic
          [w fval2] = multiplicativep(M,p,K,w0,10000,2e-4,maxiter);
          time2 = toc;
          
          feas = [feas; sum(w) min(w) fval2];
       
          tic
          [w fval3] = IPM_pK(M,p,K,1e-10,inf);
          time3 = toc;
          
          feas = [feas; sum(w) min(w) fval3];
      
          Table(repeat,:) = [time2 time3 fval2 fval3];
        
        end
      
        Table = sum(Table,1)/repeats;
      
        fid = fopen('p_sim.txt','a');
        fprintf(fid, '%3d & %3d & %4.2f &  %5.2f  &  %5.2f  &  %13.6g  &  %13.6g\n', design, m, p, Table);
        fclose(fid);
        
      end
      
    end
    
  end
   
end

