Чисельні методи у механіці

(1.67)

де y_i - наближене значення, отримане в точці із кроком h; y_2i - із кроком h/2; p - порядок методу; y(x_2i) - точний розв'язок задачі.

Формула Рунге для методу четвертого порядку

. (1.68)

1.8.2 Обчислювальна схема (алгоритм) методу Рунге-Кутта

Вибираємо початковий крок h на відрізку [a,b], задаємо точність .

Створюємо множину рівновіддалених точок (вузлів)

Знаходимо розв'язок y_i+1 за формулами при кроку h і при кроку h/2, 0 ? i ? n-1.

Перевіряємо нерівність

.

Якщо ця нерівність виконується, то беремо і продовжуємо обчислення з тим самим кроком, якщо ні, то зменшуємо початковий крок h у два рази і переходимо до пункту 3.

1.8.3 Метод прогнозу і корекції

Підправивши схему Эйлера, одержимо схему прогнозу

, (1.69)

де наближене значення . Цю формулу використовувати не можна, оскільки схема прогнозу нестійка. Тому використовуємо схему корекції

(1.70)

Оцінюючи похибки прогнозу і корекції, одержимо

- похибка корекції; (1.71)

- похибка прогнозу. (1.72)

На будь-якому кроці можна оцінити точність розв'язку. При заданому =0,0000001, наприклад, . Віднімаючи з (1.71) співвідношення (1.70), маємо

.

Уточнюємо розв'язок, виходячи з формули (1.70):

. (1.73)

Ця формула завершує схеми прогнозу і корекції.

1.8.4 Задача Коші для диференціальних рівнянь вищих порядків

Необхідно знайти функцію , що задовольняє диференціальне рівняння та додаткові умови

(1.74)

Для розв'язання цієї задачі можна застосувати таку схему зведення до системи диференціальних рівнянь першого порядку:

1) вводимо нові змінні:

(1.75)

2) розв'язуємо систему з диференціальних рівнянь першого порядку.

1.8.5 Метод скінченних різниць для розв'язання лінійних крайових задач

Маємо відрізок [a,b]. Потрібно знайти розв'язок лінійного диференціального рівняння другого порядку

, (1.76)

що задовольняє такі крайові умови:

(1.77)

Виберемо рівномірну сітку: x = a + ih, i = =0,1,2,…,n.

Нехай Апроксимуємо і у кожному внутрішньому вузлі (i = 1, 2, …, n-1) центральними різницями

,

і на кінцях відрізка - односторонніми скінченно-різницевими апроксимаціями

, .

Використовуючи ці формули, одержуємо різницеву апроксимацію вихідного крайового завдання:

(1.78)

Коефіцієнти різницевих рівнянь залежать від кроку сітки. Введемо позначення:

Перепишемо систему з урахуванням введених позначень

(1.79)

Маємо різницеву схему крайового завдання. Напишемо систему рівнянь у розгорнутій матричній формі:

Таким чином, завдання зводиться до розв'язання системи лінійних алгебраїчних рівнянь, що можна записати у вигляді Ay=d.

1.8.6 Загальна характеристика явних методів

Однокрокові і багатокрокові методи разом складають так звані явні методи чисельного інтегрування диференціальних рівнянь. Назву їх обумовлено тим, що усі вони з метою відшукання значень фазових змінних (змінних стану) на наступному кроці використовують інформацію лише про попередні точки.

Нагадаємо головні властивості кожної з двох груп явних методів.

Однокрокові методи.

1. Через те, що у методах Рунге-Кутта використовується інформація лише про останню обчислену точку інтегральної кривої, за допомогою цих методів можна починати інтегрування диференціальних рівнянь у випадку задачі Коші.

2. З тієї ж причини ці методи дозволяють у процесі інтегрування легко змінювати його крок.

3. Однак з тієї ж причини при використанні цих методів доводиться багаторазово обчислювати значення функцій і витрачати на це багато машинного часу.

4. При застосуванні цих методів важко одержати оцінку виникаючої при інтегруванні похибки обмеження (для цього потрібні додаткові обчислення приблизно того ж обсягу, що й саме інтегрування).

Методи прогнозу-корекції:

1. Через те, що у цих методах використовується інформація про декілька попередні точки інтегральної кривої, розділені за аргументом однаковими інтервалами, за їх допомогою неможливо розпочати процес інтегрування у випадку задачі Коші.

2. Оскільки у цих методах замість обчислення у попередніх точках використовуються значення, одержані при попередніх обчисленнях і збережені у пам'яті ЕОМ, вони є більш економічними у сенсі витрат машинного часу.

3. При будь-якій зміні величини кроку інтегрування доводиться тимчасово звертатися до методів Рунге-Кутта.

4. У цих методах легко одержати позитивну оцінку похибки обмеження.

Загальним для явних методів є така, ще не обговорювана обставина: явні методи при кроках інтегрування, більших за деякий крок , стають нестійкими у обчислювальному відношенні. Тобто якщо

, (1.80)

похибка обмеження методу починає катастрофічно збільшуватися з кожним кроком, внаслідок чого одержувані результати інтегрування перестають хоч якоюсь мірою відображати якісно і кількісно дійсну інтегральну криву.

Величина залежить від властивостей самої системи диференціальних рівнянь, яку потрібно чисельно інтегрувати, і тому не може бути визначена заздалегідь, якщо не відома ця система диференціальних рівнянь. Для систем лінійних диференціальних рівнянь з постійними коефіцієнтами цю величину можна наближено обчислити за формулою

, (1.81)

де означає мінімальну за величиною сталу часу системи диференціальних рівнянь, що інтегрується; - деяке число, більше за одиницю, різне для різних методів інтегрування (наприклад, для методу Рунге-Кутта 4-го порядку =5). У свою чергу, сталі часу можна визначити як величини, зворотні до значень дійсних та уявних частин (точніше - їхніх модулів) коренів характеристичного рівняння системи диференціальних рівнянь.

Таким чином, щоб організувати процес чисельного інтегрування якоїсь системи диференціальних рівнянь одним із явних методів, потрібно здійснити такі дії:

