Анализ численного решения системы дифференциальных уравнений

В настоящее время можно использовать множество программного обеспече­ния для осуществления матричных символьных вычислений [130]. Matlab (сокра­щение от англ. "Matrix Laboratory") является одним из основных инструментов более миллиона инженерных и научных работников для решения задач техниче­ских вычислений в силу его высокой интеграции с другими приложениями, а также работает на большинстве современных операционных систем, включая Linux, Mac OS и Microsoft Windows.

Реализацию программы для расчёта движения внешнего шарового слоя можно условно разделить на три этапа: загрузка начальных условий, непосредственное решение системы дифференциальных уравнений второго порядка с проверкой

74 условий отрыва шара от слоя и отрыва от барабана последнего шара и вывод ре­зультатов в виде графических фигур на экран. Исходный код программы для по­строения траекторий движения шаров внешнего слоя представлен в приложении 2.

Рисунок 2.11 Блок-схема процесса расчёта движения шаров внешнего слоя

Загрузка начальных условий (рисунок 2.11) выполняются с помощью встроен­ной функции xlsread().

После импорта данных и записи всех необходимых параметров (геометриче­ские размеры барабана, ступенек и мелющих тел, значения угловой скорости ба­рабана и коэффициента загрузки) для решения системы дифференциальных урав­нений записываются все требуемые математические формулы для приведения ис­ходной системы уравнений до пригодной для численного решения, предваритель­но определив основные матрицы модульной и угловой координат шаров и их производных как символьные с помощью служебного слова syms [98].

Далее в цикле (где S- количество шагов моделирования движения внешнего

шарового слоя) по формулам(2.44)-(2.47) производится расчёт углов, характери­зующих направление сил нормального давления на i-ый шар

Затем производится вычисление углов χ, используя формулы(1.1)-(2.54).По формулам(2.68)-(2.83) выполняется функция, выполняющая вычисление матриц расчёта A,B,C,D и Е. Решение получившейся системы на каждой итерации произ­водилось с помощью встроенного неявного метода Рунге-Кутты ode23tb,дающего удовлетворительные точность и быстродействие [100]. Полученные значения по­лярных координат шаров внешнего слоя φ. , p_iсохраняются и цикл завершается.

Заключительный этап - представление полученных результатов, предвари­тельно пересчитанных в декартовы координаты, в виде графических фигур вы­полняется с помощью функции plot() [54].

График, изображенный на рисунке2.12,показывает изменение аргумента по­лярной координаты первого шара sin(φ₁), по которому оценивалась корректность вычислений.

Рисунок 2.12 Изменение угловой полярной координаты первого шара sin(φ₁)

Метод расчёта, который был использован в программе, предполагает постоян­ную угловую скорость барабана, в силу чего на графике, представленном на ри­сунке2.12,наблюдается монотонное закономерное изменение угловых координат шара.

Траектория движения одного из шаров (рисунок2.13)позволяет определить момент отрыва шара от стенки и дальнейшее передвижение по объёму барабана.

Рисунок 2.13 Траектория движения первого шара

На рисунке2.14 показана траектория движения трех шаров внешнего шарово­го слоя.

Рисунок 2.14 Траектории движения трёх мелющих тел в барабане. Зелёный цвет - траектория первого шара, синий цвет - второго, красный - третьего шара (крайнего)

Таким образом, разработанная методика позволяет определить траектории

движения шаров внешнего слоя, их скорости, а также время отрыва и падения в зависимости от высоты и шага выступов футеровки.

2.3.