Equazioni differenziali a coefficienti non costanti

di il
473 risposte

473 Risposte - Pagina 11

  • Re: Equazioni differenziali a coefficienti non costanti

    Ho sbagliato io ad infilare i valori nella condizione dell'if, non ho considerato che prendevi traiettoria(end,1).
    Devi modificare abs(ys(:,1:3))>= r_pin1 || abs(ys(:,1:3))>= r_pin1 mettendo solo la colonna che ti serve.
  • Re: Equazioni differenziali a coefficienti non costanti

    Lo avevo pensato per un attimo ad un errore tuo, ma ho subito rimosso il pensiero...
    Vabbè, l'eccezione che conferma la regola che hai sempre ragione.

    a me serve vedere confrontare l'ultima coordinata x e y della particella con il raggio del diaframma e quelle che hanno valori più piccoli passano, quindi gli dico

    if (abs(ys(end,1)) <=r_pin1||abs(ys(end,2))<=r_pin1)

    Lui si fa i conti, non si incazza e alla fine "Particella" è una struttura 1x632 (ho detto n=1000).
    Quindi penso funzioni.
    EVVIVA!!
    ora prendo ste particelle e le faccio passare dal secondo diaframma... e si sfanculerà tutto, già lo so...


    A sproposito, ti ricordi quel casino degli assi?
    mi è rimasto che non fa il resize corretto e scrive il titolo sotto la barra degli strumenti della finestra della figura.
    Ieri una tipa mi ha chiesto quella parte di codice, lo ha copia nel suo script, NON lo ha modificato, e gli funziona che è una meraviglia.
    Ora io sono una cacca a programmare, però a questo punto non è solo colpa mia, anche l'amico (MatLab) ci mette il suo...
  • Re: Equazioni differenziali a coefficienti non costanti

    Sì sì... dai la colpa a matlab...
  • Re: Equazioni differenziali a coefficienti non costanti

    giug ha scritto:


    sì sì... dai la colpa a matlab...
    Non ho dato la colpa a lui, ma l'odio è reciproco evidentemente.
  • Re: Equazioni differenziali a coefficienti non costanti

    Invece di fargli risolvere l'equazione differenziale, visto che gli dò una volecità orientata nello spazio non gli posso fare valutare la posizione sul diaframma con una cosa diversa?

    tipo io parto dal vettore y_sorgente=[0,0,0, v_x,v_y,v_z] e gli dico di colcolarsi come è fatto questo vettore dopo una certa distanza, scrivendo qualcosa tipo:

    ys=@(distanza, y_sorgente) [bho?]

    no vero, sto delirando...
  • Re: Equazioni differenziali a coefficienti non costanti

    Devo ammettere che non ho capito cosa stai dicendo...
    Fammi un esempio con i numerini che la "velocità orientata nello spazio" con "la posizione sul diaframma" mi mandano in confusione...
  • Re: Equazioni differenziali a coefficienti non costanti

    Allora diciamo che per una particella io ho
    y_sorgente=[0,0,0, 1,2,3] %%gli zeri sono uguali a tutte le particelle gli altri numeri sono quelli estratti random

    dopo aver percorso d=1 metro lungo l'asse z quanto fa il nuovo vettore y?
    dovrebbe essere y_diaframma = [?, ?, 1, 1,2,3]
    le prime due sono quelle su cui lui deve fare la valutazione ed aventualmente ammazzarle,
    1 è la distanza percorsa
    1,2,3 è la velocità che non cambia perchè non ci sono forze, infatti per come lo abbiamo scritto lui risolve l'equazione differenziale:

    O=[0 0 0]';
    f = @(t,ys) [ys(4:6); O];

    Poi altra domanda:
    Io pensavo di poter usare k come contatore per il ciclo successivo, invece no, perchè è un numero piccolo mentre nella struttura ci sono parecchie particelle, come faccio a prendere il numero di particelle?
  • Re: Equazioni differenziali a coefficienti non costanti

    Continua a non essermi chiara la domanda. più che altro non ho capito, descrivendola come fosse una funzione, cosa vuoi che abbia in ingresso e cosa in uscita. Tu sai che si è spostata di un metro e in uscita vuoi che ti metta l'1 al terzo elemento del vettore? o deve trovare quelle che tu hai messo come punti interrogativi nell'y_diaframma?

    k è il numero di particelle sopravvissute (+1 dato che l'ultima volta che esce dal ciclo aumenta di 1 anche se ha finito). Dovrebbe essere esattamente uguale alla dimensione della struttura +1
  • Re: Equazioni differenziali a coefficienti non costanti

    giug ha scritto:


    Continua a non essermi chiara la domanda. più che altro non ho capito, descrivendola come fosse una funzione, cosa vuoi che abbia in ingresso e cosa in uscita. Tu sai che si è spostata di un metro e in uscita vuoi che ti metta l'1 al terzo elemento del vettore? o deve trovare quelle che tu hai messo come punti interrogativi nell'y_diaframma?
    In ingresso il vettore che uso come condizione iniziale, in uscita deve trovare le coordinate x e y (ossia i punti interrogativi) della particella una volta raggiunto il diaframma.
    Praticamente fargli fare quello che fa l'ode solver senza odesolver.

    giug ha scritto:


    k è il numero di particelle sopravvissute (+1 dato che l'ultima volta che esce dal ciclo aumenta di 1 anche se ha finito). Dovrebbe essere esattamente uguale alla dimensione della struttura +1
    Lo pensavo pure io ma non è così.
    o probabilmente fa qualcosa che non capisco perchè se nell'if gli faccio fare pure una matrice
    T(:,:,i)=ys;

    è di dimensioni 2x6x794 e anche la struttura è 794 ma è tutta fatta di zeri...


    Ti dò il codice:
    
    clear all
    close all
    
    uma=1.66e-27; %[Kg], unità di massa atomica
    q=1.602e-19; %[C], carica elementare
    
    m_e=5.485799e-4*uma; %[kg], massa dell'elettrone
    eV=1.6021765e-19; %[J], fattore di conversione eV-->J
    c= 2.99792458e8;  %[m/s], velocità della luce
    %KB=(8.617332*10^(-5))*eV; %[eV/K], costante di Boltzmann
    uma1=931.494028; %[MeV/c^2], fattore di conversione uma-->MeV/c^2 (massa)
    KB=8.617e-5; %[eV*K^-1], costante di Boltzmann
    KBT=1076*KB;
    
    n=1000; %numero diparticelle
    
    Ein=3;
    m_ion=1.00794*uma;
    
    q_over_m=q/m_ion;
    
    Etot_ion=zeros(n,1);
    betasquare_ion=zeros(n,1);
    v_ion=zeros(n,1);
    
    %%
    %Passaggio dal primo pinhole
    %%%%%%%%%%%%
    Pinhole1=1;
    DiametroPinhole1=0.002;
    r_pin1=DiametroPinhole1/2;
    
    %Have a for loop over a bunch of particles. For each particle, set its 
    %velocity (speed and angle) vector using the rand() function. Then use an 
    %"if" statement to decide whether the particle passed through or not and 
    %count those that do.
    
    k=1; %indice di supporto per selezionare le particelle nel ciclo interno if
    
    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);
                      
             %spread angolare 
             A = -20; B = 20;
             angolo(i) = A + (B-A) * rand(1); %random number in [A,B]
             
             teta=angolo(i); %angolo nel piano zy
             fi=angolo(i); %angolo con l'asse x
                      
             %%%%Componente z della velocità
             v_z= v_ion*cosd(teta)*cosd(fi);
              
             %calcolo la componente x              
             v_x= v_ion*sind(teta);
             
              %calcolo la componente y 
              v_y=v_ion*cosd(fi)*sind(teta);
             
              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); O]
              % 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);
                        
              %decido quali particelle passano il pimo diaframma
              
         if (abs(ys(end,1)) <=r_pin1||abs(ys(end,2))<=r_pin1)
                 Particella(i).traiettoria=ys(:,1:3);
                 Particella(i).velocita=ys(:,4:6);
                 k=k+1;
          end
    
               
    end
    
  • Re: Equazioni differenziali a coefficienti non costanti

    1keenan ha scritto:




    è di dimensioni 2x6x794 e anche la struttura è 794 ma è tutta fatta di zeri...
    ovviamente il 794 cambia ad ogni run
  • Re: Equazioni differenziali a coefficienti non costanti

    Intanto... Particella(i)->Particella(k) (dev'essere che ti funziona male il copia-incolla...
  • Re: Equazioni differenziali a coefficienti non costanti

    giug ha scritto:


    Intanto... Particella(i)->Particella(k) (dev'essere che ti funziona male il copia-incolla...

    No, ho sbagliato quando ho corretto la condizione...
  • Re: Equazioni differenziali a coefficienti non costanti

    Ma quindi così lui mi dice che su 1000 particelle ne passano 4...
  • Re: Equazioni differenziali a coefficienti non costanti

    Se k=5 sì.
    Diciamo che se lanci il codice diverse volte k varia da 2 a 5...
  • Re: Equazioni differenziali a coefficienti non costanti

    L'ho notato, ma non è possibile...
Devi accedere o registrarti per scrivere nel forum
473 risposte