Ok, lasciamo perdere la cancellazione delle righe...
se io gli faccio risolve le eq. differenziali per ogni i e glidico di salvare i dati in T lui mi fa una matrice che non contiene le soluzioni...
usando questo ciclo:
for i=1:n
K_ion=normrnd(Ein,Ein/10);
Etot_ion(i)=K_ion+((m_ion*uma1)/uma); % ho tutto in MeV
betasquare_ion(i)=1-(((m_ion*uma1)/(uma*Etot_ion(i)))^2);
v_ion(i)=sqrt(betasquare_ion(i)*c^2);
v_z=v_ion;
%calcolo la componente x
v_x= unifrnd(-0.785,0.785, n,1); %angolo compreso da -45 e +45
% returns an array of random numbers generated from the continuous uniform distributions
% with lower and upper endpoints specified by A and B, respectively.
% If A and B are arrays, R(i,j) is generated from the distribution specified
% by the corresponding elements of A and B. If either A or B is a scalar, it is expanded
% to the size of the other input. R = unifrnd(A,B,m,n,...) or R = unifrnd(A,B,[m,n,...])
%returns an m-by-n-by-... array. If A and B are scalars, all elements of R are generated from
%the same distribution. If either A or B is an array, they must be m-by-n-by-...
%calcolo la componente x
v_y= unifrnd(-0.785,0.785, n,1);
v0=[v_x, v_y, v_z]; %initial velocity
r0=zeros(n,3); %initial speed
y_sorgente=[r0,v0]; %vector of initial condition
tspan=[0,100];
O=[0 0 0]';
%%%Risolvo il moto
f = @(t,ys) [ys(4:6); O];
% the expression [y(4:6); (q/m)*cross(y(4:6),B)]
% combines two vectors of length 3 to make a vector of
% length 6(as long as y is a column vector).
options=odeset('RelTol',1e-7,'Events',@(t,ys)Event_Stop_Sorgente(t,ys,Pinhole1)); % opzione per risoluzione equazione differenziale:
%set precisione calcolo---% Sintassi per utilizzare Event_Stop con Pinhole1 come variabile in input
[t,ys] = ode23t(f,tspan,y_sorgente(i,:),options);
T(:,:,i)=ys(:,:); %Matrice tridimensionale delle traiettore di tutte le particelle
lui fa una T come quella in figura
mentre le soluzioni che calcola non hanno tutti elementi nulli a parte l'ultimo.
Allegati: