client = jcmoptimizer.Client(); 

% Definition of the search domain
design_space = {...
   struct('name', 'x1', 'type', 'continuous', 'domain', [-1.5,1.5]), ...
   struct('name', 'x2', 'type', 'continuous', 'domain', [-1.5,1.5]), ...
};

% Definition of fixed environment parameter
environment = {...     
   struct('name', 'radius', 'type', 'fixed', 'domain', 1.5) ...
};

% Definition of a constraint on the search domain
constraints = {...
   struct('name', 'circle', 'expression', 'sqrt(x1^2 + x2^2) <= radius')...
};

 % Creation of the study object with study_id 'vanilla_bayesian_optimization'
study = client.create_study( ...
    'design_space', design_space, ...
    'environment', environment, ...
    'constraints', constraints,...
    'driver','BayesianOptimization',...
    'study_name','Standard Bayesian optimization',...
    'study_id', 'vanilla_bayesian_optimization');

% Evaluation of the black-box function for specified design parameters
function observation = evaluate(study, sample)

    pause(2); % make objective expensive
    observation = study.new_observation();
    x1 = sample.x1;
    x2 = sample.x2;
    observation.add(10*2 ...
                    + (x1.^2-10*cos(2*pi*x1)) ...
                    + (x2.^2-10*cos(2*pi*x2)) ...
                   );
end  

% Run the minimization
study.set_evaluator(@evaluate);
study.run(); 


best = study.driver.best_sample;
fprintf('Best sample at x1=%0.3e, x2=%0.3e\n', best.x1, best.x2) 
% print information on all found local minima
minima = study.driver.get_minima('num_initial_samples', 20);  
array2table(cell2mat(struct2cell(minima)).','VariableNames', fieldnames(minima))

% get information on best value under uncertain inputs with standard deviation 0.1
pd = struct();
pd.distributions = {...     
    struct('type', 'normal', 'parameter', 'x1', 'mean', best.x1, 'stddev', 0.01), ...
    struct('type', 'normal', 'parameter', 'x2', 'mean', best.x2, 'stddev', 0.02) ...
    };
study.configure('parameter_distribution', pd);
 
stats = study.driver.get_statistics('quantiles', [0.5]);
fprintf("Objective mean %f, standard deviation %f, median %f\n", ...
        stats.mean(1), sqrt(stats.variance(1)), stats.quantiles(1, 1));
client.shutdown_server();