test_ei_detailed.m

Go to the documentation of this file.

function ok = test_ei_detailed
  % Tests for the generation of an empirical operator interpolant
  %
  % Test cases:
  %  - 1. Generation of one EI for one operator
  %  - 2. Generation of one EI for two operators
  %  - 3. Generation of two EIs for two operators
  
  ok = true;

  % declare minimal model
  mm                                  = [];
  mm.name                             = 'minimal_model';
  mm.mu_names                         = {'par1', 'par2', 'par3'};
  mm.mu_ranges                        = {[0 1], [0 2], [0 3]};
  mm.par1                             = 0;
  mm.par2                             = 0;
  mm.par3                             = 0;
  mm.debug                            = 0;
  mm.force_delete                     = 1;
  mm.filecache_ignore_fields_in_model = {};
  mm.detailed_simulation              = @(model, model_data) struct('u', model.get_mu(model));
  mm.L_I_local_ptr                    = @fv_implicit_space;
  mm.num_diff_flux_ptr                = @fv_num_diff_flux_gradient;
  mm.fv_impl_diff_weight              = 1.0;
  mm.fv_impl_conv_weight              = 0.0;
  mm.fv_impl_react_weight             = 0.0;

  mm.mass_matrix                      = @fv_mass_matrix;

  mm.gridtype                              = 'rectgrid';
  mm.xrange                                = [0 1];
  mm.yrange                                = [0 1];
  mm.xnumintervals                         = 30;
  mm.ynumintervals                         = 30;
  mm.bnd_rect_corner1                      = [-1 -1];
  mm.bnd_rect_corner2                      = [2 2];
  mm.bnd_rect_index                        = -2;
  mm.verbose                               = 9;
  mm.diffusivity_ptr                       = @diffusivity_exponential;
  mm.diffusivity_derivative_ptr            = @diffusivity_exponential_derivative;
  mm.diff_k0                               = 0.00;
  mm.diff_m                                = 1.0;
  mm.diff_p                                = 2.0;
  mm.decomp_mode                           = 0;
  mm.laplacian_ptr                         = @(glob, U, model) U;
  mm.laplacian_derivative_ptr              = @(glob, U, model) ones(length(U),1);
  mm.filecache_velocity_matrixfile_extract = 0;
  mm.neumann_values_ptr                    = @neumann_values_homogeneous;
  mm.c_neu                                 = 0.00;
  mm.get_inner_product_matrix              = @(detailed_data) detailed_data.W;
  mm.l2_error_sequence_algorithm           = @fv_l2_error;

  mm.new = 1;
  mm.operators_ptr = @(x) x;
  mm.init_values_algorithm = @(x) x;
  mm.error_norm = 'l2';

  bg_descr.rb_problem_type = 'LinEvol';
  bg_descr.detailed_data_constructor = @Test.DetailedData;
  bg_descr.reduced_data_constructor  = @LinEvol.ReducedData;

  dmodel = LinEvol.DetailedModel(mm);

  disp('Test 1');
  M_train      = ParameterSampling.Random(5);
  rb_generator = SnapshotsGenerator.Random(dmodel, 'random', [], true, false);
  ei_gen1      = SnapshotsGenerator.SpaceOpEvals(dmodel, 'implicit', rb_generator, mm.L_I_local_ptr, true, false);
  ei_extension1 = Greedy.Plugin.EI(ei_gen1);
  ei_extension1.stop_Mmax = 200;

  ei_greedy1   = Greedy.Algorithm(ei_extension1, M_train);
  ei_greedy1.stop_timeout = 60*60;
  ei_greedy1.stop_epsilon = 1e-9;

  bg_descr.bg_algorithm = ei_greedy1;
  red_model1     = Test.ReducedModel(dmodel, bg_descr);

  model_data     = gen_model_data(dmodel);
  try
    evalc('detailed_data1 = gen_detailed_data(red_model1, model_data);');
  catch err
    ok = false;
    disp('> failed');
  end
  if ok
    disp('> passed');
  end

  max_err_sequence = get_field_on_active_child(detailed_data1.datatree, ...
      'max_err_sequence', red_model1);
  ok = ok & max_err_sequence(end) < 1e-10;
  ok = ok & length(max_err_sequence) == 6;

  disp('Test 2');
  ei_gen2      = SnapshotsGenerator.SpaceOpEvals(dmodel, 'explicit', rb_generator, mm.L_I_local_ptr, true, false);

  ei_gen_sum    = ei_gen1 + ei_gen2;
  ei_extension2 = Greedy.Plugin.EI(ei_gen_sum);
  ei_extension2.stop_Mmax = 200;

  ei_greedy2   = Greedy.Algorithm(ei_extension2, M_train);
  ei_greedy2.stop_timeout = 60*60;
  ei_greedy2.stop_epsilon = 1e-9;

  bg_descr.bg_algorithm = ei_greedy2;
  red_model2     = Test.ReducedModel(dmodel, bg_descr);
  try
    evalc('detailed_data2 = gen_detailed_data(red_model2, model_data);');
  catch err
    ok = false;
    disp('> failed');
  end
  if ok
    disp('> passed');
  end

  disp('Test 3');
  max_err_sequence = get_field_on_active_child(detailed_data2.datatree, ...
      'max_err_sequence', red_model2);
  ok = ok & max_err_sequence(end) < 1e-10;
  ok = ok & length(max_err_sequence) == 6;

  ei_extension3 = Greedy.Plugin.EI(ei_gen2);
  ei_extension3.stop_Mmax = 300;

  ei_greedy3_exp = Greedy.Algorithm(ei_extension3, M_train);
  ei_greedy3_exp.stop_epsilon = 1e-9;
  ei_greedy3     = [ei_greedy1, ei_greedy3_exp];

  bg_descr.bg_algorithm = ei_greedy3;
  red_model3     = Test.ReducedModel(dmodel, bg_descr);
  try
    evalc('detailed_data3 = gen_detailed_data(red_model3, model_data);');
  catch err
    ok = false;
    disp('> failed');
  end
  if ok
    disp('> passed');
  end

  max_err_sequence = get_field_on_active_child(detailed_data3.datatree, ...
      'max_err_sequence', red_model3, 'explicit');
  ok = ok & max_err_sequence(end) < 1e-10;
  ok = ok & length(max_err_sequence) == 6;

  max_err_sequence = get_field_on_active_child(detailed_data3.datatree, ...
      'max_err_sequence', red_model3, 'implicit');
  ok = ok & max_err_sequence(end) < 1e-10;
  ok = ok & length(max_err_sequence) == 6;



end