| Révision | 119ee7fcf277214ee87261ad2a0aec6f113f32d7 (tree) |
|---|---|
| l'heure | 2009-07-08 01:46:51 |
| Auteur | isella |
| Commiter | isella |
I added a function which implements the ideas by Radicchi (consider three
clusters each time and connect the two smaller ones; I do this with a
probability distribution). The code is under development...be really careful!
Python-codes/kernel_time-dependent.py (diff)| @@ -5,6 +5,16 @@ | ||
| 5 | 5 | import scipy.integrate as si |
| 6 | 6 | |
| 7 | 7 | |
| 8 | + | |
| 9 | + | |
| 10 | +def calculate_probability(x): | |
| 11 | + i = s.arange(len(x))[:,s.newaxis] | |
| 12 | + j = s.arange(len(x)) | |
| 13 | + max_ij = s.maximum(i,j) | |
| 14 | + P = x[::-1].cumsum()[::-1][max_ij]/x.sum() | |
| 15 | + | |
| 16 | + return P | |
| 17 | + | |
| 8 | 18 | #optimized version of the previous function |
| 9 | 19 | # I tested that this bit does the same as the previous function but it is much faster. |
| 10 | 20 | #I modified it by removing a couple of arguments which are indeed global variables |
| @@ -14,7 +24,9 @@ | ||
| 14 | 24 | |
| 15 | 25 | #kernel=y[:,s.newaxis]*y[s.newaxis,:]/100. #this is k_ij=ij |
| 16 | 26 | |
| 17 | - kernel=y[:,s.newaxis]*y[s.newaxis,:]/100. #this is k_ij=ij | |
| 27 | + prob_mat=calculate_probability(y) | |
| 28 | + | |
| 29 | + kernel=y[:,s.newaxis]*y[s.newaxis,:]/100.*prob_mat #this is k_ij=ij | |
| 18 | 30 | |
| 19 | 31 | |
| 20 | 32 | creation=s.zeros(lensum) |
| @@ -33,6 +45,15 @@ | ||
| 33 | 45 | return out |
| 34 | 46 | |
| 35 | 47 | |
| 48 | +# test_arr=s.array([12.,15.,123.,24.,1,19.]) | |
| 49 | + | |
| 50 | +# print "test_arr is, ", test_arr | |
| 51 | + | |
| 52 | +# matrix_prob=calculate_probability(test_arr) | |
| 53 | + | |
| 54 | +# print "matrix_prob is, ", matrix_prob | |
| 55 | + | |
| 56 | + | |
| 36 | 57 | y0=s.zeros(300) |
| 37 | 58 | |
| 38 | 59 | x=s.arange(len(y0))+1. #grid in cluster size |
| @@ -71,12 +92,9 @@ | ||
| 71 | 92 | fig = p.figure() |
| 72 | 93 | axes = fig.gca() |
| 73 | 94 | |
| 74 | -#axes.plot(s.log10(n_comp_large),s.log10(r_gyr_ari_large), "bo") | |
| 95 | + | |
| 75 | 96 | axes.plot(t,concentration, "ro",label="mean number of monomers in a cluster") |
| 76 | -#axes.plot(t,number_mono,label="monodisperse approx",linewidth=2.) | |
| 77 | -# axes.loglog(time_clu,n_clus_large, "rx",label="$N_>$") | |
| 78 | -# axes.loglog(time_clu[choice_time_clu],(10.**my_conc),"k",\ | |
| 79 | -# label=r"Fit of $N_\infty $",linewidth=2.) | |
| 97 | + | |
| 80 | 98 | p.xlabel('Time') |
| 81 | 99 | p.ylabel('Concentration') |
| 82 | 100 | p.title("Evolution of number of clusters") |
Fork