Для численного интегрирования уравнения движения необходимо дисретизировать время и , используя разложение в ряд Тейлора, записать координаты расположения частиц на последующем временном шаге как функцию текущих координат, скоростей и уравнений, а в некоторых случаях – более высоких производных, следующим образом:
, (4)
.
По существу, это конечно-разностные выражения для координат и скоростей частиц.
Для интегрирования требуются два начальных условия для каждой частицы, так как уравнение движения имеет второй порядок по времени. Обычно начальным координатам и скоростям частиц приписываются определенные фиксируемые значения.
1. Далее с этими начальными значениями можно рассчитать силы, действующие на частицу в момент времени t = 0.
2. Затем, используя вышеприведенное уравнение координаты, которое обрывается после третьего слагаемого, можно получить новое значение координаты.
3. После этого по двум членам уравнения скорости можно получить новое значение скорости.
4. Ускорение, будучи пропорциональным силе, может быть непосредственно рассчитано на каждом шаге, поскольку координаты частицы известны.
Из уравнений (4) мы видим, что при таком обрыве уравнений ошибка в определении v составляет величину порядка а для .
Вследствие «обрыва» уравнений (4) слагаемые с b(t) и более высокого порядка (слагаемые третьего порядка или выше в уравнении координаты и второго порядка или выше в уравнении скорости) не не учитываются и не улулучшают результатов моделирования. Конечно, их можно рассчитать, записав конечно-разностные уравнения для ускорений, но в расчетах не требуются. Для уменьшения ошибки интегрирования по времени, сохраняя при этом шаг как можно большим, необходимы соответствующие алгоритмы.
1. Метод Верле. Сначала координата для предыдущего момента разлагается в степенной ряд в степенной ряд по координате, скорости и ускорению для текущего момента времени t, а затем складываются вместе два уравнения (для моментов и ), чтобы получить выражения для координаты и скорости частицы в момент времени :
, (5)
.
2. Метод с перешагиванием («чехарда»). Здесь скорости частиц рассчитываются в момент времени , т.е. на середине отрезка, соединяющего координаты частиц, тогда как сами координаты и ускорения определяются в моменты, кратные . Это может затруднить начало прогона моделирования. Однако, предполагая, что и , как следует ниже:
, (6)
.
Отметим, что если мы используем периодические граничные условия (как требуется при МД-моделировании макроскопическим путем) и на шаге частица возвращается в рассматривемую область с противопожарной стороны, это не приводит к каким-либо затруднением при моделировании наносистем. Следует также отметить, что в алгоритме Верле мы не получим точных значений скоростей в отличие от алгоритма с перешагиванием. Ошибки метода с перешагиванием – порядка , а СА алгоритм более устойчив, чем в методе Верле.
3. Метод скоростей Верле. Здесь, все координаты, скорости и ускорения частиц рассчитываются в моменты времени, кратные . Сначала, согласно обычным уравнениям движения, обновляются координаты частиц, затем в данных координатах вычисляется ускорение и далее из среднего значения двух ускорений рассчитывается скорость. Данный алгоритм имеет одинаковую точность для нано- и макроскопических систем, так как частицы, покидающие область с периодическими граничными условиями, не вносят каких либо ошибок в расчеты согласно данному алгоритму:
, (7)
.
В СВ- алгоритме требуется определить и сохранять все три переменные (a, v и r). Однако по стабильности этот алгоритм превосходит алгоритм Верле и алгоритм с перешагиванием, потому, что он основываться на неявном численном методе; т.е. скорость в момент времени рассчитывается из ускорения в моменты времени t и . С точки зренмия стабильности, преимущество СВ-алгоритма состоит в том, что он позволяет обновлять координаты с большим по величине шагом и, значит, можно производить более продолжительное время моделирования, используя меньшее число шагов.
4. Предикторно-корректорный метод Гира. Чтобы получать более точные траектории при том же временном шаге или иметь возможность применять большие шаги по времени, можно использовать вычислительный метод более высокого порядка, называемый предикторно-корректорным методом Гира (ПКГ). В этом методе траектории полагаются гладкими функциями времени и используется информация частиц на двух и трех предшествующих шагах времени. В принципе, сначала координаты, затем скорости и, наконец, ускорения предсказываются в момент времени , согласно записанной в конечно-разностной формуле, впрлоть до n-го порядка. Ниже представлены формулы для ПКГ третьего порядка
,
,
. (11)
Лишь после того, как положение частицы будет предсказано можно точно скорректировать ускорение aK в этот момент. Используя теорию возмущений, можно показать, что ошибка в предсказываемых координатах и скоростях пропорциональна , и, следовательно, с целью достижения максимально возможной точности мв можем внести необходимую поправку при помощи следующих уравнений.
,
,
. (12)
где с1 = 1, а остальные коэффициенты с целью определения максимальной точности можно определить методом описанным в (Gear C.W. Numerical Initial Value Prolems in Ordinary Differential Eqations. Tngltwood Cliffs. N.J.: Prentice Yfll. 1971 253 p.) В предикторно-корректорном методе третьего порядка с0 = 0, с1 = 1, с2 = 1, а в предикторно-корректорном методе четвертого порядка – с0 = 1/6, с1 = 5/6, с2 = 1/3. Кроме того, могут применяться методы более высокого порядка, но они требуют большего объема компьютерной памяти и могут потреблять больше ресурсов ЦПУ без существенного улучшения точности расчетов.
5. Выбор шага по времени . МД-моделировании важно выбирать шаг по времени ДТ при моделировании микроканонического ансамбля полная энергия системы должна сохраняться. Если величина Dt будет слишком большой, шаги станут тоже чрезвычайно большими, и частица может войти в область, запрещенную для классического движения, где потенциальная энергия является возрастающей функцией координат. Это может происходить, если сталкиваются две частицы , или частица ударяется о “барьер”,установленный внешним потенциалом. Вхождение в область, запрещенную в классической механике, означает, что новое значение потенциальной энергии стало больше ,чем максимально допустимое. В этом случае полная энергия возрастает и при больших шагах так происходит до тех пор, пока последовательность значений общей энергии не разойдется. Таким образом , в зависимости от соответствующей величины полной энергии приращение Dt следует выбирать достаточно малым, чтобы значение полной энергии всегда оставалось постоянным, и не слишком маленьким ,чтобы при моделировании не потребовалось чрезвычайно большого количества шагов. Оптимальное значение Dt обычно находится методом проб и ошибок. Одна фемтосекунда (10–15c) - хорошая пробная оценка шага , но оптимальное значение реально зависит от начальной энергии и типа рассматриваемого потенциала
, (4)
.
а для 
.
как можно большим, необходимы соответствующие алгоритмы.
разлагается в степенной ряд в степенной ряд по координате, скорости и ускорению для текущего момента времени t, а затем складываются вместе два уравнения (для моментов 
), чтобы получить выражения для координаты и скорости частицы в момент времени 
, (5)
.
, т.е. на середине отрезка, соединяющего координаты частиц, тогда как сами координаты и ускорения определяются в моменты, кратные 
и 
, как следует ниже:
, (6)
.
, (7)
.
рассчитывается из ускорения в моменты времени t и 
,
,
. (11)
, и, следовательно, с целью достижения максимально возможной точности мв можем внести необходимую поправку при помощи следующих уравнений.
,
,
. (12)