визначити (обчислити за заданими значеннями параметрів диференціальних рівнянь) сталі часу заданої системи диференціальних рівнянь;

відкинувши сталі часу, що дорівнюють нескінченності (які відповідають нульовим кореням характеристичного рівняння), знайти серед них максимальну за величиною і мінімальну ;

задатися кількістю (>1) кроків інтегрування на інтервалі часу, що дорівнює мінімальній сталій часу; при цьому необхідно, щоб було більше за сталу для заданого методу; при цьому необхідно пам'ятати також, що чим більшу за величиною ми задаємо, тим з більшою точністю здійснюватиметься процес чисельного інтегрування, але й тим більший час він займе;

розрахувати величину кроку інтегрування за формулою

; (1.82)

задатися бажаною тривалістю процесу, що виходить внаслідок інтегрування, у вигляді кількості максимальних сталих часу системи; розрахувати бажану тривалість процесу за формулою

; (1.83)

здійснити присвоювання змінним стану початкових значень;

організувати цикл за часом , який би змінювався від 0 до з кроком , усередині якого за заданими значеннями вектору у попередній точці обчислюється значення того ж вектору у наступній точці, а також на основі цих розрахунків формуються масиви, які у подальшому будуть використовуватися як результати інтегрування.

1.8.7 Жорсткі системи. Неявні методи

В інженерній практиці трапляються випадки, коли відношення максимальної сталої часу системи до мінімальної є дуже великою величиною (наприклад, більше, ніж ). Такі системи диференціальних рівнянь називаються жорсткими.

Чисельне інтегрування жорстких систем за допомогою явних методів натикається на суттєві труднощі. У цьому випадку загальна кількість кроків інтегрування, потрібна для повного інтегрування

, (1.84)

буде теж дуже великою. А через те, що похибки округлення накопичуються зі збільшенням кількості обчислень приблизно пропорційно цій кількості, на практиці вже після мільйона операцій підсумовування-множення похибки округлення можуть невпізнанно спотворити інтегральну криву. Якщо врахувати, що кожний крок інтегрування може вміщувати тисячі операцій підсумовування-множення, то вже через тисячі кроків інтегрування відносна похибка внаслідок округлення збільшується до 100 %. Подальше інтегрування стає недоцільним.

Подолати цю перешкоду можна тільки у випадку, якщо збільшити крок інтегрування до величини, у десятки разів більшої за мінімальну сталу часу системи. Для розглянутих явних методів це не є можливим, оскільки за цих умов ці методи є нестійкими.

Тому для інтегрування жорстких систем використовують так звані неявні методи.

Найпростішим шляхом розв'язування "жорстких" задач є так званий неявний метод Ейлера, за яким інтегрування на одному кроці здійснюється за неявною формулою

. (1.85)

На відміну від явних методів цей метод є стійким при величинах кроку, набагато більших за мінімальну сталу часу, хоча точність його є помірною (порядку 10 %).

2. Методи оптимізації технічних систем

2.1 Оптимальний і раціональний розв'язок

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

Розв'язок, який задовольняє всі задані обмеження, називається допустимим. З допустимих у процесі розв'язання екстремальних задач вибираються оптимальні, або раціональні розв'язки.

Під оптимальним розумітимемо такий допустимий розв'язок, який є якнайкращим з погляду вибраного критерію оптимальності. Критерій оптимальності формується на основі одного або декількох критеріїв ефективності. Критерієм ефективності може бути: маса, вартість, приведені витрати, термін виготовлення і т. п. Але не всі критерії ефективності і обмеження, що виражають вимоги до проектованих виробів, можуть бути формалізовані, наприклад, естетичність, компактність, технологічність виготовлення та ін. У цьому випадку доцільно ввести поняття раціонального розв'язку.

Під раціональним розумітимемо такий розв'язок, який виходить неформальним шляхом, тобто з урахуванням експертних (або інших неформальних) оцінок. При пошуку оптимальної конструкції за прийнятим критерієм оптимальності і обмеження можна одержати оптимальну, але нераціональну конструкцію. Наприклад, якщо вимога технологічності не знайшла кількісної оцінки і не увійшла до числа обмежень, то можна одержати проект мінімальної ваги (оптимальний за вагою), але нераціональний проект через нетехнологічність. Важча конструкція, але зручна з технологічної точки зору може виявитися раціональнішою. Стосовно механічних конструкцій, що деформуються, раціональне проектування може бути визначене як проектування, що ґрунтується засноване на принципах механіки деформованого твердого тіла і має на меті отримання оптимальної конструкції на базі вибраного проектантом критерію. Оптимізаційні задачі (і їх різні поставлення) є окремими підзадачами у процесі розроблення в цілому раціонального проекту [5,6].

2.2 Розрахункова модель

Задача проектування силової конструкції починається з вибору розрахункової моделі, яка виходить шляхом ідеалізації реальної конструкції.

При цьому необхідно розробити або прийняти ряд допоміжних моделей:

1. Моделі форми. Моделі форми елементів конструкцій є схематичним описом геометрії елемента за допомогою стандартних, типових елементів:

стрижень - тіло, поперечні розміри h (висота) і b (ширина), якого малі в порівнянні з його довжиною . Для стрижня характерне співвідношення при ; якщо , то стрижень починає працювати як пластина;

пластина - тіло, форма якого визначається серединною площиною, причому товщина набагато менша від двох інших габаритних розмірів , тобто , ;

оболонка - тіло, обмежене двома близькими криволінійними поверхнями з відстанню між ними (товщина ), значно меншого за інші габаритні розміри ;

просторове тіло - тіло довільної форми, всі габаритні розміри якого сумірні.

2. Конструктивні схеми. Класифікація конструктивних схем може бути здійснена за різними ознаками. Наприклад, залежно від використовуваних типових елементів конструктивні розрахункові схеми можуть бути: просторовими, оболонковими, пластинчастими, стрижневими і комбінованими.

До конструктивної схеми просторового тіла можна навести конструкції фундаменту, корпусної деталі машини, складного вузла з'єднання окремих елементів і т. д.

