Deleting stored intermediate values in ODE15s

I'm solving a partial differntial equation in time and space, and for various reasons related to the research, I've discretized the equation and need to use a very large number of discretizations (500,000). I am using ODE15s, and it solves very slowly due to large usage of memory. I only need a small part of the output the last discretization node as a function of time. I want ODE15s to not store entries that are not in the 500,000 node (unless they are needed for calculating the next time point). I tried using the OutputFcn option, but I am still getting the full outputs (a 500,000 x number of timepoints structure).

options = odeset('JPattern', j, RelTol=1e-7,AbsTol=1e-7, OutputFcn=@(t, y, flag) store_last_node(t, y, flag));
out = ode15s(@(t,C)PDEs(input,t,C),0:60*input.SimTime, initialize, options);


function status = store_last_node(t, y, flag)
     persistent last_node_data;
    
    if strcmp(flag, 'init')
        last_node_data = [];
    elseif isempty(flag)
        % Store only the last node value
        last_node_data = [last_node_data; [t(end), y(end)]];  % Store last node only
    elseif strcmp(flag, 'done')
        % Save the final output
        assignin('base', 'last_node_output', last_node_data);
    end
    
    status = 0;  % Continue solving
end

[Walter Roberson]

    elseif isempty(flag)
        % Store only the last node value
        if size(last_node_data, 1) < 500  %or however many
            last_node_data = [last_node_data; t(end), y(end)];  % Store last node only
        else
            last_node_data = [last_node_data(2:end,:); t(end), y(end)]

[Nick]

Thanks, when I run this, out.y is still a large matrix of 500,000 x number of time points.

[Torsten]

@Nick

If this question is still related to your former question about the advection equation, it might really be worth studying and implementing the stable 2nd order scheme from the article I supplied. The need of 500.000 grid points in a simulation to sufficiently resolve the solution is really tremendous.

[Nick]

Thank you Torsten, I'll take another look at the article and try to apply a 2nd order stable scheme. I appreciate your help.

Answer 1

[Walter Roberson]

LastN = 500;
options = odeset('JPattern', j, RelTol=1e-7,AbsTol=1e-7);
out = ode15s(@(t,C)PDEs(input,t,C), [0, 60*(input.SimTime-(LastN-1:0)), initialize, options);