Математическое моделирование УдГАУ. Каталог.
Содержание
Исходное задание и порядок выполнения 3
1 Моделирование электрической цепи 4
1.1 Исходные данные и графическая модель 4
1.2 Математическая модель 5
1.3 Решение методом Гаусса 6
1.4 Решение методами Зейделя и Якоби 7
1.5 Токи ветвей и проверка балансов 8
1.6 Алгоритмы и сравнение методов 9
2 Моделирование системы SRS на автомобиле с ППП 10
2.1 Объект и функциональная схема 10
2.2 Параметры и программные обозначения 11
2.3 Математическая модель полного привода 12
2.4 Режимы и алгоритм функционирования 13
2.5 Модель давления и аналитическое решение 14
2.6 Результаты разгона из состояния покоя 15
2.7 Аварийный режим и реакция SRS 16
2.8 Точность и выводы 17
3 Оптимизация загрузки оборудования 18
3.1 Собственная постановка для профиля ТСА 18
3.2 Максимизация прибыли 19
3.3 Доказательство оптимальности 20
3.4 Минимизация затрат 21
3.5 Графическое решение 22
3.6 Проверка затрат и рекомендации 23
Заключение 24
Список использованных источников 25
Приложение А Программа метода Гаусса 26
Приложение Б Программа итерационных методов 27
Приложение В Программа модели SRS и ППП 28
Приложение Г Программа оптимизации 30
2
Исходное задание и порядок выполнения
Контрольная работа выполнена по разделу 4 учебного пособия А. В.
Храмешина «Моделирование в агроинженерии» [1]. Рассмотрены
электрическая цепь постоянного тока, динамический процесс на автомобиле
выездной сервисной службы с постоянным полным приводом и системой
SRS, а также загрузка ремонтного оборудования.
Для варианта 99 в таблице 3 принят столбец 9 по обоим буквенным
индексам: а = 9, б = 9. В таблице 4 система выбрана по столбцу 9 — SRS; в
строке «привод б» по индивидуальному заданию принято ППП —
постоянный полный привод вместо прочерка.
Таблица 1 — Состав индивидуального задания
Задание
Исходные данные
1 Электрическая
цепь
Схема 2; R₁ = 25,8 Ом; R₂ = 18,6 Ом; R₃ = 36,8 Ом; R₄ = 42,8 Ом; E₁
= 24 В; E₂ = 32 В; E₃ = 64 В.
2 Динамический
процесс
SRS; ППП — постоянный полный привод. Режимы: тестовый,
рабочий и аварийный.
3 Оптимизация
Пункт 3.3 — загрузка оборудования. Максимум прибыли при
выборе выпуска и минимум затрат на выполнение заказа.
Данные третьего задания самостоятельно сформированы для
мастерской технического сервиса в АПК. Числа и постановки не
заимствованы из таблиц вариантов третьего задания. Производственные
нормы и экономические показатели приняты как учебные исходные данные.
Постоянный полный привод учтён в распределении тяги между осями
и проверке сцепления. SRS рассматривается как отдельная подсистема
пассивной безопасности: её выходом служит отклик удерживающего
устройства, а не тормозная сила автомобиля.
В первом задании применены прямой метод Гаусса и два
итерационных метода — Зейделя и Якоби. В работе приведены схемы,
таблицы соответствия обозначений программным переменным, расчёты в СИ
и программы для воспроизведения результатов.
3
1 Моделирование электрической цепи
1.1 Исходные данные и графическая модель
По таблицам 2 и 3 пособия выбрана схема 2. Источники ЭДС идеальны,
сопротивления постоянны. Требуется найти токи всех ветвей.
Положительные направления I₁, I₂ и I₃ приняты от узла A к узлу B; контурные
токи J₁ и J₂ направлены по часовой стрелке.
Рисунок 1 — Схема 2 с параметрами варианта 9
Таблица 2 — Параметры и программные обозначения
Параметр
Значение
Единица
В программе
R₁; R₂
25,8; 18,6
Ом
Коэффициенты A
R₃; R₄
36,8; 42,8
Ом
Коэффициенты A
E₁; E₂; E₃
24; 32; 64
В
Вектор b
J₁; J₂
Искомые
А
x[0]; x[1]
Верхняя ветвь содержит последовательно соединённые E₁, R₁ и R₄.
Между R₁ и R₄ нет разветвления, поэтому I₄ = I₁. Средняя и нижняя ветви
содержат соответственно E₂, R₂ и E₃, R₃.
4
1.2 Математическая модель
После объединения R₁ и R₄ цепь имеет два существенных узла и три
ветви. Число независимых контуров: 3 − 2 + 1 = 2. Сопротивление верхней
ветви равно:
   
Контурные токи позволяют автоматически выполнить первый закон
Кирхгофа:
   
  
Обход верхнего контура проходит по верхней ветви A → B и по
средней B → A. Нижний контур проходит по средней A → B и нижней B →
A. Второй закон Кирхгофа даёт:
󰇛󰇜 
󰇛󰇜 
Числовая система уравнений A·J = b:
 
 
87,2 -18,6
-18,6 55,4
󰇛󰇜󰢹 󰇛󰇜󰢹
Размерности согласованы: коэффициенты A выражены в омах, J — в
амперах, b — в вольтах. Определитель не равен нулю:
   