До конструктивної оболонкової схеми можуть бути віднесені конструкції купольного покриття будівлі, відсіку літака, судна і т. д.

До конструктивної пластинчастої схеми можуть бути наведені окремі підконструкції палуби судна, плоскої кришки посудини і т. д.

Прикладами конструктивних стрижневих схем є ферми мостів, зуб шестерні, корпус ракети і т. д.

До комбінованих конструктивних схем відносять схеми, що включають різні типи елементів. Наприклад, конструкція літального апарата є комбінованою (складовою) конструкцією, що складається з оболонок, стрижнів (стрингерів і шпангоутів), пластин.

3. В'язі. Всі тіла, як правило, взаємодіють із зовнішнім середовищем або іншими тілами. Якщо внаслідок яких-небудь обмежень (умов) дане тіло не може зайняти довільне положення в просторі і має довільні швидкості, то таке тіло називається невільним. У цьому випадку обмеження (умови) називають в'язями.

4. Моделі навантаження. Зовнішні сили, що діють на конструкцію, поділяються на три групи:

зосереджені сили - сили, які діють на невеликих ділянках поверхні деталі (наприклад, тиск кульки шарикопідшипника на вал);

розподілені сили - сили, прикладені до значних ділянок поверхні (наприклад, тиск рідини або газу на стінки посудини);

об'ємні або масові сили - сили, прикладені до кожної частинки матеріалу (наприклад, сили тяжіння, інерційні навантаження).

Навантаження поділяють на стаціонарні (сталі, статичні) і нестаціонарні (змінні). Вони можуть мати випадковий, або детермінований, характер.

Моделі навантаження повинні враховувати дію полів і середовищ: температурного поля, електромагнітного поля, корозійних середовищ та ін. Зокрема, температурні дії можуть викликати зміни геометричної форми тіла і зміну сил взаємодії між його точками, зміну фізико-механічних властивостей матеріалу.

5. Механічні властивості матеріалів. У розрахункову модель необхідно закладати такі механічні властивості матеріалів, які відповідають вимогам експлуатації. До них відносяться: модуль пружності, модуль зсуву, коефіцієнт Пуассона, межа текучості, межа міцності.

6. Урахування статистичної природи початкових даних у моделі. Розрахункова модель повинна бути достатньо простою, яка допускає побудову математичної моделі і в той же час достатньо адекватно відображає стан реальної конструкції. Найчастіше в моделі передбачається детермінованість форми і структури об'єкта, який проектується, геометричних розмірів елементів, зовнішніх дій і властивостей матеріалів. Насправді всі початкові дані перебувають під впливом великої кількості неврахованих факторів і тому тією чи іншою мірою мають випадковий характер. Наприклад, змінність геометричних розмірів і форми елементів конструкції виконує істотну роль при втраті стійкості форм рівноваги конструкції. Майже всі зовнішні навантаження є випадковими. Випадковий характер мають і механічні характеристики матеріалу. Випадковий характер початкових даних у розрахунковій моделі може бути врахований при виборі коефіцієнтів запасу згідно з діючими нормами проектування або застосування статистичних методів аналізу стану конструкції з урахуванням статистичних початкових даних, що входять до розрахункової моделі.

2.3 Поставлення задач параметричної оптимізації

Кінцевою метою проектування є створення раціональних виробів, конструкцій і т. п., виходячи з наявних ресурсів і можливостей. Щоб добитися якнайкращого результату, необхідно оптимізувати на всіх етапах проектування.

У процесі проектування виникає різноманітність різних поставлень задач оптимізації:

- вибір оптимальних геометричних форм;

- вибір оптимальних структур конструкцій;

- оптимальний розподіл внутрішніх зусиль за рахунок попередніх напружень;

- підбір матеріалів;

- створення конструкцій мінімальної ваги із заданою надійністю;

- оптимальне армування конструкцій і т. д.

При проектуванні розрізняють задачі визначення структури і визначення значень внутрішніх параметрів (параметрів елемента). Якщо серед варіантів структури відшукується найкращий, то таку задачу називають структурною оптимізацією. Розрахунок внутрішніх параметрів, оптимальних з позиції деякого критерію при заданій структурі об'єкта, називають параметричною оптимізацією.

2.4 Критерій ефективності. Цільова функція

Для оцінки проектних рішень у процесі їх пошуку використовуються різні показники. Ними можуть бути показники ваги, вартості, надійності, технологічності, естетичності і т. п. Показники можна розбити на два класи: кількісні і якісні.

Кількісні показники дозволяють виконати оцінку проектного рішення кількісно на основі аналізу математичної моделі ефективності. Наприклад, наявність функції ваги виробу дозволяє кількісно порівнювати варіанти, які розробляються.

Якісні показники не дозволяють виконати оцінку проектного рішення кількісно, оскільки немає математичній моделі ефективності. Наприклад, оцінка механічної системи без розробленої математичної моделі технології виготовлення здійснюється якісно експертами -технологами.

Під критерієм ефективності розумітимемо показник, який дозволяє кількісно виконувати оцінку ефективності проектних рішень.

2.5 Обмеження

Працездатною (допустимою) вважається така система, для якої виконуватимуться всі умови, записані у формі обмежень типу рівностей або нерівностей:

. (2.1)

Наведемо основні види обмежень на поведінку пружних технічних систем.

1. Обмеження на міцність. Вимогу виконання умови міцності можна сформулювати в такій розмірній формі:

, (2.2)

де - число підконструкцій у заданій конструкції; - номер варіанта дії на систему; - максимальне еквівалентне напруження в об'ємі i-ї підконструкції, яке визначається за прийнятою гіпотезою або теорією міцності при j-му варіанті дії; - допустиме напруження для матеріалу i-ї підконструкції при j-му варіанті дії. Якщо ввести функцію обмежень

,

то обмеження (2.2) можна привести до канонічної форми (2.1):

. (2.3)

2. Обмеження на жорсткість. Обмеження на пружні переміщення можна записати у формі

, (2.4)

