Занятие 4. Вычисление параметров для молекулярной механики

Главная

1.

Дана оптимизированная структура этана в виде z-matrix.
Вместо значений длин и углов связей стоят переменные.
Цель - создать около 20 разных файлов для расчёта энергии в Gamess с разными значениями по длине одной из связей.
Для этого в представленных координатах одна связь или угол отличается названием переменной, но не её значением.

2.

Автоматизируем процесс с помощью скрипта на bash.
Составила файл-заготовку для размножения (et.inp).

3.

Tекстовый файл скрипта (make_b.bash).
Запустила скрипт:
bash ./make_b.bash
Получился 21 inp файл, и в каждом разное значение для переменной сс.

4.

Чтобы запустить расчёт для этих файлов, для этого перед строчкой с, нужно записать "done".
Вставка запуска Gamess:
gms b_${i}.inp 1 > b_${i}.log
Запустила скрипт. Для извлечения значений энергии из log файла воспользовалась awk.
awk '/TOTAL ENERGY =/{print $4}' b_${i}.log
Запустила скрипт. На экране должно появилось 21 значение энергии.
Чтобы выводить и значение длины связи, добавила перед вызовом awk распечатку переменой nb.
Распечатала переменную и несколько пробелов без переноса строки:
1 echo -n "$nb "
Запустила скрипт и перенаправила поток в файл :
bash ./make_b.bash > bond

5.

Есть зависимость энергии молекулы от длины одной связи.
зависимость нужно построить в gnuplot.
Построила зависимость энергии от длины связи : plot "bond".
Появился график с точками похожими на параболу.
Теперь надо было найти коэффициенты в функции f(x)=a+k(x-b)^2, которые бы позволили наиболее близко описать наблюдаемую зависимость.
Для этого воспользовалась возможностями Gnuplot.
Сначала задала функцию в развернутом виде, в строке gnuplot ввела: f(x)=a + k*x*x - 2*k*x*b + k*b*b
И задала стартовые значения коэффициентов:
a=-80
k=1
b=1.5
Провела подгонку коэффициентов под имеющиеся точки в файле bond: fit f(x) "bond" via a,k,b
Сохранила значения коэффициентов.
Построила графики функции и значений энергии из Gamess.

Рис.1.Зависимость энергии молекулы от длины С-С связи

Параметры f(x):
a = -79.7652
k = 0.563608
b = 1.55432
Функция не совсем совпадает с данными значениями. Возможно, потому что на самом деле зависимость энергии связи от ее длины задается более сложной функцией.

6.

Проделала аналогичные операции для валентного угла HCH, его значения должны изменяться от 109.2 до 113.2. (Значения).

Рис.2.Зависимость энергии молекулы от величины ССН угла

Параметры f(x):
a = -79.7647
k = 3.56076e-05
b = 111.38
Здесь апроксимация прошла успешно.

7.

Проделала аналогичные операции для торсионного угла d3, его значения должны изменяться от -180 до 180 c шагом 12. (Значения).

Рис.3.Зависимость энергии молекулы от величины торсионного угла d3

Это график без подгонки.
Количество минимумов у функции: 3.