KerMor  0.9
Model order reduction for nonlinear dynamical systems and nonlinear approximation
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Speed.m
Go to the documentation of this file.
1 namespace testing{
2 
3 
4 /* (Autoinserted by mtoc++)
5  * This source code has been filtered by the mtoc++ executable,
6  * which generates code that can be processed by the doxygen documentation tool.
7  *
8  * On the other hand, it can neither be interpreted by MATLAB, nor can it be compiled with a C++ compiler.
9  * Except for the comments, the function bodies of your M-file functions are untouched.
10  * Consequently, the FILTER_SOURCE_FILES doxygen switch (default in our Doxyfile.template) will produce
11  * attached source files that are highly readable by humans.
12  *
13  * Additionally, links in the doxygen generated documentation to the source code of functions and class members refer to
14  * the correct locations in the source code browser.
15  * However, the line numbers most likely do not correspond to the line numbers in the original MATLAB source files.
16  */
17 
18 class Speed {
35  public:
36 
37 
38  public: /* ( Static ) */
39 
40 
41  static function [colvec<double>times , doublee ] = KernelExpCustomBaseEval(kernels.KernelExpansion kexp,integer numpts) {
42  if nargin < 2
43  numpts = 1000;
44  if nargin < 1
45  kexp = kernels.KernelExpansion;
46  kexp.Centers.xi= rand(400,40);
47  kexp.Ma= rand(30,40);
48  o = general.Orthonormalizer;
49  A = o.orthonormalize(rand(40,40));
50  L = chol(A^t*diag(randi(10,1,40))*A);
51  kexp.Base= L;
52  end
53  end
54  if ~kexp.HasCustomBase
55  error(" No custom base set to compare to. ");
56  end
57  x = rand(size(kexp.Centers.xi,1),numpts);
58  dbase = kexp.toTranslateBase;
59  runs = 10;
60  t = zeros(2,runs);
61  e = zeros(2,runs);
62  for r = 1:runs
63  tic;
64  fx1 = kexp.evaluate(x);
65  t(1,r) = toc;
66  tic;
67  fx2 = dbase.evaluate(x);
68  t(2,r) = toc;
69  er = Norm.L2(fx1-fx2);
70  e(1,r) = max(er);
71  fxin = Norm.L2(fx1);
72  fxin(fxin == 0) = 1;
73  e(2,r) = max(er./fxin);
74  end
75  times = mean(t,2);
76  e = max(e,[],2);
77  }
101  static function pt = TryCatch(loopsize) {
102 
103  t = [];
104 
105  if nargin < 1
106  loopsize = 1000000;
107  end
108 
109  tic
110  for ii = 1:loopsize
111  A = ii;
112  end
113  t(end+1) = toc;
114 
115  tic
116  try
117  for ii = 1:loopsize
118  A = ii.^2;
119  end
120  catch ME
121  end
122  t(end+1) = toc;
123 
124  tic
125  for ii = 1:loopsize
126  try
127  A = ii.^2;
128  catch ME
129  end
130  end
131  t(end+1) = toc;
132 
133  tic
134  for ii = 1:loopsize
135  if ~ii
136  /* Never here.
137  * elseif 0
138  * % Never here.
139  * elseif mod(ii,2)
140  * B = ii; */
141  else
142  C = ii.^2;
143  end
144  end
145  t(end+1) = toc;
146 
147  pt = PrintTable;
148  pt.HasHeader= true;
149  pt.HasRowHeader= true;
150  pt.addRow(," Plain for "," try(for) "," for(try) "," for(if) ");
151  tc = num2cell(t);
152  pt.addRow(" Time ",tc[:],[" %6gs "]);
153  reltimes = num2cell(t./min(t));
154  pt.addRow(" Rel factor ",reltimes[:],[" %3.3g "]);
155  }
172  static function pt = BinaryvsMatSave(double fill,integer n,integer m) {
173 
174  if nargin < 3
175  m = 1000;
176  if nargin < 2
177  n = 2000;
178  if nargin < 1
179  fill = .5;
180  end
181  end
182  end
183 
184  a = KerMor.App;
185 
186  runs = 2;
187  pt = PrintTable;
188  pt.addRow(," Direct "," -v4 "," -v6 "," binary ");
189  pt.HasHeader= true;
190  pt.HasRowHeader= true;
191  for k = 1:runs
192  A = full(sprand(n,m,fill));
193 
194  f = fullfile(a.DataDirectory," tmp.mat ");
195  ti = tic;
196  save(f," A ");
197  t[1] = toc(ti);
198  d = dir(f);
199  s[1] = d.bytes;
200  delete(f);
201 
202  f = fullfile(a.DataDirectory," tmp2.mat ");
203  ti = tic;
204  save(f," A "," -v4 ");
205  t[2] = toc(ti);
206  d = dir(f);
207  s[2] = d.bytes;
208  delete(f);
209 
210  f = fullfile(a.DataDirectory," tmp3.mat ");
211  ti = tic;
212  save(f," A "," -v6 ");
213  t[3] = toc(ti);
214  d = dir(f);
215  s[3] = d.bytes;
216  delete(f);
217 
218  f = fullfile(a.DataDirectory," tmp4.mat ");
219  ti = tic;
220  fh = fopen(f," w+ ");
221  fwrite(fh,A," double ");
222  fclose(fh);
223  t[4] = toc(ti);
224  d = dir(f);
225  s[4] = d.bytes;
226  delete(f);
227  pt.addRow(" Time ",t[:],[" %6gs "]);
228  pt.addRow(" Size ",s[:],[@(v)sprintf(" %5gk ",v/1024)]);
229  end
230  }
256  static function FindVecInMatrix(integer n,integer m) {
257 
258  if nargin < 2
259  m = 12000;
260  if nargin < 1
261  n = 5000;
262  end
263  end
264  A = rand(n,m);
265  pos = [1 round(m/2) m];
266  for kpos = 1:3
267  b = A(:,pos(kpos));
268 
269  time = tic;
270  hlp = bsxfun(@eq,A,b);
271  i = find(sum(hlp,1) == n,1);
272  time1 = toc(time);
273 
274  time = tic;
275  i = find(sum(A == repmat(b,1,m)) == n,1);
276  time2 = toc(time);
277 
278  time = tic;
279  ind = strfind(reshape(A,1,[]),b^t);
280  round((ind+n-1)/n);
281  time3 = toc(time);
282 
283  fprintf(" Dims: %dx%d, match pos: %d, bsxfun time: %f, repmat time: %f, strfind time: %f\n ",...
284  n,m,pos(kpos), time1,time2,time3);
285  end
286  }
300  static function AffParamMatrix(integer n) {
301 
302  if nargin < 1
303  n = 5000;
304  end
305 
306  d = 50;
307 
308  str = [" 4*t "," 4*t + mu(1) "," exp(sum(mu)) "," sin(t) "," exp(-t)*mu(2) "," t*mu(2) "," cos(t) "," exp(-mu(2)) "];
309  theta[1] = @(t,mu)4*t;
310  theta[2] = @(t,mu)4*t + mu(1);
311  theta[3] = @(t,mu)exp(sum(mu));
312  theta[4] = @(t,mu)sin(t);
313  theta[5] = @(t,mu)exp(-t)*mu(2);
314  theta[6] = @(t,mu)t*mu(2);
315  theta[7] = @(t,mu)cos(t);
316  theta[8] = @(t,mu)exp(-mu(2));
317  nc = length(theta);
318 
319  am = general.AffParamMatrix;
320 
321  for i=1:nc
322  mat[i] = rand(d,d);/* #ok */
323 
324  am.addMatrix(str[i], mat[i]);
325  end
326 
327  v = rand(d,1);
328  params = rand(3,n);
329 
330  /* % Linear */
331  fprintf(" Testing linear case...\n ");
332  t = tic;
333  for i=1:n
334  coeff = zeros(1,5);
335  for k=1:nc
336  fun = theta[k];
337  coeff(k) = fun(params(1,i),params(2:3,i));
338  end
339  dummy = lincomb_sequence(mat,coeff)*v;/* #ok */
340 
341  end
342  t1 = toc(t);
343 
344  t = tic;
345  for i=1:n
346  dummy = am.compose(params(1,i),params(2:3,i))*v;/* #ok */
347 
348  end
349  t2 = toc(t);
350 
351  fprintf(" Direct: %gs, AffParamMatrix: %gs\n ",t1,t2);
352 
353  /* % Quadratic */
354  fprintf(" Testing quadratic case...\n ");
355  t = tic;
356  for i = 1:nc
357  for j = 1:nc
358  mat2[i,j] = mat[i]*mat[j];/* #ok */
359 
360  end
361  end
362  for i=1:n
363  coeff = zeros(1,5);
364  for k=1:nc
365  fun = theta[k];
366  coeff(k) = fun(params(1,i),params(2:3,i));
367  end
368  dummy = lincomb_sequence2(mat2,coeff,coeff)*v;/* #ok */
369 
370  end
371  t1 = toc(t);
372 
373  t = tic;
374  am = am*am;
375  for i=1:n
376  dummy = am.compose(params(1,i),params(2:3,i))*v;/* #ok
377  *dummy = am(params(1,i),params(2:3,i))*v;%#ok */
378 
379  end
380  t2 = toc(t);
381 
382  fprintf(" Direct: %gs, AffParamMatrix: %gs\n ",t1,t2);
383 
384  function res = lincomb_sequence(seq,sigma)
385  /* function res = lincomb_sequence(seq,sigma)
386  *
387  * function performing a linear combination of the elements in the
388  * cell array seq with coefficients in sigma result is a
389  * vector/matrix the same size as the entries of seq.
390  * if sigma = 0, the component is not touched, i.e. the component
391  * may also be an empty matrix. */
392 
393  /* Bernard Haasdonk 15.5.2007 */
394 
395  Q = length(sigma);
396  res = seq[1]* sigma(1);
397  for q=2:Q
398  if sigma(q)~=0
399  res = res + sigma(q)*seq[q];
400  end;
401  end;
402  end
403 
404  function res = lincomb_sequence2(seq,sigma1,sigma2)
405  /* function res = lincomb_sequence2(seq,sigma1,sigma2)
406  *
407  * function performing a linear combination of the elements in the
408  * 2d cell array seq with coefficients in sigma1 and sigma2 result is a
409  * vector/matrix the same size as the entries of seq.
410  * size of seq is length(sigma1) x length(sigma2)
411  * if some sigma = 0, the component is not touched, i.e. the component
412  * may also be an empty matrix. */
413 
414  /* Bernard Haasdonk 15.5.2007 */
415 
416  Q1 = length(sigma1);
417  Q2 = length(sigma2);
418  res = zeros(size(seq[1,1]));
419  for q1=1:Q1
420  if sigma1(q1)~=0
421  for q2=1:Q2
422  if sigma2(q2)~=0
423  res = res + sigma1(q1)*sigma2(q2)*seq[q1,q2];
424  end;
425  end;
426  end;
427  end;
428  end
429 
430  }
445  static function res = GaussMexSpeedTest1Arg(sx,sy,iter) {
446  if nargin < 3
447  iter = 50;
448  if nargin < 2
449  sy = 100;
450  if nargin < 1
451  sx = 5000;
452  end
453  end
454  end
455  k = kernels.GaussKernel(1);
456  x = rand(sx,sy);
457 
458  fprintf(" One argument speed test with sx=%d, sy=%d and %d iterations\n ",sx,sy,iter);
459  tmex = zeros(1,iter); tmex2 = zeros(1,iter);
460  tmexp = zeros(1,iter); tmat = zeros(1,iter);
461  pi = ProcessIndicator(" Iterating.. ",iter);
462  for i=1:iter
463  t = tic;
464  Kmex = k.dontuse_evaluate(x);
465  /* Kmex = k.evaluateIntel(x); */
466  tmex(i) = toc(t);
467 
468  t = tic;
469  Kmex2 = k.dontuse_evaluateDirect(x);
470  tmex2(i) = toc(t);
471 
472  t = tic;
473  KmexP = k.evaluateMex(x);
474  tmexp(i) = toc(t);
475 
476  t = tic;
477  K = k.evaluate(x,[]);
478  tmat(i) = toc(t);
479  pi.step;
480  end
481  pi.stop;
482  fprintf([" 1: %1.5fs - Mex straight\n2: %1.5fs - Mex time opt\n "...
483  " 3: %1.5fs - Mex time opt openmp\n4: %1.5fs - Matlab time\n "...
484  " Difference norms: 1-4=%1.5f, 2-4=%1.5f, 3-4=%1.5f\n "],...
485  mean(tmex2),mean(tmex),mean(tmexp),mean(tmat),...
486  norm(Kmex2-K),norm(Kmex-K),norm(KmexP-K));
487 
488  res = true;
489  }
504  static function res = GaussMexSpeedTest2Arg(sx,sy1,sy2,iter) {
505  if nargin < 4
506  iter = 40;
507  if nargin < 3
508  sy2 = 100;
509  if nargin < 2
510  sy1 = 100;
511  if nargin < 1
512  sx = 500;
513  end
514  end
515  end
516  end
517 
518  k = kernels.GaussKernel(1);
519  x = rand(sx,sy1);
520  y = rand(sx,sy2);
521 
522  fprintf(" Two argument speed test with sx=%d, sy1=%d, sy2=%d and %d iterations\n ",sx,sy1,sy2,iter);
523  tmex = zeros(1,iter); tmex2 = zeros(1,iter);
524  tmexp = zeros(1,iter); tmat = zeros(1,iter);
525  fprintf(" Iteration ");
526  for i=1:iter
527  fprintf(" %d ",i);
528 
529  t = tic;
530  Kmex = k.dontuse_evaluate(x,y);
531  /* Kmex = k.evaluateIntel(x,y); */
532  tmex(i) = toc(t);
533 
534  t = tic;
535  Kmex2 = k.dontuse_evaluateDirect(x,y);
536  tmex2(i) = toc(t);
537 
538  t = tic;
539  KmexP = k.evaluateMex(x,y);
540  tmexp(i) = toc(t);
541 
542  t = tic;
543  K = k.evaluate(x,y);
544  tmat(i) = toc(t);
545 
546  end
547  fprintf(" done!\n ");
548  fprintf([" 1: %1.5fs - Mex straight\n2: %1.5fs - Mex time opt\n "...
549  " 3: %1.5fs - Mex time opt openmp\n4: %1.5fs - Matlab time\n "...
550  " Difference norms: 1-4=%1.5f, 2-4=%1.5f, 3-4=%1.5f\n "],...
551  mean(tmex2),mean(tmex),mean(tmexp),mean(tmat),...
552  norm(Kmex2-K),norm(Kmex-K),norm(KmexP-K));
553 
554  res = true;
555  }
570 };
571 }
572 
logical HasHeader
Flag that determines if the first row should be used as table header.
Definition: PrintTable.m:266
Speed: Collects tests regarding speed of different methods and strategies.
Definition: Speed.m:18
static function [ colvec< double > times , double e ] = KernelExpCustomBaseEval(kernels.KernelExpansion kexp,integer numpts)
Tests the evaluation speed and determines the evaluation error of a kernel expansion and this expansi...
Definition: Speed.m:41
static function res = GaussMexSpeedTest1Arg(sx, sy, iter)
Tests the speed of the c implementations of evaluate for gaussians.
Definition: Speed.m:445
An integer value.
static function res = GaussMexSpeedTest2Arg(sx, sy1, sy2, iter)
Tests the speed of the c implementations of evaluate for gaussians.
Definition: Speed.m:504
static function pt = BinaryvsMatSave(double fill,integer n,integer m)
BinaryvsMatSaveSpeed: Tests the speeds and storage size required for storing double matrices...
Definition: Speed.m:172
KernelExpansion()
Default constructor.
static function FindVecInMatrix(integer n,integer m)
Created for test purposes of finding a vector in a matrix.
Definition: Speed.m:256
static function AffParamMatrix(integer n)
Tests the evaluation speed of the general.AffParamMatrix compared to loop-type evaluations of affine ...
Definition: Speed.m:300
PrintTable: Class that allows table-like output spaced by tabs for multiple rows. ...
Definition: PrintTable.m:17
function addRow(varargin)
Adds a row to the current table.
Definition: PrintTable.m:611
static function rowvec< double > n = L2(matrix< double > x)
Returns the discrete norm for each column vector in x.
Definition: Norm.m:39
Global configuration class for all KerMor run-time settings.
Definition: KerMor.m:17
static function KerMor theinstance = App()
The singleton KerMor instance.
Definition: KerMor.m:910
Norm: Static class for commonly used norms on sets of vectors.
Definition: Norm.m:17
static function pt = TryCatch(loopsize)
TryCatch: Demonstrate how slow try-catch blocks are.
Definition: Speed.m:101
ProcessIndicator: A simple class that indicates process either via waitbar or text output...
KernelExpansion: Base class for state-space kernel expansions.