де - максимальне узагальнене переміщення поверхні i-ї підконструкції при j-му варіанті дії; - допустиме узагальнене переміщення. Під узагальненим переміщенням розуміємо різні переміщення (лінійні, кутові, відносні, абсолютні).

Обмеження (2.4) приводиться до канонічного вигляду

, (2.5)

якщо взяти функцію обмежень за жорсткістю у вигляді

.

3. Умова місцевої стійкості. Під місцевою втратою стійкості розуміємо явище втрати стійкості на локальній ділянці системи, що деформується, типу випинання. Умову місцевої стійкості можна записати в розмірній формі

, (2.6)

де - параметр зовнішньої j-ї дії на i-ту підконструкцію; - мінімальне значення критичного параметра j-ї зовнішньої дії на i-ту підконструкцію.

Обмеження (2.6) можна записати в канонічній формі

, (2.7)

якщо взяти функцію обмежень на місцеву стійкість у вигляді

.

4. Обмеження на частоти власних коливань. Частоти власних коливань системи є одними з найважливіших динамічних характеристик системи. Обмеження на частоти власних коливань мають вигляд

, (2.8)

де - нижча частота власних j-х коливань пружної системи; - допустиме найменше значення частоти власних j-х коливань системи, що призначається як розрахункове значення частоти вимушених j-х коливань.

Умова (2.8) в канонічній формі має вигляд

, (2.9)

.

5. Обмеження на параметри. При пошуку оптимальних рішень в рамках прийнятих фізичних і математичних моделей повинне виконуватися певне співвідношення керованих (ті, що змінні у процесі оптимізації) і некерованих (ті, що незмінні у процесі оптимізації) параметрів у формі ( - задане число)

. (2.10)

Обмеження на геометричні параметри (2.10) можна записати в канонічній формі

, (2.11)

.

2.6 Задачі оптимізації

У процесі розв'язання задачі оптимізації звичайно необхідно знайти оптимальні значення деяких параметрів, що визначають дану задачу [5,6]. При розв'язанні інженерних задач їх прийнято називати проектними параметрами. Як проектні параметри можуть бути, зокрема, значення лінійних розмірів об'єкта, маси, температури і т. п. Число проектних параметрів характеризує розмірність (і ступінь складності) задачі оптимізації.

Вибір оптимального рішення або порівняння двох альтернативних рішень проводиться за допомогою деякої залежної величини (функції), яка визначається проектними параметрами. Ця величина називається цільовою функцією (або критерієм якості). У процесі розв'язання задачі оптимізації повинні бути знайдені такі значення проектних параметрів, при яких цільова функція має мінімум (або максимум). Таким чином, цільова функція - це глобальний критерій оптимальності в математичних моделях, за допомогою яких описуються інженерні задачі.

Цільову функцію можна записати у вигляді

(2.12)

Прикладами цільової функції, що трапляються в інженерних розрахунках, є міцність або маса конструкції, потужність установки і т. п.

У разі одного проектного параметра () цільова функція (2.12) є функцією однієї змінної, і її графік - деяка крива на площині. При цільова функція є функцією двох змінних, і її графіком є поверхня.

Необхідно зазначити, що цільова функція не завжди може бути представлена у вигляді формули. Іноді вона може набувати тільки деякі дискретні значення, задаватися у вигляді таблиці і т. п. У всіх випадках вона повинна бути однозначною функцією проектних параметрів.

Цільових функцій може бути декілька. Наприклад, при проектуванні виробів машинобудування одночасно вимагається забезпечити максимальну надійність, мінімальну матеріаломісткість, максимальний корисний об'єм (або вантажопідйомність). Деякі цільові функції можуть виявитися несумісними. У таких випадках необхідно вводити пріоритет тієї або іншої цільової функції.

Виділяють два типи задач оптимізації: безумовні і умовні. Безумовна задача оптимізації полягає у відшуканні максимуму або мінімуму дійсної функції (2.12) від дійсних змінних і визначенні відповідних значень аргументів на деякій множині -вимірного простору.

Умовні задачі оптимізації, або задачі з обмеженнями, - це такі задачі, при формулюванні яких задаються деякі умови (обмеження) на множині . Ці обмеження задаються сукупністю деяких функцій, що задовольняють рівняння або нерівності.

Обмеження-рівності виражають залежність між проектними параметрами, яка повинна враховуватися при знаходженні рішення. Ці обмеження відображають закони природи, наявність ресурсів і т. п.

У результаті обмежень область проектування , яка визначається всіма проектними параметрами, може бути істотно зменшена відповідно до фізичної суті задачі. Число обмежень-рівностей може бути довільним. Їх можна записати у вигляді

. (2.13)

У ряді випадків із цих співвідношень можна виразити одні проектні параметри через інші. Це дозволяє виключити деякі параметри з процесу оптимізації, що приводить до зменшення розмірності задачі і полегшує її рішення. Аналогічно можуть вводитися також обмеження-нерівності, що мають вигляд

(2.14)

Необхідно зазначити особливість у відшуканні рішення за наявності обмежень. Оптимальне рішення тут може відповідати або локальному екстремуму (максимуму або мінімуму) усередині області проектування, або значенню цільової функції на межі області. Якщо ж обмеження відсутні, то шукається оптимальне рішення на всій області проектування, тобто знаходиться глобальний екстремум.

2.7 Чисельний пошук екстремуму функції однієї змінної

Більшість задач оптимізації зводиться до пошуку найбільшого (або найменшого) значення деякої функції. Методи математичного аналізу зручні для розв'язання цієї проблеми, коли функція задається явно і є при цьому диференційованою. Коли ж функція задається табличним способом або аналітично громіздкою формулою, ефективними є чисельні методи розв'язання.

Функція однієї змінної називається унімодальною на відрізку , якщо на ньому знаходиться єдина точка , в якій функція набуває мінімального значення.

2.7.1 Метод золотого перетину

Золотим перетином відрізка називається розподіл його точкою на дві нерівні частини так, щоб відношення усього відрізка до більшої частини дорівнювало відношенню більшої частини до меншої

Рисунок 2.1 - Ілюстрація золотого відношення

