2014-01-25
684
NEXT
NEXT
NEXT
NEXT
NEXT
NEXT
DIM p(n)
p(1) =.9 'вероятность первого состояния
FOR i = 2 TO n 'цикл заполнения остальных вероятностей
p(i) = 1 / 10 / n
s = s - p(i) * LOG(p(i)) 'подсчет энтропии
delt =.1 / n 'шаг изменения вероятностей
nstep = 1000000 'число шагов цепи Маркова
FOR j = 1 TO nstep 'цикл блуждания
i = INT(RND * n) + 1 'случайный выбор номера вероятн ости
dp = delt * (2 * RND - 1) 'изменение вероятности
p(i) = p(i) + dp
'страховка на отрицательность
IF p(i) <= 0 THEN p(i) = p(i) - dp: GOTO 1
ds = -dp * (s + LOG(p(i))) 'вариация энтропии
IF ds < 0 THEN p(i) = p(i) - dp: GOTO 1 'неудачный шаг
s = 0 'действия в случае удачного шага
summa = 0
FOR k = 1 TO n
p(k) = p(k) / (1 + dp)
s = s - p(k) * LOG(p(k))
summa = summa + p(k) 'сумма вероятностей для контроля
PRINT "энтропия s = ", s 'печать результатов
PRINT "энтропия равновесия s0 = ", LOG(n)
FOR i = 1 TO n
PRINT USING "#.####"; p(i)
PRINT "сумма вероятностей summa = ", summa
Пример экспериментальных данных.
число состояний n = 10
энтропия s = 2.301872
энтропия равновесия s0 = 2.302585
0.0999
0.1000
0.1001
0.0999
0.0999
0.1000
0.0999
0.0999
0.0999
0.1000
сумма вероятностей summa =.9994525
После релаксации система пришла в состояние с равными вероятностями.
 |
Задача 4. Моделирование канонического ансамбля системы независимых спинов.
'Программа: спины в термостате
INPUT "число спинов n=", n: DIM s(n) 'ввод числа спинов и создание массива
e = 1 'энергия спина
INPUT "e/kT=", ekT 'ввод относительной температуры
FOR i = 1 TO n 'заполнение массива состояний (все вверх)
s(i) = 1
nstep = 100000 'число МК шагов
n1 = n 'сколько спинов вверх
FOR j = 1 TO nstep 'цикл блуждания
i = INT(n * RND) + 1 'случайный выбор номера спина
dE = -2 * s(i) 'изменение энергии системы при повороте спина
'если энергия уменьшилась, то переворот остается и
'уменьшаем число спинов вверху, заканчиваем шаг
IF dE < 0 THEN s(i) = -s(i): n1 = n1 + dE / 2: GOTO 1
P12 = EXP(-dE * ekT) 'вероятность перехода при увеличении E
'разыгрывание вероятности перехода
IF RND <= P12 THEN s(i) = -s(i): n1 = n1 + dE / 2: GOTO 1
'подсчет среднего числа спинов вверху
1 n1cp = n1cp + n1 / nstep
'печать результатов
PRINT "число спинов в верхнем состоянии n1cp=", n1cp
PRINT "в последний момент n=", n1
PRINT "среднее значение энергии системы Ecp=", 2 * n1cp - n
'теорет. значение энергии
E0 = -n * e * (EXP(ekT) - EXP(-ekT)) / (EXP(ekT) + EXP(-ekT))
PRINT "теоретическое значение энергии Ecp0=", E0
число спинов n = 10
e/kT=.01 высокая температура
число спинов в верхнем состоянии n1cp = 4.95226
в последний момент n1 = 6
среднее значение энергии системы Ecp = - 9,547901E-02
теоретическое значение энергии Ecp0 = - 9,999669E-02
спины разделились примерно поровну, энергия примерно равна нулю
число спинов n = 10
e/kT = 1 средняя температура
число спинов в верхнем состоянии n1cp = 1.197173
в последний момент n1 = 2
среднее значение энергии системы Ecp = -7.605654
теоретическое значение энергии Ecp0 = -7.615942
теоретические и экспериментальные значения практически совпадают
число спинов n=10
e/kT=10 низкая температура
число спинов в верхнем состоянии n1cp =.00103
в последний момент n1 = 0
среднее значение энергии системы Ecp = -9.99794
теоретическое значение энергии Ecp0 = -10
вверху нет спинов, они «сконденсировались» вниз, энергия минимальна
