avm99963 | 802583e | 2022-05-12 10:56:56 +0200 | [diff] [blame] | 1 | program main |
| 2 | implicit none |
| 3 | integer*4, PARAMETER :: L = 32 |
| 4 | integer*2 :: S(1:L, 1:L) |
| 5 | real*8 :: MAGNE, ENERG, TEMP, E, DIFE, DELTA, M, SUMA, W(-8:8), genrand_real2 |
| 6 | integer*4 :: PBC(0:L+1), I, J, IMC, MCTOT, IPAS, N, SEED |
| 7 | character :: TEMPSTRING*100, str*20 |
| 8 | |
| 9 | if (command_argument_count() < 1) then |
| 10 | print *, "Error: please supply a command-line argument with the" |
| 11 | print *, "temperature, via './mc1 {TEMPERATURE}', where {TEMPERATURE}" |
| 12 | print *, "is the temperature." |
| 13 | stop |
| 14 | endif |
| 15 | |
| 16 | ! Inicialitzem algunes variables sobre el problema |
| 17 | call get_command_argument(1, TEMPSTRING) ! Ex: 2.0d0 |
| 18 | TEMPSTRING = trim(adjustl(TEMPSTRING)) |
| 19 | read (TEMPSTRING, *) TEMP |
| 20 | SEED = 234567 |
| 21 | MCTOT = 10000 |
| 22 | N = L*L |
| 23 | |
| 24 | ! Obrim el fitxer on escriurem els resultats a cada iteració |
| 25 | open(unit = 12, file = "out_data/" // "SIM-L-" // trim(str(L)) // "-TEMP-" // & |
| 26 | trim(str(int(1000*TEMP))) // "-MCTOT-" // trim(str(MCTOT)) //".out") |
| 27 | |
| 28 | ! Inicialitzem variable que ens fa latent la periodicitat |
| 29 | PBC(0) = L |
| 30 | PBC(L + 1) = 1 |
| 31 | do I = 1, L |
| 32 | PBC(I) = I |
| 33 | enddo |
| 34 | |
| 35 | ! Cache dels valors de l'exponencial |
| 36 | do I = -8, 8 |
| 37 | W(I) = exp(-float(I)/TEMP) |
| 38 | enddo |
| 39 | |
| 40 | ! Inicialitzem la matriu d'spins aleatòriament |
| 41 | call WRITECONFIG(S, L) |
| 42 | |
| 43 | E = ENERG(S, L, PBC) |
| 44 | M = MAGNE(S, L) |
| 45 | print *, "Energia inicial:", E, "Magnet. inicial:", M |
| 46 | write (12, *) "0", E, M |
| 47 | |
| 48 | ! Iterem amb el mètode de Montecarlo |
| 49 | do IMC = 1, MCTOT |
| 50 | do IPAS = 1, N |
| 51 | ! LOS LOOPS HACEN COSAS |
| 52 | I = INT(genrand_real2()*L) + 1 |
| 53 | J = INT(genrand_real2()*L) + 1 |
| 54 | SUMA = S(PBC(I - 1), J) + S(PBC(I + 1), J) + S(I, PBC(J - 1)) + S(I, PBC(J + 1)) |
| 55 | DIFE = 2*SUMA*S(I, J) |
| 56 | |
| 57 | if (DIFE > 0) then |
| 58 | DELTA = genrand_real2() |
| 59 | if (DELTA >= W(int(DIFE))) then |
| 60 | ! NO ho acceptem |
| 61 | cycle |
| 62 | endif |
| 63 | endif |
| 64 | |
| 65 | ! Sí ho acceptem: |
| 66 | S(I, J) = -S(I, J) |
| 67 | E = E + DIFE |
| 68 | M = M + 2*S(I, J) |
| 69 | enddo |
| 70 | if (mod(IMC, 1000) == 0) then |
| 71 | print *, "Iter:", IMC, "Energia:", E, "Magn:", INT(M) |
| 72 | endif |
| 73 | write (12, *) IMC, E, M |
| 74 | enddo |
| 75 | |
| 76 | ! Tanquem el fitxer de sortida |
| 77 | close(12) |
| 78 | endprogram main |