Решение единственно. Выполняется строгое диагональное
преобладание: 87,2 > 18,6 и 55,4 > 18,6. Оно обеспечивает сходимость
используемых итерационных методов.
5
1.3 Решение методом Гаусса
Составляем расширенную матрицу. В программе используется выбор
главного элемента по текущему столбцу; для данной системы перестановка
строк не требуется.
87,2 -18,6 -8
-18,6 55,4 -32
Для исключения J₁ из второй строки находим множитель:
  
Выполняем преобразование L₂ := L₂ − m₂₁L₁. Новый диагональный
элемент и правая часть:
  
 󰇛󰇜
87,2 -18,6 -8
0 51,432568807339 -33,706422018349
Обратный ход: сначала находим второй контурный ток, затем первый.
При вычислениях коэффициенты до вывода не округляются:
 
 
 󰇛󰇜 
Подстановка в исходную систему с неокруглёнными значениями:
 
 
Невязка в выполненном расчёте не превышает 10⁻¹³ В. Отрицательные
значения контурных токов означают направление, противоположное
выбранному обходу. Это допустимый физический результат, а не ошибка
метода.
6
1.4 Решение методами Зейделя и Якоби
Начальное приближение J⁽⁰⁾ = (0; 0)ᵀ А. В методе Зейделя новое
значение J₁ сразу используется при расчёте J₂:
󰣈 󰇛󰣈󰇜
󰣈 󰇛󰣈󰇜
Критерии остановки проверяются одновременно: δ = ‖J⁽ᵏ⁺¹⁾ − J⁽ᵏ⁾‖∞ ≤
10⁻¹⁰ А и ρ = ‖AJ − b‖∞ ≤ 10⁻⁹ В. Предельное число итераций — 10 000.
Таблица 3 — Последовательность приближений Зейделя
k
J₁, А
J₂, А
0
0,0000000000
0,0000000000
1
-0,0917431193
-0,6084191700
2
-0,2215206028
-0,6519906717
3
-0,2308145240
-0,6551110135
5
-0,2315277664
-0,6553504775
8
-0,2315314418
-0,6553517115
10
-0,2315314431
-0,6553517120
Условия выполнены на 10-й итерации: δ = 8,978·10⁻¹¹ А, ρ = 5,607·10⁻¹⁰
В. Максимальное отличие от метода Гаусса равно 6,926·10⁻¹² А.
В методе Якоби обе новые компоненты вычисляются только по
старому вектору; первая формула сохраняется, во второй используется J₁⁽ᵏ⁾:
󰣈 󰇛󰣈󰇜
Якоби достигает тех же допусков за 19 итераций: ρ = 5,323·10⁻¹⁰ В.
Таким образом, использованы два разных итерационных способа, а метод
Гаусса служит дополнительной прямой проверкой.
 󰇛󰇜 
   
Числа qЗ и qЯ — спектральные радиусы матриц итераций. Меньшее
значение для Зейделя объясняет более быстрое уменьшение ошибки.
7
1.5 Токи ветвей и проверка балансов
Таблица 4 — Результаты расчёта токов
Ветвь
Ток, А
Фактическое направление
E₁, R₁, R₄
-0,231531443
От B к A
E₂, R₂
-0,423820269
От B к A
E₃, R₃
0,655351712
От A к B
Первый закон Кирхгофа:
 
Разность потенциалов U = φB − φA одинакова при расчёте по любой
ветви:
󰇛󰇜 
 
 
Мощность сопротивления всегда положительна: P = RI². Отдаваемая
мощность источника равна EI, если положительный ток совпадает с
направлением его ЭДС.
Таблица 5 — Баланс мощностей
Элемент
Мощность, Вт
Характер
R1
1,383056
Потребление
R2
3,340999
Потребление
R3
15,805080
Потребление
R4
2,294371
Потребление
E1
-5,556755
Поглощение
E2
-13,562249
Поглощение
E3
41,942510
Отдача
  
