Directional sensitivities

tests/TestStateJumps.m

@@ -195,7 +195,7 @@ function testBounceball(testCase)
             sensTimes = [t0 sensTimes tEnd];
             fdStep = generateFDstep(numel(x0), numel(p), 'hy', 1e-6, 'hp', 1e-6);
             sensFun = generateSensitivityFunction(datahandle, sol, 'method', 'VDE', 'FDstep', fdStep, ...
-                'CalcGy', true, 'CalcGp', true, 'Gmatrices_intermediate', true);
+                'CalcGy', true, 'CalcGp', true, 'Gmatrices_intermediate', false);
             sens = sensFun(sensTimes);
             Gy = {sens.Gy};
             Gp = {sens.Gp};
@@ -257,7 +257,7 @@ function testSensitivitiesSimpleVDE(testCase)
             rtol2 = 1e-5;
 
             sensFun = generateSensitivityFunction(datahandle, sol, 'method', 'VDE', ...
-                'CalcGy', true, 'CalcGp', true, 'Gmatrices_intermediate', true);
+                'CalcGy', true, 'CalcGp', true, 'Gmatrices_intermediate', false);
 
             t1 = sol.switches(1);
             t1Minus = leftLimit(t1);
@@ -293,7 +293,7 @@ function testSensitivitiesSimpleEND_piecewise(testCase)
             rtol2 = 1e-5;
 
             sensFun = generateSensitivityFunction(datahandle, sol, 'method', 'END_piecewise', ...
-                'CalcGy', true, 'CalcGp', true, 'Gmatrices_intermediate', true);
+                'CalcGy', true, 'CalcGp', true, 'Gmatrices_intermediate', false);
 
             t1 = sol.switches(1);
             t1Minus = leftLimit(t1);

toolbox/examples/canonicalExample/analyticalSensitivities/analyticalSensitivities_test.m

@@ -7,7 +7,7 @@
 tspan         = [0 20];
 initialvalues = [1;0];
 parameters    = 5.437;
-sol           = solveODE(datahandle, tspan, initialvalues, parameters); 
+sol           = solveODE(datahandle, tspan, initialvalues, parameters);
 
 %% Precalculations for sensitivities
 dim_y = size(sol.y, 1);
@@ -16,8 +16,10 @@
 
 %% Sensitivities VDE
 integrator_options = odeset( 'AbsTol', 1e-14,'RelTol', 1e-12);
-sensitivities_function_VDE = generateSensitivityFunction(datahandle, sol, 'FDstep', FDstep, 'integrator_options', integrator_options, 'method', 'VDE', 'Gmatrices_intermediate', true);
+sensitivities_function_VDE = generateSensitivityFunction(datahandle, sol, 'FDstep', FDstep, ...
+    'integrator_options', integrator_options, 'method', 'VDE', 'Gmatrices_intermediate', true, 'legacy', true);
 sensitivities_VDE = sensitivities_function_VDE(20);
+
 %% Analytical solution
 yA_1 = @(t) (1/300)*t.^3 + 1;
 yA_2 = 0;
@@ -89,4 +91,3 @@
 sensitivities_VDE.Up{1} - U1_p;
 sensitivities_VDE.Up{2} - U2_p;
 sensitivities_VDE.Gp - Gp_tf_t0;
-

toolbox/examples/canonicalExample/canonicalExampleAnalyticSolution.m

@@ -123,36 +123,5 @@
     yFull{idxModel} = @(t) [y{1, idxModel}(t); y{2, idxModel}(t)];
 end
 
-funcSolPiecewise = @(t) evalPiecewise(t, yFull, length(y0), tSwitch);
-end
-
-function y = evalPiecewise(t, solPiece, dimy, sw)
-% Make sure eval points and switching times are sorted in ascending order.
-[t, tSortIdx] = sort(t);
-sw = sort(sw);
-y = zeros(dimy, length(t));
-
-done = false;
-idxModelStart = 1;
-for idxSw=1:length(sw)
-    % Since eval points and switching times are sorted, we look for the first eval point that lies after the next switch.
-    idxModelEnd = (idxModelStart - 1) + find(t(idxModelStart:end) >= sw(idxSw), 1);
-    % If no such point exists, then we stay in the same model until the end.
-    if isempty(idxModelEnd)
-        y(:, idxModelStart:end) = solPiece{idxSw}(t(idxModelStart:end));
-        done = true;
-        break
-    end
-    % Otherwise, evaluate points before switch in the current model and update model for next iteration.
-    y(:, idxModelStart:idxModelEnd-1) = solPiece{idxSw}(t(idxModelStart:idxModelEnd-1));
-    idxModelStart = idxModelEnd;
-end
-
-% If we actually reach the last model or stay in first model, then we have to evaluate the remaining points in that model.
-if ~done
-    y(:, idxModelStart:end) = solPiece{end}(t(idxModelStart:end));
-end
-
-% Reverse the sorting
-y(:, tSortIdx) = y(:, 1:end);
+funcSolPiecewise = @(t) evalPiecewiseFunc(t, yFull, length(y0), tSwitch);
 end

