porsche_script.m

function res=porsche_script(step,varargin)
%function res=porsche_script(step,varargin)
%
%
% GRID FUNCTIONS
% step 1: plot of the macrogrid
% step 2: plot of the microgrid
% step 3: plot microgrid elementwise 
% step 4: reshape a grid and plot it
% DETAILED SIMULATIONS
% step 10: FEM simulation of the model with the reference grid
% step 11: FEM simulation with the reshaped grid, plot of the original grid
% RB SIMULATIONS
% step 20: RB simulation of the model with the reference grid
% step 21: RB simulation with the reshaped grid, plot of the original grid
% COMPARISON OF DETAILED AND RB SIMULATION
% step 30: compare detailed and RB solution of reference car
% step 31: compare detailed and RB solution of original car
% OPTIMIZATION
% step 40: Detailed Optimization
% step 41: RB Optimization N=9
% step 42: RB Optimization N=16
% step 43: RB Optimization N=25
% step 44: Comparison of detailed and RB Optimization
%
%
% Oliver Zeeb, 01.09.11

if nargin == 0
  step = 1;
end;

if nargin < 2
  params = [];
end;

%set some parameters
params = [];
params.dimrange = 1;
pdeg = 1;
qdeg = 2 * pdeg; % quadrature_degree
params.pdeg = pdeg;
params.qdeg = qdeg;
params.no_plot = 0; %plotting active

res = [];