Источник E₃ отдаёт 41,942510 Вт, а источники E₁ и E₂ поглощают
соответственно 5,556755 и 13,562249 Вт. Оставшаяся мощность рассеивается
на резисторах. Невязка энергетического баланса меньше 10⁻¹² Вт.
8
1.6 Алгоритмы и сравнение методов
Рисунок 2 — Прямой и итерационные алгоритмы
Рисунок 3 — Уменьшение невязки
Таблица 6 — Сравнение методов
Показатель
Гаусс
Зейдель
Якоби
J₁, А
-0,231531443
-0,231531443
-0,231531443
J₂, А
-0,655351712
-0,655351712
-0,655351712
Итерации
Прямой
10
19
Все методы дают одинаковые токи в пределах заданной точности. Проверки по
законам Кирхгофа и балансу мощности подтверждают правильность знаков и соединений.
Полные программы приведены в приложениях А и Б.
9
2 Моделирование системы SRS на автомобиле с ППП
2.1 Объект и функциональная схема
Объект — автомобиль выездной службы технического сервиса АПК с
постоянным полным приводом. По таблице 4 выбрана система SRS. Блок
пассивной безопасности анализирует датчики и активирует удерживающие
устройства при соответствующем событии; его калибровка зависит от
автомобиля [2].
Рисунок 4 — ППП и SRS в общей модели автомобиля
Модель состоит из двух связанных по входному сигналу частей:
разгона полноприводного автомобиля и реакции SRS на заданный удар.
Начальное состояние — покой. В тестовом и рабочем режимах выполняется
разгон, в аварийном после разгона вводится фронтальный удар.
Тяга передаётся обеим осям постоянно в долях 50:50. SRS не управляет
этим распределением. Её расчётный выход — избыточное давление p(t).
Взаимодействие подушки с человеком и критерии травмирования не
моделируются.
Горизонтальная дорога, неизменные нагрузки на оси и линейное
сопротивление движению приняты как упрощения. Пробуксовка отдельных
колёс, дифференциалы и упругость шин не рассчитываются. Численные
параметры ниже относятся только к учебной модели.
10
2.2 Параметры и программные обозначения
Таблица 7 — Параметры и интервалы исследования
Величина
Принято
Интервал
Имя
Масса m, кг
2000 ТР; 2400 РР
2000–2600
m
Сцепление μ
0,65
0,30–0,80
MU
Качение f
0,015
0,010–0,025
F
Линейное сопротивление c, Н·с/м
40
20–60
C
Радиус колеса r, м
0,34
0,32–0,36
R
Доля нагрузки спереди β
0,55
0,50–0,60
BETA
Доля тяги спереди α
0,50
Постоянна
ALPHA
Суммарный момент T*, Н·м
900 ТР; 1800 РР
0–1800
demand
Постоянная привода τ, с
0,25
0,15–0,40
TAU
Начало удара tс, с
6
Задано
TC
Замедление удара aуд, м/с²
120
80–160
DEC
Порог aп, м/с²; Δvп, м/с
30; 2,4
20–40; 1,8–3,0
AP; DV
Задержка tд, с
0,003
0,002–0,005
DELAY
Масштаб P₀, Па
80 000
60 000–100 000
P0
τн; τс, с
0,008; 0,070
0,006–0,012; 0,05–0,09
TR; TD
Пороги, задержка и параметры давления назначены для
вычислительного эксперимента; они не являются калибровками серийной
SRS. g = 9,81 м/с². Основной шаг разгона h = 0,01 с; давления — 0,0005 с.
Таблица 8 — Переменные состояния и результаты
Обозначение
Смысл и единица
В программе
s; v
Путь, м; скорость, м/с
y[0]; y[1]
T
Суммарный момент, Н·м
R·force·(1−exp(−t/TAU))
w; D
Изменение скорости, м/с; исправность
DEC·(t−TC); healthy
p; u
Давление, Па; время после команды, с
p; u
11
2.3 Математическая модель полного привода
Статические нормальные нагрузки Nп = βmg, Nз = (1 − β)mg. При
равных радиусах колёс силы тяги равны Fп = αT/r и Fз = (1 − α)T/r. Для
сохранения заданного деления момента оба ограничения сцепления должны
выполняться одновременно:
  
 󰇛󰇛󰇜󰇜
 󰇛󰇜 
В рабочем режиме Nп = 12 949,2 Н, Nз = 10 594,8 Н. Поэтому Fпред =
min(16 833,96; 13 773,24) = 13 773,24 Н. Задание 1800 / 0,34 = 5294,118 Н
ниже этого предела; на каждую ось приходится 2647,059 Н.
Инерционность передачи момента зададим апериодическим звеном:
 󰇛󰇜󰇛󰇜 
󰇛󰇜  󰇟󰇛󰇜󰇠
До удара и после начала движения уравнения имеют вид:
  󰇛󰇜
При v = 0 скорость удерживается нулевой, пока тяга не превысит fmg.
Начало движения при F∞ > fmg определяется точно:
 󰇛󰇜
Обозначим u = t − tп, k = c/m, λ = 1/τ, A₀ = (F∞ − fmg)/m. Для u ≥ 0
аналитическое решение:
󰇛󰇜 󰇛󰇜
󰇟󰇛󰇜󰇛󰇜󰇛󰇜󰇠
 󰇟󰇛󰇜󰇛󰇜󰇠
