rbmatlab  1.16.09
 All Classes Namespaces Files Functions Variables Modules Pages
Output.m
1 classdef Output < handle
2  % result class for computations executed by an
4  %
5  % This class provides methods to transform and visualize the gathered data
6  % for publication in LaTeX/Tikz documents.
7 
8  properties
9 
10  % values computed by Postprocess.StochasticAssessment.Assessment
11  %
12  % Usually this is a cell array of structure arrays. For cells correspond to
13  % the @ref Assessment.rsamples "rsamples" and the structure array
14  % corresponds to the test samples given by @ref Assessment.M_test "M_test".
15  values;
16 
17  % cell array of all condition numbers of reduced simulaton system matrices.
18  rd_conds;
19 
20  rtime_est = {}; % time for the computation of the error estimator. By default this is not set.
21 
22  % the object of type Assessment with which the values are generated
23  %
24  % @todo Maybe we should cache this, too...
25  assessment;
26 
27  % maxmimum error values to be plotted. Above this values error display is
28  % cropped.
29  %
30  % @default = 1;
31  error_bound = 1;
32 
33  % maxmimum estimator values to be plotted. Above this values error display
34  % is cropped.
35  %
36  % @default = 100;
37  estimator_bound = 100;
38 
39  % define stability-region as error being smaller than sqrt(diffmax * area),
40  % e.g. diffmax = 4, area = 2e-7
41 
42  % stability region
43  %
44  % If the averaged error exceeds this value, the visualization is cropped
45  % here. This happens especially for small reduced basis sizes. (Default =
46  % '1e-2')
47  stab_limit = 1e-2;
48 
49  end
50 
51  methods
52 
53  function output = Output(values, assessment)
54  % function output = Output(values, assessment)
55  % Constructor
56  %
57  % Parameters:
58  % values: the structure with the gathered data by the
60  % assessment: the generating object of type Postprocess.StochasticAssessment.Assessment
61 
62  output.values = values;
63  output.assessment = assessment;
64 
65  end
66 
67  function plot_landscape(this)
68  %function plot_landscape(this)
69  % plots a two dimensional error landscape in case of two dimensional
70  % @ref Assessment.plot_fields "plot_fields"
71  %
72  % Used attributes:
73  % - #error_bound and
74  % - #stab_limit
75 
76  errs = this.get_slice('max_err');
77  i = logical(errs>this.error_bound);
78  errs(i) = this.error_bound;
79  i = logical(isnan(errs));
80  errs(i) = this.error_bound;
81  reshape(errs, size(this.values));
82  C = ones(size(errs));
83  i = logical(errs<this.stab_limit);
84  C(i) = 2;
85 
86  figure;
87  if length(this.assessment.samples) == 2
88  if isempty(find(C==1,1)); % i.e. all stable
89  surf(this.assessment.samples{1},this.assessment.tsamples{2},errs');
90  else % some not stable
91  surf(this.assessment.samples{1},this.assessment.tsamples{2},errs',C');
92  end;
93  shading interp;
94  % figure, pcolor(Ms, Ns,C);
95  set(gca,'Zscale','log');
96  xlabel(this.assessment.pf_descr{1});
97  ylabel(this.assessment.pf_descr{2});
98  zlabel(this.assessment.error_label);
99  else
100  plot(this.assessment.rsamples, errs');
101  xlabel(this.assessment.pf_descr{1});
102  ylabel(this.assessment.error_label);
103  set(gca,'Yscale','log');
104  end
105  run_name = strrep(this.assessment.run_name,'_',' ');
106  title([run_name,' L-infty([0,T],L2) error']);
107 
108  figure;
109  if length(this.assessment.samples) == 2
110  surf(this.assessment.samples{1}, this.assessment.samples{2}, real(this.get_slice('max_err_inds')));
111  xlabel(this.assessment.pf_descr{1});
112  ylabel(this.assessment.pf_descr{2});
113  zlabel('time index');
114  else
115  plot(this.assessment.samples, this.get_slice('max_err_inds'));
116  xlabel(this.assessment.pf_descr{1});
117  zlabel('time index');
118  end
119 
120  title('time index of maximum l2-error');
121  %set(gca,'Zscale','log');
122  end
123 
124  function print_table(this, datafile)
125  % function print_table(this, datafile)
126  % prints out a CSV table with gathered data.
127  %
128  % An example for a pgfplotstable graphic that can interpret the table output
129  % is shown in the following code snippet:
130  % @code
131  % \pgfplotstabletypeset[
132  % font={\footnotesize},
133  % columns={dim1,dim2,{time_avg},{err_max},{offtime_total}},
134  % columns/{dim1}/.style={
135  % column name={N},
136  % postproc cell content/.append code={%
137  % \ifnum\pgfplotstablerow=0
138  % \pgfkeyssetvalue{/pgfplots/table/@cell content}{$H={}$#2}%
139  % \fi
140  % }%
141  % },
142  % columns/{dim2}/.style={
143  % column name={M},
144  % postproc cell content/.append code={%
145  % \ifnum\pgfplotstablerow=0
146  % \pgfkeyssetvalue{/pgfplots/table/@cell content}{\ensuremath{-}}%
147  % \fi
148  % },
149  % },
150  % columns/{time_avg}/.style={
151  % column name={\o-run--time[s]}
152  % },
153  % columns/{err_max}/.style={
154  % column name={max. error}
155  % },
156  % columns/{offtime_total}/.style={
157  % column name={offline time[h]},
158  % divide by=3600
159  % },
160  % every even row/.style={
161  % before row={\rowcolor[gray]{0.9}}},
162  % every head row/.style={
163  % before row=\toprule,after row=\midrule},
164  % every last row/.style={
165  % after row=\bottomrule},
166  % display columns/2/.style={dec sep align},
167  % display columns/3/.style={zerofill,precision=2},
168  % empty cells with={\ensuremath{-}}
169  % ]{table.dat};
170  % @endcode
171  %
172  % Parameters:
173  % datafile: a string specifying the path and filename of the data table to be stored
174  %
175  % This function prints out a table with data for the average, min and max
176  % values of
177  % - reduced simulation time,
178  % - reduced simulation reconstruction time,
179  % - reduced simulation error estimator computation time,
180  % - maximum error estimator and the
181  % - maximum "true" error
182  % Furthermore, the offline time is printed out split up into
183  % - the total time spent in the offline phase,
184  % - the time spent for (pre-)computation of detailed simulations
185  % - the time spent in the EI-greedy algorithm
186  % - the time spent in the reduced basis generation algorithm
187  ratime = this.reduce_slice(@mean, 'rtime');
188  rmtime = this.reduce_slice(@min, 'rtime');
189  rMtime = this.reduce_slice(@max, 'rtime');
190  rratime = this.reduce_slice(@mean, 'rrtime');
191  rrmtime = this.reduce_slice(@min, 'rrtime');
192  rrMtime = this.reduce_slice(@max, 'rrtime');
193  restatime = this.reduce_slice(@mean, 'rtime_est');
194  restmtime = this.reduce_slice(@min, 'rtime_est');
195  restMtime = this.reduce_slice(@max, 'rtime_est');
196  etaa = this.reduce_slice(@mean, 'max_Delta');
197  etam = this.reduce_slice(@min, 'max_Delta');
198  etaM = this.reduce_slice(@max, 'max_Delta');
199  erra = this.reduce_slice(@mean, 'max_err');
200  errm = this.reduce_slice(@min, 'max_err');
201  errM = this.reduce_slice(@max, 'max_err');
202 
203  dimension = ones(length(ratime),2);
204  ratio = ones(length(ratime),1);
205  simulation = cell(length(ratime),1);
206  offt = ones(length(ratime), 4);
207  asm = this.assessment;
208  rbdd = get_by_description(asm.detailed_data, 'rb');
209 
210  for comb = 1:length(asm.rsamples)
211  pf_val = asm.rsamples(:,comb);
212  for i = 1:length(pf_val)
213  asm.update_rbmodel(asm.plot_fields{i}, asm.rsamples(i, comb));
214  end
215  N = asm.rbmodel.N;
216  M = asm.rbmodel.M;
217  simulation{comb} = 'Reduced';
218  dimension(comb,1) = N;
219  dimension(comb, 2) = M;
220  if length(pf_val) == 1
221  ratio(comb) = asm.rsamples(comb);
222  else
223  ratio(comb) = comb;
224  end
225  [offt(comb,1), offt(comb,2), offt(comb,3)] = offtime(asm.detailed_data, asm.rbmodel);
226  end
227  if offt(1,2) < 50
228  offt(:,1) = 90;
229  elseif offt(1,2) < 65
230  offt(:,1) = 117;
231  elseif offt(1,2) < 1300
232  offt(:,1) = 887;
233  else
234  offt(:,1) = 1394;
235  end
236  offt(:,4) = sum(offt, 2);
237  datime = this.reduce_slice(@mean, 'dtime', [], 1);
238  dmtime = this.reduce_slice(@min, 'dtime', [], 1);
239  dMtime = this.reduce_slice(@max, 'dtime', [], 1);
240  Hsize = asm.detailed_data.model_data.grid.nelements;
241  ddimension = Hsize;
242  dratio = -1;
243  print_datatable(datafile,...
244  {...
245  'ratio', 'dim1', 'dim2', ...
246  'time_avg', 'time_max', 'time_min', 'rec_time_avg', 'rec_time_min', 'rec_time_max', ...
247  'est_time_avg', 'est_time_min', 'est_time_max', 'eta_avg', 'eta_min', 'eta_max', ...
248  'err_avg', 'err_min', 'err_max', ...,
249  'offtime_dsims', 'offtime_ei', 'offtime_rb', 'offtime_total' ...
250  }, ...
251  [...
252  dratio, 0, ddimension, ...
253  datime, dmtime, dMtime, -1, -1 , -1, ...
254  -1, -1, -1, 0, 0, 0, ...
255  0, 0, 0, ...
256  0, 0, 0, 0; ...
257  ratio(:), dimension, ...
258  ratime(:), rmtime(:), rMtime(:), rratime(:), rrmtime(:), rrMtime(:), ...
259  restatime(:), restmtime(:), restMtime(:), etaa(:), etam(:), etaM(:), ...
260  erra(:), errm(:), errM(:), ...
261  offt...
262  ]');
263 
264  end
265 
266  function print_landscape(this, datafile, fieldname)
267  % function print_landscape(this, datafile, fieldname)
268  % equivalent to the plot_landscape() function put prints out a table
269  % which can be interpreted e.g. by pgfplots/Tikz.
270  %
271  % An example for a pgfplots graphics that can interpret the table output
272  % is shown in the following code snippet:
273  % @code
274  % \begin{axis}[zmode=log,log basis z=10,axis on top,
275  % view={-50}{-30}, axis lines=box,
276  % xlabel={$N$}, ylabel={$M$},
277  % zlabel={$\eta_{N,M,M'}$}]
278  % \addplot3[surf]
279  % table[x=N, y=M, z=error] {error_landscape.dat};
280  % @endcode
281  %
282  % Parameters:
283  % datafile: a string specifying the path and filename of the data table to be stored
284  % fieldname: the fieldname in the #values attribute that shall be
285  % plotted @default 'max_err'
286  if nargin <=2
287  fieldname = 'max_err';
288  end
289  errs = this.reduce_slice(@max, fieldname);
290  i = logical(errs>this.error_bound);
291  errs(i) = this.error_bound;
292  i = logical(isnan(errs));
293  errs(i) = this.error_bound;
294 
295  rsamples = this.assessment.rsamples;
296  blocklength = length(this.assessment.samples{1});
297 
298  print_datatable(datafile, [this.assessment.plot_fields, {'error'}], ...
299  [rsamples; errs(:)'], blocklength);
300 % C = ones(size(errs));
301 % i = logical(errs<this.stab_limit);
302 % C(i) = 2;
303  end
304 
305  function print_3d_curve(this, datafile, fieldname)
306  % function print_landscape(this, datafile, fieldname)
307  % equivalent to the print_landscape() function put prints out a table
308  % which can be interpreted for a 3d-curve plot e.g. by pgfplots/Tikz.
309  %
310  % This function can be used to display certain curves on the error
311  % landscape generated by print_landscape().
312  %
313  % An example for a pgfplots graphics that can interpret the table output
314  % is shown in the following code snippet:
315  % @code
316  % \begin{axis}[zmode=log,log basis z=10,axis on top,
317  % view={-50}{-30}, axis lines=box,
318  % xlabel={$N$}, ylabel={$M$},
319  % zlabel={$\eta_{N,M,M'}$}]
320  % \addplot3[mark=0]
321  % table[x=N, y=M, z=error] {error_curve.dat};
322  % @endcode
323  %
324  % Parameters:
325  % datafile: a string specifying the path and filename of the data table to be stored
326  % fieldname: the fieldname in the #values attribute that shall be
327  % plotted @default 'max_err'
328  if nargin <=2
329  fieldname = 'max_err';
330  end
331  errs = this.reduce_slice(@max, fieldname);
332  i = logical(errs>this.error_bound);
333  errs(i) = this.error_bound;
334  i = logical(isnan(errs));
335  errs(i) = this.error_bound;
336 
337  Ns = ones(length(errs), 1);
338  Ms = ones(length(errs), 1);
339  rsamples = this.assessment.rsamples;
340  for i=1:length(Ns)
341  c = rsamples(i);
342  this.assessment.update_rbmodel([], c);
343  Ns(i) = this.assessment.rbmodel.N;
344  Ms(i) = this.assessment.rbmodel.M;
345  end
346 
347  print_datatable(datafile, {'c', 'N', 'M', 'error'}, ...
348  [rsamples(:), Ns, Ms, errs(:)]');
349  end
350 
351  function slic = get_slice(this, fieldname, cell_range, struct_range)
352  % function slic = get_slice(this, fieldname, cell_range, out_range)
353  % helper function slicing a field in the #values attribute
354  %
355  % Parameters:
356  % fieldname: name of the field in the #values attributed to be sliced
357  % cell_range: cell indices to be extracted from the field.
358  % If empty all cell entries are returned
359  % @default []
360  % struct_range: Usually the #values attribute is a cell of structure arrays.
361  % arrays. This range defines a slice in the structure
362  % array. If empty all array entries are returned.
363  % @default []
364  %
365  % Return values:
366  % slic: the sliced #values attribute
367  if nargin <= 2 || isempty(cell_range)
368  cell_range = 1:length(this.values(:));
369  end
370  if nargin <= 3 || isempty(struct_range)
371  struct_range = 1:length(this.values{1}(:));
372  end
373  if isequal(fieldname, 'rtime_est')
374  slic = cellfun(@(x) x(struct_range), this.rtime_est(cell_range), 'UniformOutput', false);
375  else
376  slic = cellfun(@(x) [x(struct_range).(fieldname)], this.values(cell_range), 'UniformOutput', false);
377  end
378  slic = cell2mat(slic');
379  end
380 
381  function rslice = reduce_slice(this, red_func, fieldname, cell_range, struct_range, direction)
382  % function rslice = reduce_slice(this, red_func, fieldname, cell_range, struct_range, direction)
383  % extends the get_slice() function by also applying a reduction function
384  % on each matrix entry.
385  %
386  % Parameters:
387  % fieldname: name of the field in the #values attributed to be sliced
388  % red_func: function pointer who is applied to each matrix entry in
389  % the slice, e.g. '@@max' or '@@mean'
390  % cell_range: cell indices to be extracted from the field. If empty
391  % all cell entries are returned.
392  % @default []
393  % struct_range: Usually the #values attribute is a cell of structure arrays.
394  % arrays. This range defines a slice in the structure
395  % array. If empty all array entries are returned.
396  % @default []
397  % direction: integer specifying in which direction a matrix shall be
398  % reduced. In case it is equal to '2', the matrix is
399  % transposed before 'red_func' is applied.
400  % @default 1
401  %
402  % Return values:
403  % rslice: the sliced and reduced #values attribute
404  if nargin <=3
405  cell_range = [];
406  end
407  if nargin <= 4
408  struct_range = [];
409  end
410  if nargin <= 5
411  direction = 1;
412  end
413  slic = get_slice(this, fieldname, cell_range, struct_range);
414  if direction == 2
415  slic = slic';
416  rslice = zeros(1,size(slic,1));
417  else
418  rslice = zeros(size(slic,1),1);
419  end
420  for i = 1:size(slic,1)
421  rslice(i) = red_func(slic(i,:));
422  end
423  end
424 
425 
426  function merge(this, other)
427  % function merge(this, other)
428  % helper function merging this Output object with another one.
429  %
430  % The following rules apply during merging:
431  % -# If 'values.max_Delta' is all 'Inf', copy 'max_Delta' and
432  % 'max_Delta_ind' from 'other'.
433  % -# If 'values.max_err' is all 'Inf', copy 'max_err' and
434  % 'max_err_ind' from 'other'.
435  % -# If one of 'rconds', 'dconds' or 'rtime_est' is empty, copy it
436  % from 'other'
437  % -# If one of 'rtime', 'rrtime' or 'dtime' is set to all-'NaN', copy it
438  % from 'other'
439  % -# If 'this' and 'other' have valid 'rtime' (i.e. not-'Nan') and only
440  % one of both has computed estimators, compute 'rtime_est' by
441  % subtracting one from the other.
442  %
443  % Parameters:
444  % other: another object of type Postprocess.StochasticAssessment.Output to be
445  % merged into this one.
446  cpfields = {};
447  rtime_exists = false;
448  % check what we ar missing:
449  if all([this.values{end}(:).max_Delta] == Inf)
450  cpfields = [ cpfields, {'max_Delta', 'max_Delta_ind'}];
451  end
452  if all([this.values{end}(:).max_err] == Inf)
453  cpfields = [ cpfields, {'max_err', 'max_err_ind'}];
454  end
455  if isempty([this.values{end}(:).rconds])
456  cpfields = [ cpfields, {'rconds'}];
457  end
458  if isempty([this.values{end}(:).dconds])
459  cpfields = [ cpfields, {'dconds'}];
460  end
461  if all(isnan([this.values{end}(:).rtime]))
462  cpfields = [ cpfields, {'rtime'}];
463  else
464  rtime_exists = true;
465  end
466  if all(isnan([this.values{end}(:).rrtime]))
467  cpfields = [ cpfields, {'rrtime'}];
468  end
469  if all(isnan([this.values{end}(:).dtime]))
470  cpfields = [ cpfields, {'dtime'}];
471  end
472 
473  for i = 1:length(cpfields)
474  for j = 1:length(this.values(:))
475  [this.values{j}(:).(cpfields{i})] = other.values{j}(:).(cpfields{i});
476  end
477  end
478 
479  if rtime_exists && any(~isnan([other.values{end}(:).rtime]))
480  if this.assessment.compute_estimates ~= other.assessment.compute_estimates
481  this.rtime_est = cell(size(this.values));
482  for i = 1:length(this.values(:))
483  this.rtime_est{i} = abs([this.values{i}(:).rtime] - [other.values{i}(:).rtime]);
484  end
485  end
486  end
487 
488  if isempty(this.rtime_est)
489  this.rtime_est = other.rtime_est;
490  end
491 
492  if isempty(this.rd_conds)
493  this.rd_conds = other.rd_conds;
494  end
495  end
496  end
497 end
Tools for post-processing data, i.e. data extraction and visual enhancements for publication.
Definition: Assessment.m:1
result class for computations executed by an Postprocess.StochasticAssessment.Assessment object...
Definition: Output.m:19
a class used to compute reduced several reduced simulations over a huge parameter sample extracting u...
Definition: Assessment.m:19
Tools for gathering and storing data from a huge set of randomly generated reduced simulations...
Definition: Assessment.m:2
rsamples
the "real" collection of vectors used for the plot_fields variation
Definition: Assessment.m:175
plot_fields
is a vector of field names with at most 2 elements over which the error landscape is computed or a ce...
Definition: Assessment.m:36