Mercurial > repos > public > sbplib
changeset 167:15baeb35f74e feature/grids
Merge in changes from default.
| author | Jonatan Werpers <jonatan@werpers.com> |
|---|---|
| date | Tue, 23 Feb 2016 13:25:43 +0100 |
| parents | 7cb97c1988d9 (current diff) 344bde2f9d9b (diff) |
| children | ba1ae5b2c45e |
| files | |
| diffstat | 15 files changed, 236 insertions(+), 44 deletions(-) [+] |
line wrap: on
line diff
diff -r 7cb97c1988d9 -r 15baeb35f74e +anim/animate.m --- a/+anim/animate.m Tue Feb 23 10:55:46 2016 +0100 +++ b/+anim/animate.m Tue Feb 23 13:25:43 2016 +0100 @@ -5,28 +5,67 @@ % it will be called for increasnig t. %Todo: make it catch up and produce warnings if it lags behind? Instead of just requesting the next target time -function animate(F, tstart, tend, time_modifier , frame_rate) - if ~exist('time_modifier') - time_modifier = 1; - end + + +% If adapt is true time_modifier is treated as an upper bound +function animate(F, tstart, tend, time_modifier, target_frame_rate) + default_arg('time_modifier', 1); + default_arg('target_frame_rate',30); + + % t is simulation time + % tau is real time + + time_modifier_bound = time_modifier; + dTau_target = 1/target_frame_rate; % Real time between frames + + rs = util.ReplaceableString(); + rs.appendFormat(' tau: %d\n'); + rs.appendFormat(' target tau: %d\n'); + rs.appendFormat(' Target fps: %.2f\n'); + rs.appendFormat(' Actual fps: %.2f\n'); + rs.appendFormat('Target time_modifier: %d\n'); + rs.appendFormat('actual time_modifier: %d\n'); + + animation_start = tic(); + prevTau = 0; + targetTau = 0; + tauFrameStart = -dTau_target; + t = F(tstart); - if ~exist('frame_rate') - frame_rate = 30; + while t < tend + % Sleep until the frame should start + pause(targetTau-toc(animation_start)); + tau = toc(animation_start); + dTau = tau - tauFrameStart; + + % Calculate error in tau + e_Tau = tau - targetTau; + + % Regulate time_modifier based on e_Tau + % time_modifier = min(time_modifier_bound, max(0.5, abs(1-e_Tau/dTau)) * time_modifier); + + % Mark the start of the frame + tauFrameStart = tau; + + dt_target = dTau_target*time_modifier; % Targeted simulation time between frames + + t_prev = t; + t = F(t + dt_target); % Run simulation + + % Calculate when this frame should end and the next start. (this depends on what simulation time we ended up on) + dt = t-t_prev; + % targetTau = targetTau + dt/time_modifier; + targetTau = targetTau + dTau_target; + + % Update information about this frame + tau = toc(animation_start); + rs.updateParam(tau, targetTau, 1/dTau_target, 1/dTau, time_modifier_bound, time_modifier); end - frame_time = 1/frame_rate; - dt = frame_time*time_modifier; - animation_start = tic(); - t = F(tstart); - while t < tend - t = F(t + dt); - t_left = (t-tstart)/time_modifier-toc(animation_start); - pause(t_left) - end + % Final time reporting time_to_animate = toc(animation_start); expected_time = tend/time_modifier; - fprintf('\n'); fprintf('Time to animate: %.3f\n', time_to_animate) fprintf('Expected time : %.3f\n', expected_time)
diff -r 7cb97c1988d9 -r 15baeb35f74e +noname/animate.m --- a/+noname/animate.m Tue Feb 23 10:55:46 2016 +0100 +++ b/+noname/animate.m Tue Feb 23 13:25:43 2016 +0100 @@ -1,8 +1,8 @@ -% hand = animate(discretization, time_modifier, Tend, dirname, opt) +% hand = noname.animate(discretization, time_modifier, Tend, dirname, opt) % % Example: -% animate(discr,timemodifier,tend) -% animate(discr,1, [tstart tend],'my_mov', opt) +% noname.animate(discr,timemodifier,tend) +% noname.animate(discr,1, [tstart tend],'my_mov', opt) function hand = animate(discretization, time_modifier, Tend, dirname, opt) default_arg('time_modifier',1); @@ -67,8 +67,6 @@ if makemovies save_frame(); end - % pause(0.1) - str = util.replace_string(str,'t = %.5f',ts.t); if do_pause pause
diff -r 7cb97c1988d9 -r 15baeb35f74e +noname/calculateSolution.m --- a/+noname/calculateSolution.m Tue Feb 23 10:55:46 2016 +0100 +++ b/+noname/calculateSolution.m Tue Feb 23 13:25:43 2016 +0100 @@ -1,12 +1,13 @@ % Calculates the solution of discretization for a given set of ms ts and orders. % discrHand -- function handle to a Discretization constructor -% method -- time stepping method % m -- grid parameter % order -- order of accuracy of the approximtion % T -- time to calculate solution for +% tsOpt -- options for the time stepper creation. % input paramters m, t, order may all be vectors. -function [] = calculateSolution(filename, discrHand, method, m, T_in, order, force_flag) +function [] = calculateSolution(filename, name, discrHand, m, T_in, order, tsOpt, force_flag) default_arg('force_flag',false); + default_arg('tsOpt', []); if exist(filename,'file') && ~force_flag fprintf('File ''%s'' already exist.',filename); @@ -19,8 +20,6 @@ sf = SolutionFile(filename); - - orderWidth = findFieldWidth('%d',order); mWidth = findFieldWidth('%d',m); TWidth = findFieldWidth('%d',T_in); @@ -30,12 +29,12 @@ T = sort(T_in); % Make sure times are sorted discr = discrHand(m(j),order(i)); - k_max = discr.getTimestep(method); + k_max = discr.getTimestep(tsOpt); % Do we want to to save the initial conditions? if T(1) == 0 snapshot = discr.getTimeSnapshot(0); - saveToFile(sf, method, order(i), m(j),T(1), snapshot, NaN, NaN, discr); + saveToFile(sf, name, order(i), m(j),T(1), snapshot, NaN, NaN, discr); T(1) = []; end @@ -47,7 +46,8 @@ if is_int_multiples fprintf('Calculating time series in increments\n'); else - fprintf('Restarting for each time in timeseries\n'); + fprintf('RESTARTING for each time in timeseries\n'); + fprintf('If this is not what you want try giving T in integer multiples.\n'); end % T now contains all the times we need to step to, @@ -56,7 +56,8 @@ if is_int_multiples % Times are integer multiples, we can save time [k,N] = alignedTimestep(k_max,T(1)); - ts = discr.getTimestepper(method,k); + tsOpt.k = k; + ts = discr.getTimestepper(tsOpt); runtime = 0; for l = 1:length(T) end_step = N * time_multiples(l); @@ -65,20 +66,21 @@ ts.stepN(end_step-ts.n,true); runtime = runtime + toc(clock_start); snapshot = discr.getTimeSnapshot(ts); - saveToFile(sf, method, order(i), m(j),T(l), snapshot, runtime, k, discr); + saveToFile(sf, name, order(i), m(j),T(l), snapshot, runtime, k, discr); fprintf('Done! (%.3fs)\n',runtime); end else % Times are not interger multiples, we have to start from 0 every time. for l = 1:length(T) [k,N] = alignedTimestep(k_max,T(l)); - ts = discr.getTimestepper(method,k); + tsOpt.k = k; + ts = discr.getTimestepper(tsOpt); fprintf('[order = %-*d, m = %-*d, T = %-*d]: ',orderWidth,order(i),mWidth,m(j),TWidth,T(l)); clock_start = tic(); [v,t] = ts.stepN(N-ts.n,true); runtime = toc(clock_start); snapshot = discr.getTimeSnapshot(ts); - saveToFile(sf, method, order(i), m(j),T(l), snapshot, runtime, k, discr); + saveToFile(sf, name, order(i), m(j),T(l), snapshot, runtime, k, discr); fprintf('Done! (%.3fs)\n',runtime); end @@ -89,11 +91,11 @@ end -function saveToFile(sf, method, order, m, T, snapshot, runtime, k, discr) - key.method = method; - key.order = order; - key.m = m; - key.T = T; +function saveToFile(sf, name, order, m, T, snapshot, runtime, k, discr) + key.name = name; + key.order = order; + key.m = m; + key.T = T; entry.repr = snapshot; entry.runtime = runtime;
diff -r 7cb97c1988d9 -r 15baeb35f74e +scheme/Beam2d.m --- a/+scheme/Beam2d.m Tue Feb 23 10:55:46 2016 +0100 +++ b/+scheme/Beam2d.m Tue Feb 23 13:25:43 2016 +0100 @@ -1,4 +1,4 @@ -classdef SchmBeam2d < noname.Scheme +classdef Beam2d < noname.Scheme properties m % Number of points in each direction, possibly a vector N % Number of points total @@ -26,7 +26,7 @@ end methods - function obj = SchmBeam2d(m,lim,order,alpha,opsGen) + function obj = Beam2d(m,lim,order,alpha,opsGen) default_arg('opsGen',@sbp.Higher); default_arg('a',1);
diff -r 7cb97c1988d9 -r 15baeb35f74e +scheme/Wave.m --- a/+scheme/Wave.m Tue Feb 23 10:55:46 2016 +0100 +++ b/+scheme/Wave.m Tue Feb 23 13:25:43 2016 +0100 @@ -1,4 +1,4 @@ -classdef SchmWave < noname.Scheme +classdef Wave < noname.Scheme properties m % Number of points in each direction, possibly a vector h % Grid spacing @@ -20,7 +20,7 @@ end methods - function obj = SchmWave(m,xlim,order,alpha) + function obj = Wave(m,xlim,order,alpha) default_arg('a',1); [x, h] = util.get_grid(xlim{:},m);
diff -r 7cb97c1988d9 -r 15baeb35f74e +time/Timestepper.m --- a/+time/Timestepper.m Tue Feb 23 10:55:46 2016 +0100 +++ b/+time/Timestepper.m Tue Feb 23 13:25:43 2016 +0100 @@ -72,7 +72,7 @@ end function evolve_without_progress(obj, tend) - while obj.t < tend - obj.k/100 + while obj.t < tend - obj.k/2 obj.step(); end end @@ -84,7 +84,7 @@ steps_since_update = 0; last_update = tic(); s = util.replace_string('',' %d %%',0); - while obj.t < tend - obj.k/100 + while obj.t < tend - obj.k/2 obj.step(); steps_since_update = steps_since_update + 1;
diff -r 7cb97c1988d9 -r 15baeb35f74e +util/ReplaceableString.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/+util/ReplaceableString.m Tue Feb 23 13:25:43 2016 +0100 @@ -0,0 +1,50 @@ +classdef ReplaceableString < handle + properties + n + fmt + param + end + + methods + function obj = ReplaceableString(fmt, varargin) + default_arg('fmt', ''); + obj.n = 0; + + obj.setFormat(fmt); + obj.param = varargin; + end + + function setFormat(obj, fmt) + obj.fmt = fmt; + end + + function appendFormat(obj, fmt) + obj.fmt = [obj.fmt, fmt]; + end + + function update(obj, fmt, varargin) + obj.setFormat(fmt); + obj.param = varargin; + + obj.display(); + end + + function updateParam(obj, varargin) + obj.param = varargin; + obj.display(); + end + + function display(obj) + reverseStr = repmat(sprintf('\b'), 1, obj.n); + newStr = padStr(sprintf(obj.fmt, obj.param{:}),obj.n); + fprintf([reverseStr, newStr]); + + obj.n = length(newStr); + end + end + +end + +function b = padStr(a, n) + b = sprintf('%-*s', n, a); +end
diff -r 7cb97c1988d9 -r 15baeb35f74e Dictionary.m --- a/Dictionary.m Tue Feb 23 10:55:46 2016 +0100 +++ b/Dictionary.m Tue Feb 23 13:25:43 2016 +0100 @@ -72,7 +72,8 @@ end end otherwise - error('Unsupported indexing operator: %s',S.type); + B = builtin('subsref', obj, S); + % error('Unsupported indexing operator: %s',S.type); end end
diff -r 7cb97c1988d9 -r 15baeb35f74e SolutionFile.m --- a/SolutionFile.m Tue Feb 23 10:55:46 2016 +0100 +++ b/SolutionFile.m Tue Feb 23 13:25:43 2016 +0100 @@ -34,15 +34,18 @@ end function stupidSave(obj) + % Read file contents matObj = matfile(obj.filename,'Writable',true); keys = obj.matfile.keys; entries = obj.matfile.entries; + % Delete the file if exist(obj.filename,'file') delete(obj.filename); end + % Rewrite the file matObj = matfile(obj.filename,'Writable',true); matObj.keys = keys; matObj.entries = entries;
diff -r 7cb97c1988d9 -r 15baeb35f74e axPos.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/axPos.m Tue Feb 23 13:25:43 2016 +0100 @@ -0,0 +1,26 @@ +function axPos(n,m) + a = {}; + aPos = []; + for i = 1:n*m + a{i} = subplot(n, m, i); + aPos(i,:) = a{i}.Position; + end + + dx = aPos(1,3); + dy = aPos(1,4); + + x = unique(aPos(:,1)); + y = unique(aPos(:,2)); + + + fprintf('dx: %f\n', dx); + fprintf('dy: %f\n', dy); + + fprintf(' x: '); + fprintf('%f ', x); + fprintf('\n'); + + fprintf(' y: '); + fprintf('%f ', y); + fprintf('\n'); +end
diff -r 7cb97c1988d9 -r 15baeb35f74e minors.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/minors.m Tue Feb 23 13:25:43 2016 +0100 @@ -0,0 +1,31 @@ +function [minor, sub] = minors(A) + [n, m] = size(A); + + if n ~= m + error('A must be square'); + end + + sub = {}; + ks = {}; + + ind = 1:n; + for k = 1:n + C = nchoosek(ind,k); + for i = 1:size(C,1) + ks{end + 1} = k; + sub{end + 1} = A(C(i,:),C(i,:)); + end + end + + for i = 1:length(sub) + fprintf('%d:\n', ks{i}); + disp(sub{i}) + + minor(i) = det(sub{i}); + end +end + + +% A is positive semidefinite if all minors are non-negative +% A is negative semidefinite if all odd minors are non-positive and all even minors are non-negative +
diff -r 7cb97c1988d9 -r 15baeb35f74e runtestsAll.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/runtestsAll.m Tue Feb 23 13:25:43 2016 +0100 @@ -0,0 +1,22 @@ +function res = runtestsAll() + l = dir(); + + dirNames = {l([l.isdir]).name}; + + packages = {}; + for i = 1:length(dirNames) + if dirNames{i}(1) == '+' + packages{end+1} = dirNames{i}(2:end); + end + end + + rootSuite = matlab.unittest.TestSuite.fromFolder(pwd); + packageSuites = {}; + for i = 1:length(packages) + packageSuites{i} = matlab.unittest.TestSuite.fromPackage(packages{i}); + end + + ts = [rootSuite, packageSuites{:}]; + + res = ts.run(); +end
diff -r 7cb97c1988d9 -r 15baeb35f74e savepng.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/savepng.m Tue Feb 23 13:25:43 2016 +0100 @@ -0,0 +1,5 @@ +function savepng(h, filename, dpi) + print(h,filename,'-dpng',sprintf('-r%d',dpi)); + % Let print size in inches as input parameter + % Smaller boundingbox +end
diff -r 7cb97c1988d9 -r 15baeb35f74e semiDefIneq.m --- /dev/null Thu Jan 01 00:00:00 1970 +0000 +++ b/semiDefIneq.m Tue Feb 23 13:25:43 2016 +0100 @@ -0,0 +1,13 @@ +function ineq = semiDefIneq(A) + [m, sub] = minors(A); + + ineqsys = true; + for i = 1:length(m) + ineqsys = ineqsys & m(i) >= 0; + end + + ineq = simplify(ineqsys); + + str = toString(ineq); + fprintf('%s\n',strjoin(strsplit(str,' & '), '\n')); +end \ No newline at end of file
diff -r 7cb97c1988d9 -r 15baeb35f74e toString.m --- a/toString.m Tue Feb 23 10:55:46 2016 +0100 +++ b/toString.m Tue Feb 23 13:25:43 2016 +0100 @@ -19,6 +19,8 @@ str = cell2string(value); elseif isa(value,'function_handle') str = func2str(value); + elseif isa(value,'sym') + str = char(value); else warning('No string representation for class ''%s''', class(value)) str = 'NO_STR_REP';