Движение нейтрона в твердом теле

[Xumukys]

Здравствуйте, уважаемые форумчане!
Помогите, пожалуйста, составить прогу на VBA Excel по движению электрона в твердом теле.
Объем тела задается координатами некоторой прямоугольной области xmin, xmax; ymin, ymax. Задается также средняя длина свободного пробега нейтрона Lamda.
1.Разыгрывается случайная начальная координата x0 внутри заданного интервала [xmin, xmax] с помощью выражения q = a + r(b-a), в котором следует положить a = xmin, b = xmax.
2.Аналогично разыгрывается случайная начальная координата y0 внутри заданного интервала [ymin, ymax].
3.Разыгрывается угол φ случайного направления движения нейтрона на плоскости xy внутри интервала углов [0, 2pi] с помощью выражения q = a + r(b-a), в котором следует положить a = 0, b = 2pi.
4.Разыгрывается случайное значение длины свободного пробега L в соответствии с выражением L= - Lamda*ln(1-r).
5.Вычисляются значения координат точки конца свободного пробега по формулам:
x = x0 + L cosφ, y = y0 + L sinφ
6.Если частица оказывается за пределами хотя бы одного из интервалов [xmin, xmax] или [ymin, ymax], то рассмотрение ее движения заканчивается. Если обе координаты остаются внутри объема, в точке с этими координатами происходит столкновение нейтрона с ядром атома вещества. Рассмотрим здесь три возможных случайных варианта последствий такого столкновения. Для моделирования случайности каждого варианта разделим интервал случайного числа r, т.е. интервал [0, 1] на три подинтервала: [0, r1], (r1, r2], (r2, 1] (круглая скобка означает, что данный конец интервала в сам интервал не включен). Разыгрываем случайное число r и если оно попадает в первый подинтервал, реализуется первый вариант, если во второй, то второй, если в третий, то третий.
7.Первый вариант – в результате столкновения происходит упругое изотропное рассеяние нейтрона ядром атома. В этом случае координата точки столкновения принимается за новую начальную координату:
x0 = x, y0 = y
и моделирование движения этого нейтрона повторяется, начиная с п. 3 данного алгоритма.
8. Второй вариант – нейтрон поглощается ядром атома вещества и его существование в качестве свободной частицы прекращается.
9.Третий вариант – столкновение нейтрона с ядром приводит к делению этого ядра и к появлению еще одного нейтрона. В этом случае оба нейтрона (старый и новый) начинают движение с координаты, определяемой x0 = x, y0 = y , движение каждого из них воспроизводится, начиная с п. 3.

Моделируемая цепная ядерная реакция происходит во времени, поэтому необходимо предусмотреть в программе счетчик времени в некоторых условных единицах. При запуске программы счетчик должен принимать нулевое значение, а при каждом обращении к п. 3 он должен увеличиваться на единицу.
Вводится зависимость величины границы между вторым и третьим подинтервалами r2 от относительного количества нейтронов:
r2=r02(1+kv* ((n-nv)/nv)
Здесь r02 – первоначально задаваемая величина границы между вторым и третьим интервалами, n – текущее количество нейтронов, nv – задаваемое количество нейтронов, вблизи которого должны происходить флуктуации текущего количества n, kv – задаваемый коэффициент пропорциональности, определяющий эффективность системы автоматического регулирования.
Исходными данными являются средняя длина свободного пробега Lamda, параметры интервалов r1 и r02, максимальное модельное время tmax, начальное количество нейтронов n0, задаваемое среднее значение количества нейтронов nv (рабочий режим реактора), а также коэффициент kv, определяющий эффективность системы автоматического регулирования режима работы реактора.
Границы моделируемой области можно в данном случае задать фиксированными в тексте программы (xmin = 0, xmax = 1, ymin = 0, ymax = 1)
r = Rnd(1)

Помогите, пожалуйста, очень нужно.