󰇛󰇜󰇛󰇜
До tп имеем s = v = 0. Размерности: T/r, fmg и cv — силы в ньютонах;
правая часть уравнения скорости — м/с². Аналитические формулы
используются для контроля численного расчёта.
12
2.4 Режимы и алгоритм функционирования
Таблица 9 — Постановка вычислительных опытов
Режим
Условия движения
Холостой
ход
Tзад = 0; автомобиль в покое
ТР
Без груза: m = 2000 кг, Tзад = 900 Н·м
РР
С грузом: m = 2400 кг, Tзад = 1800 Н·м
АР
До 6 с — РР; затем заданный удар
Рисунок 5 — Алгоритм функционирования учебной модели
В АР с момента tс = 6 с тяга снимается, а движение на коротком ударном участке
задаётся замедлением aуд = 120 м/с² до v = 0. Это внешняя сила столкновения, не
ограниченная сцеплением колёс; уравнение обычного разгона на этом участке заменяется
ударной моделью.
В течение события w = ∫aуд dt. В модели команда фиксируется при D = 1, aуд ≥ 30
м/с² и w ≥ 2,4 м/с. Для выбранного постоянного импульса это момент t* = 6 + 2,4/120 =
6,020 с. После задержки tд = 0,003 с начинается отклик при t₀ = 6,023 с.
При D = 0 команда в учебной логике блокируется и отмечается неисправность; это
отдельная проверка, а не универсальное правило реальных блоков SRS. В ТР и РР
ударного события нет, поэтому p = 0.
13
2.5 Модель давления и аналитическое решение
После однократной команды введём локальное время u = t − t₀. До t₀
принимаем p = 0. После команды сосредоточенная модель источника
давления и выхода газа имеет вид:
 󰇛󰇜󰇟󰇛󰇜󰇠
󰇛󰇜 
Каждое слагаемое имеет размерность Па/с. Модель описывает форму
переходного процесса, не разделяя изменение объёма, температуры газа и
механическое взаимодействие оболочки с человеком.
Умножим уравнение на интегрирующий множитель exp(u/τс):
󰇟󰇛󰇜󰇠 󰇛󰇜󰇛󰇜
После интегрирования от 0 до u получаем:
󰇛󰇜 󰇟󰇛󰇜󰇠󰇛󰇜
Условие максимума dp/du = 0 приводит к exp(−u/τн) = τн/(τн + τс).
Следовательно:
 󰇛󰇜 
  
  
Коэффициент P₀ = 80 000 Па является масштабом модели и не равен
максимальному давлению. После максимума давление убывает и при u → ∞
стремится к нулю.
Численная модель
Для y′ = F(t,y) используется метод Рунге—Кутты четвёртого порядка:
 󰇛󰇜 󰇛󰇜
 󰇛󰇜 󰇛󰇜
 󰇛󰇜
Для разгона y = (s,v), для давления y = p. Начало движения, начало
удара и команда SRS заданы точно; через переключения интегрирование
одним неразделённым шагом не выполняется.
14
2.6 Результаты разгона из состояния покоя
Расчёт разгона проведён на интервале 0–6 с. В ТР тяга ниже и
автомобиль легче; в РР перевозится сервисное оборудование. Начало
движения: tп = 0,029465 с в ТР и 0,017259 с в РР.
Рисунок 6 — Скорость и путь в тестовом и рабочем режимах
Таблица 10 — Контрольные значения движения
t, с
v ТР, м/с
s ТР, м
v РР, м/с
s РР, м
0
0,0000
0,0000
0,0000
0,0000
1
0,8469
0,3387
1,5084
0,6116
2
1,9890
1,7570
3,5154
3,1237
4
4,2172
7,9779
7,4495
14,1103
6
6,3582
18,5676
11,2549
32,8359
За 6 с автомобиль в ТР достигает 6,358168 м/с и проходит 18,567580 м; в РР —
11,254914 м/с и 32,835885 м. Увеличение момента в РР перевешивает увеличение массы.
При холостом ходе s = v = 0 во всём интервале.
15
2.7 Аварийный режим и реакция SRS
До tс = 6 с АР совпадает с РР. Для θ = t − 6 до остановки v = v₆ − 120θ и
s = s₆ + v₆θ − 60θ². Затем v = 0, s сохраняется. Продолжительность ударного
импульса определяется скоростью перед ударом:
  
  
Рисунок 7 — Ударный импульс и отклик SRS
Перемещение за ударный участок равно v₆²/(2·120) = 0,527805 м. Это величина
упрощённой ударной модели, а не расчёт тормозного пути или деформации реального
автомобиля. Распознавание происходит через 20 мс, начало отклика — через 23 мс,
максимум давления — через 41,218 мс после начала удара.
Предел сцепления ППП относится к разгону. Он не ограничивает заданную силу
столкновения. SRS в этой модели не меняет траекторию кузова; поэтому наличие
давления само по себе не доказывает безопасность для человека.
16
2.8 Точность и выводы
Таблица 11 — Ошибки RK4 относительно аналитического разгона
Режим
h, с
Ошибка s, м
Ошибка v, м/с
TR
0,020
2,630e-08
3,656e-09
TR
0,010
1,644e-09
2,284e-10
TR
0,005
1,022e-10
1,417e-11
RR
0,020
4,662e-08
6,541e-09
RR
0,010
2,913e-09
4,088e-10
RR
0,005
1,808e-10
2,536e-11
Таблица 12 — Максимальная ошибка давления за 0,2 с
Шаг, с
Максимальная ошибка, Па
0,00100
4,508e-03
0,00050
2,821e-04
0,00025
1,764e-05
При уменьшении шага вдвое ошибка давления уменьшается примерно
в 16 раз, что согласуется с четвёртым порядком RK4. При h = 0,0005 с
максимальная абсолютная ошибка составляет 0,000282 Па, или около
5,10·10⁻⁹ от pмакс.
Проверка ППП при меньшем сцеплении μ = 0,30 в РР даёт Fпред =
6356,88 Н > 5294,118 Н. Тяга не ограничивается и на нижней границе
выбранного диапазона. Предельное значение для данного момента μкр =
5294,118 / (0,9·2400·9,81) = 0,249845. При меньшем μ модель уменьшает
доступное задание тяги.
При D = 0 моделируется диагностическая неисправность без команды.
В холостом, тестовом и рабочем режимах давление остаётся нулевым. В АР
давление имеет один максимум и затем убывает; повторной команды в
пределах события нет.
Для ТСА модель показывает различие между проверкой трансмиссии и
диагностикой удерживающей системы. ППП определяет тяговые
возможности автомобиля, а SRS — реакцию пассивной защиты на событие.
Для применения к конкретной машине нужны её характеристики и заводская
документация; учебные пороги не используются как настройки.
17
3 Оптимизация загрузки оборудования
3.1 Собственная постановка для профиля ТСА
Выбран пункт 3.3 методички. Ремонтная мастерская восстанавливает
валы, пальцы шарнирных соединений и втулки сельскохозяйственной
техники. Составим две самостоятельные постановки: выбор объёма выпуска
на максимум условной прибыли и распределение заданного заказа на
минимум затрат.
Рисунок 8 — Две постановки задачи загрузки оборудования
В постановке А каждая деталь последовательно проходит токарную
операцию, шлифование и контроль. Фонды времени уже уменьшены на
заранее предусмотренные простои; переменные x₁, x₂, x₃ обозначают число
восстановленных валов, пальцев и втулок.
В постановке Б рассматривается отдельный заказ на 40 валов и 30
втулок. Одна и та же операция может выполняться на универсальном станке
либо на ЧПУ. Разделение между станками допускается по деталям, но не
внутри одной детали.
Все исходные числа сформированы самостоятельно. Постановки А и Б
относятся к разным плановым задачам; их фонды времени и финансовые
результаты не складываются. Прибыль и затраты заданы в рублях, времена
— в секундах, объёмы — в штуках.
18
3.2 Максимизация прибыли
Таблица 13 — Исходные данные по выпуску деталей
Показатель
Вал x₁
Палец x₂
Втулка x₃
Фонд
Токарная операция, с/шт.
1200
900
600
57 600 с
Шлифование, с/шт.
600
1080
480
44 640 с
Контроль, с/шт.
300
240
180
15 240 с
Прибыль, руб./шт.
1380
1494
852
—
Предельный спрос, шт.
30
26
24
—
Условная прибыль на деталь задана после вычета переменных затрат;
постоянные затраты считаются неизменными. Целевая функция:
 
