Out of memory issue

Hi guys!

For research I'm running a ode45 simulation. The script is running fine however for more accurace I would like to decrease the length of the time interations. However, when I try to run the script with decreased time steps (Fs= 100000) my operating system seems to have some issues with the memory allocation. My script is posted below.

%% Script written by M.F.Duijn
driven_fre=@(t,y,tau,Delta,eta,g,J)[ ...
    Delta/(tau*pi)+2*y(1,:).*y(2,:)-g*y(1,:); ...
    y(2,:).^2+eta(t)-(pi*tau*y(1,:)).^2+tau*J.*y(1,:) ...
    ]/tau;

% parameters used in Fig. 4a
tau=10; % in ms!!
Delta=1;
J=0;
g=3;
eta=1;
frequency=[0:1 :(4.5*freq0)];% driving frequency in Hz [note that 30 Hz is not exactly the osc. frequency for eta=1]
drive=@(eta,d_eta,f,t) eta+d_eta*sin(2*pi*f/1000*t); % recall that we use ms
d_eta=0.5;   % driving amplitude

% for i=1:numel(eta)

    for j=1:numel(frequency)

        % define the integration time dependent on the driving frequency
         if frequency(j)
           T=[0 200000]/frequency(j); % roughly 200 periods of the driver
         else
            T=[0 200000]/freq0;
        end
      

        [t,y]=ode45(@(t,y)driven_fre(t,y,tau,Delta,...
            @(t)drive(eta,d_eta,frequency(j),t),g,J),T(1):1/100000:T(end),[0.1,1]);

        % remove the transient by selecting the last half
        ind=t>=t(end)-t(fix(end/2));
        t=t(ind); t=t-t(1);
        y=y(ind,:);
        Fs= 1/(1/100);
     
        
        % convert into the 'real' r and v variables
        r=y(:,1)*1000;      % this is to generate 'proper' Hz because time is in ms;
        
        [a_mean(j),freq(j)]= peaks(r,t)
         
         [pospks,posidx] = findpeaks(r,t); %find positive peaks
         [negpks,negidx] = findpeaks(-r,t);%find negative peaks 
         negpks=-negpks;                   %Convert negative peaks to right value
         
          for k= 1:(length(posidx)-1)
             dist(k)= (posidx(k+1)-posidx(k));
             dist(k)= dist(k)/1000;                  %Convert to seconds to compute frequency
          end          
         mean_dist= mean(dist);
         freq(j)= 1/mean_dist ;   
         a_mean(j)=mean(pospks)-mean(negpks);        % clear dist posidx negidx pospks negpks
        
    end
% end

There is probably a more subtle way to write this code however, I'm not that experienced with efficient coding. I would love to have more knowledge about this but I would really appriciate a starting point.

Im running this script on a 64-bit operating system with 16GB RAM, in my opinion this most be enough to run these kind of scripts.

Any helps is really appreciated!

Answer 1

[Jan]

What is freq0?

For j=2 the time span is:

T = [0 200000] / 1   % This is [0, 2e5]

In

[t,y]=ode45(@(t,y)driven_fre(t,y,tau,Delta, @(t)drive(eta,d_eta,frequency(j),t),g,J), ...
    T(1):1/100000:T(end), ...  % <- are you sure?
    [0.1,1]);

You dive this time span in steps of the size 1/1e5. This mean, that you ask ODE45 for 2e5*1e5 = 2e10 points for the output. For the time vector this needs 160 GB (8 bytes per double). For the trajectory y the double size is required. With 640 GB RAM one iteration might run successfully, but with 16 GB there is no chance.

Why do you want to output the trajectory with such a huge number of time steps?

[Maarten Duijn]

Freq0 is the natural frequency of the oscillator with the specific parameters. I calculate this before running this script with a ode45 without driver frequency. Freq0 is normally is about 30.7 or 35.2 Hz depending on the parameters used.

T = [0 200000] / frequency(j) what I try to do here is to make sure that I have around 200 periods of the driver frequency. The reason I am simulating with 1/1e5 is that I have noticed that the results are really inconsistent with different time steps.

[Jan]

Do you understand, that requesting a trajectory which needs 380 GB of RAM, must exhaust your memory?

If the results are not consistent, the problem might be stiff. Then forcing the ODE45 integrator to use tiny step sizes is a bad idea. You need another solver in this case.

Answer 2

[Anshika Chaurasia]

Hi Maarten,

It is my understanding that you are facing out of memory issues. Yo can refer to following link to resolve 'out of memory' error:

Why do I get 'Out of Memory' errors in MATLAB? - MATLAB Answers - MATLAB Central (mathworks.com)