rbmatlab  1.16.09
 All Classes Namespaces Files Functions Variables Modules Pages
advection_fv_output.m
Go to the documentation of this file.
1 function advection_fv_output(step)
2 % small script implementing a simple advection example
3 % for producing matrices to be used in the RB-DS framework
4 % Discretization with FV Functions.
5 %
6 % mu_names = {'vx_weight','vy_weight'};
7 % mu_ranges = [0,1]^2
8 %
9 % possible steps:
10 %
11 % step = 1; % initialization of model and plot of model data
12 % step = 2; % detailed simulation of lin_evol
13 % step = 2.5; % conversion to DS primal model, lin_ds detailed simulation
14 % step = 3; % conversion to DS primal model, DS detailed simulation
15 % and comparison with lin_evol
16 % step = 4 % compute several trajectories.
17 % step = 4.5 % compute reduced basis based on trajecories.
18 % step = 4.75 % comparison of detailed and reduced simulation
19 % step = 5 % figures of detailed simulations based on step 4
20 % step = 6 % figures of reduced simuations based on step 4
21 % step = 7 % time measurements based on step 4
22 % step = 8; % compare original and accelerated ds_model
23 % step = 9; % generate greedy basis on subset of parameter domain
24 % step = 10; % reduced simulations and error plots
25 % step = 11; % test of matlab ode solver for resulting sparse system
26 % step = 12; % further plots for MCMDS paper
27 
28 % Bernard Haasdonk 25.8.2009
29 
30 % cartesian grid
31 
32 if nargin < 1
33  step = 1;
34 end;
35 
36 params = [];
37 model = [];
38 %params.coarse_factor = 2; % fact 2 => ndofs = 16384, nt = 1024, 15s detailed
39 %params.coarse_factor = 1; % fact 1 => ndofs = 32768, nt = 2048, 60s detailed
40 params.coarse_factor = 2;
41 %verbose(1);
42 
43 switch step
44 
45  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
46  %step = 1; % initialization of model and plot of model data
47  case 1
48  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
49 
50 % verbose(1,'initializaton of lin_evol_model and plot model_data');
51  disp('initializaton of lin_evol_model and plot model_data');
52  model = advection_fv_output_model(params);
53  model_data = gen_model_data(model);
54  params.axis_equal = 1;
55  params.axis_tight = 1;
56  plot(model_data.grid,params);
57  % inspect(model_data.grid)
58 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
59  %step = 2; % detailed simulation of lin_evol
60 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
61  case 2
62 % verbose(1,'detailed simulation of lin_evol_model and plot');
63  disp('detailed simulation of lin_evol_model and plot');
64  model = advection_fv_output_model(params);
65  model_data = gen_model_data(model);
66  % model.debug = 1;
67  % disp('model debugging turned on.')
68  model.cone_weight = 0.0;
69  model.vx_weight = 0.0;
70  model.vy_weight = 0.0;
71  model.vx_weight = 1;
72  model.vy_weight = 1;
73  model.cone_weight = 1.0;
74 % model.vx_weight = 1;
75 % model.vy_weight = 0;
76  sim_data = detailed_simulation(model,model_data);
77  plot_sim_data(model,model_data,sim_data,[]);
78 
79 
80 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
81 %step = 2.5; % conversion to DS primal model, DS detailed simulation
82 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
83  case 2.5
84  disp('initialization of lin_ds model from lin_evol model');
85  % verbose(1,'initialization of lin_ds model from lin_evol model and ...
86  % plot');
87 
88  lin_evol_model = advection_fv_output_model(params);
89  %lin_evol_model.debug = 1;
90  lin_evol_model_data = gen_model_data(lin_evol_model);
91 
92  % set some additional fields in model which will be copied into
93  % ds model depending on whether constants are known or not:
94  estimate_constants = 0;
95  if estimate_constants
96  %%% if constant estimaton has to be performed:
97  lin_evol_model.estimate_lin_ds_nmu = 4;
98  lin_evol_model.estimate_lin_ds_nX = 10;
99  lin_evol_model.estimate_lin_ds_nt = 4;
100  else
101  %%% otherwise:
102  % the following values are rough bounds, so could be
103  % non-rigorous, then larger search must be performed.
104  lin_evol_model.state_bound_constant_C1 = 1;
105  lin_evol_model.output_bound_constant_C2 = 1.2916;
106  end;
107 
108  ds_model = lin_ds_from_lin_evol_model(lin_evol_model);
109  %ds_model.u_function_ptr = @(params) 0.1; % define new
110 % ds_model.theta = 1;
111  ds_model.theta = 0;
112  ds_model_data = gen_model_data(ds_model);
113  mu = [1,1];
114  %mu = [0.5,1,1];
115  % set mu in model and base_model!!!
116  ds_model = ds_model.set_mu(ds_model,mu);
117 
118  ds_sim_data = detailed_simulation(ds_model,ds_model_data);
119 
120  params.no_lines = 1;
121  params.clim = [0,1];
122  params.plot_title='detailed implicit lin_ds';
123  lin_evol_from_ds_sim_data.U = ds_sim_data.X;
124  lin_evol_from_ds_sim_data.y = ds_sim_data.Y;
125  plot_sim_data(lin_evol_model,lin_evol_model_data,lin_evol_from_ds_sim_data,params);
126 
127  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
128 %step = 3; % conversion to DS primal model, DS detailed simulation
129  % and comparison with lin_evol
130 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
131 
132  case 3 % initialization of DS model and lin_ds detailed simulation
133  disp('initialization of lin_ds model from lin_evol model and plot');
134 % verbose(1,'initialization of lin_ds model from lin_evol model and ...
135 % plot');
136 
137  lin_evol_model = advection_fv_output_model(params);
138  %lin_evol_model.debug = 1;
139  lin_evol_model_data = gen_model_data(lin_evol_model);
140 
141  % set some additional fields in model which will be copied into
142  % ds model depending on whether constants are known or not:
143  estimate_constants = 0;
144  if estimate_constants
145  %%% if constant estimaton has to be performed:
146  lin_evol_model.estimate_lin_ds_nmu = 4;
147  lin_evol_model.estimate_lin_ds_nX = 10;
148  lin_evol_model.estimate_lin_ds_nt = 4;
149  else
150  %%% otherwise:
151  % the following values are rough bounds, so could be
152  % non-rigorous, then larger search must be performed.
153  lin_evol_model.state_bound_constant_C1 = 1;
154  lin_evol_model.output_bound_constant_C2 = 1.2916;
155  end;
156 
157  ds_model = lin_ds_from_lin_evol_model(lin_evol_model);
158 % ds_model.theta = 0.0;
159  %ds_model.u_function_ptr = @(params) 0.1; % define new
160  ds_model_data = gen_model_data(ds_model);
161  %mu = [0.5,1,1];
162  mu = [1,1];
163  % set mu in model and base_model!!!
164  ds_model = ds_model.set_mu(ds_model,mu);
165 
166  ds_sim_data = detailed_simulation(ds_model,ds_model_data);
167  lin_evol_model = lin_evol_model.set_mu(lin_evol_model,mu);
168  lin_evol_sim_data = detailed_simulation(lin_evol_model,...
169  lin_evol_model_data);
170  err = lin_evol_sim_data.U - ds_sim_data.X;
171  disp(['max error in DOFs of lin_evol and ds simulation: ',...
172  num2str(max(abs(err(:))))]);
173 
174  params.plot_title='detailed lin_evol';
175  params.no_lines = 1;
176  plot_sim_data(lin_evol_model,lin_evol_model_data,lin_evol_sim_data,params);
177  params.plot_title='detailed lin_ds';
178  lin_evol_from_ds_sim_data.U = ds_sim_data.X;
179  lin_evol_from_ds_sim_data.y = ds_sim_data.Y;
180  plot_sim_data(lin_evol_model,lin_evol_model_data,lin_evol_from_ds_sim_data,params);
181  params.plot_title='error explicit lin_ds and lin_evol';
182  error_sim_data = lin_evol_sim_data;
183  error_sim_data.U = error_sim_data.U-ds_sim_data.X;
184  plot_sim_data(lin_evol_model,lin_evol_model_data,error_sim_data,params);
185 
186  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
187  %step = 4 % reduced Basis generation by single trajectories and
188  %comparison of detailed and reduced simulation
189  %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
190 
191  case 4 % generate some trajectories
192 
193  ds_model = fast_model(params);
194  ds_model.enable_error_estimator = 0;
195 % lin_evol_model = advection_fv_output_model(params);
196  %lin_evol_model.debug = 1;
197 % lin_evol_model_data = gen_model_data(lin_evol_model);
198 
199  % set some additional fields in model which will be copied into
200  % ds model depending on whether constants are known or not:
201 % estimate_constants = 0;
202 % if estimate_constants
203 % %%% if constant estimaton has to be performed:
204 % lin_evol_model.estimate_lin_ds_nmu = 4;
205 % lin_evol_model.estimate_lin_ds_nX = 10;
206 % lin_evol_model.estimate_lin_ds_nt = 4;
207 % else
208 % %%% otherwise:
209 % % the following values are rough bounds, so could be
210 % % non-rigorous, then larger search must be performed.
211 % lin_evol_model.state_bound_constant_C1 = 1;
212 % lin_evol_model.output_bound_constant_C2 = 1.2916;
213 % end;
214 %
215 % ds_model = lin_ds_from_lin_evol_model(lin_evol_model);
216  %ds_model.u_function_ptr = @(params) 0.1; % define new
217  ds_model_data = gen_model_data(ds_model);
218 % mu = [1,1,0.5];
219 % mu = [0.5,0.5,0.5];
220 
221  %mus = {[0,0,0],[1,1,1],[1,1,0.5],[0,0.5,1]};
222  mus = {[0,0],[1,1],[1,0.5],[0.5,1]};
223  for m = 1:length(mus)
224  disp(['computing trajectory number ',num2str(m)]);
225  mu = mus{m};
226  ds_model = ds_model.set_mu(ds_model,mu);
227  %ds_model.theta = 0;
228  %ds_model.theta = 1;
229  %disp('set theta to 1!!!');
230  tic;
231  ds_sim_data = detailed_simulation(ds_model,ds_model_data);
232  t_det = toc;
233  disp(['time for full simulation: ',num2str(t_det)]);
234  save(fullfile(rbmatlabresult,['trajectory',num2str(m)]),...
235  'ds_model','ds_model_data','ds_sim_data');
236  end;
237 
238  case 4.5 % compute reduced basis.
239  % make sure that the previous trajectories are computed!!
240 
241  disp('PCA of snapshots takes a while...');
242  load(fullfile(rbmatlabresult,'trajectory1'),...
243  'ds_model','ds_model_data','ds_sim_data');
244  % smaller epsilon does result in more than 1 vector for mu=0,0,0
245  epsilon = 1e-4;
246  RB = PCA_fixspace(ds_sim_data.X,[],ds_model_data.G,[],[],epsilon);
247  %RB = orthonormalize_qr(ds_sim_data.X,ds_model_data.G);
248  % select correct orthogonal set
249  K = RB' * ds_model_data.G * RB;
250  i = find(abs(diag(K)-1)>1e-2);
251  i = min(i);
252  disp(['found ',num2str(size(RB,2)),' basis vectors after traj. 1']);
253  if isempty(i)
254  i = size(RB,2)+1;
255  end;
256 % keyboard;
257  RB = RB(:,1:(i-1));
258  disp(['reduced to ',num2str(size(RB,2)),' basis vectors.']);
259 % keyboard;
260 
261  load(fullfile(rbmatlabresult,'trajectory2'),...
262  'ds_model','ds_model_data','ds_sim_data');
263  RB2 = PCA_fixspace(ds_sim_data.X,RB,ds_model_data.G,[],[],epsilon);
264  %RB = orthonormalize_qr(ds_sim_data.X,ds_model_data.G);
265  % select correct orthogonal set
266  K = RB2' * ds_model_data.G * RB2;
267  i = find(abs(diag(K)-1)>1e-2);
268  i = min(i);
269  if isempty(i)
270  i = size(RB2,2)+1;
271  end;
272  disp(['found ',num2str(size(RB2,2)),' basis vectors after traj. 2']);
273 % keyboard;
274  RB2 = RB2(:,1:(i-1));
275  disp(['reduced to ',num2str(size(RB2,2)),' basis vectors after traj. 2']);
276  RB = [RB, RB2];
277 % RB = orthonormalize(RB);
278 
279  V = RB;
280  W = ds_model_data.G * V;
281  rbfname = fullfile(rbmatlabresult,'advection_fv_output_basis');
282  save(rbfname,'RB');
283 
284  ds_model.RB_basis_filename = rbfname;
285  ds_model.RB_generation_mode = 'file_load';
286 
287  detailed_data = gen_detailed_data(ds_model,ds_model_data);
288 
289  model = ds_model;
290  save(fullfile(rbmatlabresult,'advection_fv_output_detailed_data_2traj'),...
291  'model','detailed_data')
292 
293  case 4.75 % compare reduced and detailed simulation
294 
295  load(fullfile(rbmatlabresult,'advection_fv_output_detailed_data_2traj'),...
296  'model','detailed_data')
297  load(fullfile(rbmatlabresult,'trajectory2'),...
298  'ds_sim_data');
299 
300  ds_model = model;
301 
302  %perform single reduced simulation for given parameter
303  reduced_data = gen_reduced_data(ds_model,detailed_data);
304  ds_model.N = reduced_data.N;
305  tic,
306 % ds_model.nt = ds_model.nt/4;
307 % disp('set nt to quarter!!');
308  rb_sim_data = rb_simulation(ds_model,reduced_data);
309  t_red = toc;
310  rb_sim_data = rb_reconstruction(ds_model,detailed_data,rb_sim_data)
311  disp(['time for reduced simulation: ',num2str(t_red)]);
312 
313  % plot of reduced trajectory
314  lin_evol_model = ds_model.base_model;
315  lin_evol_model_data = gen_model_data(lin_evol_model);
316 
317  params.no_lines = 1;
318  params.show_colorbar = 1;
319  params.axis_equal = 1;
320  params.plot = lin_evol_model.plot;
321  % plot single snapshot
322  %lin_evol_model.plot(rb_sim_data.X(:,1),lin_evol_model_data.grid, ...
323  % params);
324 
325  % plot sequence:
326  params.title = 'reduced solution';
327  plot_sequence(rb_sim_data.X,lin_evol_model_data.grid,params)
328 
329  % plot sequence:
330  params.title = 'error';
331  plot_sequence((rb_sim_data.X-ds_sim_data.X),...
332  lin_evol_model_data.grid,params)
333 
334  % plot of reduced output
335  figure,plot(rb_sim_data.time,[rb_sim_data.Y;ds_sim_data.Y]);
336  legend('reduced output','detailed output');
337 
338  keyboard;
339 
340 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
341  case 5 % generate plots of detailed data for paper:
342  % step = 5 % figures of detailed simuations based on step 4
343 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
344 
345  disp('generate plots of detailed_data for morepas poster');
346 
347  load(fullfile(rbmatlabresult,'advection_fv_output_detailed_data_2traj'),...
348  'model','detailed_data')
349 
350 % model = fast_model(params);
351 % rbfname = 'advection_fv_output_basis';
352 % load(fullfile(rbmatlabresult,rbfname),'ds_model','ds_model_data');
353 
354  % make sure, that these trajectories exist.
355  trajectories = [3,4];
356 
357  for tr = 1:length(trajectories);
358  trname = ['trajectory',num2str(trajectories(tr))];
359  % plots of trajectory
360  load(fullfile(rbmatlabresult,trname));
361  lin_evol_model = ds_model.base_model;
362  lin_evol_model_data = gen_model_data(lin_evol_model);
363 
364  params.no_lines = 1;
365  params.show_colorbar = 0;
366  params.axis_equal = 1;
367  params.plot = lin_evol_model.plot;
368  lin_evol_model.plot(lin_evol_model_data.grid,ds_sim_data.X(:,1), ...
369  params);
370  set(gca,'Xtick',[])
371  set(gca,'Ytick',[])
372  saveas(gcf,fullfile(rbmatlabresult,[trname,'_t1.eps']),'epsc');
373  close(gcf);
374 
375  lin_evol_model.plot(lin_evol_model_data.grid,ds_sim_data.X(:,64), ...
376  params);
377  set(gca,'Xtick',[])
378  set(gca,'Ytick',[])
379  saveas(gcf,fullfile(rbmatlabresult,[trname,'_tmiddle.eps']),'epsc');
380  close(gcf);
381 
382  lin_evol_model.plot(lin_evol_model_data.grid,ds_sim_data.X(:,end), ...
383  params);
384  set(gca,'Xtick',[])
385  set(gca,'Ytick',[])
386  saveas(gcf,fullfile(rbmatlabresult,[trname,'_tend.eps']),'epsc');
387  close(gcf);
388 
389  p = ds_model.plot_sim_data(ds_model, ds_model_data, ds_sim_data,params);
390  % Ylim = get(gca,'Ylim');
391  % Ylim = [Ylim(1),Ylim(2)*1.1];
392  set(gca,'Ylim',[0,0.09]);
393  set(p(2),'Color',[0,0.8,0]);
394  saveas(gcf,fullfile(rbmatlabresult,[trname,'_output.eps']),'epsc');
395  close(gcf);
396  close(gcf);
397  end;
398 
399 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
400 % step = 6 % figures of reduced simuations based on step 4
401 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
402  case 6 % generate plots of reduced model for poster/paper
403 
404  disp('generate plots of reduced simulations output estimation');
405 
406  load(fullfile(rbmatlabresult,'advection_fv_output_detailed_data_2traj'),...
407  'model','detailed_data')
408  ds_model = model;
409 
410  % as reference
411 % load(fullfile(rbmatlabresult,'trajectory1'));
412  load(fullfile(rbmatlabresult,'trajectory2'));
413  ds_model.enable_error_estimator = 1;
414 % ds_model.error_estimation = 1;
415 % ds_model.RB_generation_mode = 'file_load'
416 % ds_model.RB_basis_filename = fullfile(...
417 % rbmatlabresult,'advection_fv_output_basis');
418 % detailed_data = gen_detailed_data(ds_model,ds_model_data);
419  reduced_data = gen_reduced_data(ds_model,detailed_data);
420 
421  Ns = [150:-10:10];
422 
423  Ns = Ns(find(reduced_data.N>=Ns));
424 
425  for nind = 1:length(Ns)
426  % perform reduced simulation
427  ds_model.N = Ns(nind);
428  rb_sim_data = rb_simulation(ds_model,reduced_data);
429  rb_sim_data = rb_reconstruction(ds_model,detailed_data,rb_sim_data)
430 
431  params.no_lines = 1;
432  params.show_colorbar = 1;
433  params.axis_equal = 1;
434 
435  % plot of reduced output
436  figure;
437  p = plot(rb_sim_data.time,...
438  [rb_sim_data.Y;...
439  rb_sim_data.Y+rb_sim_data.DeltaY; ...
440  rb_sim_data.Y-rb_sim_data.DeltaY; ...
441  ds_sim_data.Y]);
442  colors = {[1,0,0],[1,0,0],[1,0,0],[0,1,0]};
443  linestyles = {'-','-.','-.','-'};
444  widths = [2,2,2,2];
445  for i = 1:length(p)
446  set(p(i),'Color',colors{i},'Linestyle',linestyles{i},...
447  'Linewidth',widths(i));
448  end;
449  % set marker on rb_sim_data:
450  Xdata = get(p(1),'XData');
451  Ydata = get(p(1),'YData');
452  Zdata = get(p(1),'ZData');
453  set(p(1),'XData',Xdata(1:4:end));
454  set(p(1),'YData',Ydata(1:4:end));
455  set(p(1),'ZData',Zdata(1:4:end));
456  set(p(1),'Marker','o');
457 
458  legend(['y\^(t): reduced output'],...
459  'y\^(t)+\Delta_y(t)',...
460  'y\^(t)-\Delta_y(t)',...
461  ['y(t): detailed output'],'Location','NorthWest');
462  title(['k = ',num2str(ds_model.N)],'Fontsize',15);
463  saveas(gcf,...
464  fullfile(rbmatlabresult,...
465  ['output_estimation',num2str(ds_model.N),'.eps']),'epsc');
466  close(gcf);
467  end;
468 
469 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
470 % step = 7 % time measurements
471 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
472  case 7
473  disp('time measurements for morepas poster');
474 
475 % clear all;
476  nreps = 10;
477  t_detailed = zeros(1,nreps);
478 
479  % initialize detailed model and simulate
480 
481  load(fullfile(rbmatlabresult,'advection_fv_output_detailed_data_2traj'),...
482  'model','detailed_data')
483  ds_model = model;
484 
485 % lin_evol_model = advection_fv_output_model(params);
486  %lin_evol_model.debug = 1;
487  %lin_evol_model_data = gen_model_data(lin_evol_model);
488 
489  % set some additional fields in model which will be copied into
490  % ds model depending on whether constants are known or not:
491 % estimate_constants = 0;
492 % if estimate_constants
493 % %%% if constant estimaton has to be performed:
494 % lin_evol_model.estimate_lin_ds_nmu = 4;
495 % lin_evol_model.estimate_lin_ds_nX = 10;
496 % lin_evol_model.estimate_lin_ds_nt = 4;
497 % else
498 % %%% otherwise:
499 % % the following values are rough bounds, so could be
500 % % non-rigorous, then larger search must be performed.
501 % lin_evol_model.state_bound_constant_C1 = 1;
502 % lin_evol_model.output_bound_constant_C2 = 1.2916;
503 % end;
504 %
505 % ds_model = lin_ds_from_lin_evol_model(lin_evol_model);
506 % % set accelerated data functions from below
507 % ds_model.A_function_ptr = @(model,model_data) ...
508 % eval_affine_decomp(@A_components,...
509 % @A_coefficients,...
510 % model,model_data);
511 %
512 % ds_model.B_function_ptr = @(model,model_data) ...
513 % eval_affine_decomp(@B_components,...
514 % @B_coefficients,...
515 % model,model_data);
516 % %ds_model.u_function_ptr = @(params) 0.1; % define new
517  ds_model_data = gen_model_data(ds_model);
518 
519  disp('skipping detailed simulations...')
520 % if 0
521 
522  for r = 1:nreps
523  disp(['repetition = ',num2str(r)])
524  %mu = [1,1,0.5]-rand(1,3)*0.01;
525  mu = [1,0.5]-rand(1,2)*0.01;
526  ds_model = ds_model.set_mu(ds_model,mu);
527  tic,
528  ds_sim_data = detailed_simulation(ds_model,ds_model_data);
529  t_detailed(1,r) = toc;
530  end;
531 
532 % end;
533 
534  disp('--------------------------------------------------')
535  disp('detailed simulation')
536  disp(t_detailed);
537  disp(['mean t_detailed=',num2str(mean(t_detailed))]);
538  disp('--------------------------------------------------')
539 
540  % reduced simulation with and without error estimation
541 
542  clear('ds_sim_data');
543  ds_model.RB_generation_mode = 'file_load';
544  ds_model.RB_basis_filename = ...
545  fullfile(rbmatlabresult,'advection_fv_output_basis');
546  detailed_data = gen_detailed_data(ds_model,ds_model_data);
547  reduced_data = gen_reduced_data(ds_model,detailed_data);
548 
549  Ns = [150:-10:10]; %,30,20,10];
550  Ns = Ns(find(reduced_data.N>=Ns));
551 
552  t_reduced_with_err_est = zeros(length(Ns),nreps);
553  t_reduced_without_err_est = zeros(length(Ns),nreps);
554 
555  % reduced simulation with and without error estimation
556 
557  for err_est = 0:1;
558  ds_model.error_estimation = err_est;
559  disp(['err_est= ',num2str(err_est)]);
560  for n_ind = 1:length(Ns)
561  disp(['M = ',num2str(Ns(n_ind))]);
562  ds_model.N = Ns(n_ind);
563  for r = 1:nreps
564  disp(['repetition = ',num2str(r)])
565  %mu = [1,1,1]-rand(1,3)*0.01;
566  mu = [1,1]-rand(1,2)*0.01;
567  ds_model = ds_model.set_mu(ds_model,mu);
568  tic;
569  rb_sim_data = rb_simulation(ds_model,reduced_data);
570  if size(rb_sim_data.Xr,1)~=Ns(n_ind)
571  disp('N dimension not set correctly!');
572  keyboard;
573  end;
574  if err_est == 1
575  t_reduced_with_err_est(n_ind,r) = toc;
576  else
577  t_reduced_without_err_est(n_ind,r) = toc;
578  end;
579  end;
580  end;
581  end;
582 
583  % output measurements
584 
585  disp('--------------------------------------------------')
586  disp('detailed simulation')
587  disp(t_detailed);
588  disp(['mean t_detailed=',num2str(mean(t_detailed))]);
589 
590  disp('--------------------------------------------------')
591  disp('reduced simulation with error estimation')
592  disp(t_reduced_with_err_est);
593  disp(['mean t_reduced_with_err_est=',...
594  num2str(mean(t_reduced_with_err_est,2)')]);
595  disp(['Ns =',num2str(Ns)]);
596 
597  disp('--------------------------------------------------')
598  disp('reduced simulation without error estimation')
599  disp(t_reduced_without_err_est);
600  disp(['mean t_reduced_without_err_est=',...
601  num2str(mean(t_reduced_without_err_est,2)')]);
602  disp(['Ns =',num2str(Ns)]);
603  disp('--------------------------------------------------')
604 
605  save(fullfile(rbmatlabresult,'time_measurements'));
606 
607 % dim_x = 65536, nt = 2048
608 %
609 %mean t_detailed=64.3806
610 %--------------------------------------------------
611 %reduced simulation with error estimation
612 % 3.4542 3.5969 3.5659 3.4216 3.6170
613 % 3.3648 3.3816 3.4218 3.4391 3.4552
614 % 3.1643 3.1312 3.1244 3.1433 3.1720
615 % 3.0594 3.0496 3.0628 3.0486 3.0387
616 % 2.9993 2.9953 3.0078 2.9413 2.9603
617 % 2.9555 2.9093 2.9091 2.9510 2.9256%%%
618 %
619 %mean t_reduced_with_err_est=3.5311 3.4125 3.147 3.0518 2.9808 2.9301
620 %Ns =60 50 40 30 20 10
621 %--------------------------------------------------
622 %reduced simulation without error estimation
623 % 2.3975 2.4220 2.2924 2.4490 2.3302
624 % 2.3322 2.3198 2.3425 2.3056 2.3512
625 % 2.2227 2.2734 2.2210 2.2582 2.2504
626 % 2.2531 2.2180 2.2429 2.2212 2.2612
627 % 2.2283 2.1979 2.2581 2.2241 2.1957
628 % 2.2274 2.2014 2.2577 2.1856 2.2268
629 %
630 %mean t_reduced_without_err_est=2.3782 2.3303 2.2451 2.2393 2.2208 2.2198
631 %Ns =60 50 40 30 20 10
632 %--------------------------------------------------
633 
634  case 8 % compare original and accelerated detailed model
635 
636  disp('compare original and accelerated lin_ds models');
637 
638  % clear all;
639  % initialize detailed model and simulate
640  lin_evol_model = advection_fv_output_model(params);
641  estimate_constants = 0;
642  if estimate_constants
643  %%% if constant estimaton has to be performed:
644  lin_evol_model.estimate_lin_ds_nmu = 4;
645  lin_evol_model.estimate_lin_ds_nX = 10;
646  lin_evol_model.estimate_lin_ds_nt = 4;
647  else
648  %%% otherwise:
649  % the following values are rough bounds, so could be
650  % non-rigorous, then larger search must be performed.
651  lin_evol_model.state_bound_constant_C1 = 1;
652  lin_evol_model.output_bound_constant_C2 = 1.2916;
653  end;
654 
655  ds_model = lin_ds_from_lin_evol_model(lin_evol_model);
656  ds_model.cone_weight = lin_evol_model.cone_weight;
657  %ds_model.u_function_ptr = @(params) 0.1; % define new
658  ds_model_data = gen_model_data(ds_model);
659 
660  % set accelerated data functions
661  ds_model_fast = lin_ds_from_lin_evol_model(lin_evol_model);
662  ds_model_fast.cone_weight = lin_evol_model.cone_weight;
663  ds_model_fast.A_function_ptr = @(model,model_data) ...
664  eval_affine_decomp(@A_components,...
665  @A_coefficients,...
666  model,model_data);
667 
668  ds_model_fast.B_function_ptr = @(model,model_data) ...
669  eval_affine_decomp(@B_components,...
670  @B_coefficients,...
671  model,model_data);
672 
673  % ds_model_fast.C_function_ptr = @(model,model_data) ...
674  % eval_affine_decomp(@C_components,...
675  % @C_coefficients,...
676  % model,model_data);
677 
678  % disp('skipping detailed simulations...')
679  % if 0
680 
681  %mu = [1,1,0.5]-rand(1,3)*0.01;
682  mu = [1,0.5]-rand(1,2)*0.01;
683  ds_model = ds_model.set_mu(ds_model,mu);
684  ds_model_fast = ds_model_fast.set_mu(ds_model_fast,mu);
685  tic,
686  ds_sim_data = detailed_simulation(ds_model,ds_model_data);
687  t_detailed = toc;
688  tic,
689  ds_sim_data_fast = detailed_simulation(ds_model_fast,ds_model_data);
690  t_detailed_fast = toc;
691 
692  disp(['original t = ',num2str(t_detailed),...
693  ', accelerated t = ',num2str(t_detailed_fast)]);
694  if ~isequal(ds_sim_data,ds_sim_data_fast)
695  disp('simulation results of original and accelerated model differ!')
696  else
697  disp('simulation results of original and accelerated model equal!')
698  end;
699 
700  keyboard;
701 
702  case 9 % generate basis from POD-Greedy
703 
704  filecache_clear;
705  ds_model = fast_model(params);
706  ds_model.cone_range = [0.4,0.6];
707  ds_model.base_model.cone_range = ds_model.cone_range;
708  % ds_model.
709  % shrink range to 'single point'
710 % r = [0.4,0.6]; % => 2 -> 0.8 error nach 10 basisvektoren
711 % r = [0.0,1.0]; % => 4.6 -> 1.9 error nach 10 basisvektoren
712 % r = [0.45,0.55]; % => 2.38 -> 0.51 error nach 10 basisvektoren
713  %ds_model.mu_ranges = {[0,1],[0.4,0.6],[0.4,0.6]};
714  ds_model.mu_ranges = {[0.4,0.6],[0.4,0.6]};
715  % ds_model.rb_basis_generation_ptr = ...
716  % @basisgen_POD_greedy_uniform_fixed;
717 
718  % ds_model.RB_generation_mode = 'uniform_fixed';
719  % ds_model.RB_num_intervals = [2,2,2];
720 
721  ds_model_data = gen_model_data(ds_model);
722 
723  disp('generating basis from POD-Greedy');
724 
725  % for basis generation, set further fields:
726 
727  ds_model.verbose = 5; % => only dots in greedy loop
728  ds_model.RB_generation_mode = 'greedy_uniform_fixed';
729  ds_model.RB_numintervals = [2,1,1]; % => 0.5 in cone_pos obtained
730  ds_model.RB_numintervals = [1,1];
731 % ds_model.RB_numintervals = [4,4,4];
732  ds_model.RB_stop_epsilon = 1e-5;
733 % ds_model.RB_stop_timeout = 3600*6;
734 % ds_model.RB_stop_Nmax =150;
735  ds_model.RB_stop_timeout = 600;
736  ds_model.RB_stop_Nmax =10;
737  ds_model.RB_extension_algorithm = @RB_extension_PCA_fixspace;
738  ds_model.RB_error_indicator = 'estimator';
739  ds_model.enable_error_estimator = 1;
740  % default:
741  % ds_model.get_estimator_from_sim_data = @(sim_data) sim_data.DeltaX(end);
742 
743  detailed_data = gen_detailed_data(ds_model,ds_model_data);
744 
745  model = ds_model;
746  dfname = fullfile(rbmatlabresult,'advection_fv_output_lin_ds_detailed');
747  save(dfname,'model','detailed_data');
748 
749  disp('successfully generated and saved detailed_data');
750 
751  params.plot_title = 'reduced basis';
752  params.plot = model.base_model.plot;
753  lin_evol_model_data = gen_model_data(model.base_model);
754  plot_sequence(detailed_data.V,lin_evol_model_data.grid,params)
755 
756  keyboard;
757 
758  case 10 % reduced simulations and error plots
759 
760  %basis_from_step = 9; % step 9
761  basis_from_step = 4; % step 4
762 
763  if basis_from_step == 9
764  dfname = fullfile(rbmatlabresult,'advection_fv_output_lin_ds_detailed');
765  load(dfname);
766  disp('loaded basis from step 9')
767  else
768  load(fullfile(rbmatlabresult,'advection_fv_output_detailed_data_2traj'),...
769  'model','detailed_data')
770 
771 % model = fast_model(params);
772  model.enable_error_estimator = 1;
773  model_data = gen_model_data(model);
774  % dfname = fullfile(rbmatlabresult,'advection_fv_output_detailed_data_146');
775 % model.RB_basis_filename = ...
776 % fullfile(rbmatlabresult,'advection_fv_output_basis');
777 % model.RB_generation_mode = 'file_load';
778 % detailed_data = gen_detailed_data(model,model_data);
779  % dfname = fullfile(rbmatlabresult,'advection_fv_output_detailed');
780  % load(dfname);
781  tmp = load(fullfile(rbmatlabresult,'trajectory2'));
782  ds_sim_data = tmp.ds_sim_data;
783  mu_from_file = get_mu(tmp.ds_model);
784  disp('loaded basis from step 4')
785  model = model.set_mu(model,mu_from_file);
786  end;
787  % trajectory2 above is made with this
788  % model = model.set_mu(model,[1,1,1]);
789  % model = model.set_mu(model,[1,1,1]);
790  % trajectory3 above is made with this:
791  %model = model.set_mu(model,[1,1,0.5]);
792  %model = model.set_mu(model,[0.5,0.6,0.6]);
793  %model = model.set_mu(model,[0.5,0.4,0.4]);
794  %model = model.set_mu(model,[0.5,0.5,0.5]);
795  model.N = size(detailed_data.V,2);
796  %model.N = 150;
797  reduced_data = gen_reduced_data(model,detailed_data);
798  rb_sim_data = rb_simulation(model,reduced_data);
799 
800  params.no_lines = 1;
801  params.show_colorbar = 1;
802  params.axis_equal = 1;
803 
804  % plot of reduced output
805  figure;
806  p = plot(rb_sim_data.time,...
807  [rb_sim_data.Y;...
808  rb_sim_data.Y+rb_sim_data.DeltaY; ...
809  rb_sim_data.Y-rb_sim_data.DeltaY]);
810  colors = {[1,0,0],[1,0,0],[1,0,0],[0,1,0]};
811  linestyles = {'-','-.','-.','-'};
812  widths = [2,2,2,2];
813  for i = 1:length(p)
814  set(p(i),'Color',colors{i},'Linestyle',linestyles{i},...
815  'Linewidth',widths(i));
816  end;
817  % set marker on rb_sim_data:
818  Xdata = get(p(1),'XData');
819  Ydata = get(p(1),'YData');
820  Zdata = get(p(1),'ZData');
821  set(p(1),'XData',Xdata(1:4:end));
822  set(p(1),'YData',Ydata(1:4:end));
823  set(p(1),'ZData',Zdata(1:4:end));
824  set(p(1),'Marker','o');
825 
826  legend(['y\^(t): reduced output'],...
827  'y\^(t)+\Delta_y(t)',...
828  'y\^(t)-\Delta_y(t)','Location','NorthWest');
829  title(['k = ',num2str(model.N)],'Fontsize',15);
830 
831  % table of error, error estimator and effectivity
832 
833  reduced_data = gen_reduced_data(model,detailed_data);
834 
835  Ns = 150:-10:10;
836  Ns = Ns(find(reduced_data.N>=Ns));
837 
838  err = zeros(1,length(Ns));
839  err_estim = zeros(1,length(Ns));
840  effectivity = zeros(1,length(Ns));
841 
842  for ni = 1:length(Ns)
843  model.N = Ns(ni);
844  rb_sim_data = rb_simulation(model, reduced_data);
845  err(ni) = abs(ds_sim_data.Y(end) - rb_sim_data.Y(end));
846  err_estim(ni) = rb_sim_data.DeltaY(end);
847  effectivity(ni) = err_estim(ni)/err(ni);
848  end;
849  Ns
850  err
851  err_estim
852  effectivity
853  disp('finished computing effectivities')
854  keyboard;
855 
856  case 11 % test of matlab ode solver for resulting sparse system
857 
858  params.coarse_factor = 8; % nake model smaller
859  model = fast_model(params);
860  model.RB_generation_mode = 'from_detailed_data';
861 
862  % detailed simulation:
863  model_data = gen_model_data(model);
864  model_data.RB = zeros(size(model_data.G,1),0);
865  detailed_data = gen_detailed_data(model,model_data);
866 
867  model.decomp_mode = 0; % complete
868  model = model.set_time(model,0)
869  %model = model.set_mu(model,[0.5,1,1]);
870  model = model.set_mu(model,[1,1]);
871 
872  A = @(model,model_data,t) ...
873  model.A_function_ptr(model.set_time(model,t),model_data);
874 
875  B = @(model,model_data,t) ...
876  model.B_function_ptr(model.set_time(model,t),model_data);
877 
878  x0 = model.x0_function_ptr(model,model_data);
879  tic;
880  [t,x] = ode15s(@(t,x) A(model,model_data,t)*x + B(model,model_data,t), ...
881  [0 1], x0);
882  time_ode15 = toc;
883  disp(['ode15 takes ',num2str(time_ode15),' sec']);
884 
885  params.plot_title = 'ode15 solution';
886  params.plot = model.base_model.plot;
887  params.axis_equal = 1;
888  lin_evol_model_data = gen_model_data(model.base_model);
889  plot_sequence((x'),lin_evol_model_data.grid,params)
890 
891  tic;
892  sim_data = detailed_simulation(model,model_data);
893  time_rbmatlab = toc;
894  disp(['rbmatlab takes ',num2str(time_rbmatlab),' sec']);
895  params.plot_title = 'rbmatlab solution';
896  params.plot = model.base_model.plot;
897  params.axis_equal = 1;
898  lin_evol_model_data = gen_model_data(model.base_model);
899  plot_sequence(sim_data.X,lin_evol_model_data.grid,params)
900 
901  % further plots for MCMCS paper
902  case 12
903 
904  load(fullfile(rbmatlabresult,'advection_fv_output_detailed_data_2traj'),...
905  'model','detailed_data')
906 
907  % plot of error estimator over parameter variation
908 % model = fast_model(params);
909 % model_data = gen_model_data(model);
910  % dfname = fullfile(rbmatlabresult,'advection_fv_output_detailed_data_146');
911 % model.RB_basis_filename = ...
912 % fullfile(rbmatlabresult,'advection_fv_output_basis');
913 % model.RB_generation_mode = 'file_load';
914 % detailed_data = gen_detailed_data(model,model_data);
915 % save(fullfile(rbmatlabresult,...
916 % 'advection_fv_output_detailed_2trajectories.mat'),...
917 % 'model','detailed_data')
918 % keyboard;
919  % dfname = fullfile(rbmatlabresult,'advection_fv_output_detailed');
920  % load(dfname);
921  tmp = load(fullfile(rbmatlabresult,'trajectory2'));
922  ds_sim_data = tmp.ds_sim_data;
923  disp('loaded basis from step 4')
924  model.enable_error_estimator = 1;
925  % factors = 0:0.05:1;
926  reduced_data = gen_reduced_data(model,detailed_data);
927  Ns = [1,8,16,32];
928  Ns = Ns(find(reduced_data.N>=Ns));
929  factors = 0:0.05:1;
930  err_estim = zeros(length(factors),length(Ns));
931  for ni = 1:length(Ns)
932  model.N = Ns(ni);
933  for i =1:length(factors)
934  %model = model.set_mu(model, factors(i)*[1,1,1]);
935  model = model.set_mu(model, factors(i)*[1,1]);
936  rb_sim_data = rb_simulation(model, reduced_data);
937  % err(ni) = abs(ds_sim_data.Y(end) - rb_sim_data.Y(end));
938  err_estim(i,ni) = rb_sim_data.DeltaY(end);
939  % effectivity(ni) = err_estim(ni)/err(ni);
940  end;
941  end;
942 % keyboard;
943  p = plot(factors,err_estim');
944  colors = {[1,0,0],[0,0.6,0],[0,0,1],[0,0,0]};
945  linestyles = {'-','-.','-','-.'};
946  for i = 1:min(length(Ns),4);
947  set(p(i),'Linewidth',2,'Color',colors{i},'Linestyle',linestyles{i});
948  end;
949 
950  xlabel('parameter \mu_1=\mu_2=\mu_3');
951  ylabel('\Delta_y(\mu)');
952  title('error estimator from parameter sweep')
953  legend({'k=1','k=8','k=16','k=32'},'Location','NorthWest');
954  saveas(gcf,fullfile(rbmatlabresult,'parameter_sweep.eps'),'epsc');
955  keyboard;
956  otherwise
957  disp('step number unknown');
958 end;
959 
960 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
961 % generation of ds_model from lin_evol_model and replace
962 % by fast vectors
963 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
964 
965 function ds_model = fast_model(params);
966  % fast_model
967 lin_evol_model = advection_fv_output_model(params);
968 %lin_evol_model.debug = 1;
969 %lin_evol_model_data = gen_model_data(lin_evol_model);
970 
971 % set some additional fields in model which will be copied into
972 % ds model depending on whether constants are known or not:
973 estimate_constants = 0;
974 if estimate_constants
975  %%% if constant estimaton has to be performed:
976  lin_evol_model.estimate_lin_ds_nmu = 4;
977  lin_evol_model.estimate_lin_ds_nX = 10;
978  lin_evol_model.estimate_lin_ds_nt = 4;
979 else
980  %%% otherwise:
981  % the following values are rough bounds, so could be
982  % non-rigorous, then larger search must be performed.
983  lin_evol_model.state_bound_constant_C1 = 1;
984  lin_evol_model.output_bound_constant_C2 = 1.2916;
985 end;
986 
987 ds_model = lin_ds_from_lin_evol_model(lin_evol_model);
988 ds_model.cone_weight = lin_evol_model.cone_weight;
989 % set accelerated data functions from below
990 ds_model.A_function_ptr = @(model,model_data) ...
991  eval_affine_decomp(@A_components,...
992  @A_coefficients,...
993  model,model_data);
994 
995 ds_model.B_function_ptr = @(model,model_data) ...
996  eval_affine_decomp(@B_components,...
997  @B_coefficients,...
998  model,model_data);
999 %ds_model.u_function_ptr = @(params) 0.1; % define new
1000 
1001 %ds_model.theta = 1;
1002 %disp('set theta to 1!!!');
1003 %keyboard;
1004 
1005 
1006 
1007 
1008 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1009 % auxialiary coefficient functions for acceleration of
1010 % advection model: explicit implementation of coefficients
1011 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
1012 
1013 function Acomp = A_components(model,model_data)
1014  % A_components
1015 Acomp = model_data.A_components;
1016 
1017 function Acoeff = A_coefficients(model);
1018  % A_coefficients
1019 
1020 %model.base_model.decomp_mode = 2;
1021 %model.base_model.t = model.t;
1022 %mu = get_mu(model);
1023 %model.base_model = model.set_mu(model.base_model,mu);
1024 %[L_I_coeff, L_E_coeff, b_coeff] = ...
1025 % model.base_model.operators_ptr(model.base_model, );
1026 %% L_E = Id + Delta t * A
1027 %Acoeff = L_E_coeff(2:end)/model.base_model.dt;
1028 %keyboard;
1029 %Acoeff = -model.base_model.dt * ...
1030 % model.base_model.velocity_coefficients_ptr(model.base_model);
1031 %if model.t>0
1032  Acoeff = - [model.vx_weight, model.vy_weight]*(1-model.t);
1033 % keyboard;
1034 %end;
1035 
1036 function Bcomp = B_components(model,model_data)
1037  % B_components
1038 Bcomp = model_data.B_components;
1039 
1040 function Bcoeff = B_coefficients(model);
1041  % B_coefficients
1042 
1043 % hmmm the coefficients are now called twice in A and B
1044 % should somehow be cached later...
1045 %model.base_model.decomp_mode = 2;
1046 %model.base_model.t = model.t;
1047 %mu = get_mu(model);
1048 %model.base_model = model.set_mu(model.base_model,mu);
1049 %[L_I_coeff, L_E_coeff, b_coeff] = ...
1050 % model.base_model.operators_ptr(model.base_model, );
1051 %% L_E = Id + Delta t * A
1052 %Bcoeff2 = b_coeff/model.base_model.dt;
1053 vcoeff = - [model.vx_weight, model.vy_weight]*(1-model.t);
1054 Q_0 = model.base_model.cone_number;
1055 dcoeff = zeros(1,Q_0);
1056 max_pos = model.cone_weight * (Q_0-1)+1;
1057 t = model.t;
1058 for q = 1:Q_0
1059  dcoeff(q) = (1-(max_pos-q))*(1-t) * ((max_pos>=q) && (max_pos < q+1)) ...
1060  + (1+(max_pos-q))*(1-t) * ((max_pos>=q-1) && (max_pos < q));
1061 end;
1062 v = -(vcoeff(:)*(dcoeff(:)'));
1063 Bcoeff = v(:);
1064 %if model.t>0
1065 %keyboard;
1066 %end;
1067 
1068 %function Ccomp = C_components(model,model_data)
1069 %Ccomp = model_data.C_components;
1070 %
1071 %function Ccoeff = C_coefficients(model);
1072 %model.base_model.decomp_mode = 2;
1073 %model.base_model.t = model.t;
1074 %mu = get_mu(model);
1075 %model.base_model = model.set_mu(model.base_model,mu);
1076 %Ccoeff = ...
1077 % model.base_model.operators_output(model.base_model);
1078 
1079 
function r = verbose(level, message, messageId)
This function displays messages depending on a message-id and/or a level. Aditionally you can set/res...
Definition: verbose.m:17
function [ RBext , dummy ] = RB_extension_PCA_fixspace(model, detailed_data)
function computing a RB basis extension for given parameters by the POD-Greedy algorithm.
function advection_fv_output(step)
small script implementing a simple advection example for producing matrices to be used in the RB-DS f...
function p = plot_sequence(varargin)
plotting a sequence of data slices on polygonal 2d grid (constructed from params if empty) and provid...
Definition: plot_sequence.m:17
function p = plot_sim_data(model, model_data, sim_data, plot_params)
function performing the plot of the simulation results as specified in model.
Definition: plot_sim_data.m:17
Customizable implementation of an abstract greedy algorithm.
Definition: DuneRBLeafNode.m:1