How to take 100 samples of a plot.

I have a task to take 100 samples of a plot at a fixed value of parameters a2, omega1, and omega2 and I will be changing a1 and take 100 samples. It consumes so much time so I am wondering is their a way to take a 100 samples of a plot.

clear all; close all; clc;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%                       Simulation                          %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%avgnum=10;
%for avg=1:avgnum
OMEGA1= 1:0.2:6.4;
OMEGA2= 6.6:0.2:12;
%Erroravg_nonzero_x1_y1=zeros([1,length(OMEGA1)]);
%Erroravg_zero_x1_y1=zeros([1,length(OMEGA1)]);
%Erroravg_nonzero_x2_y2=zeros([1,length(OMEGA2)]);
%Erroravg_zero_x2_y2=zeros([1,length(OMEGA2)]);
 for num1=1:length(OMEGA1)
     for num2=1:length(OMEGA2)
%value of constants
a1=0.2;a2=0.3;
%omega1=5;omega2=4;
omega1=OMEGA1(num1);
omega2=OMEGA2(num2);
%omega1= 1:0.5:10;
%omega2= 11:0.5:20;
G=1;C12=0.01;C21=0.02;
dt=0.01;   %step size
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
x1(1)=0.5;
y1(1)=0.5;
x2(1)=1;
y2(1)=1;
for i=2:1000
    x1(i)=x1(i-1)+((a1-x1(i-1)^2-y1(i-1)^2)*x1(i-1)-omega1*y1(i-1)+G*C12*(x2(i-1)-x1(i-1)))*dt;
    y1(i)=y1(i-1)+((a1-x1(i-1)^2-y1(i-1)^2)*y1(i-1)+omega1*x1(i-1)+G*C12*(y2(i-1)-y1(i-1)))*dt;
    x2(i)=x2(i-1)+((a2-x2(i-1)^2-y2(i-1)^2)*x2(i-1)-omega2*y2(i-1)+G*C21*(x1(i-1)-x2(i-1)))*dt;
    y2(i)=y2(i-1)+((a2-x2(i-1)^2-y2(i-1)^2)*y2(i-1)+omega2*x2(i-1)+G*C21*(y1(i-1)-y2(i-1)))*dt;
end
     end
 end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%                       Observation                        %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
