Создание модели возникновения Солнечной системы из межзвездного газа на базе численного моделирования с учетом гравитационного взаимодействия частиц

Курсовой проект - Авиация, Астрономия, Космонавтика

Другие курсовые по предмету Авиация, Астрономия, Космонавтика

? элемент поверхности dА и предположим, что среднее количество движения, пересекающее в единицу времени нашу поверхность слева направо, равно S+, а S--среднее количество движения пересекающее поверхность в единицу времени справа налево. Тогда средняя сила F равна

 

F = S+ - S- (5.5)

 

и среднее давление определяется выражением

 

P = (5.6)

 

где Fn обозначает компоненту силы, нормальную к элементу поверхности. (В двумерном случае давление равно плотности потока импульса через единичный отрезок, а не через единичную площадку.) Что служит соответствующей мерой давления для потенциала Леннарда - Джонса?

Один из способов реализации этого метода определения среднего давления состоит в следующем. Представим себе, что на ребрах элементарной ячейки установлены четыре воображаемые поверхности. Тогда мы можем модифицировать подпрограмму реriоdiс и вычислять поток импульса за один шаг по времени.

 

procedure periodic(var xtemp,ytemp,xflux,yflux:real;,py,Lx,Ly:real);xtemp Ly then:=ytemp-Ly;:=yflux+py;

end;;

 

Среднее давление можно найти, добавив в основную программу следующие инструкции:

:=0.0;:=0.0;

Другой способ вычисления давлении вытекает из теоремы вириала:

 

, (5.7)

 

где ri- координата i-й частицы, Fi-полная сила, действующая на частицу i со стороны всех остальных частиц, и сумма берется по всем N частицам. Вывод формулы (5.7) дан в приложении 6А. Чтобы с помощью вириала вычислить давление, добавим в подпрограмму Verlet следующую инструкцию:

 

virial:=virial+x[i]*ax[i]+y[i]*ay[i];

 

Для вычисления среднего давления добавим в подпрограмму оutput следующие инструкции:

Итак, мы встретили уже два примера-(6.4) и (6.7), содержащие связь макроскопических величин со средними по времени от функций координат и скоростей частиц системы. В равновесном состоянии эти средние не зависят от времени. В гл. 15 и 16 мы будем рассматривать второй вид усреднения - среднее по статистическим ансамблям. Связь этих двух методов усреднения кратко обсуждается в гл. 15.

При моделировании методом молекулярной динамики процесс установления равновесия зачастую может занимать существенную часть от общего времени счета. Как правило, в качестве начальных условий удобнее всего выбирать равновесную конфигурацию из какого-либо прежнего расчета, которая отвечает температуре и плотности, близким к требуемым. Следующая подпрограмма запоминает последнюю конфигурацию пропускаемого варианта и может быть использована в качестве последней подпрограммы, вызываемой программой md.

Следующие инструкции можно добавить в подпрограмму initiаl, с тем чтобы можно было воспользоваться какой-то прежней равновесной конфигурацией в качестве начальной конфигурации нового расчета.

Необходимо только один раз провести моделирование молекулярной динамики, отправляясь от системы с произвольной конфигурацией. Последующие расчеты для различных состояний можно начинать, модифицируя условия, сформированные на момент завершения соответствующего предыдущего расчета. Вопрос о том, как изменять температуру системы при неизменной плотности, рассмотрен в задаче.

ЗАДАЧА 1. Качественные свойства жидкости и газа

а. Во всех пунктах данной задачи для получения своей начальной конфигурации используйте следующие инструкции dаtа и rеаd. Инструкция RЕАD читает элементы данных из инструкций DАТА и присваивает их значения соответствующим переменным. Каждая инструкция READ читает очередной элемент из инструкции DАТА. Напишите короткую программу, в которой используются инструкции RЕАD и DАТА, и удостоверьтесь, что вы правильно понимаете работу этих инструкций.

Выберите N = 16, Lх = Lу = 6, rsса1е = 1 и ?t = 0.02. Чему равна приведенная плотность? Выберите nsnар = 20 и проведите моделирование на протяжении по крайней мере 200 шагов по времени. Чему равна начальная температура системы? Каков характер ваших снимков? Как долго система достигает равновесия? Каков ваш критерий равновесия? Вычислите полную энергию системы, температуру и давление. Не учитывайте при усреднении по времени неравновесные конфигурации.

б. Как было найдено в п. а, полная энергия отрицательна. Какой смысл имеет знак полной энергии? Повторите моделирование с теми же начальными условиями (rscа1е = 1), но с Lх =Ly = 30. Чему равна полная энергия в этом случае? Тому же, что и в п. а? Если нет, то почему? Опишите характер полученных снимков. Заполняют ли частицы ящик равномерно за время порядка 200 шагов или у них наблюдается тенденция к образованию капелек?

в. Увеличьте начальное расстояние между частицами. Используйте Lx = Lу = 30 и rsсаlе = 2. Чему равна приведенная плотность системы? Оцените начальную плотность капельки. Какая получилась полная энергия -отрицательная или положительная? Если величину rsсаlе сделать достаточно большой, капелька должна в конце концов испариться, даже если полная энергия отрицательна. В реальной системе капелька находилась бы в равновесии со своим паром и обе фазы существовали бы совместно. По мере уменьшения плотности начальной капельки все больше частиц внутри капельки будет испаряться и размер капельки будет сокращаться. Разумеется, нужно с осторожностью относиться к любым заключениям о структуре системы, основанным на моделировании только с 16 частицами.

г. Положите Lх = Lу = 30 и rsса1е = 0.8. Чему равна приведенная плотность системы? Оцените начальную плотность капельки. Поскольку вначале частицы ближе друг к другу, не