lin_evol_opt_rb_fd_derivative.m

function U_der = lin_evol_opt_rb_fd_derivative(model, detailed_data,reduced_data, mu)
%function U_der = lin_evol_opt_rb_fd_derivative(model, reduced_data, mu)
%
% function calculating the derivative solution via a finite difference
% procedure for the parameter mu.
%
% Markus Dihlmann 03.02.2011

if isempty(mu)
    mu = get_mu(model);
end

old_mu = get_mu(model);



compute_derivative_old = model.compute_derivative_info;

model.compute_derivative_info =0;

sim_data = rb_simulation(model, reduced_data);
sim_data = rb_reconstruction_derivative(model,detailed_data, sim_data);

U_fix=sim_data.U;

h=0.0001; %stepsize

for i=1:length(mu)
    %is mu+h inside mu_range?
    if (mu(i)+h>=model.mu_ranges{i}(1))&&(mu(i)+h<=model.mu_ranges{i}(2))
        mu(i) = mu(i) + h;
        model = set_mu(model,mu);
        sim_data = rb_simulation(model, reduced_data);
        sim_data = rb_reconstruction_derivative(model, detailed_data,sim_data);
        mu(i) = mu(i) -h;
        model = set_mu(model,mu);
        
        U_der{i} = (sim_data.U-U_fix)./h;
    else
        mu(i) = mu(i) - h;
        model = set_mu(model,mu);
        sim_data = rb_simulation(model, reduced_data);
        sim_data = rb_reconstruction_derivative(model, detailed_data,sim_data);
        mu(i) = mu(i) +h;
        model = set_mu(model,mu);
        
        U_der{i} = (-sim_data.U+U_fix)./h;
    end
    
end

model = set_mu(model, old_mu);
model.compute_derivative_info = compute_derivative_old;