Система ограничений по времени и спросу:
 
 
 
      

Размерности: (с/шт.)·шт. = с в ограничениях и (руб./шт.)·шт. = руб. в
функции прибыли. Сначала решается непрерывная задача линейного
программирования, затем проверяется целочисленность. В данном случае
оптимальная вершина целая.
Разделив активные ограничения на 300, 120 и 60, получаем 4x₁ + 3x₂ +
2x₃ = 192; 5x₁ + 9x₂ + 4x₃ = 372; 5x₁ + 4x₂ + 3x₃ = 254. Вычитание уравнений
даёт:
  
   
   
 
Объёмы не превышают спрос. Независимый целочисленный перебор
всех 31·27·25 = 20 925 планов даёт тот же результат. Программа перебора
использует только стандартную библиотеку Python.
19
3.3 Доказательство оптимальности
Обозначим загрузки трёх постов через L₁, L₂ и L₃. Коэффициенты
прибыли точно представляются положительной линейной комбинацией норм
времени:
 
 
 
Следовательно, для любого допустимого плана:



План (24; 20; 18) полностью использует все три фонда и достигает этой
границы. Поэтому он глобально оптимален как для непрерывной, так и для
целочисленной задачи. Коэффициенты 0,5; 0,7; 1,2 имеют размерность руб./с.
Рисунок 9 — Загрузка при двух допустимых планах
Таблица 14 — Сравнение выпуска
Показатель
Исходный план
Оптимальный
Валы; пальцы; втулки
20; 16; 14
24; 20; 18
Токарная загрузка, с
46 800
57 600
Шлифование, с
36 000
44 640
Контроль, с
12 360
15 240
Прибыль, руб.
63 432
78 336
Прирост прибыли равен 14 904 руб., или 23,50 %. Все три ресурса становятся
полностью загруженными. Увеличение только одного фонда не гарантирует такого же
прироста без проверки остальных ограничений.
20
3.4 Минимизация затрат
Для отдельного заказа требуется восстановить 40 валов и 30 втулок.
Примем u₁, u₂ — количества на универсальном станке; z₁, z₂ — количества на
ЧПУ.
Таблица 15 — Данные по альтернативным станкам
Показатель
Универсальный
ЧПУ
Время одного вала, с
1500
900
Время одной втулки, с
1020
660
Затраты на вал, руб.
680
440
Затраты на втулку, руб.
460
280
Доступный фонд, с
72 000
45 000
Целевая функция и ограничения:
 
  
 
 
 