switch step
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    %%%%%%%%%%%%%%%%% GRID FUNCTIONS %%%%%%%%%%%%%%%%%
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    case 1 %plot of the macrogrid
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
        disp('=====================')
        disp('plot of the macrogrid')
        disp('=====================')
        model = porsche_model;
        model_data = gen_model_data(model);
        plot_params.color=[1 0 0];
        disp('plotting the macrogrid')
        plot(model_data.macrogrid,plot_params);
        axis equal;
   
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    case 2 %plot of the microgrid
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
        disp('=====================')
        disp('plot of the microgrid')
        disp('=====================')
        model = porsche_model;
        model_data = gen_model_data(model);
        disp('plotting the microgrid')
        plot(model_data.grid);
        axis equal;
        
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    case 3 % plot microgrid elementwise 
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
        disp('==============================')
        disp('plot the microgrid elementwise')
        disp('==============================')
        clear
        model = porsche_model;
        model_data = gen_model_data(model);  

        microgrid = model_data.grid;
        macrogrid = model_data.macrogrid;
        mic2mac=micro2macro_map(microgrid,macrogrid);

        for k=1:macrogrid.nelements
            elids=find(mic2mac==k);
               
            figure(1)
            axis equal
            gridp=gridpart(model_data.grid,elids);
            plot(gridp)
            axis([0 750 -50 350])
        
            disp(['Macrotriangle Nr', num2str(k)])

            disp('press enter to continue')
            pause;
            
            %close(1)
            hold on
        end
        
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%     
    case 4 %reshape a grid and plot it
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
        disp('==============================')
        disp('  reshape a grid and plot it  ')
        disp('==============================')
        model = porsche_model;
        model_data = gen_model_data(model);
        
        %reset the mus
        mu_mod = [0,30];
        model = change_mu(model, mu_mod);
        
        %reshape the grid
        grid_reshaped = grid_reshape(model, model_data.grid);
        
        %plot reshaped grid
        disp('plotting the original grid')
        figure(1)
        plot(grid_reshaped);
        axis equal
        title('original grid')
        
        %plot reference grid
        disp('plotting the reference grid')
        figure(2)
        plot(model_data.grid)
        axis equal
        title('reference grid')
        
        
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    %%%%%%%%%%%%% DETAILED SIMULATIONS %%%%%%%%%%%%%%%
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%    
        
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    case 10 % finite element assembly, solution and plot 
            % with the reference car
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
        disp('==============================================')
        disp('  FEM solution and plot of the reference car  ')
        disp('==============================================')

          %%%%%%%%%%%%%%% initialize model %%%%%%%%%%%%%%%%%%%
          model = porsche_model(params);
          model_data = gen_model_data(model);

          % do an FEM-simulation
          disp('calculating detailed simulation...')
          sim_data = detailed_simulation(model,model_data);

          
          %%%%%%%%%%%%%%% plotting %%%%%%%%%%%%%%%%%%%
          disp('plotting the PDE solution...')
          plot_sim_data(model,model_data,sim_data,1)
          % calculate and plot pressure and velocity
          disp('calculating gradient (velocity) and pressure')
          [sim_data.gradx, sim_data.grady]=disc_gradient([sim_data.uh.grid.X,sim_data.uh.grid.Y],sim_data.uh);
          sim_data.pressure=compute_pressure(model,sim_data.gradx, sim_data.grady);
          disp('plotting the velocity and pressure...')
          plot_sim_data(model,model_data,sim_data,3)
    
          res = sim_data;
          disp('porsche_script step 10 complete! Inspect workspace...')
          keyboard
          
   
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    case 11 %do a FEM Simulation with the reshaped grid 
            %and plot the solution on the reshaped grid
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    
        disp('=============================================')
        disp('  FEM solution and plot of the original car  ')
        disp('=============================================')

        %%%%%%%%%%%%%%% initialize model %%%%%%%%%%%%%%%%%%%
        model = porsche_model;
        model_data = gen_model_data(model);
        
        %reset the mus
        mu_mod = [0,30];
        model=change_mu(model, mu_mod);
        
        disp('calculating detailed simulation...')
        sim_data = detailed_simulation(model,model_data);
  
        %%%%%%%%%%%%%%% plotting %%%%%%%%%%%%%%%%%%%
        disp('plotting the PDE solution...')
        plot_sim_data(model,model_data,sim_data,2);
        
        disp('calculating gradient (velocity) and pressure')
        [sim_data.gradx,sim_data.grady]=disc_gradient([sim_data.uh.grid.X, sim_data.uh.grid.Y],sim_data.uh); % OLI: 14.02.11 ONLY VALID FOR PDEG = 1 !
        sim_data.pressure=compute_pressure(model, sim_data.gradx, sim_data.grady);
        disp('plotting the velocity and pressure...')
        plot_sim_data(model,model_data,sim_data,4)
        
        res=sim_data;
        disp('porsche_script step 11 complete! Inspect workspace...')
        keyboard
        
       

    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    %%%%%%%%%%%%%%%%% RB SIMULATIONS %%%%%%%%%%%%%%%%%
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
    
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    case 20  % RB-Simulation of the reference car
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    
        disp('=============================================')
        disp('  RB-Simulation of the reference car  ')
        disp('=============================================')

          %%%%%%%%%%%%%%% initialize model %%%%%%%%%%%%%%%%%%%
          disp('Setting the model...');
          model=porsche_model;
          disp('generating model_data');
          model_data=gen_model_data(model);
        
          %%%%%%%%%%%% RB-Simulation %%%%%%%%%%%%%%
          disp('generating detailed_data...');
          detailed_data=gen_detailed_data(model,model_data);
          disp('generating reduced_data...');
          reduced_data=gen_reduced_data(model,detailed_data);
          disp('doing a reduced simulation');
          rb_sim_data=rb_simulation(model,reduced_data);
          disp('rb_reconstruction...');
          rb_sim_data=rb_reconstruction(model,detailed_data,rb_sim_data);
          

          
          %%%%%%%%%%%%%%% plotting %%%%%%%%%%%%%%%%%%%
          disp('plotting the PDE solution...')
          plot_sim_data(model,model_data,rb_sim_data,1);
          
          % calculate and plot pressure and velocity
          disp('calculating gradient (velocity) and pressure')
          [rb_sim_data.gradx, rb_sim_data.grady]=disc_gradient([rb_sim_data.uh.grid.X, rb_sim_data.uh.grid.Y],rb_sim_data.uh);

          rb_sim_data.pressure=compute_pressure(model,rb_sim_data.gradx, rb_sim_data.grady);
          disp('plotting the velocity and pressure...')
          plot_sim_data(model, model_data, rb_sim_data,3)
    
          res = rb_sim_data;
          
          disp('porsche_script step 20 complete! Inspect workspace...')
          keyboard
    
          
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    case 21 %do a RB-Simulation with the original grid 
            %and plot the solution on the original grid
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    
        disp('=============================================')
        disp('      RB-Simulation of the original car      ')
        disp('=============================================')

        disp('Setting the model...');
        model = porsche_model;
        disp('generating model_data');
        model_data = gen_model_data(model);
        disp('generating detailed_data...');
        detailed_data=gen_detailed_data(model,model_data);
        disp('generating reduced_data...');
        reduced_data=gen_reduced_data(model,detailed_data);
        
        %reset the mus
        mu_mod=[0,30];
        model=change_mu(model, mu_mod);
        
        %RB-Simulation
        disp('doing a reduced simulation');
        rb_sim_data=rb_simulation(model,reduced_data);
        disp('rb_reconstruction...');
        rb_sim_data=rb_reconstruction(model,detailed_data,rb_sim_data);
        
        %%%%%%%%%%%%%%% plotting %%%%%%%%%%%%%%%%%%%
        disp('plotting the PDE solution...')
        plot_sim_data(model,model_data,rb_sim_data,2);
        
        
        disp('calculating gradient (velocity) and pressure')
        [rb_sim_data.gradx,rb_sim_data.grady]=disc_gradient([rb_sim_data.uh.grid.X, rb_sim_data.uh.grid.Y],rb_sim_data.uh); % OLI: 14.02.11 ONLY VALID FOR PDEG = 1 !
        rb_sim_data.pressure=compute_pressure(model, rb_sim_data.gradx, rb_sim_data.grady);
        disp('plotting the velocity and pressure...')
        plot_sim_data(model,model_data,rb_sim_data,4)
        
        res=rb_sim_data;
        
        disp('porsche_script step 21 complete! Inspect workspace...')
        keyboard
    

        
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    %%%% COMPARISON OF DETAILED AND RB SIMULATION %%%%
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
    
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    case 30  % Comparison with the reference car
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    
        disp('=====================================================================')
        disp(' Comparison between detailed an RB-simulation with the reference car ')
        disp('=====================================================================')

          %%%%%%%%%%%%%%% initialize model %%%%%%%%%%%%%%%%%%%
          disp('Setting the model...');
          model=porsche_model;
          disp('generating model_data');
          model_data=gen_model_data(model);
          
          %%%%%%%%% Detailed Simulation %%%%%%%%%%%
          disp('doing a detailed_simulation...');
          sim_data = detailed_simulation(model,model_data);
        
          %%%%%%%%%%%% RB-Simulation %%%%%%%%%%%%%%
          disp('generating detailed_data...');
          detailed_data=gen_detailed_data(model,model_data);
          disp('generating reduced_data...');
          reduced_data=gen_reduced_data(model,detailed_data);
          disp('doing a reduced simulation');
          rb_sim_data=rb_simulation(model,reduced_data);
          disp('rb_reconstruction...');
          rb_sim_data=rb_reconstruction(model,detailed_data,rb_sim_data);
          
          %comparison
          comparison_sim_data = sim_data;
          comparison_sim_data.uh=copy(sim_data.uh); %need to copy the fdf-class uh seperatly
          comparison_sim_data.uh.dofs = sim_data.uh.dofs-rb_sim_data.uh.dofs;
          
          %%%%%%%%%%%%%%% plotting %%%%%%%%%%%%%%%%%%%
          figure(1);
          disp('plotting the detailed solution...')
          plot_sim_data(model,model_data,sim_data,1);
          title('FEM-Solution on reference domain');
         
          figure(2);
          disp('plotting the reduced solution...')
          
          plot_sim_data(model,model_data,rb_sim_data,1);
          title('RB-Solution on reference domain');
          
          figure(3);
          disp('plotting comparison...')
          plot_sim_data(model,model_data,comparison_sim_data,1);
          title('RB-Comparison (FEM-RB)');

          res = {rb_sim_data, sim_data};
          disp('porsche_script step 30 complete! Inspect workspace...')
          keyboard
          
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    case 31  % Comparison of detailed and RB-solution of the original car
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
        disp('====================================================================')
        disp(' Comparison between detailed an RB-simulation with the original car ')
        disp('====================================================================')

          %%%%%%%%%%%%%%% initialize model %%%%%%%%%%%%%%%%%%%
          disp('Setting the model...');
          model=porsche_model;
          disp('generating model_data');
          model_data=gen_model_data(model);
          

          %%%%%%%%%%%% calculate RB-data %%%%%%%%%%%%%%
          disp('generating detailed_data...');
          detailed_data=gen_detailed_data(model,model_data);
          disp('generating reduced_data...');
          reduced_data=gen_reduced_data(model,detailed_data);
          
          %reset the mus
          mu_mod=[0,30];
          model=change_mu(model, mu_mod);
          
          %%%%%%%%% Detailed Simulation %%%%%%%%%%%
          disp('doing a detailed_simulation...');
          sim_data = detailed_simulation(model,model_data);
          %%%%%%%%%%%%% RB-Simulation %%%%%%%%%%%%%
          disp('doing a reduced simulation');
          rb_sim_data=rb_simulation(model,reduced_data);
          disp('rb_reconstruction...');
          rb_sim_data=rb_reconstruction(model,detailed_data,rb_sim_data);
          
          %comparison
          comparison_sim_data = sim_data;
          comparison_sim_data.uh=copy(sim_data.uh); %need to copy the fdf-class uh seperatly
          comparison_sim_data.uh.dofs = sim_data.uh.dofs-rb_sim_data.uh.dofs;
          
          %%%%%%%%%%%%%%% plotting %%%%%%%%%%%%%%%%%%%
          figure(1);
          disp('plotting the detailed solution...')
          plot_sim_data(model,model_data,sim_data,2);
          title('FEM-Solution on original domain');
          
          figure(2);
          disp('plotting the reduced solution...')
          plot_sim_data(model,model_data,rb_sim_data,2);
          title('RB-Solution on original domain');
          
          figure(3);
          disp('plotting comparison...')
          plot_sim_data(model,model_data,comparison_sim_data,2);
          title('RB-Comparison (FEM-RB)');
    
          res = {rb_sim_data, sim_data};
          disp('porsche_script step 31 complete! Inspect workspace...')
          keyboard
          
          
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    %%%%             OPTIMIZATION                 %%%%
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% 
    
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    case 40  % Detailed Optimization
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%     
    disp('==================================')
    disp('      Detailed Optimization       ')
    disp('==================================')
    
    disp('---')
    disp('setting the model... (hit RETURN to continue)');
    pause
    model=porsche_model;
    disp('...done');
    disp('---')
    disp('generating model_data...(hit RETURN to continue)');
    pause
    model_data=gen_model_data(model);
    disp('...done');
    disp('---')
    % optimization parameters
    model.optimization.constraint_circle_midpoint = [176.1103;161.1864];
    model.optimization.constraint_circle_radius = 10;
    mu_start=[178.06;151.38];
    model=set_mu(model,mu_start);
        
    disp('optimization with SQP... (CAREFUL: This can take a while!!! Hit RETURN to continue)')
    pause
    model.optimization.opt_mode='detailed';
    J_mu_0 = model.optimization.objective_function(model,model_data);
    opt_data=SQP_algorithm(model,model_data);
    res=opt_data;
    disp('...done');
    disp('---')
    disp('plotting the iterates... (hit RETURN to continue)')
    pause
    figure(1)
    plot(opt_data.mu_k(1,:),opt_data.mu_k(2,:),'--bs','LineWidth',2,...
                'MarkerEdgeColor','r',...
                'MarkerFaceColor','w',...
                'MarkerSize',5)
    title('Detailed Optimization, N=9 mu_0=[178.06;151.38]')
    disp('...done');
    disp('---')
    
    
    disp('-----------------------------------------')
    disp(' some information about the optimization ')
    disp('-----------------------------------------')
    disp(['mu_0          = (',num2str(mu_start(1)),', ', num2str(mu_start(2)),')'])
    disp(['mu_opt        = (', num2str(opt_data.mu_opt(1)), ', ', num2str(opt_data.mu_opt(2)),')'])
    disp(['J(0)          = ', num2str(J_mu_0, '%10.5e\n')])
    disp(['J(mu_opt)     = ', num2str(opt_data.obj_function, '%10.5e\n')])
    disp(['nr_iterations : ', num2str(opt_data.nr_iterations)])
    disp(['time online   : ', num2str(opt_data.time_elapsed),' sec'])
    disp('-----------------------------------------')
    
    disp('-----------------------------------------')
    disp('porsche_script step 40 complete! Inspect workspace...')
    keyboard
    
    
    
    
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    case 41  % RB Optimization  N=9
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  
    disp('====================================')
    disp('      Reduced Optimization  N=9     ')
    disp('====================================')
    
    disp('---')
    disp('setting the model... (hit RETURN to continue)');
    pause
    model=porsche_model;
    % optmization parameters
    model.optimization.constraint_circle_midpoint = [176.1103;161.1864];
    model.optimization.constraint_circle_radius = 10;
    mu_start=[178.06;151.38];
    model=set_mu(model,mu_start);
    %setting 3 points in x- and y-direction for the basis --> N=9
    midpoint_x=model.optimization.constraint_circle_midpoint(1);
    midpoint_y=model.optimization.constraint_circle_midpoint(2);
    snapshot_x=[midpoint_x-20, midpoint_x, midpoint_x+20];
    snapshot_y=[midpoint_y-20, midpoint_y, midpoint_y+20];
    snapshot_combination=createCombinations(snapshot_x, snapshot_y);
    for k=1:size(snapshot_combination,2)
        RB_mu_list{k} = snapshot_combination(:,k)';
    end
    
    model.RB_mu_list=RB_mu_list;    
    disp('...done');
    disp('---')
    disp('generating model_data...(hit RETURN to continue)');
    pause
    model_data=gen_model_data(model);
    disp('...done');
    disp('---')
    disp('generation of offline_data (hit RETURN to continue)')
    pause
    tic;
    disp('...generating detailed_data...');
    detailed_data=gen_detailed_data(model,model_data);
    disp('...generating reduced_data...');
    reduced_data=gen_reduced_data(model,detailed_data);
    time_offline = toc;
    disp('...done');
    disp('---')
    
    disp('optimization with SQP... (hit RETURN to continue)')
    pause
    model.optimization.opt_mode='reduced';
    J_mu_0 = model.optimization.objective_function(model,reduced_data);
    opt_data=SQP_algorithm(model,reduced_data);
    res=opt_data;
    disp('...done');
    disp('---')
    disp('plotting the iterates... (hit RETURN to continue)')
    pause
    figure(2)
    plot(opt_data.mu_k(1,:),opt_data.mu_k(2,:),'--bs','LineWidth',2,...
                'MarkerEdgeColor','r',...
                'MarkerFaceColor','w',...
                'MarkerSize',5)
    title('RB-Optimization, N=9 mu_0=[178.06;151.38]')
    disp('...done');
    disp('---')
    
    
    disp('-----------------------------------------')
    disp(' some information about the optimization ')
    disp('-----------------------------------------')
    disp(['mu_0          = (',num2str(mu_start(1)),', ', num2str(mu_start(2)),')'])
    disp(['mu_opt        = (', num2str(opt_data.mu_opt(1)), ', ', num2str(opt_data.mu_opt(2)),')'])
    disp(['J(0)          = ', num2str(J_mu_0, '%10.5e\n')])
    disp(['J(mu_opt)     = ', num2str(opt_data.obj_function, '%10.5e\n')])
    disp(['nr_iterations : ', num2str(opt_data.nr_iterations)])
    disp(['time online   : ', num2str(opt_data.time_elapsed),' sec'])
    disp(['time offline  : ', num2str(time_offline),' sec'])
    disp(['N             : ', num2str(length(model.RB_mu_list))])
    disp('-----------------------------------------')
    
    disp('-----------------------------------------')
    disp('porsche_script step 41 complete! Inspect workspace...')
    keyboard
    
    
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    case 42  % RB Optimization  N=16
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  
    disp('====================================')
    disp('      Reduced Optimization  N=16     ')
    disp('====================================')
    
    disp('---')
    disp('setting the model... (hit RETURN to continue)');
    pause
    model=porsche_model;
    % optimization parameters
    model.optimization.constraint_circle_midpoint = [176.1103;161.1864];
    model.optimization.constraint_circle_radius = 10;
    mu_start=[178.06;151.38];
    model=set_mu(model,mu_start);
    %setting 4 points in x- and y-direction for the basis --> N=16
    x_min=model.optimization.constraint_circle_midpoint(1)-20;
    x_max=model.optimization.constraint_circle_midpoint(1)+20;
    y_min=model.optimization.constraint_circle_midpoint(2)-20;
    y_max=model.optimization.constraint_circle_midpoint(2)+20;
    snapshot_x=linspace(x_min,x_max,4);
    snapshot_y=linspace(y_min,y_max,4);
    snapshot_combination=createCombinations(snapshot_x, snapshot_y);
    for k=1:size(snapshot_combination,2)
        RB_mu_list{k} = snapshot_combination(:,k)';
    end
    
    model.RB_mu_list=RB_mu_list;    
    disp('...done');
    disp('---')
    disp('generating model_data...(hit RETURN to continue)');
    pause
    model_data=gen_model_data(model);
    disp('...done');
    disp('---')
    disp('generation of offline_data (hit RETURN to continue)')
    pause
    tic;
    disp('...generating detailed_data...');
    detailed_data=gen_detailed_data(model,model_data);
    disp('...generating reduced_data...');
    reduced_data=gen_reduced_data(model,detailed_data);
    time_offline = toc;
    disp('...done');
    disp('---')
            
    disp('optimization with SQP... (hit RETURN to continue)')
    pause
    model.optimization.opt_mode='reduced';
    J_mu_0 = model.optimization.objective_function(model,reduced_data);
    opt_data=SQP_algorithm(model,reduced_data);
    res=opt_data;
    disp('...done');
    disp('---')
    disp('plotting the iterates... (hit RETURN to continue)')
    pause
    figure(2)
    plot(opt_data.mu_k(1,:),opt_data.mu_k(2,:),'--bs','LineWidth',2,...
                'MarkerEdgeColor','r',...
                'MarkerFaceColor','w',...
                'MarkerSize',5)
    title('RB-Optimization, N=16 mu_0=[178.06;151.38]')
    disp('...done');
    disp('---')
    
    
    disp('-----------------------------------------')
    disp(' some information about the optimization ')
    disp('-----------------------------------------')
    disp(['mu_0          = (',num2str(mu_start(1)),', ', num2str(mu_start(2)),')'])
    disp(['mu_opt        = (', num2str(opt_data.mu_opt(1)), ', ', num2str(opt_data.mu_opt(2)),')'])
    disp(['J(0)          = ', num2str(J_mu_0, '%10.5e\n')])
    disp(['J(mu_opt)     = ', num2str(opt_data.obj_function, '%10.5e\n')])
    disp(['nr_iterations : ', num2str(opt_data.nr_iterations)])
    disp(['time online   : ', num2str(opt_data.time_elapsed),' sec'])
    disp(['time offline  : ', num2str(time_offline),' sec'])
    disp(['N             : ', num2str(length(model.RB_mu_list))])
    disp('-----------------------------------------')
    
    disp('-----------------------------------------')
    disp('porsche_script step 42 complete! Inspect workspace...')
    keyboard
    
    
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    case 43  % RB Optimization  N=25
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%  
    disp('=====================================')
    disp('      Reduced Optimization  N=25     ')
    disp('=====================================')
    
    disp('---')
    disp('setting the model... (hit RETURN to continue)');
    pause
    model=porsche_model;
    % optimization parameters
    model.optimization.constraint_circle_midpoint = [176.1103;161.1864];
    model.optimization.constraint_circle_radius = 10;
    mu_start=[178.06;151.38];
    model=set_mu(model,mu_start);
    %setting 5 points in x- and y-direction for the basis --> N=25
    midpoint_x=model.optimization.constraint_circle_midpoint(1);
    midpoint_y=model.optimization.constraint_circle_midpoint(2);
    snapshot_x=[midpoint_x-40, midpoint_x-20, midpoint_x, midpoint_x+20, midpoint_x+40];
    snapshot_y=[midpoint_y-40, midpoint_y-20, midpoint_y, midpoint_y+20, midpoint_y+40];
    snapshot_combination=createCombinations(snapshot_x, snapshot_y);
    for k=1:size(snapshot_combination,2)
        RB_mu_list{k} = snapshot_combination(:,k)';
    end
    model.RB_mu_list=RB_mu_list;    
    disp('...done');
    disp('---')
    
    disp('generating model_data...(hit RETURN to continue)');
    pause
    model_data=gen_model_data(model);
    disp('...done');
    disp('---')
    disp('generation of offline_data (hit RETURN to continue)')
    pause
    tic;
    disp('...generating detailed_data...');
    detailed_data=gen_detailed_data(model,model_data);
    disp('...generating reduced_data...');
    reduced_data=gen_reduced_data(model,detailed_data);
    time_offline = toc;
    disp('...done');
    disp('---')
        
    disp('optimization with SQP... (hit RETURN to continue)')
    pause
    model.optimization.opt_mode='reduced';
    J_mu_0 = model.optimization.objective_function(model,reduced_data);
    opt_data=SQP_algorithm(model,reduced_data);
    res=opt_data;
    disp('...done');
    disp('---')
    disp('plotting the iterates... (hit RETURN to continue)')
    pause
    figure(2)
    plot(opt_data.mu_k(1,:),opt_data.mu_k(2,:),'--bs','LineWidth',2,...
                'MarkerEdgeColor','r',...
                'MarkerFaceColor','w',...
                'MarkerSize',5)
    title('RB-Optimization, N=25 mu_0=[178.06;151.38]')
    disp('...done');
    disp('---')
    
    
    disp('-----------------------------------------')
    disp(' some information about the optimization ')
    disp('-----------------------------------------')
    disp(['mu_0          = (',num2str(mu_start(1)),', ', num2str(mu_start(2)),')'])
    disp(['mu_opt        = (', num2str(opt_data.mu_opt(1)), ', ', num2str(opt_data.mu_opt(2)),')'])
    disp(['J(0)          = ', num2str(J_mu_0, '%10.5e\n')])
    disp(['J(mu_opt)     = ', num2str(opt_data.obj_function, '%10.5e\n')])
    disp(['nr_iterations : ', num2str(opt_data.nr_iterations)])
    disp(['time online   : ', num2str(opt_data.time_elapsed),' sec'])
    disp(['time offline  : ', num2str(time_offline),' sec'])
    disp(['N             : ', num2str(length(model.RB_mu_list))])
    disp('-----------------------------------------')
    
    disp('-----------------------------------------')
    disp('porsche_script step 43 complete! Inspect workspace...')
    keyboard
    
    
    
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    case 44  % Comparison of detailed and RB Optimization
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    disp('==============================================')
    disp('  Comparison of detailed and RB optimization  ')
    disp('==============================================')
    
    disp('---')
    disp('setting the model... (hit RETURN to continue)');
    pause
    model=porsche_model;
    disp('...done');
    disp('---')
    disp('generating model_data...(hit RETURN to continue)');
    pause
    model_data=gen_model_data(model);
    %%% optimization parameters %%%
    model.optimization.constraint_circle_midpoint = [176.1103;161.1864];
    model.optimization.constraint_circle_radius = 10;
    mu_start=[178.06;151.38];


    midpoint_x=model.optimization.constraint_circle_midpoint(1);
    midpoint_y=model.optimization.constraint_circle_midpoint(2);
    model=set_mu(model,mu_start);
    disp('...done');
    disp('---')
        
    
    %---------------------
    %------ detailed -----
    %---------------------    
    disp('detailed optimization with SQP... (CAREFUL! This can take a while!!! Hit RETURN to continue)')
    pause
    model.optimization.opt_mode='detailed';
    J_mu_0_det = model.optimization.objective_function(model,model_data);
    opt_data_detailed=SQP_algorithm(model,model_data);
    disp('...done');
    disp('---')
    
    %----------------------
    %------- N9 -----------
    %----------------------
    disp('N=9: generation of offline_data (hit RETURN to continue)')
    pause
    
    snapshot_x_N9=[midpoint_x-20, midpoint_x, midpoint_x+20];
    snapshot_y_N9=[midpoint_y-20, midpoint_y, midpoint_y+20];
    snapshot_combination_N9=createCombinations(snapshot_x_N9, snapshot_y_N9);
    for k=1:size(snapshot_combination_N9,2)
        RB_mu_list_N9{k} = snapshot_combination_N9(:,k)';
    end
    model.RB_mu_list=RB_mu_list_N9;    
    tic;
    disp('N=9: ...generating detailed_data...');
    detailed_data_N9=gen_detailed_data(model,model_data);
    disp('N=9: ...generating reduced_data...');
    reduced_data_N9=gen_reduced_data(model,detailed_data_N9);
    time_offline_N9 = toc;
    disp('...done');
    disp('---')
    disp('N=9: RB optimization with SQP... (hit RETURN to continue)')
    pause
    model.optimization.opt_mode='reduced';
    J_mu_0_red_N9 = model.optimization.objective_function(model,reduced_data_N9);
    opt_data_red_N9 = SQP_algorithm(model,reduced_data_N9);
    disp('...done');
    disp('---')
    
    
    %----------------------
    %------- N16 ----------
    %----------------------
    disp('N=16: generation of offline_data (hit RETURN to continue)')
    pause
    x_min_N16=model.optimization.constraint_circle_midpoint(1)-20;
    x_max_N16=model.optimization.constraint_circle_midpoint(1)+20;
    y_min_N16=model.optimization.constraint_circle_midpoint(2)-20;
    y_max_N16=model.optimization.constraint_circle_midpoint(2)+20;
    snapshot_x_N16=linspace(x_min_N16,x_max_N16,4);
    snapshot_y_N16=linspace(y_min_N16,y_max_N16,4);
    snapshot_combination_N16=createCombinations(snapshot_x_N16, snapshot_y_N16);
    for k=1:size(snapshot_combination_N16,2)
        RB_mu_list_N16{k} = snapshot_combination_N16(:,k)';
    end
    model.RB_mu_list=RB_mu_list_N16;    
    tic;
    disp('N=16: ...generating detailed_data...');
    detailed_data_N16=gen_detailed_data(model,model_data);
    disp('N=16: ...generating reduced_data...');
    reduced_data_N16=gen_reduced_data(model,detailed_data_N16);
    time_offline_N16 = toc;
    disp('...done');
    disp('---')
    disp('N=16: RB optimization with SQP... (hit RETURN to continue)')
    pause
    model.optimization.opt_mode='reduced';
    J_mu_0_red_N16 = model.optimization.objective_function(model,reduced_data_N16);
    opt_data_red_N16 = SQP_algorithm(model,reduced_data_N16);
    disp('...done');
    disp('---')
    %----------------------
    %------- N25 ----------
    %----------------------
    disp('N=25: generation of offline_data (hit RETURN to continue)')
    pause
    
    snapshot_x_N25=[midpoint_x-40, midpoint_x-20, midpoint_x, midpoint_x+20, midpoint_x+40];
    snapshot_y_N25=[midpoint_y-40, midpoint_y-20, midpoint_y, midpoint_y+20, midpoint_y+40];
    snapshot_combination_N25=createCombinations(snapshot_x_N25, snapshot_y_N25);
    for k=1:size(snapshot_combination_N25,2)
        RB_mu_list_N25{k} = snapshot_combination_N25(:,k)';
    end
    model.RB_mu_list=RB_mu_list_N25;    
    
    tic;
    disp('N=25: ...generating detailed_data...');
    detailed_data_N25=gen_detailed_data(model,model_data);
    disp('N=25: ...generating reduced_data...');
    reduced_data_N25=gen_reduced_data(model,detailed_data_N25);
    time_offline_N25 = toc;
    disp('...done');
    disp('---')
    
    disp('N=25: RB optimization with SQP... (hit RETURN to continue)')
    pause
    model.optimization.opt_mode='reduced';
    J_mu_0_red_N25 = model.optimization.objective_function(model,reduced_data_N25);
    opt_data_red_N25 = SQP_algorithm(model,reduced_data_N25);
    disp('...done');
    disp('---')
    
    %-------------------------
    
    res={opt_data_detailed,opt_data_red_N9,opt_data_red_N16,opt_data_red_N25}; %setting the returning variable
    
    %plotting the iterates
    disp('plotting the iterates... (hit RETURN to continue)')
    pause
    figure(1)
    plot(opt_data_detailed.mu_k(1,:),opt_data_detailed.mu_k(2,:),'-rs','LineWidth',2,...
                'MarkerEdgeColor','r',...
                'MarkerFaceColor','w',...
                'MarkerSize',5)
    hold on
    plot(opt_data_red_N9.mu_k(1,:),opt_data_red_N9.mu_k(2,:),'--cs','LineWidth',2,...
                'MarkerEdgeColor','c',...
                'MarkerFaceColor','w',...
                'MarkerSize',5)
    plot(opt_data_red_N16.mu_k(1,:),opt_data_red_N16.mu_k(2,:),'-.gs','LineWidth',2,...
                'MarkerEdgeColor','g',...
                'MarkerFaceColor','w',...
                'MarkerSize',5)            
    plot(opt_data_red_N25.mu_k(1,:),opt_data_red_N25.mu_k(2,:),'--bs','LineWidth',2,...
                'MarkerEdgeColor','b',...
                'MarkerFaceColor','w',...
                'MarkerSize',5)            
    legend('detailed', 'reduced, N=9', 'reduced, N=16','reduced, N=25','Location','NorthWest')
    disp('...done');
    disp('---')
    
    %calculating and plotting the error
    disp('plotting the error... (hit RETURN to continue)')
    pause
    k_det = opt_data_detailed.nr_iterations;
    k_red_N9 = opt_data_red_N9.nr_iterations;
    k_red_N16 = opt_data_red_N16.nr_iterations;
    k_red_N25 = opt_data_red_N25.nr_iterations;
    k_min_N9= min(k_det,k_red_N9);
    k_min_N16= min(k_det,k_red_N16);
    k_min_N25= min(k_det,k_red_N25);
    % calculating error N9
    differenz_mu_N9(:,1:k_min_N9) = opt_data_red_N9.mu_k(:,1:k_min_N9) - opt_data_detailed.mu_k(:,1:k_min_N9);
    for k=1:k_min_N9
        norm_diff_N9(k)=norm(differenz_mu_N9(:,k));
    end
    % calculating error N16
    differenz_mu_N16(:,1:k_min_N16) = opt_data_red_N16.mu_k(:,1:k_min_N16) - opt_data_detailed.mu_k(:,1:k_min_N16);
    for k=1:k_min_N16
        norm_diff_N16(k)=norm(differenz_mu_N16(:,k));
    end
    % calculating error N25
    differenz_mu_N25(:,1:k_min_N25) = opt_data_red_N25.mu_k(:,1:k_min_N25) - opt_data_detailed.mu_k(:,1:k_min_N25);
    for k=1:k_min_N25
        norm_diff_N25(k)=norm(differenz_mu_N25(:,k));
    end

    % plot of the error
    figure(2)
    plot_achse_N9=0:1:k_min_N9-1;
    plot_achse_N16=0:1:k_min_N16-1;
    plot_achse_N25=0:1:k_min_N25-1;
    
    plot(plot_achse_N9,norm_diff_N9,'c')
    hold on
    plot(plot_achse_N16,norm_diff_N16,'g')
    plot(plot_achse_N25,norm_diff_N25,'b')
    legend('reduced, N=16', 'reduced, N=9','reduced, N=25')
    title('Error');
    xlabel('k');
    disp('...done');
    disp('---')
    
    % displaying information about the optimization
    disp('-----------------------------------------')
    disp(' some information about the optimization ')
    disp('-----------------------------------------')
    disp(['mu_0               = (',num2str(mu_start(1)),'  , ', num2str(mu_start(2)),')'])
    
    disp(['mu_opt_detailed    = (', num2str(opt_data_detailed.mu_opt(1)), ', ', num2str(opt_data_detailed.mu_opt(2)),')'])
    disp(['mu_opt_red_N9      = (', num2str(opt_data_red_N9.mu_opt(1)), ', ', num2str(opt_data_red_N9.mu_opt(2)),')'])
    disp(['mu_opt_red_N16     = (', num2str(opt_data_red_N16.mu_opt(1)), ', ', num2str(opt_data_red_N16.mu_opt(2)),')'])
    disp(['mu_opt_red_N25     = (', num2str(opt_data_red_N25.mu_opt(1)), ', ', num2str(opt_data_red_N25.mu_opt(2)),')'])
    disp('-----------------------------------------')
    disp(['J_det(mu_0)        = ', num2str(J_mu_0_det, '%10.5e\n')])
    disp(['J_red_N9(mu_0)     = ', num2str(J_mu_0_red_N9, '%10.5e\n')])
    disp(['J_red_N16(mu_0)    = ', num2str(J_mu_0_red_N16, '%10.5e\n')])
    disp(['J_red_N25(mu_0)    = ', num2str(J_mu_0_red_N25, '%10.5e\n')])
    disp('-----------------------------------------')
    disp(['J_det(mu_opt)      = ', num2str(opt_data_detailed.obj_function, '%10.5e\n')])
    disp(['J_red_N9(mu_opt)   = ', num2str(opt_data_red_N9.obj_function, '%10.5e\n')])
    disp(['J_red_N16(mu_opt)  = ', num2str(opt_data_red_N16.obj_function, '%10.5e\n')])
    disp(['J_red_N25(mu_opt)  = ', num2str(opt_data_red_N25.obj_function, '%10.5e\n')])
    disp('-----------------------------------------')
    disp(['nr_iterations_det     : ', num2str(opt_data_detailed.nr_iterations)])
    disp(['nr_iterations_red_N9  : ', num2str(opt_data_red_N9.nr_iterations)])
    disp(['nr_iterations_red_N16 : ', num2str(opt_data_red_N16.nr_iterations)])
    disp(['nr_iterations_red_N25 : ', num2str(opt_data_red_N25.nr_iterations)])
    disp('-----------------------------------------')
    disp(['time_online_det       : ', num2str(opt_data_detailed.time_elapsed),' sec'])
    disp(['time_online_red_N9    : ', num2str(opt_data_red_N9.time_elapsed),' sec'])
    disp(['time_online_red_N16   : ', num2str(opt_data_red_N16.time_elapsed),' sec'])
    disp(['time_online_red_N25   : ', num2str(opt_data_red_N25.time_elapsed),' sec'])
    disp('-----------------------------------------')
    disp(['time offline N9   : ', num2str(time_offline_N9),' sec'])
    disp(['time offline N16  : ', num2str(time_offline_N16),' sec'])
    disp(['time offline N25  : ', num2str(time_offline_N25),' sec'])
    disp('-----------------------------------------')
    disp('-----------------------------------------')

    
    disp('porsche_script step 44 complete! Inspect workspace...')
    keyboard
        
    otherwise
        error('unknown step number');
end;