. (2.15)

Число називають золотим відношенням. Його значення відоме: (рис. 2.1). На відрізку можна визначити дві симетрично розміщені точки і , що реалізують золотий перетин. Їх знаходимо за формулами

,

.

Опишемо алгоритм мінімізації функції однієї змінної , вважаючи її унімодальною на відрізку . Початковий відрізок ділимо точками (перша точка) і (друга точка) за правилом золотого перетину і у цих точках обчислюємо значення функцій і . Порівняння цих значень дозволяє відкинути або інтервал , якщо , або інтервал , якщо . Довжина відрізка, що залишився, зменшиться у разів:

.

Після цього процес повторюємо.

На інтервалі, що залишився, уже є одна точка, що робить його золотий перетин: є друга точка золотого перетину відрізка ,а - перша точка золотого перетину відрізка . Знаючи одну з точок золотого перетину, іншу можна знайти за однією із вищезгаданих формул та обчислити значення у знову знайденій точці (значення в іншій точці вже обчислено на попередньому кроці).

Таким чином, на кожному кроці, починаючи з другого, потрібно лише одне обчислення функції , і інтервал невизначеності зменшується в разів:

Точка мінімуму

Після кроків маємо довжину інтервалу невизначеності

(2.16)

Процес золотого перетину продовжуємо до того часу, поки інтервал стане менше деякого заданого числа , названого точністю.

З формули випливає, що при довжина відрізка, який залишився, наближається до нуля як геометрична прогресія зі знаменником , тобто метод золотого перетину завжди збігається, причому лінійно. Число кроків , що забезпечують задану точність знаходження точки мінімуму, повинно задовольняти нерівність

(2.17)

Використовуючи цю формулу, можна знайти необхідне число кроків для забезпечення необхідної точності . Однак на практиці часто роблять інакше: визначають межу відрізка за формулою (2.16) і потім порівнюють із заданою точністю .

Метод золотого перетину гарантує знаходження мінімуму у найнесприятливіших умовах, коли функція є не тільки недиференційованою, але і навіть має розриви першого роду.

2.7.2 Алгоритм мінімізації функції методом золотого перетину

Знаходимо точки і за формулами

де

Обчислюємо значення функції та і порівнюємо їх:

якщо то за інтервал невизначеності беремо відрізок і знаходимо його золотий перетин, вважаючи, що , та обчислюємо

якщо за інтервал невизначеності беремо відрізок і знаходимо його золотий перетин, вважаючи, що

,

та обчислюємо

Перевіряємо умову . Якщо ця нерівність виконується, то за точку мінімуму беремо величину

якщо ні, то повертаємося до пункту 1.

Процес повторюємо до того часу, поки для деякого не справдиться .

2.7.3 Метод Ньютона

Для функцій однієї змінної класичний підхід при пошуку значень x в точках перегину функції полягає в розв'язанні рівняння

. (2.18)

Робота алгоритму починається у точці , яка являє собою початкове наближення координати стаціонарної точки або корінь рівняння (2.18). Потім будується лінійна апроксимація функції в точці , і точка, в якій апроксимуюча лінійна функція обертається в нуль, береться як наступне наближення. Якщо точка взята як поточне наближення до стаціонарної точки, то лінійна функція, яка апроксимує функцію у точці , записується у вигляді

. (2.19)

Прирівнявши праву частину рівняння (2.19) до нуля, одержимо таке наближення:

.

На жаль, залежно від вибору початкової точки і виду функції алгоритм може як збігатися до істиної стаціонарної точки, так і розходитися.

2.7.4 Апроксимація кривими. Кубічна апроксимація

Раніше була зроблена спроба знайти малий інтервал, в якому знаходиться мінімум функції. Тут розглянемо інший підхід, згідно з яким використовується декілька значень функції у певних точках для апроксимації функції звичним поліномом, принаймні, в невеликій області значень. Потім положення мінімуму функції апроксимується положенням мінімуму полінома, оскільки останній обчислити простіше. Найбільшого поширення набули методи оцінювання з використанням квадратичної і кубічної апроксимацій, оскільки побудова апроксимуючого полінома вище третього порядку стає дуже складною процедурою. У даному розділі розглядається метод оцінювання з використанням кубічної апроксимації, що забезпечує велику точність [5]. Для кубічної апроксимації в цьому методі використовуються значення функції і її похідної, обчислені у двох точках (p, q). Розглянемо задачу мінімізації на прямій, тобто мінімізацію функції

; (2.20)

,

де - задана точка; - заданий напрям; - крок. Отже,

. (2.21)

Припускаємо, що відомі такі значення:

. (2.22)

Цю інформацію можна використовувати для побудови кубічного полінома

, (2.23)

який апроксимуватиме функцію . Якщо так, то рівняння, які визначають а,b,c,d, мають такий вигляд:

(2.24)

Ці рівняння мають таке розв'язання:

(2.25)

Точки перегину кубічного полінома є розв'язком рівняння

Отже, якщо r є точкою мінімуму кубічного полінома, то

(2.26)

Одне із значень (2.26) відповідає мінімуму. Друга похідна дорівнює

. (2.27)

Якщо ми вибираємо додатний знак, то при

друга похідна буде . Тоді:

. (2.28)

Кращі чисельні результати можуть бути одержані при використанні такої еквівалентної формули:

Якщо , то необхідно вибрати значення q додатним, тобто зробити крок у напрямку спадання функції , інакше значення q необхідно вибрати від'ємним. Значення q повинно бути таким, щоб інтервал містив мінімум. Це буде справедливо, якщо або якщо .

Якщо жодна з цих умов не виконана, то значення q подвоюється, це повторюється до того часу, поки даний інтервал не міститиме мінімум. Для визначення початкового значення q Давідон, Флетчер і Пауелл запропонували вибирати q таким чином:

(2.29)

де - оцінка найменшого значення істинного мінімуму ; - константа, значення якої вибирається таким, що дорівнює 2 або 1.

Ця ітераційна процедура має такі кроки:

1. Знайти

і

