function Y = nerv_p(P,lambda, epsilon)
% compute NeRV layout
%
% lambda controls the tradeoff; lambda=1 maximizes recall, lambda=0
% maximizes precision, default lambda=0.9;
%
% epsilon is a remedy for zeros in P (avoid numerical problems in
% precision, default=1e-10
% 
% Copyright (c) 2014, Zhirong Yang (Aalto University)
% All rights reserved.

if ~exist('lambda', 'var') || isempty(lambda)
    lambda = 0.9;
end

if ~exist('epsilon', 'var') || isempty(epsilon)
    epsilon = 1e-10;
end

method = 'nerv';
P = preprocess_input_similarities((P+P')/2, method);

theta = 2;
rng('default');
Y0 = randn(size(P,1),2)*1e-4;
check_step = 1;
tol = 1e-4;
max_time = inf;
verbose = true;
optimizer = 'sd';
recordYs = false;


fprintf('===============================================\n');
fprintf('Initialization stage\n');
fprintf('===============================================\n');
attr = 4;
max_iter = 30;
Y1 = ne_wrapper(P, method, Y0, attr, theta, max_iter, check_step, tol, max_time, verbose, optimizer, recordYs, lambda, epsilon);

fprintf('===============================================\n');
fprintf('Long-run stage\n');
fprintf('===============================================\n');
attr = 1;
max_iter = 100;
Y = ne_wrapper(P, method, Y1, attr, theta, max_iter, check_step, tol, max_time, verbose, optimizer, recordYs, lambda, epsilon);