toolbox/examples/canonicalExample/sensitivities_CanonicalExample.m

@@ -17,9 +17,9 @@
 %% Generation of sensitivity functions
 tic;
 integrator_options = odeset( 'AbsTol', 1e-14,'RelTol', 1e-12);
-sensitivities_function_ENDfull = generateSensitivityFunction(datahandle, sol, 'FDstep', FDstep, 'method', 'END_full', 'calcGy', true,'calcGp',true, 'Gmatrices_intermediate', true, 'save_intermediates', true);
-sensitivities_function_ENDpiece = generateSensitivityFunction(datahandle, sol, 'FDstep', FDstep, 'method', 'END_piecewise', 'calcGy', true,'calcGp',true, 'Gmatrices_intermediate', true, 'save_intermediates', false);
-sensitivities_function_VDE = generateSensitivityFunction(datahandle, sol, 'FDstep', FDstep, 'integrator_options', integrator_options, 'method', 'VDE', 'calcGy', true,'calcGp',false, 'Gmatrices_intermediate', true, 'save_intermediates', false);
+sensitivities_function_ENDfull = generateSensitivityFunction(datahandle, sol, 'FDstep', FDstep, 'method', 'END_full', 'calcGy', true,'calcGp',true, 'Gmatrices_intermediate', false, 'save_intermediates', true);
+sensitivities_function_ENDpiece = generateSensitivityFunction(datahandle, sol, 'FDstep', FDstep, 'method', 'END_piecewise', 'calcGy', true,'calcGp',true, 'Gmatrices_intermediate', false, 'save_intermediates', false);
+sensitivities_function_VDE = generateSensitivityFunction(datahandle, sol, 'FDstep', FDstep, 'integrator_options', integrator_options, 'method', 'VDE', 'calcGy', true,'calcGp',false, 'Gmatrices_intermediate', false, 'save_intermediates', false);
 toc;
 
 %% Sensitivty calculations
@@ -35,8 +35,8 @@
 toc;
 
 %% Plot sensitivities initial values
-sensitivities_function_ENDpiecewise = generateSensitivityFunction(datahandle, sol, 'FDstep', FDstep, 'method', 'END_piecewise', 'Gy', true,'Gp',false, 'Gmatrices_intermediate', true, 'save_intermediates', true);
-sensitivities_function_VDE = generateSensitivityFunction(datahandle, sol, 'FDstep', FDstep, 'method', 'VDE', 'calcGy', true,'calcGp',false, 'Gmatrices_intermediate', true, 'save_intermediates', true);
+sensitivities_function_ENDpiecewise = generateSensitivityFunction(datahandle, sol, 'FDstep', FDstep, 'method', 'END_piecewise', 'Gy', true,'Gp',false, 'Gmatrices_intermediate', false, 'save_intermediates', true);
+sensitivities_function_VDE = generateSensitivityFunction(datahandle, sol, 'FDstep', FDstep, 'method', 'VDE', 'calcGy', true,'calcGp',false, 'Gmatrices_intermediate', false, 'save_intermediates', true);
 
 t_plot = 0:0.01:20;
 %sensitivities_END_plot = sensitivities_function_ENDpiecewise(t_plot);
@@ -97,8 +97,8 @@
 
 %% Plot sensitivities parameters
 
-%sensitivities_function_ENDpiecewise = generateSensitivityFunction(datahandle, sol, 'FDstep', FDstep, 'method', 'END_piecewise', 'Gy', false,'Gp',true, 'Gmatrices_intermediate', true, 'save_intermediates', true);
-sensitivities_function_VDE = generateSensitivityFunction(datahandle, sol, 'FDstep', FDstep, 'method', 'VDE', 'calcGy', false,'calcGp',true, 'Gmatrices_intermediate', true, 'save_intermediates', true);
+%sensitivities_function_ENDpiecewise = generateSensitivityFunction(datahandle, sol, 'FDstep', FDstep, 'method', 'END_piecewise', 'Gy', false,'Gp',true, 'Gmatrices_intermediate', false, 'save_intermediates', true);
+sensitivities_function_VDE = generateSensitivityFunction(datahandle, sol, 'FDstep', FDstep, 'method', 'VDE', 'calcGy', false,'calcGp',true, 'Gmatrices_intermediate', false, 'save_intermediates', true);
 
 t_plot = 0:0.01:20;
 %sensitivities_END_plot = sensitivities_function_ENDpiecewise(t_plot);