2. Перевірити, чи виконується умова , і якщо вона не виконується, то необхідно здійснювати пошук уздовж напрямку d. Вибрати q з виразу (2.29). При цьому необхідно "вгадати" .

3. Обчислити

і .

4. Якщо або , то інтервал, що містить мінімум, знайдений. Інакше необхідно замінити q на 2q і повернутися до кроку 3.

5. Використання рівняння (2.28) для апроксимації точки мінімуму на інтервалі (0, q) значенням r.

6. Якщо

,

де - задана точність, то необхідно зупинитися.

7. Повернутися на крок 5, використовуючи інтервал (0, r), якщо , або використовуючи інтервал (r, q), якщо .

На кроці 6 здійснюється перевірка значення похідної. Попередні перевірки приводять до зупинення тоді, коли положення мінімуму не змінюється.

2.8 Чисельні методи пошуку екстремуму функції декількох змінних

Методи, орієнтовані на розв'язання задач безумовної оптимізації, можна розділити на три широкі класи відповідно до типу інформації, яка використовується при реалізації того або іншого методу [5,6]:

1. Методи прямого пошуку, які гуртуються на обчисленні тільки значень цільової функції.

2. Градієнтні методи, в яких використовуються значення перших похідних цільовій функції.

3. Методи другого порядку, в яких використовуються значення перших і других похідних цільовій функції.

Із методів прямого пошуку звичайно відзначають три методи, які є ефективними і можуть буди використані для широкого числа практичних застосувань:

1) пошук за симплексом, або S²-метод;

2) метод пошуку Хука-Джівса;

3) метод зв'язаних напрямів Пауелла.

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

У процесі пошуку за S²-методом послідовно оперують регулярними симплексами (множина (n+1)-ї рівновіддаленої точки в n-вимірному просторі) у просторі керованих змінних. При реалізації методу Хука-Джівса використовується фіксована множина (координатних) напрямів, які вибираються рекурсивним способом. Метод Пауелла орієнтований на розв'язання задач з квадратичними цільовими функціями і збігається за кінцеве число ітерацій. До загальних особливостей всіх трьох методів необхідно віднести відносну простоту відповідних обчислювальних процедур, які легко реалізуються і швидко коригуються. З іншого боку, реалізація вказаних методів часто вимагає дуже великої кількості обчислень значень функції. Ця обставина призводить до необхідності розгляду методів, що гуртуються на використанні градієнта цільової функції.

До градієнтних методів належать такі:

1) метод найшвидшого спуску; 2) метод Флетчера-Рівса (метод зв'язаних градієнтів); 3) методи змінної метрики (метод Давідона-Флетчера-Пауелла (ДФП), метод Бройдена-Флетчера-Шенно (БФШ) та ін.

У градієнтних методах для визначення найкращого напрямку пошуку використовується відома властивість напрямку градієнта - напрямок градієнта є найшвидшим зростанням функції. Отже, протилежний напрямок є напрямком найшвидшого спадання функції. Напрямок градієнта перпендикулярний в будь-якій точці лінії постійного рівня, оскільки уздовж цієї лінії функція постійна. Тому якщо ми знаходимося в точці на деякому кроці процесу оптимізації, то пошук мінімуму функції здійснюється уздовж напрямку . На кроці i точка мінімуму апроксимується точкою . Наступною апроксимацією є точка

, (2.30)

де - поточне наближення до розв'язку ; л_i - параметр, який характеризує дпвжину кроку; - напрямок пошуку. Значення л_i може бути знайдене за допомогою одного з методів одновимірного пошуку.

До методів другого порядку відносять:

1) метод Ньютона; 2) модифікований метод Ньютона; 3) метод Марквардта, що є комбінацією методів Коші і Ньютона.

Зважаючи на великий обсяг інформації, у цьому розділі розглянемо два градієнтних методи: 1) метод найскорішого спуску, на прикладі якого можна продемонструвати окремі прийоми, що використовуються при реалізації різних градієнтних алгоритмів; 2) метод ДФП, який відрізняється найширшим використанням при розв'язанні різноманітніших задач.

2.8.1 Метод найшвидшого спуску

Припустимо, що в деякій точці простору керованих змінних вимагається визначити напрямок найшвидшого локального спуску, тобто найбільшого локального зменшення цільової функції. Розкладемо цільову функцію навколо точки в ряд Тейлора:

і відкинемо члени другого порядку і вище. Очевидно, що локальне зменшення цільової функції визначається другим доданком, оскільки значення фіксоване. Найбільше зменшення F асоціюється з вибором такого напрямку в (4.1), якому відповідає найбільша від'ємна величина скалярного добутку другої складової розкладу. З властивості скалярного добутку випливає, що названий вибір забезпечується при

і другий доданок набуває вигляду

.

Розглянутий випадок відповідає найшвидшому локальному спуску. Визначення л доцільно проводити на кожній ітерації

.

Блок-схема методу найшвидшого спуску наведена на рис. 2.2. Для пошуку мінімуму функції

у напрямку з точки використовується метод квадратичної апроксимації.

У точці , і ми вибираємо довжину кроку такою, щоб крок "перекрив" мінімум функції . Похідна

.

У алгоритмі перевіряється умова "перекриття" мінімуму, яка виконується, якщо або , або . Якщо мінімум не потрапив у відрізок , то подвоюється, і це повторюється стільки разів, скільки необхідно для виконання умови "перекриття".

Упевнившись, що відрізок містить мінімум, як третю точку візьмемо точку , і потім проводиться квадратична апроксимація з метою пошуку мінімальної точки. У процесі пошуку припускається збіжність до екстремуму, тому для ефективності процедури зменшується довжина кроку. При цьому поділ кроку саме вибирається довільно.

Рисунок 2.2 - Блок-схема методу найшвидшого спуску

Критерій завершення кожної ітерації повинен лише призводити до зменшення значення функції, в порівнянні із значенням це дозволить скоротити об'єм обчислень порівняно з об'ємом обчислень при точному пошуку мінімуму, який може виявитися дуже значним.

