scripts/convdiff.m

Go to the documentation of this file.

% small script demonstrating the convdiff example from the M2AN Paper, that is
% also implemented in Dune. Later it will be possible to use the Dune
% implementation through a mex interface.

% Martin Drohmann 06.05.2009
% based on burgers_fv.m by
% Bernard Haasdonk 14.8.2007

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%% Select here, what is to be performed                        %%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
step = 0 % initialize the model_data structure. This step needs to be executed
          % once per Matlab session only
%step = 1 % single detailed simulation with given data and plot. Run
          % this with varying parameters mu until sure that scheme
          % is stable. Modify dt or the data-functions accordingly,
          % until a nice parameter-domain with uniformly stable
          % detailed scheme is obtained.
%step = 2 % generate dummy reduced basis from single trajectory and check, if
          % ei_interpolation with projection on this space maintains
          % result. A simple reduced simulation can also be
          % performed. All results should be visually identical
%step = 3 % generate reduced basis
%step = 6 % time measurement of reduced simulation and
          % use reduced basis in rb_demo_gui
%step = 7 % generate error-landscape over varying N and M
          % can take several hours!!!
%step = 8 % do runtime comparisons between detailed and reduced simulations

%steps = {2,5,7,8}
steps = {0,3,4};

for si=1:length(steps)

  step = steps{si};

  % output-filenames in rbmatlabresult
  detailedfname = 'convdiff_detailed_interpol.mat';

  %% parameters for visualization
  plot_params.show_colorbar = 1;
  plot_params.colorbar_mode = 'EastOutside';
  plot_params.plot          = @fv_plot;

  switch step
   case 0 % initialize model data
    descr             = convdiff_descr;
    [dmodel, rmodel] = gen_models(descr);

    model_data = gen_model_data(dmodel);
   case 1 % single detailed simulation and plot
    disp('performing single detailed simulation')
    tic
    sim_data  = detailed_simulation(dmodel, model_data);
    toc
    plot_sim_data(dmodel, model_data, sim_data, plot_params);

   case 2 % construct dummy reduced basis by single trajectory and simulate
    disp('detailed interpolated simulation for basis construction:')
    % mu_test      = cellfun(@mean, dmodel.mu_ranges);
    %  mu_test2  = cellfun(@min, params.mu_ranges);
    simple_bg_descr.rb_problem_type           = 'LinEvol';
    simple_bg_descr.detailed_data_constructor = @SimpleDetailedData;
    simple_bg_descr.do_orthonormalize         = true;
    simple_bg_descr.mu_list                   = {sim_data.mu};
    rsimple_model = gen_reduced_model(dmodel, simple_bg_descr);
    detailed_data = gen_detailed_data(rsimple_model, model_data);

    disp('reduced simulation:')

    tic;
    reduced_data = gen_reduced_data(rsimple_model, detailed_data);
    disp(['gen_reduced_data timining: ', num2str(toc)]);
    rsimple_model.N = get_rb_size(detailed_data);

    reduced_data = extract_reduced_data_subset(rsimple_model, reduced_data);
    set_mu(rsimple_model, sim_data.mu);

    tic;
    rb_sim_data = rb_simulation(rsimple_model, reduced_data);
    disp(['rb_simulation timing:      ', num2str(toc)]);
    tic;
    rb_sim_data = rb_reconstruction(rsimple_model, detailed_data, rb_sim_data);
    disp(['rb_reconstruction timing:  ', num2str(toc)]);

    plot_params.plot_title = 'reduced simulation result';
    plot_sim_data(dmodel, model_data, rb_sim_data, plot_params);
    plot_params.plot_title = 'detailed simulation result';
    plot_sim_data(dmodel, model_data, sim_data, plot_params);

   case 3 % reduced basis
    disp('constructing reduced basis')
    detailed_data = gen_detailed_data(rmodel, model_data);
    %  detailed_data = rb_basis_generation(detailed_data, ...
    %                  params);
    save(fullfile(rbmatlabresult,detailedfname),...
         'detailed_data','rmodel');
    leaf_dd = get_leaf(detailed_data, rmodel);
    max_err_sequence = get_field_on_active_child(detailed_data.datatree, 'max_err_sequence', rmodel);
    plot(max_err_sequence);
    set(gca,'Yscale','log');
    title('RB-generation error convergence');
   case 4
    mu_test   = cellfun(@mean, rmodel.mu_ranges);
    load(fullfile(rbmatlabresult,detailedfname));
    disp('reduced simulation:')
    rmodel.N = get_rb_size(detailed_data, rmodel);
    reduced_data = gen_reduced_data(rmodel, detailed_data);

    reduced_data = extract_reduced_data_subset(rmodel, reduced_data);
    tic;

    set_mu(rmodel, mu_test);
    rb_sim_data = rb_simulation(rmodel,reduced_data);
    t = toc;
    disp(['time for online phase: t = ',num2str(t)]);

    disp('full simulation:')
    set_mu(dmodel, get_mu(rmodel));
    tic;
    sim_data = detailed_simulation(dmodel,model_data);
    t = toc;
    disp(['time for detailed simulation: t = ',num2str(t)]);

    demo_rb_gui(rmodel,detailed_data,[],plot_params);

   case 7 % training-error landscape
    disp('warning: takes a few hours!');
    load(fullfile(rbmatlabresult,detailedfname));

    model.N = params.RB_stop_Nmax;

    offline_data = gen_reduced_data(model,detailed_data);
    range_params.plot_fields = { 'N', 'M' };
    range_params.max         = [size(detailed_data.RB,2), ...
                                size(detailed_data.BM{1},2) ];
    range_params.sample_size = [ 7, 12 ];
    range_params.mu_set_size = 100;

    testdir = ['convdiff_test_', range_params.mu_set_size ];

    params.run_name    = [testdir, '_step7'];
    params.tictoctable = true;

    output = stochastic_error_estimation(model, detailed_data, offline_data,...
                                         testdir, range_params, params);

  otherwise
    error('step-number is unknown!');

  end; % switch


end; % for si
%| \docupdate