Ndatamin=4;Ndatamax=31;
Data=Ndatamin:1:Ndatamax;
Data=Data/Ndatamax;
Error_nonzero_x1=zeros([1,Ndatamax-Ndatamin+1]);%initialization
Error_zero_x1=zeros([1,Ndatamax-Ndatamin+1]);%initialization
Error_zero_y1=zeros([1,Ndatamax-Ndatamin+1]);%initialization
Error_nonzero_x2=zeros([1,Ndatamax-Ndatamin+1]);%initialization
Error_zero_x2=zeros([1,Ndatamax-Ndatamin+1]);%initialization
Error_nonzero_y2=zeros([1,Ndatamax-Ndatamin+1]);%initialization
Error_zero_y2=zeros([1,Ndatamax-Ndatamin+1]);%initialization
count=0;
for N_m=Ndatamin:Ndatamax
    count=count+1;
    N_measurements=N_m;
    N_basis=31;
    %index=randi([20000,30000],1,N_m);
    index=randi([100,499],1,N_m);
    Xdot1=zeros([1,N_measurements]);
    Ydot1=zeros([1,N_measurements]);
    Xdot2=zeros([1,N_measurements]);
    Ydot2=zeros([1,N_measurements]);
    
    for ni=1:N_measurements
        Xdot1(ni)=(x1(index(ni)+1)-x1(index(ni)))/dt;
        Ydot1(ni)=(y1(index(ni)+1)-y1(index(ni)))/dt;
        Xdot2(ni)=(x2(index(ni)+1)-x2(index(ni)))/dt;
        Ydot2(ni)=(y2(index(ni)+1)-y2(index(ni)))/dt;
    end
    M=zeros([N_measurements,N_basis]);
    for i=1:N_measurements
        for j=1:N_basis
            if j==1
               M(i,j)=1;
            elseif j==2
               M(i,j)=x1(index(i));
            elseif j==3
               M(i,j)=y1(index(i));
            elseif j==4
               M(i,j)=x2(index(i));
            elseif j==5
               M(i,j)=y2(index(i));
            elseif j==6
               M(i,j)=x1(index(i))^2;
            elseif j==7
                M(i,j)=x2(index(i))^2;
            elseif j==8
                M(i,j)=y1(index(i))^2;
            elseif j==9
                M(i,j)=y2(index(i))^2;
            elseif j==10
                M(i,j)=x1(index(i))*x2(index(i));
            elseif j==11
                M(i,j)=x1(index(i))*y1(index(i));   
            elseif j==12
                M(i,j)=x1(index(i))*y2(index(i));  
            elseif j==13
                M(i,j)=x2(index(i))*y1(index(i));  
            elseif j==14
                M(i,j)=x2(index(i))*y2(index(i)); 
            elseif j==15
                M(i,j)=y1(index(i))*y2(index(i)); 
            elseif j==16
                M(i,j)=x1(index(i))^3;  
            elseif j==17
                M(i,j)=y1(index(i))^3;
            elseif j==18
                M(i,j)=x2(index(i))^3; 
            elseif j==19
                M(i,j)=y2(index(i))^3; 
            elseif j==20
                M(i,j)=x1(index(i))^2*x2(index(i));
            elseif j==21
                M(i,j)=x1(index(i))^2*y1(index(i));
            elseif j==22
                M(i,j)=x1(index(i))^2*y2(index(i));
            elseif j==23
                M(i,j)=x2(index(i))^2*x1(index(i));
            elseif j==24
                M(i,j)=x2(index(i))^2*y1(index(i));
            elseif j==25
                M(i,j)=x2(index(i))^2*y2(index(i));
            elseif j==26
                M(i,j)=y1(index(i))^2*x1(index(i));
            elseif j==27
                M(i,j)=y1(index(i))^2*x2(index(i));
            elseif j==28
                M(i,j)=y1(index(i))^2*y2(index(i));
            elseif j==29
                M(i,j)=y2(index(i))^2*x1(index(i));
            elseif j==30
                M(i,j)=y2(index(i))^2*x2(index(i));
            else   
                M(i,j)=y2(index(i))^2*y1(index(i));             
            end
        end
    end
end
    figure
    hold on
    plot(x1,'r')
    plot(x2,'g')

[Image Analyst]

Define "take".

[Heya :)]

It means that I have to run the code 100 times.

Answer 1

[Sindar]

First, what are you actually doing with M, Xdot1, etc.? Currently, the code overwrites them each N_m loop, without doing anything with the results

Try this: define sampled vectors, so you can later vectorize things:

index=randi([100,499],1,N_m);
x1_samp0 = x1(index);
x1_samp1 = x1(index+1);
...

Xdot1=(x1_samp1-x1_samp0)/dt;
...


M(:,1)=1;
M(:,2)=x1_samp0;
...
M(:,6)=x1_samp0.^2;

edit: changed my mind on other option; don't do it. Here are some general tips, though:

• if you have a set of conditions that can be expressed as variable=first, variable=second, etc. -- use a switch instead of if-elseif
• If you have a loop running over a set of conditions and operations which have nothing else in common -- just run each line individually, then fill in the "else" portion after

for j=1:N_basis
    if j==1
       M(i,j)=1;
    elseif j==2
       M(i,j)=x1(index(i));
    elseif j==3
       M(i,j)=y1(index(i));
    end
end

does the same as the below, but cleaner and faster:

M(i,1)=1;
M(i,2)=x1(index(i));
M(i,3)=y1(index(i));

• many, many, many for loops can be replaced with vector operations. It's kinda the point of Matlab