SCM.m

function SCM_data = SCM(model,params,operator,innerproduct_matrix)
%function SCM_data = SCM(model,params,operator,innerproduct_matrix)

% Offline programm of the successive constraint method (scm) for stability
% factors.
%
% Requires 
%  model : standard model
%  params: parameters for the scm algorithm (see below)
%  operator:  (function handle or model_data) returning the system matrix 
%              function handle could be (operator =@ (model) model.operators(model,operator);
%             the funtion operator should depend on model.decomp_mode        
%  innerproduct_matrix: if operator==model_data, the innerpoduct_matrix is optional
% 
%
% For simple cases, the call SCM_data=SCM(model,params,model_data) is also
% possible.
%
% Returns SCM_data, which is necessary for scm online computations 
% via SCMonline  (not SCM_online).
% 
% See scm_coercive_demo() or scm_infsup_demo() for demonstrations
%
%
% scm parameters (default values are set in scm_infsup_greedy
% or scm_coercive_demo):
%
%  params.type = 'coercive'  :    Decides which algorithm is used and 
%              = 'infsup'         which constant should be computed 
%    
%  params.Ma   = int           : number of parameters with an already
%                                computed constant near to mu
%  params.Mp   = int           : number of parameters without computed
%                                constant near to mu
%  params.NumberOfParameters = int  : size of the test set used in the SCM
%                                 Greedy procedure
%  params.eps  = double         : error boundary of the Greedy procedure
%
%  params.max_iterations = int  : maximal number of Greedy iterations
%
%  params.dirichlet_gids        : Indices to remove the dirichlet boundary
%                               conditions
%
%  params.parallel    = 0        non parallel version
%                     = 1        parallel (pool stays active after
%                                                              computation)
%                     = 2        parallel (closing pool afterwards)
%
%
%
% David Kreplin 05.01.2015


% Allows the call SCM_data=SCM(model,params,model_data);
if ~isa(operator, 'function_handle')
    
    if nargin == 3
        innerproduct_matrix = model.get_inner_product_matrix(operator);
    end
    operator =@ (model) model.operators(model,operator);   
end
  
if ~isfield(params, 'type')
    params.type = 'infsup';
    warning('No scm type defined! Setting type = "infsup". For coercive constants type="coercive" will be faster!');
end


if ~isfield(params, 'dirichlet_gids')
    params.dirichlet_gids = [];
end

%coercive case
if strcmp(params.type,'coercive');
    SCM_data = scm_coercive_greedy(model, params, operator, innerproduct_matrix);
end

%infsup case
if strcmp(params.type,'infsup');
    SCM_data = scm_infsup_greedy(model, params, operator, innerproduct_matrix);
end

disp('-- SCM: finished!');


end