Метод найшвидшого спуску не рекомендується як "серйозна" оптимізаційна процедура, оскільки для практичного застосування "працює" дуже повільно. Це пояснюється тим, що властивість найшвидшого спуску є лише локальною властивістю і тому необхідна часта зміна напряму, що і приводить у результаті до неефективної обчислювальної процедури.

2.8.2 Метод Давідона-Флетчера-Пауелла

У цьому методі напрямок пошуку на кроці i є напрямком , де - додатно певна симетрична матриця, яка відновлюється на кожному кроці.

Почнемо пошук з початкової точки , взявши як початкову матрицю (як правило, одиничну матрицю). Ітераційна процедура може бути представлена таким чином (замість далі пишеться ):

1. На кроці i є точка і додатно визначена симетрична матриця H_i.

2. Як напрямок пошуку береться напрямок d_i=-H_ig_i.

3. Щоб знайти функцію , що мінімізує функцію F(x_i+d_i), необхідно виконати одновимірний пошук уздовж прямої x_i+d_i.

4. Припустити, що v_i= d_i.

5. Припустити, що x_i+1=x_i+v_i.

6. Знайти F(x_i+1) і g_i+1. Завершити процедуру, якщо величини або достатньо малі. Інакше необхідно продовжити.

7. Припустити, що u_i=g_i+1-g_i.

8. Обновити матрицю H таким чином:

H_i+1=H_i+A_i+B_i,

; .

9. Збільшити i на одиницю і повернутися на крок 2.

Одновимірний пошук необхідно здійснювати кубічною апроксимацією, але пошук закінчується, не досягнувши граничної збіжності, а після того, як буде одержано покращення значення функції.

2.9 Методи оптимізації за наявності обмежень

2.9.1 Поняття штрафної функції

Розглянемо задачу нелінійного програмування такого вигляду:

мінімізувати (2.31) при обмеженнях . (2.32)

Початкова задача (2.31), (2.32) умовної оптимізації перетвориться в задачу пошуку мінімуму без обмежень (у задачу безумовної оптимізації) функції

(2.33)

де - штрафна функція; - штрафний параметр.

Необхідно, щоб при порушенні обмежень "штрафувала" функцію Z, тобто збільшувала її значення. У цьому випадку мінімум Z знаходитиметься усередині області обмежень. Функція , що задовольняє цю умову, може бути не єдиною. Розглянемо основні типи штрафів.

1. Нескінченний бар'єр. Відповідний вираз набуває нескінченно великих значень в недопустимих точках і нульові значення в допустимих точках. На практиці у якості "нескінченність" використовується велике додатне число, що допускає запис в ЕОМ, наприклад:

. (2.34)

2. Логарифмічний штраф

. (2.35)

Цей штраф додатний при всіх x, таких, що , і від'ємний при .

3. Штраф, заданий зворотною функцією

(2.36)

не має від'ємних значень в допустимій області.

4. Штраф типу квадрата зрізання

, . (2.37)

Розглянемо розв'язок задачі на основі, наприклад, штрафу, заданого зворотною функцією (2.36). Тоді функція (2.33) набуде вигляду

. (2.38)

Якщо x набуває допустимих значень, тобто значення, для яких , то Z набуває значення, які більші за відповідні значення (істинної цільової функції задачі), і різницю можна зменшити за рахунок того, що r може бути дуже малою величиною. Але якщо x набуває значень, які хоча і є допустимими, але близькі до границі області обмежень і, принаймні, одна з функцій близька до нуля, тоді значення функції , і, отже, значення функції Z стануть дуже великими. Таким чином, вплив функції полягає у створенні "гребню з крутими краями" уздовж кожної границі області обмежень. Отже, якщо пошук починається з допустимої точки, і здійснюється пошук мінімуму функції без обмежень, то мінімум буде досягнутий усередині допустимої області для задачі з обмеженнями. Вважаючи r достатньо малою величиною, для того, щоб вплив був малим у точці мінімуму, можна зробити точку мінімуму функції без обмежень такою, що збігається з точкою мінімуму функції з обмеженнями.

2.9.2 Метод SUMT

Результати попереднього розділу показують, що можна розв'язати задачу мінімізації з обмеженнями, розв'язуючи, для послідовності значень r, яка прагне до нуля, задачу без обмежень наступного вигляду: мінімізувати

.

Метод SUMT (sequential unconstrained minimisation technique) вперше запропонований Керролом в 1961 р. Його ідеї були досліджені Фіакко і Маккорміком, які розглянули теоретичні питання і збіжність методу, а також створили практичну систему для його реалізації. На практиці необхідно побудувати обчислювальний метод, що використовує теоретичну властивість збіжності, яка розглянута в попередньому розділі. Для заданих функції і обмежень необхідно вибрати початкове значення , щоб сформувати функцію, яка мінімізується без обмежень методом ДФП, який був наведений раніше. Визначивши мінімум функції , необхідно зменшити значення r шляхом , де константа . Потім необхідно мінімізувати функцію , знову використовуючи метод ДФП. Таким чином, буде розроблена ітераційна процедура. На к-му кроці мінімізується функція , мінімум якої знаходиться у точці . Важливо, що її можна використовувати надалі як першу точку в ітераційній процедурі мінімізації функції , де . Тепер зрозуміло, що послідовність спадає і прямує до нуля, отже, послідовність точок мінімуму сходитиметься до розв'язання задачі з обмеженнями. Блок-схема методу наведена на рис. 2.3.

Рисунок 2.3 - Блок-схема методу SUMT

Передбачається, що на початку процедури є допустима точка. Важливо, щоб у процесі подальших обчислень одержані точки належали допустимій області. Метод ДФП є градієнтним методом мінімізації, який використовує при одномірному пошуку кубічну апроксимацію. Тоді у міру наближення точки x до границі усередині допустимої області , а у міру наближення точки x до границі зовні допустимої області . Таким чином, якщо пошук здійснюється уздовж прямої, що з'єднує дві точки, одна з яких лежить всередині, а інша поза областю обмежень, то кубічна апроксимація виявляється неприйнятною, оскільки функція розривна уздовж даної прямої.

