blob: 554890de907ec22208eba12b7b39ef7294dc087d [file] [log] [blame]
program main
use, intrinsic :: iso_c_binding, only: sp=>c_float, dp=>c_double
use iso_fortran_env
implicit none
integer*4 :: L
real*8 :: magne, energ, TEMP, E, DIFE, DELTA, M, SUMA, W(-8:8)
integer*4 :: I, J, IMC, MCTOT, MCINI, MCD, IPAS, N, SEED, NSEED, SEED0, pbc
integer :: genrand_int31
character(100) :: NOM
logical :: SAVEFINALCONF
integer :: SUMI
real(dp) :: SUME, SUME2, SUMM, SUMAM, SUMM2, VARE, VARM, CV, CHI
! This should be replaced by L but it doesn't work dynamically, so we're
! setting 64 as the upper limit which should be higher than L.
integer*2 :: S(1:64, 1:64)
! Inicialitzem algunes variables sobre el problema
namelist /DADES/ L, NOM, TEMP, NSEED, SEED0, MCTOT, MCINI, MCD, SAVEFINALCONF
read(nml = DADES, unit = 5)
N = L*L
! Cache dels valors de l'exponencial
do I = 0, 8
W(I) = (2**30 - 1 + 2**30)*exp(-float(I)/TEMP)
enddo
SUMI = 0
SUME = 0
SUME2 = 0
SUMM = 0
SUMAM = 0
SUMM2 = 0
do SEED = SEED0, SEED0 + NSEED - 1
print *, "Seed: ", SEED
! Inicialitzem la matriu d'spins aleatòriament
call generateSpinMatrix(S, L, SEED)
E = energ(S, L)
M = magne(S, L)
! Iterem amb el mètode de Montecarlo
do IMC = 1, MCTOT
do IPAS = 1, N
! LOS LOOPS HACEN COSAS
I = mod(genrand_int31(), L) + 1
J = mod(genrand_int31(), L) + 1
SUMA = S(PBC(I - 1, L), J) + &
S(PBC(I + 1, L), J) + &
S(I, PBC(J - 1, L)) + &
S(I, PBC(J + 1, L))
DIFE = 2*SUMA*S(I, J)
if (DIFE > 0) then
DELTA = genrand_int31()
if (DELTA >= W(int(DIFE))) then
! NO ho acceptem
cycle
endif
endif
! Sí ho acceptem:
S(I, J) = -S(I, J)
E = E + DIFE
M = M + 2*S(I, J)
enddo
if (IMC > MCINI .and. mod(IMC, MCD) == 0) then
SUMI = SUMI + 1
SUME = SUME + E
SUME2 = SUME2 + E*E
SUMM = SUMM + M
SUMAM = SUMAM + abs(M)
SUMM2 = SUMM2 + M*M
endif
enddo
enddo
! Normalitzem els promitjos
SUME = SUME/SUMI
SUME2 = SUME2/SUMI
SUMM = SUMM/SUMI
SUMAM = SUMAM/SUMI
SUMM2 = SUMM2/SUMI
VARE = SUME2 - SUME**2
VARM = SUMM2 - SUMM**2
CV = (SUME2 - SUME**2)/real(N*TEMP**2)
CHI = (SUMM2 - SUMAM**2)/real(N*TEMP)
! Guardem a un fiter els promitjos
open(unit = 13, file = "data_out/" // trim(NOM) // ".res")
! Valors: L, T, <E>, <E**2>, Var(E), <M>, <|M|>, <M**2>, Var(M), C_V, CHI
write(13, *) L, TEMP, SUME, SUME2, VARE, SUMM, SUMAM, SUMM2, VARM, CV, CHI
if (SAVEFINALCONF) then
call writeConfig(S, L, "data_out/" // trim(NOM) // ".conf")
endif
! Tanquem els fitxers de sortida
close(13)
endprogram main