Исключим u₁ = 40 − z₁ и u₂ = 30 − z₂. Получим эквивалентную задачу
двух переменных:
 
 
 
    
Минимизация затрат равносильна максимизации экономии S = 240z₁ +
180z₂. Значение 41 000 руб. — алгебраическая стоимость обработки всего
заказа на универсальном станке; сам такой план недопустим по его фонду
времени.
На одну секунду ЧПУ вал даёт экономию 240/900 = 0,266667 руб./с,
втулка — 180/660 = 0,272727 руб./с. В данной задаче выгодно сначала
загрузить ЧПУ втулками, затем использовать оставшееся время для валов.
21
3.5 Графическое решение
Построим допустимую область в координатах z₁, z₂. Прямые
ограничений определяют выпуклый многоугольник. Линия цели S = const
перемещается в направлении увеличения экономии.
Рисунок 10 — Допустимая область и оптимальное распределение
Таблица 16 — Затраты в вершинах непрерывной области
z₁, шт.
z₂, шт.
12,4000
0,0000
40,0000
0,0000
40,0000
13,6364
28,0000
30,0000
0,0000
30,0000
0,0000
18,2353
Минимум достигается на пересечении границы спроса по втулкам z₂ =
30 и полного фонда ЧПУ:
 
 󰇛󰇜 
Оптимальная вершина (28; 30) целая. Дробные вершины таблицы
нужны для анализа непрерывной задачи и не предлагаются как физические
производственные планы.
22
3.6 Проверка затрат и рекомендации
На универсальном станке остаются u₁ = 40 − 28 = 12 валов и u₂ = 30 −
30 = 0 втулок. Проверим выполнение заказа, загрузку и затраты:
  
   
  


Для доказательства глобального минимума оценим максимальную
экономию:
󰇛󰇜󰇛󰇜
󰇛󰇜 
 