Необхідно ретельно досліджувати такі питання при використанні методу ДФП в даній задачі. Цілком ефективним є наступний прийом. Нехай є точка p і напрям пошуку . Наступна точка необхідна для здійснення кубічної апроксимації. Почнемо із значення = 2 (подвоєний крок у методі Ньютона) і перевіримо, чи є точка q допустимою, тобто чи виконується нерівність для всіх j. Якщо вона виконується, то не змінюється, але якщо нерівність не виконується, то замінюється на , знаходиться нова точка q і знову виконується перевірка, поки допустима точка q не буде знайдена. Вибір значення а не цілком очевидний. Вибір може бути успішним, при довжина кроку стає близькою до відстані до найближчої границі області обмежень і тому є безпечною для інтерполяційної процедури. Вибір початкового значення r важливий з погляду скорочення числа ітерацій при мінімізації функції . Якщо спочатку r вибрано дуже малим, то метод збігатиметься дуже швидко. Проте такий вибір може призвести до серйозних ускладнень при обчисленнях. Для малих r функція швидко змінюватиметься довкола мінімуму, що може викликати ускладнення при використанні градієнтного методу. Дуже велике значення r може привести до того, що штрафна функція в рівнянні (6.8) стане домінуючою. Тому "розумний" вибір початкового значення r дуже важливий. Для багатьох задач можна припустити, що . Більш раціональніший підхід ґрунтується на тому, щоб зрозуміти, що якщо початкова точка x лежатиме поблизу мінімуму функції

,

то градієнт функції буде малий:

.

Квадрат норми цього вектора

і мінімум буде досягнутий при

.

Таке початкове значення r, яке припускають Фіакко і Маккормік, повинне давати позитивні результати в загальному випадку. Функція мінімізується до того часу, поки два послідовні значення F₁ і F₂ не стануть такими, що . Ця умова, звичайно, може бути змінена. Програма закінчує роботу, коли

.

3. Метод скінченних елементів у задачах механіки деформованого твердого тіла

3.1 Загальні поняття і класифікація задач обчислювальної механіки

Механіка як фундаментальна наукова дисципліна може бути поділена на три основні частини [7].

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

Обчислювальна механіка. Сфери застосування обчислювальної механіки можуть бути виділені таким чином.

Механіка суцільного середовища вивчає континуальні моделі тіл. Дві традиційні сфери застосування механіки суцільного середовища - механіка рідин і твердих тіл; остання включає механіку конструкцій, які, очевидно, складаються з твердих тіл. Обчислювальна механіка твердого тіла використовує науково-прикладний підхід, тоді як обчислювальна механіка конструкцій акцентує увагу на технологічних застосуваннях до аналізу і дизайну конструкцій.

Обчислювальна механіка рідин має справу із задачами, що описують рівновагу і рух рідин і газів, включаючи такі добре розроблені галузі, як гідродинаміку, аеродинаміку, фізику атмосфери, фізику горіння і вибуху.

Механіка зв'язаних задач - одна з останніх сфер застосування чисельних методів. Мається на увазі, що ця сфера включає фізичні системи, що виходять за класичні границі механіки твердих тіл і рідин, розглядає взаємодію рідин і конструкцій.

Нарешті, механіка систем має справу з механічними об'єктами, природними або штучними. Приклад штучних систем - літаки, мости, двигуни, машини, мікрочіпи, роботи і т. д.; біологічних систем - тварини або дерева, що вивчаються з позицій біомеханіки.

Статика і динаміка. Задачі механіки суцільного середовища можуть бути класифіковані відповідно до того, береться до уваги ефект інерції чи ні, тобто можуть бути розглянуті задачі статики і динаміки.

У динаміці явно розглядається залежність від часу, оскільки обчислення сил інерції (а так само, можливо, демпфірувальних сил) вимагає обчислення похідних за часом.

Задачі статики так само можуть бути тимчасово залежними, проте силами інерції нехтують через їх малу величину. Статичні задачі можуть бути так само поділені на строго статичні і квазістатичні. У математичне формулювання останніх час явно не входить; замість нього у разі потреби використовуються тимчасово подібний параметр, що характеризує історію процесу. У інших квазістатичних задачах, таких, як розрахунок деформації повзучості, циклічна утома, пластичність, залежна від швидкості навантаження, потрібна реалістичніша оцінка часу, проте і в цьому випадку інерційними силами можна нехтувати.

Лінійність і нелінійність. Задачі механіки можуть бути так само поділені на лінійні і нелінійні:

Лінійна механіка має справу з проблемами, в яких відгук системи лінійний в причинно-наслідковому значенні. Наприклад, якщо прикладені зовнішні сили подвоюються, то результуючі переміщення і внутрішні напруження так само подвоюються. У протилежному разі задачі класифікуються як нелінійні.

Методи дискретизації. Остання класифікація задач обчислювальної механіки базується на методі дискретизації, за допомогою якого безперервна математична модель перетвориться в дискретну модель, що складається з кінцевого числа ступенів вільності. Найбільш поширеними методами просторової дискретизації є:

метод скінченних елементів;

метод граничних елементів;

метод скінченних різниць;

метод скінченних об'ємів.

У даний час в лінійній механіці твердого тіла метод кінцевих елементів найбільш поширений, тоді як застосування методу граничних елементів для розв'язання даних задач перебуває на другому місці. Для нелінійних задач метод скінченних елементів є найефективнішим і домінуючим.

Класичний метод скінченних різниць майже повністю втратив своє значення при вирішенні практичних задач механіки твердого тіла. Це твердження, проте, неправильне для механіки рідини і газів, де різницеві методи до цього часу значно поширені. Метод скінченних об'ємів, що базується на законах збереження, застосовується для розв'язання сильно нелінійних задач механіки рідини і газів.

Варіанти МСЕ. Термін "метод скінченних елементів" насправді визначає широкий спектр обчислювальних технологій відповідно до деяких загальних властивостей, які будуть розглянуті надалі [8,9]. Класифікація МСЕ щодо механіки конструкцій залежно від типу шуканої функції може бути представлена таким чином. Формулювання МСЕ: