Innanzi tutto, grazie ad entrambi per gli interventi e gli spunti di riflessione e risoluzione.
Approfitto della pausa di Pasqua per vedere se e cosa cambia codificando i vostri suggerimenti.
Detto questo...
Quanto al tuo problema, dovrei analizzare il tuo codice e ricavarne la stabilità e l'approssimazione dei metodi numerici che usi, ma ora non ne avrei neppure il tempo.
cerco di riassumere brevemente, perché non credo ci voglia molto ad una persona preparata come te per capire dove possa essere il problema... se il problema gli viene spiegato da una persona chiara come me
Il codice riguarda una simulazione che ho fatto un paio di anni fa per la tesi e a cui ogni tanto rimetto mano per fare calcoli su regimi temporali differenti.
La simulazione riguarda l'evoluzione temporale di una "matrice densità".
Per evoluzione temporale si intende che
1) al "tempo zero" la matrice viene "popolata" da valori (di tipo double) inferiori ad 1 ma maggiori di zero in tutte le sue posizioni, tranne che nell'origine, in cui si mette un altro valore (compreso tra zero e il valore uniformemente distribuito sul resto della "griglia)
2) per i tempi successivi al primo, l'evoluzione temporale comporta che il valore di densità nel sito generico al tempo (t+1) è dato dalla media delle densità che ci sono nei suoi quattro primi vicini al tempo (t), dove il "conteggio" dei primi vicini è fatto con condizioni periodiche al contorno (è quindi un "toro" e non un rettangolo); questa evoluzione inoltre vale per tutti i siti, eccetto il sito nell'origine, il cui valore rimane fisso nel tempo.
L'evoluzione, quindi, consiste nella "propagazione" del valore inserito nell'origine all'interno del reticolo su cui invece sì è distribuito un valore uniforme. Con le ovvie analogie all'equazione del calore e tutto il resto e con tutta una serie di pippe matematiche di analogie con i random walk su reticolo e compagni bella, su cui non mi dilungo in quanto inutili alla fine della comprensione del problema.
Facendo notare a margine che la griglia ha un lato intero e dispari (per via del fatto che un asse, in entrambe le direzioni, deve essere inteso come asse di simmetria), il problema presenta una "enorme" simmetria, al punto che si può concentrare l'analisi su un ottavo del reticolo... ovvero una degli otto spicchi in cui è suddiviso dai suoi assi di simmetria. Concentrando l'analisi su questo singolo spicchio, si dovrà fare attenzione alle "nuove condizioni al contorno"
Detto questo, io faccio un tot di evoluzioni della matrice densità e poi faccio il conteggio di queste densità.
Dal conteggio di questi valori double che (in pratica) si stanno man mano avviando ad essere nulli, io ricavo il numero N di particelle contenute nella griglia. La densità del singolo sito è infatti interpretabile (visto che il sito è puntiforme) con il numero di particelle che lo occupano (numero per modo di dire... in quanto frazione di intero).
Il numero di particelle, come le singole densità, sarà un double che si avvia a diventare nullo.
Facendo la differenza tra il numero di particelle iniziali (un intero dell'ordine di qualche migliaio) e questo numero di particelle libere (double e piccolo), io mi ricavo il numero di particelle "virtualmente" intrappolate (dove virtualmente va inteso nel senso che, come spiegato sopra, la densità nell'origine non varia...).
Il problema, credo, non è dato dal fatto che io sommi un numero molto grande di singoli valori double. Almeno questo lo credo perché altre osservabili che dipendono da questa somma non mi danno problemi e sono vicine ai valori previsti teoricamente.
Io credo che il problema sia legato al fatto che alla fine io faccio la differenza appena spiegata tra un numero intero elevato e un numero double piccolo... perdendomi l'approssimazione, cosa molto grave perché poi di questi numeri approssimati devo fare logaritmi e altri "cazzi e mazzi" per i quali la precisione è fondamentale.
Comprendendo che (cosa che sto facendo prima di vedere di modificare il tipo di dato) io potrei semplicemente stamparmi il numero di particelle libere e usare questo per il datafit (modificando poco i conti teorici per ottenere quantità che dipendano dal numero di particelle libere e non dal numero di particelle intrappolate).
Pertanto, da grande ignorante in materia vi lascio due domande:
1) la somma di una elevata quantità di numeri double piccoli e dello stesso ordine può dare problemi?
2) la differenza tra un numero intero elevato e un numero double piccolo, allo stesso modo, può dare problemi?
Sto approfondendo le mie conoscenze di queste problematiche anche in virtù del fatto che questo tipo di problematiche sta pian piano diventando il mio lavoro
Grazie in anticipo per i vostri consigli di cui farò tesoro