Найденный план достигает нижней границы. Перебор 41·31 = 1271
целого распределения и независимое решение линейной программы
подтверждают этот результат.
Таблица 17 — Сравнение двух допустимых распределений
Показатель
Исходный план
Оптимум
Универсальный: валы; втулки
20; 10
12; 0
ЧПУ: валы; втулки
20; 20
28; 30
Универсальный: загрузка, с
40 200
18 000
ЧПУ: загрузка, с
31 200
45 000
Затраты, руб.
32 600
28 880
Резерв универсального, с
31 800
54 000
Экономия составляет 3720 руб., или 11,41 %. Рекомендуемый план: все
30 втулок и 28 валов на ЧПУ, остальные 12 валов — на универсальном
станке. ЧПУ полностью загружен; резерв универсального станка равен 15 ч.
В реальном планировании следует уточнить нормы времени, стоимость
инструмента и переналадок. При изменении исходных данных задачу
требуется решить повторно.
23
Заключение
Выполнены три задания варианта 9 для профиля «Технический сервис
в АПК». В задании 2 по индивидуальному условию принят постоянный
полный привод ППП; выбранная по варианту система SRS сохранена.
В цепи по схеме 2 получены токи I₁ = −0,231531443 А, I₂ =
−0,423820269 А, I₃ = 0,655351712 А, I₄ = I₁. Токи верхней и средней ветвей
направлены от B к A, нижней — от A к B. Решения Гаусса, Зейделя и Якоби
согласуются; методы Зейделя и Якоби потребовали 10 и 19 итераций
соответственно. Баланс мощности составляет 22,823506 Вт.
В динамической модели тяга передаётся обеим осям в долях 50:50 с
проверкой сцепления. В рабочем режиме автомобиль начинает движение из
покоя и за 6 с достигает 11,254914 м/с, пройдя 32,835885 м. Полный привод
отражён в распределении момента и ограничении доступной силы тяги.
Для заданного аварийного события модель SRS формирует команду
при t = 6,023 с. Давление достигает 55,343 кПа при t = 6,041218 с и затем
убывает. Ударный участок и давление рассчитаны раздельно: удерживающая
система не используется как источник тормозной силы кузова.
Аналитическое решение и сгущение сетки подтверждают корректность
численного расчёта.
В задаче максимизации прибыли получен выпуск 24 валов, 20 пальцев
и 18 втулок. Прибыль равна 78 336 руб.; все три фонда использованы
полностью. По сравнению с исходным допустимым планом прирост
составляет 14 904 руб.
В задаче минимизации затрат заказ распределён следующим образом:
на ЧПУ — 28 валов и 30 втулок, на универсальный станок — 12 валов.
Минимум затрат равен 28 880 руб.; экономия относительно сравниваемого
плана — 3720 руб. Оптимальность обоих планов доказана аналитическими
оценками и проверена перебором.
Подпись студента: __________________________
Дата выполнения: ___________________________
Затраченное время: _______________________ ч
24
Список использованных источников
1. Храмешин А. В. Моделирование в агроинженерии: учебное пособие.
Ижевск: Ижевская ГСХА, 2019. 50 с. Раздел 4, таблицы 2–4; пункт 3.3.
Предоставленный файл Modelirovanie_v_agroinzhenerii.doc.
2. Bosch Engineering. Passive safety systems. Официальное описание
удерживающих систем, датчиков и индивидуальной калибровки.
Электронный ресурс:
https://www.bosch-engineering.com/services/mobility-solutions/body-and-comfort/passive-safety-systems/
Дата обращения: 30.09.2026.
Пояснения к исходным данным
По пособию определены электрическая схема, численные параметры
цепи и система SRS. Тип привода ППП задан индивидуально. Источник 2
использован для описания назначения SRS; численные пороги и
характеристики учебной динамической модели из него не заимствованы.
Данные оптимизации составлены самостоятельно для технического
сервиса сельскохозяйственной техники. Постановка А описывает
последовательные операции, постановка Б — альтернативные способы
выполнения отдельного заказа.
Программы расчёта
Таблица 18 — Назначение полных программ
Файл
Вход
gauss.py
A, b
seidel.py
A, b, допуски, jacobi
srs.py
Режим, шаг, параметры
equipment.py
Нормы, фонды, затраты
В приложениях А–Г приведены полные исполняемые тексты на Python.
Программа srs.py занимает две страницы; обе части объединяются в один
файл без добавления текста заголовков. Для запуска достаточно стандартной
библиотеки Python 3.10 или новее.
Числа в программах записаны с десятичной точкой. В расчётах
сохраняются неокруглённые значения. Программы не зависят от файла Word
и воспроизводят представленные результаты.
25
Приложение А Программа метода Гаусса
Сохранить текст в gauss.py. Матрица и правая часть соответствуют
схеме 2 варианта 9. На каждом шаге выбирается максимальный по модулю
опорный элемент.
A = [[87.2, -18.6], [-18.6, 55.4]]
b = [-8., -32.]
def solve(A, b):
n = len(b)
a = [list(row) + [b[i]] for i, row in enumerate(A)]
for k in range(n):
q = max(range(k, n), key=lambda i: abs(a[i][k]))
if abs(a[q][k]) < 1e-12:
raise ValueError('Singular matrix')
a[k], a[q] = a[q], a[k]
for i in range(k+1, n):
w = a[i][k] / a[k][k]
for j in range(k, n+1):
a[i][j] -= w*a[k][j]
x = [0.]*n
for i in range(n-1, -1, -1):
x[i] = (a[i][n] - sum(a[i][j]*x[j]
for j in range(i+1, n))) / a[i][i]
return x
if __name__ == '__main__':
x = solve(A, b)
print('J =', x)
print('residual =', max(abs(sum(c*v for c,v in zip(r,x))-d)
for r,d in zip(A,b)))
Контрольный результат:
J = [-0.23153144314725793, -0.6553517119591877]
residual < 1e-13
Токи ветвей после запуска: I₁ = J[0], I₂ = J[1] − J[0], I₃ = −J[1], I₄ = I₁.
Копия расширенной матрицы позволяет сохранить исходные A и b для
проверки невязки.
Метод Гаусса является прямым. В этой работе он применяется
совместно с итерационными алгоритмами следующего приложения.
26
Приложение Б Программа итерационных методов
Файл seidel.py. По умолчанию используется Зейдель. Для Якоби вызвать solve(A, b,
jacobi=True). history хранит k, J₁, J₂, невязку и изменение решения.
A = [[87.2, -18.6], [-18.6, 55.4]]
b = [-8., -32.]
def solve(A, b, ex=1e-10, er=1e-9, limit=10000, jacobi=False):
x = [0.]*len(b)
history = [[0, *x, max(map(abs,b)), 0.]]
for k in range(1, limit+1):
old = x.copy()
for i in range(len(b)):
if abs(A[i][i]) < 1e-15:
raise ValueError('Zero diagonal')
left = old if jacobi else x
s = sum(A[i][j]*left[j] for j in range(i))
s += sum(A[i][j]*old[j] for j in range(i+1,len(b)))
x[i] = (b[i]-s) / A[i][i]
delta = max(abs(u-v) for u,v in zip(x,old))
rho = max(abs(sum(c*v for c,v in zip(row,x))-d)
for row,d in zip(A,b))
history.append([k,*x,rho,delta])
if delta <= ex and rho <= er:
return x, history
raise RuntimeError('No convergence')
if __name__ == '__main__':
x, hist = solve(A,b)
print('J =',x)
print('k, J1,J2, residual, delta =', hist[-1])
Контрольные показатели:
Seidel: k = 10; rho = 5.606804e-10 V
Jacobi: k = 19; rho = 5.322924e-10 V
В методе Якоби left содержит старый вектор; в методе Зейделя — уже
обновляемый. Это единственное различие формул шага; критерии остановки одинаковы.
27
Приложение В Программа модели SRS и ППП
Файл srs.py. Часть 1: параметры, пределы тяги обеих осей, шаг RK4 и разгон. TR —
тестовый режим, RR — рабочий, IDLE — неподвижный холостой ход.
from math import exp, log
G, F, C, R, TAU = 9.81, .015, 40., .34, .25
BETA, ALPHA, MU = .55, .50, .65
TC, DEC, DV, DELAY = 6., 120., 2.4, .003
AP = 30.
P0, TR, TD = 80000., .008, .070
def params(mode):
m, demand = (2000., 900.) if mode == 'TR' else (2400., 1800.)
if mode == 'IDLE': demand = 0.
cap = min(MU*BETA*m*G/ALPHA,
MU*(1-BETA)*m*G/(1-ALPHA))
force = min(demand/R, cap)
return m, force
def rk4(fun, t, y, h):
a = fun(t, y)
b = fun(t+h/2, [v+h*w/2 for v,w in zip(y,a)])
c = fun(t+h/2, [v+h*w/2 for v,w in zip(y,b)])
d = fun(t+h, [v+h*w for v,w in zip(y,c)])
return [v+h*(a0+2*b0+2*c0+d0)/6
for v,a0,b0,c0,d0 in zip(y,a,b,c,d)]
def drive(mode, h=.01, end=TC):
m, force = params(mode)
f0 = F*m*G
if force <= f0: return [[0.,0.,0.],[end,0.,0.]]
t0 = -TAU*log(1-f0/force)
def rhs(t,y):
return [y[1], (force*(1-exp(-t/TAU))-f0-C*y[1])/m]
t, y = t0, [0.,0.]
rows = [[0.,0.,0.], [t,*y]]
while t < end-1e-12:
dt = min(h,end-t)
y = rk4(rhs,t,y,dt); t += dt
rows.append([t,*y])
return rows
28
Приложение В Продолжение программы
Часть 2 того же файла srs.py: аналитический контроль, давление и аварийное
событие. Возвращаемое fire равно None при блокировке команды.
def exact_drive(mode,t=TC):
m, force = params(mode); f0 = F*m*G
if force <= f0: return 0.,0.
t0 = -TAU*log(1-f0/force)
if t <= t0: return 0.,0.
u, k, lam = t-t0, C/m, 1/TAU
a = (force-f0)/m
e, q = exp(-k*u), exp(-lam*u)
v = a*((1-e)/k-(e-q)/(lam-k))
s = a*(u/k-(1-e)/k**2
-((1-e)/k-(1-q)/lam)/(lam-k))
return s,v
def pressure(u):
return P0*(1-exp(-u/TR))*exp(-u/TD) if u>=0 else 0.
def pressure_numeric(h=.0005,end=.2):
def rhs(u,y):
return [P0/TR*exp(-(1/TR+1/TD)*u)-y[0]/TD]
u,y=0.,[0.]; rows=[[u,y[0]]]
while u<end-1e-12:
dt=min(h,end-u); y=rk4(rhs,u,y,dt); u+=dt
rows.append([u,y[0]])
return rows
def crash(h=.01,healthy=True):
_,s0,v0=drive('RR',h)[-1]
stop=TC+v0/DEC
ok = healthy and DEC>=AP and v0>=DV
fire=TC+DV/DEC+DELAY if ok else None
return s0,v0,stop,s0+v0*v0/(2*DEC),fire
if __name__=='__main__':
for mode in ['IDLE','TR','RR']:
z=drive(mode)[-1]; ref=exact_drive(mode)
print(mode,z,'exact:',ref)
print('AR: s6,v6,tstop,sstop,tfire =',crash())
peak=TR*log(1+TD/TR)
print('peak after fire, pressure:',peak,pressure(peak))
for h in [.001,.0005,.00025]:
err=max(abs(p-pressure(u)) for u,p in pressure_numeric(h))
print('pressure h,error:',h,err)
29
Приложение Г Программа оптимизации
Файл equipment.py. Первая функция перебирает планы выпуска, вторая — способы
распределить фиксированный заказ между двумя станками.
from itertools import product
TIME = [[1200,900,600],[600,1080,480],[300,240,180]]
CAP = [57600,44640,15240]
PROFIT = [1380,1494,852]
DEMAND = [30,26,24]
def maximize():
best = (-1,None)
for x in product(*(range(n+1) for n in DEMAND)):
used=[sum(a*b for a,b in zip(row,x)) for row in TIME]
if all(u<=c for u,c in zip(used,CAP)):
value=sum(p*n for p,n in zip(PROFIT,x))
if value>best[0]: best=(value,x)
return best
def minimize():
best=(float('inf'),None)
for z1 in range(41):
for z2 in range(31):
u1,u2=40-z1,30-z2
if 1500*u1+1020*u2>72000: continue
if 900*z1+660*z2>45000: continue
cost=680*u1+460*u2+440*z1+280*z2
if cost<best[0]: best=(cost,(u1,u2,z1,z2))
return best
if __name__=='__main__':
print('Maximum:',maximize())
print('Minimum:',minimize())
Контрольный результат:
Maximum: (78336, (24, 20, 18))
Minimum: (28880, (12, 0, 28, 30))
В первом результате компоненты — валы, пальцы и втулки. Во втором — валы и
втулки на универсальном станке, затем валы и втулки на ЧПУ.
Команды запуска:
python gauss.py
python seidel.py
python srs.py
python equipment.py