Решение одномерного уравнения теплопроводности в mathcad

Видео:8.1 Решение уравнения теплопроводности на отрезкеСкачать

Решение одномерного уравнения теплопроводности в mathcad

11.3.1. Параболические и гиперболические уравнения

Разработчики впервые применили дополнительные встроенные функции для решения параболических и гиперболических уравнений в частных производных в версии Mathcad 11, отлично осознавая значимость этих задач для современного исследователя и инженера. Предусмотрены два варианта решения: при помощи вычислительного блока Given/pdeso l ve , а также при помощи встроенной функции numo l . Первый путь проще в применении и нагляднее, зато второй позволяет автоматизировать процесс решения уравнений в частных производных, например, если нужно включить его в качестве составного шага в более сложную Mathcad-программу.

Вычислительный блок Given /pdesolve

Встроенная функция pdesolve применяется в рамках вычислительного блока, начинающегося ключевым словом Given , и пригодна для решения различных гиперболических и параболических уравнений. Она предназначена для решения одномерного уравнения (или системы уравнений) в частных производных (того, которое определит пользователь в рамках вычислительного блока Given ), зависящего от времени t и пространственной координаты х , имеет целый набор различных аргументов и работает следующим образом:

• pdesolve(u,x,xrange,t,trange,[xpts],[tpts])) — Возвращает скалярную (для единственного исходного уравнения) или векторную (для системы уравнений) функцию двух аргументов (x,t) , являющуюся решением дифференциального уравнения (или системы уравнений) в частных производных. Результирующая функция получается интерполяцией сеточной функции, вычисляемой согласно разностной схеме:

• u — явно заданный вектор имен функций (без указания имен аргументов), подлежащих вычислению. Эти функции, а также граничные условия (в форме Дирихле или Неймана) должны быть определены пользователем перед применением функции pdeso l ve в вычислительном блоке после ключевого слова Given . Если решается не система уравнений в частных производных, а единственное уравнение, то, соответственно, вектор и должен содержать только одно имя функции и вырождается в скаляр;
• х — пространственная координата (имя аргумента неизвестной функции);
• xrange — пространственный интервал, т. е. вектор значений аргумента х для граничных условий. Этот вектор должен состоять из двух действительных чисел (представляющих левую и правую границу расчетного интервала);
• t — время (имя аргумента неизвестной функции);
• t range — расчетная временная область: вектор значений аргумента t, который должен состоять из двух действительных чисел (представляющих левую и правую границу расчетного интервала по времени);
• xpts — количество пространственных точек дискретизации (может не указываться явно, в таком случае будет подобрано программой автоматически);
• tpts — количество временных слоев, т. е. интервалов дискретизации по времени (также может не указываться пользователем явно).

В качестве примера использования функции pdeso l ve (листинг 11.4) используем то же самое одномерное уравнение теплопроводности (11.5) с граничными и начальными условиями (11.6) и (11.7).

Листинг 11.4. Решение одномерного уравнения теплопроводности

Для корректного использования функции pdeso l ve предварительно, после ключевого слова Given , следует записать само уравнение и граничные условия при помощи логических операторов (для их ввода в Mathcad существует специальная панель). Обратите внимание, что уравнение должно содержать имя неизвестной функции u(x,t) вместе с именами аргументов (а не так, как она записывается в пределах встроенной функции pdeso l ve ). Для идентификации частных производных в пределах вычислительного блока следует использовать нижние индексы, например, uxx(,t) , для обозначения второй производной функции и по пространственной координате х .

Как видно из рис. 11.14, на котором изображены результаты расчетов по листингу 11.4, встроенная функция с успехом справляется с уравнением диффузии, отыскивая уже хорошо знакомое нам решение. Заметим, что использование встроенной функции pdeso l ve связано с довольно громоздкими вычислениями, которые могут отнимать существенное время.

Как вы можете заметить, выбирать величину шага по пространственной и временной переменным может как сам алгоритм, так и пользователь (неявным образом, через число узлов сетки). Читателю предлагается повторить вычисления листинга 11.4 для различных комбинаций параметров (главным образом, числа узлов сетки), чтобы проверить, в каких случаях алгоритм встроенной функции справляется с задачей, выдавая верное решение, а в каких дает сбой.

Рис. 11.14. Решение уравнения диффузии тепла при помощи встроенной функции pdesoдve (листинг 11.4)

Пример: волновое уравнение

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

Здесь неизвестная функция u(x,t) описывает динамику смещения профиля струны относительно невозмущенного (прямолинейного) положения, а параметр с характеризует материал, из которого изготовлена струна.

Как вы видите, уравнение (11.11) содержит производные второго порядка, как по пространственной координате, так и по времени. Для того чтобы можно было использовать встроенную функцию pdesolve , необходимо переписать волновое уравнение в виде системы двух уравнений в частных производных, введя вторую неизвестную функцию v=ut . Программа для решения волнового уравнения приведена в листинге 11.5, а результат— на рис. 11.15.

Листинг 11.5. Решение волнового уравнения

Рис. 11.15. Решение волнового уравнения (продолжение листинга 11.5)

Встроенная функция numol

Альтернативный вариант решения дифференциальных уравнений в частных производных заключается в применении еще одной встроенной функции numo |, которая реализует тот же самый алгоритм сеток, позволяя вручную задать большинство его параметров:

• numol(xrange,xpts,trange,tpts,Npde,Nae,rhs,init,bc) — Возвращает матрицу решения дифференциального уравнения в частных производных, представляющую искомую сеточную функцию в каждой точке по пространственной (по строкам) и временной координате (по столбцам). Если решается не одно уравнение, а система уравнений, то результатом является составная матрица, образованная путем слияния (слева-направо) матриц со значениями каждой искомой сеточной функции:

• Npde — общее количество дифференциальных уравнений в частных производных в системе;
• Nae — общее количество дополнительных алгебраических уравнений, которые также могут входить в систему;
• rhs — векторная функция, определяющая систему дифференциальных и алгебраических уравнений (формат этого и двух следующих матричных параметров объяснен в листинге 11.9);
• init — векторная функция, определяющая начальные условия для каждой неизвестной функции;
• be — функциональная матрица граничных условий;
• xrange — пространственный интервал, т. е. вектор значений аргумента х для граничных условий. Этот вектор должен состоять из двух действительных чисел (представляющих левую и правую границу расчетного интервала);
• xpts — количество пространственных точек дискретизации (может не указываться явно, в таком случае будет подобрано программой автоматически);
• trange — расчетная временная область: вектор значений аргумента t, который должен состоять из двух действительных чисел (представляющих левую и правую границу расчетного интервала по времени);
• tpts — количество временных слоев, т. е. интервалов дискретизации по времени (также может не указываться пользователем явно);

Пример решения волнового уравнения при помощи функции numol приведен в листинге 11.6, особое внимание в котором мы призываем уделить формату представления векторов rhs , init и be , а также принципу извлечения отдельных сеточных решений из матрицы-результата. График решения, показанный на рис. 11.16, полезно сравнить с результатом применения вычислительного блока из предыдущего раздела (см. листинг 11.5 и рис. 11.15).

Листинг 11.6 . Решение волнового уравнения при помощи функции numol

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

Рис. 11.16. Решение волнового уравнения (продолжение листинга 11.6)

Именно в целях визуализации решения параболических и гиперболических уравнений в частных производных использование функции numol наиболее полезно. График решения динамических уравнений (зависящих от времени t) выглядит намного эффектнее и воспринимается несравненно лучше, если он оформлен в виде анимации. Для создания анимационных роликов расчетное время следует выразить через константу FRAME и затем применить команду View / Animate (Вид / Анимация) (см. разд. 13.3.2).

Видео:Mathcad-09. Пример: уравненияСкачать

Решение одномерного уравнения теплопроводности в mathcad

Глава 5. Решение дифференциальных уравнений

5.8 Функции решения параболических и гиперболических уравнений

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

Уравнения в частных производных можно разделить на три типа:

1) параболические, содержащие первую производную по одной переменной и вторую – по другой, причем все производные входят в уравнение с одинаковым знаком;

2) гиперболические, содержащие первую производную по одной переменной и вторую – по другой, входящие в уравнение с разными знаками;

3) эллиптические, содержащие вторые производные, причем одного знака.

Функции решения параболических и гиперболических уравнений

MathCAD 11 включает в себя две функции для решения параболических и гиперболических уравнений, pdesolve и numol .

Функция pdesolve используется в составе вычислительного блока Given – pdesolve для решения параболических и гиперболических уравнений (или систем уравнений) в частных производных, имеющих в качестве аргументов, как правило, время t и пространственную координату x .

Обращение к этой функции:

возвращает скалярную (для одного уравнения) или векторную (для системы уравнений) функцию, являющуюся решением уравнения (или системы уравнений). Здесь u –явно заданный вектор имен функций (без указания имен аргументов), подлежащих вычислению. Эти функции, а также граничные условия должны быть определены внутри вычислительного блока Given – pdesolve ; х – пространственная координата; x range – вектор значений аргумента х для граничных условий. Он должен состоять из двух чисел, представляющих две границы расчетного интервала; t – время (имя второго аргумента неизвестной функции); t range – вектор значений аргумента t для граничных условий, состоящих из двух чисел, представляющих две границы расчетного интервала; x pts – количество пространственных точек дискретизации (может не указываться); t pts – количество временных слоев (может не указываться).

Пример использования функции pdesolve приведен на рис. 5.19 запись вычислительного блока с функцией pdesolve аналогична записи блока с функцией Odesolve . Результаты расчета показаны на рис. 5.20.

Решение одномерного волнового уравнения

Здесь w -перемещение, v -скорость перемещения

где тогда

Представим первое уравнение как систему двух

уравнений первого порядка

Given

граничные условия

Рис. 5.1 9 Пример использования функции pdesolve

Единичное решение волнового уравнения

Сетка решений волнового уравнения на временном и пространственном интервалах

Рис. 5. 20 Результаты решения волнового уравнения

Обратите внимание на то, что уравнения внутри вычислительного блока должны записываться с аргументами. Для идентификации частных производных следует использовать нижний индекс, например, – вторая производная функции u по пространственной координате х.

Недостатком функции pdesolve (как и функции Odesolve ) является невозможность ее использования в составе выражения – программы для многократного решения дифференциального уравнения. При необходимости многократного решения обыкновенных дифференциальных уравнений в состав программного модуля можно включать функции Rkadapt или Bulstoer.

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

Функция numol предназначена для решения тех же уравнений, что и функция pdesolve .

Обращение к этой функции:

возвращает матрицу решения дифференциального уравнения в частных производных в каждой точке по пространственной (по строкам) и временной (по столбцам) координате. Если решается не одно уравнение, а система уравнений, то результатом решения является составная матрица, образованная путем слияния (слева направо) со значениями каждой искомой сеточной функции. Здесь x range – вектор значений аргумента х для граничных условий. Он должен состоять из двух чисел, представляющих две границы расчетного интервала; t range – вектор значений аргумента t для граничных условий, состоящих из двух чисел, педставляющих две границы расчетного интервала; x pts – количество пространственных точек дискретизации (может не указываться); t pts – количество временных слоев (может не указываться); N pde – количество дифференциальных уравнений в частных производных в системе; N ae – количество дополнительных алгебраических уравнений, входящих в систему; rhs – вектор правых частей уравнений; init – векторная функция, определяющая начальные условия для каждой неизвестной функции; bc – функциональная матрица граничных условий.

Вектор граничных условий может иметь значения трех типов:

– rhs содержит вторые пространственные производные: граничные условия (или Дирихле « D », или Неймана « N ») требуются по одному с каждой стороны интервала интегрирования;

– rhs содржит первые пространственные производные: граничные условия Дирихле на левой или правой границе интервала, на другой стороне NA ;

– нет пространственных производных – граничные условия не требуются.

Функциональная матрица bc содержит три столбца, имеющих ледующий вид:

– ( init _ left ( t ) init _ right ( t ) » D «) – для граничных условий Дирихле;

– ( init _ left ( t ) init _ right ( t ) » N «) – для граничных условий Неймана.

Пользоваться функцией numol намного сложнее, чем функцией pdesolve .

Граничное условие называется условием Дирихле, если задано значение функции, или Неймана, если задана первая производная функции.

Видео:Лекция №1.1 Явная и неявная схемы для уравнения теплопроводностиСкачать

Интернет-версия справочника «Теплоэнергетика и теплотехника»

А.П.Солодов, В.Ф.Очков « Mathcad / Дифференциальные модели»

Москва, Издательство МЭИ ( publish @ admin . mpei . ac . ru тел. (095) 361-16-81), 2002 г., – 239 с.: илл.

A . Solodov , V . Ochkov « Differential Models »

I n English, Springer Publishing House ( Springer-Verlag, spring , 200 4)

Глава 1. Дифференциальные математические модели

1.2. Законы в дифференциальной форме

1.3. Модели роста

1.4. Законы сохранения

1.5. Одномерные стационарные модели. ТВЭЛ ядерного реактора.

Глава 2. Интегрируемые в квадратурах дифференциальные уравнения

2.1. Перечень интегрируемых уравнений

2.2. Линейное уравнение первого порядка

2.3. Линейные однородные уравнения с постоянными коэффициентами

2.4. Линейные неоднородные уравнения с постоянными коэффициентами

2.5. Уравнения с разделяющимися переменными

2.6. Однородные дифференциальные уравнения

2.7. Уравнения, допускающие понижение порядка

Глава 3. Динамическая модель системы с тепловыделением

3.2. Математическая модель

3.3. Фазовый портрет системы. Устойчивые и неустойчивые состояния равновесия

3.4. Представление множества равновесных состояний

3.5. Построение бифуркационной диаграммы

3.6. Трехмерное представление равновесных состояний в форме «катастрофы сборки»

3.7. Скачки состояний при плавных изменениях параметров

3.8. Временная эволюция системы с тепловыделением

Глава 4. Жесткие дифференциальные уравнения

4.2. Метод rkfixed. Неустойчивость решения

4.3. Метод rkadapt

4.4. M етод stiffr. Решение модельной жесткой задачи

4.5. Метод stiffr . Решение уравнений химической кинетики.

4.6. Явные и неявные методы

Глава 5. Теплопередача в окрестности критической точки поперечно обтекаемой трубы

5.2. Интегральное уравнение теплового пограничного слоя

5.3. Математическое описание задачи

5.4. Распределение скорости внешнего потока по окружности трубы

5.5. Вывод соотношений для критической точки

5.6. Приведение математического описания к безразмерному виду

5.7. Представление правой части дифференциального уравнения в форме алгоритма оптимизации

5.8. Численное интегрирование с помощью встроенной функции Odesolve

Глава 6. Уравнение Фолкнера-Скэн. Трение и теплообмен в пограничном слое

6.2. Математическая формулировка

6.3. Сведение краевой задачи к начальной задаче методом sbval

6.4. Решение начальной задачи методом rkfixed

6.5. Построение поля течения

6.6. Пограничный слой на проницаемой поверхности.

6.7. Уравнение теплового пограничного слоя

6.8. Закон теплообмена

6.9. Неприятности при решении краевой задачи посредством функции Odesolve

Глава 7. Уравнение Рэлея. Гидродинамическая неустойчивость

7.2. Уравнения гидродинамики для свободного сдвигового потока

7.3. Метод малых возмущений. Линеаризация.

7.4. Переход в комплексную область

7.5. Численное интегрирование в комплексной области. Программа Euler

7.6. Интегрирование и поиск собственных значений

7.7. Возвращение в действительную область

Глава 8. Кинематические волны концентрации в ионитовом фильтре

8.2. Уравнение сохранения концентрации

8.3. Волновое уравнение для концентрации

8.4. Уравнение сохранения в безразмерной форме

8.5. Уравнение изотермы адсорбции

8.6. Решение волнового уравнения методом характеристик (Пример 1 анимации, Пример 2 анимации)

Глава 9. Кинематические ударные волны

9.2. Уравнение сохранения в конечно-разностной форме.

9.3. Разрывные решения. Ударные волны

9.4. Метод Мак-Кормака. Вычислительная программа McCrm

9.5. Ударные волны концентрации примесей в фильтре

9.6. Ударные волны на автостраде

9.7. Гравитационные пузырьковые течения. Ударные волны паросодержания

Глава 10. Численное моделирование температурного поля платы компьютера

10.2. Встроенные функции для решения уравнений в частных производных

10.3. Представление дифференциального уравнения теплопроводности в конечно-разностной форме

10.4. Метод Гаусса-Зайделя. Программа Plate

10.5. Тепловая модель платы компьютера

10.6. Задача об орбитальной платформе. Функция bvalfit

11.2. Математическая формулировка задачи

11.3. Дискретное представление

11.4. Метод прогонки. Вычислительные программы Coef и SysTRD

11.5. Компьютерное моделирование периодических тепловых воздействий

Приложение. Встроенные функции для решения обыкновенных дифференциальных уравнений

Основу книги составляют инженерно-физические проекты, выполненные в среде Mathcad по примерно одинаковой схеме:

q Постановка задачи и формулировка физической модели

q Формулировка дифференциальной математической модели, т.е. модели в форме дифференциальных уравнений

q Интегрирование дифференциальных уравнений

q Визуализация результатов.

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

Хотя проекты, включенные в книгу, являются учебными, они содержат исследовательскую интригу, имеют важные практические приложения и не являются тривиальными задачами в вычислительном аспекте. Мы надеемся, что удалось найти подходящее сочетание инженерной актуальности, с одной стороны, и относительно высокой математической (вычислительной) сложности, с другой, чтобы обеспечить достаточную мотивацию для обращения к какому-либо из современных математических программных пакетов.

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

Книга состоит из одиннадцати глав и приложения.

В главе первой «Дифференциальные математические модели» обсуждается происхождение дифференциальных уравнений. Они не появляются сразу на рабочем столе инженера в готовом для интегрирования виде, а чаще всего должны быть составлены, выведены самим исследователем – как математические модели устройств, технологических процессов или аппаратов. Возможно, именно этот этап будет самым сложным в инженерном проекте, и мы считали полезным привести примеры того, как из нечеткого словесного описания рождаются дифференциальные математические модели исследуемых объектов.

Контрольный объем для двухмерной задачи

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

Логистическая модель с учетом случайных колебаний фактора роста

Глава вторая «Интегрируемые в квадратурах дифференциальные уравнения» содержит перечень и краткий рецептурный справочник для уравнений указанного типа. Показано, как можно применить символьный процессор Mathcad ’а при математических выкладках в этом случае.

Глава третья «Динамическая модель системы с тепловыделением» посвящена системам, которые могут развиваться по катастрофическим сценариям – с воспламенением и взрывом.

Поверхность равновесных состояний (вид сбоку и вид сверху)

Эволюция системы в случае трех равновесных состояний

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

Быстрая эволюция отклонений, возникающих в различных точках X медленно осциллирующего решения

В главе четвертой «Жесткие дифференциальные уравнения» подробно рассмотрена проблема численного интегрирования этого особого типа уравнений. Несмотря на практическую важность вопроса, студенту или инженеру трудно найти доступное описание сущности феномена. Mathcad имеет эффективную встроенную функцию для интегрирования жестких уравнений, но его справочная система также почти ничего не содержит по обсуждаемой проблеме. Материалы данной главы восполняют в некоторой степени этот пробел.

Интегрирование жестких уравнений химической кинетики.

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

Глава пятая «Теплопередача в окрестности критической точки поперечно обтекаемой трубы» моделирует коллизию, перед которой может оказаться начинающий инженер-исследователь при решении реальной задачи. В учебных примерах правые части дифференциальных уравнений всегда задаются простыми аналитическими выражениями. В реальных задачах правая часть почти всегда вообще не может быть записана как некоторое, пусть и сложное, аналитическое выражение, а представляется в форме сложного алгоритма. Поэтому мы считали полезным показать в действии структуры Mathcad’а (программные модули, встроенные функции оптимизации и т.п.), которые приходится использовать при решении реальных задач с дифференциальными уравнениями.

Идеальное обтекание цилиндра

Распределение по окружности трубы коэффициента теплоотдачи, теплового потока и температурного напора

Векторное поле скорости и функция тока в пограничном слое при вдуве

Глава шестая «Уравнение Фолкнера-Скэн. Трение и теплообмен в пограничном слое» посвящена численному анализу фундаментальной задачи гидродинамики и тепломассообмена. Совместно с разделом 1.4 «Законы сохранения» шестая глава образует краткий теоретический курс конвективного теплообмена. С точки зрения овладения техникой работы в Mathcad ’е, главной темой является численное решение краевой двухточечной задачи для системы обыкновенных дифференциальных уравнений с применением встроенной функции, позволяющей свести краевую проблему к задаче с начальными условиями.

В главе седьмой «Уравнение Рэлея. Гидродинамическая неустойчивость» методом малых возмущений (методом линеаризации) анализируется устойчивость течений в слое сдвига. Обнаруживаемые при анализе прогрессирующие возмущения в потоке являются предшественниками турбулентности. В методическом отношении эта глава может служить элементарным введением в общую теорию устойчивости. Рассматривается техника работы с дифференциальными уравнениями в комплексной области. Показано, как интегрировать уравнения с комплексными коэффициентами численным методом и как интерпретировать получившиеся комплекснозначные решения.

Распределение комплекснозначной амплитуды возмущений скорости в слое сдвига.

Глава восьмая «Кинематические волны концентрации в ионитовом фильтре» развивает тему интегрирования уравнений в частных производных, открытую в первой главе задачей об ударных волнах на автостраде. Сформулирована физическая и математическая модель ионитовых фильтров – устройств, обеспечивающих получение очень чистой воды для паросиловых установок тепловых и атомных электростанций. Методом характеристик построено решение уравнения в частных производных, описывающее пространственно-временную эволюцию поля концентраций. Показано, как нелинейность дифференциальной модели приводит к появлению разрывных неоднозначных решений.

Опрокидывание волны концентрации при расчете ионитового фильтра методом характеристик (случай изотермы адсорции легмюровского типа).

Глава девятая «Кинематические ударные волны» продолжает тему волновых уравнений. Разработана Mathcad -программа для численного интегрирования дифференциальных уравнений в частных производных, обеспечивающая эффективное воспроизведение ударных волн. Благодаря этому доводится до конца исследование задач об ударных волнах на автострадах и ударных волнах в фильтрах.

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

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

В главе десятой «Численное моделирование температурного поля платы компьютера» на примере теплопроводности рассматриваются численные методы для уравнений в частных производных второго порядка. Стандартные средства Mathcad ’а дополняются программой, позволяющей моделировать стационарные двухмерные температурные поля для существенно расширенного круга задач и корректно описывать разнообразные условия теплового взаимодействия объекта с окружающей средой.

Температурное поле компьютерной платы

Температурное поле в грунте вокруг трубопровода

Глава одиннадцатая «Температурные волны» посвящена нестационарным температурным полям. Инженерные приложения этого раздела теплофизики очень разнообразны. Например, в энергетике для оптимального управления пусковыми режимами необходимо рассчитывать нестационарные температурные поля в элементах машин и оборудования, с целью исключить нарушения зазоров в движущихся элементах из-за неодинакового расширения или предотвратить возникновение разрушающих термических напряжений в массивных деталях.

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

Авторы стремились показать в этой книге на реальных примерах, как можно эффективно использовать Mathcad на всех этапах разработки инженерного проекта:

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

· при получении аналитических решений на этапах, где это возможно (или где это позволяют возможности символьного процессора Mathcad)

· при численном решении, когда аналитические решения невозможны или неэффективны

· при презентации результатов с помощью многообразных и простых в использовании средств визуализации.

Мы следовали убеждению, что основная современная тенденция в инженерных исследованиях и проектировании, как и в инженерном образовании — все более широкое использование компьютерного моделирования, с тем чтобы максимально учесть принципиально важные эффекты и работать с полными математическими моделями объектов. Инженер-исследователь не должен находиться под прессом опасений, что получившаяся модель слишком сложна в формальном, вычислительном плане и что в связи с этим придется существенно «упрощать» задачу. Может быть, до такой степени, что модель совсем перестанет быть похожей на реальный объект.

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

Главное соображение в пользу этого – замечательная общность математических моделей. В центре книги – серия проектов о феномене кинематических ударных волн, когда однотипные уравнения моделируют и образование пробок на автострадах, и формирование фронта концентрации, грозящего прорывом через фильтр, и вспенивание газированной двухфазной среды в стакане или в корпусе энергетической установки. Мы начинаем книгу с обсуждения экспоненциальной и логистической кривых – элементарных моделей, лежащих в основе математического описания экономической конъюнктуры или распространения новой технологии. Мы подробно анализируем варианты эволюции системы с тепловыделением по катастрофическим сценариям, хотя такого же типа модели применяются и в экономике и социологии. Возможно, что будет полезной некоторая «диффузия» математических методов и моделей из инженерных дисциплин в экономические и социальные.

В конце книги приведен список литературы. Кроме прямых ссылок на источники, он включает перечень общей литературы по численному анализу и компьютерному моделированию, по-видимому, далеко не полный. Он содержит те относящиеся к нашей теме книги, которые были прочитаны авторами, показались им интересными и, следовательно, в явной или неявной форме были использованы при подготовке пособия. «Цель расчетов – не числа, а понимание» — этот девиз книги Р.Хемминга [ Ошибка! Источник ссылки не найден. ] мы хотели бы сделать главной идеей и нашей книги.

Видео:Решение первой начально-краевой задачи для одномерного уравнения теплопроводности.Скачать

Глава 11. Температурные волны

Видео:Регулярный режим. Численное решение дифференциальное уравнение теплопроводностиСкачать

11.1. Введение

Для оптимального управления пусковыми и переходными (т.е. нестационарными) режимами энергетических или теплотехнологических установок необходимо рассчитывать нестационарные температурные поля в элементах машин и оборудования. Прогнозирование температурного поля позволяет избежать недопустимого повышения температуры или возникновения слишком больших перепадов температуры. Характерным примером служит управление пуском мощной паровой турбины на тепловой электростанции. Естественное стремление оперативно ввести резервную энергетическую мощность наталкивается на ограничения – могут произойти недопустимые изменения осевых зазоров в проточной части турбины из-за неодинакового расширения или возникнуть недопустимые температурные напряжения в массивных деталях ротора и статора турбины.

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

Например, в проблеме Управляемого Термоядерного Синтеза плотность теплового потока на тепловоспринимающих твердых поверхностях может достигать 10 8 Вт/м 2 .

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

Большие потоки и высокие температуры возникают при лазерной или электронно-лучевой обработке деталей с целью упрочнения их поверхности. Похожие процессы имеют место при изготовлении микросхем.

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

Мы рассмотрим далее в качестве модели описанных выше процессов одномерную нестационарную задачу теплопроводности с внутренними источниками теплоты ( REF _Ref510596708 h Рис. 11 . 1 ).

Рис. 11 . 1 . Одномерная нестационарная задача теплопроводности

Предполагается, что пространственные изменения происходят только вдоль оси координат x . Боковая поверхность считается адиабатически изолированной. Если же требуется учесть теплоообмен на боковой поверхности, то теплосъем следует имитировать внутренним стоком, точно так же как в предыдущей главе.

Чтобы обеспечить универсальность модели, воспользуемся численным методом и разработаем для этой цели специальную Mathcad -программу. В ее основе будет лежать знаменитый в численном анализе метод прогонки – сверхбыстрый алгоритм для решения больших систем уравнений со специальной диагональной структурой.

Видео:Решение задачи Коши для уравнения теплопроводности (Часть 1)Скачать

11.2. Математическая формулировка задачи

В качестве исходной формулировки можно принять уравнение энергии, в котором следует исключить оператор конвективного переноса и учесть одномерность задачи, t = t ( x , τ ). В результате получится:

На левом и правом торце объекта (см. REF _Ref510596708 h Рис. 11 . 1 ) происходит тепловое взаимодействие со средой, и здесь необходимо задать граничные условия. Универсальным способом описать самые различные ситуации будет применение граничных условий третьего рода на левом ( x = 0) и правом ( x = L ) торцах объекта:

В этих соотношениях приравнены значения плотности теплового потока,

q поступающего из окружающей среды и вычисленного по уравнению Ньютона-Рихмана (правые части)

q и отводимого внутрь тела посредством теплопроводности и вычисленного по закону Фурье (левые части).

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

Видео:Mathcad-10. Пример: дифференциальные уравненияСкачать

11.3. Дискретное представление

Для численного интегрирования дифференциального уравнения в частных производных REF _Ref510601305 h ( 11 . 1 ) необходимо подготовить его конечно-разностную аппроксимацию. Мы получим нужный результат, записывая закон сохранения энергии для малого, но конечного контрольного объема. Здесь было бы полезно просмотреть еще раз аналогичные вычисления, которые уже проводились ранее (см.главы 1 и 9).

Рис. 11 . 2 . Контрольный объем XE «Контрольный объем» и потоки энергии

Итак, мы собираемся записать баланс энергии для выделенного контрольного объема ( REF _Ref510602514 h Рис. 11 . 2 ), сформированного в окрестности узла P . Последний находится сейчас в фокусе нашего внимания, но выведенные далее соотношения будут годиться для любого внутреннего узла.

От соседних, западного W и восточного E узлов поступают потоки теплоты за счет теплопроводности. Внутри контрольного объема выделяется теплота при действии внутреннего источника. В результате произойдет увеличение тепловой энергии в рассматриваемом объеме, что будет обнаружено по повышению температуры – от T 0 _P до T_P за время Δτ :

Величина баланса bal должна быть нулевой. Большими греческими «лямбда» Λ обозначены относительные значения коэффициента теплопроводности, вычисленные для контрольных поверхностей, проходящих между узлами. Например, для восточной поверхности:

,

где величина λ без индекса означает некоторое характерное значение теплопроводности.

Если коэффициент теплопроводности постоянен, то достаточно указать значение λ, а величинам Λ приписать единичные значения. Если же коэффициент теплопроводности сильно зависит от температуры или даже испытывает разрывы из-за слоистой структуры материала, среднегармонический способ осреднения REF _Ref510624678 h ( 11 . 4 ) обеспечивает вычисление потоков с хорошей точностью.

Уравнение REF _Ref510609786 h ( 11 . 3 ) записано для внутренних узлов. Ниже будет получено аналогичное уравнение для граничных узлов REF _Ref510673605 h ( 11 . 8 ) .

Неявная схема . Подчеркнем, что температуры в других узлах, т.е. T_W и T_E – неизвестные величины из «будущего», такие же как T_P . Поэтому REF _Ref510609786 h ( 11 . 3 ) представляет собой уравнение с тремя неизвестными. Всего таких уравнений можно записать столько же, сколько имеется неизвестных температур в узлах сетки. Для нахождения температур необходимо решить систему уравнений.

Можно повторить то же другими словами: не существует явного выражения для неизвестных, а требуется решать уравнения для их определения. Такие схемы в численном анализе называются неявными. Они обладают важным свойством устойчивости при счете, хотя и сложны из-за необходимости решать систему уравнений.

Явная схема. Если в уравнении REF _Ref510609786 h ( 11 . 3 ) разности температур во втором и третьем операторах (теплопроводности) записать для «прошлого» и снабдить обозначением «0», то получится явная схема: в каждом уравнении содержится одно единственное неизвестное значение T_P (в первом операторе, описывающем изменение энергии во времени). Программа и вычисления будут очень простыми, но если временной шаг превысит некоторое значение, при счете возникнут прогрессирующие паразитные осцилляции. Ограничения на шаг довольно обременительны, поэтому в настоящее время предпочитают неявные схемы.

Подробное обсуждение явной и неявной схем и демонстрационный вычислительный пример приведены в главе 4.

Вернемся к анализу уравнения баланса REF _Ref510609786 h ( 11 . 3 ) . Далее приводятся два варианта записи этого уравнения. Первое – в мнемонической форме, второе – в индексной, необходимой для составления программы. Неизвестные величины в этих уравнениях – температуры в узлах для следующего момента времени. Остальные, известные величины входят в коэффициенты уравнения. Среди этих известных величин – значение температуры в рассматриваемом узле T 0 _P , взятое с предыдущего шага по времени, либо из начальных условий.

Сравнивая выражения REF _Ref510609786 h ( 11 . 3 ) и REF _Ref510622313 h ( 11 . 5 ) , получают значения коэффициентов A , B , C , D . Эти операции можно выполнить в Mathcad ’е с помощью оператора collect , как показано ниже. Предварительно вводится обозначение для безразмерного шага по времени, т.н. сеточного числа Фурье:

Структура получившейся линейной системы уравнений с коэффициентами A , B , C , D оказывается трехдиагональной, как в демонстрационном примере для сетки из пяти узлов:

При большом числе узлов матрица коэффициентов окажется почти пустой. Например, для ста улов всего три процента ячеек будет занято числами, остальные останутся нулевыми. Ясно, что при использовании классического Гауссова исключения в основном будут перемножаться и складываться нули. Однако существует специальный вариант метода исключения, учитывающий трехдиагональную структуру матрицы и называемый методом прогонки [ Ошибка! Источник ссылки не найден. , Ошибка! Источник ссылки не найден. , Ошибка! Источник ссылки не найден. , Ошибка! Источник ссылки не найден. , Ошибка! Источник ссылки не найден. ]. Mathcad -программа этого метода приведена на REF _Ref510673830 h Рис. 11 . 5 .

Рис. 11 . 3 . Баланс энергии для поверхностного узла

Для граничных узлов контрольные объемы оказываются половинного размера, как показано на REF _Ref510625054 h Рис. 11 . 3 . Кроме того, для контрольной грани, попадающей на поверхность объекта, необходимо специальным образом записать тепловой поток:

Это уравнение является дискретным аналогом граничных условий, заданных выше соотношениями REF _Ref510626179 h ( 11 . 2 ) . Индекс « inner » означает ближайший внутренний узел: ( n -1) правого торца, 2 – для левого торца. Мы не будем повторять вычисления для коэффициентов в уравнениях для поверхностных узлов. Они полностью аналогичны тем, что выполнены для внутренних узлов. Окончательные выражения можно прочитать в тексте программы Coef ( REF _Ref510628079 h Рис. 11 . 4 ).

Видео:6-2. Метод сетокСкачать

11.4. Метод прогонки. Вычислительные программы Coef и SysTRD

Mathcad -программа Coef

Эта подпрограмма ( REF _Ref510628079 h Рис. 11 . 4 ) вычисляет коэффициенты A , B , C , D для системы уравнений, структура которой иллюстрируется соотношением REF _Ref510627923 h ( 11 . 7 ) .

Величина T 0 в списке формальных параметров означает вектор начальных значений температуры. Индексация начинается с единицы, т.е. значение переменной Origin среды Mathcad следует установить в единицу. Размерность вектора T 0 равна числу узлов сетки.

Величины Tf и Bi – векторы из двух элементов, задающие соответственно температуры среды и числа Био на левом и правом торце объекта ( REF _Ref510596708 h Рис. 11 . 1 ). Число Био – это безразмерный коэффициент теплоотдачи: Bi = α·Δx /λ.

Цикл for обеспечивает вычисления во внутренних узлах по формулам REF _Ref513089696 h ( 11 . 6 ) . Отдельно рассчитываются коэффициенты для граничных узлов, исходя из соотношения REF _Ref510627923 h ( 11 . 7 ) .

Параметр iTime сообщает подпрограмме, на каком шаге по времени находится процесс вычислений. Это важно для вычисления температуры среды, контактирующей с правым торцом. Задавая в векторе Pulse ненулевую амплитуду Ampl и некоторое значение частоты ν, мы сможем имитировать периодические тепловые воздействия на объект.

Программа Coef возвращает массив, собранный посредством функции augment из векторов A , B , C , D .

Рис. 11 . 4 . Подпрограмма для вычисления матрицы коэффициентов

Mathcad -программа SysTRD

Рис. 11 . 5 . Метод прогонки XE «Прогонки метод»

Коэффициенты A , B , C , D для трехдиагональной системы уравнений подготовлены подпрограммой Coef . В теле процедуры значения коэффициента D замещаются вычисленными значениями температуры в узлах сетки. Итак, процедура SysTRD возвращает вектор решений системы.

Mathcad -программа TimeHistory

Программа TimeHistory организует вычисления по временным шагам. Параметр nTime задает число шагов по времени. Функция TimeHistory возвращает распределения температуры по координате x для последовательных моментов времени iTime как столбцы матрицы F .

Рис. 11 . 6 . Главная программа

На каждом шаге по времени вызывается процедура Coef для вычисления коэффициентов A , B , C , D . Затем происходит обращение к решателю системы уравнений методом прогонки SysTRD . Обновляется вектор значений температуры T 0 на предыдущем временном шаге и добавляется столбец к матрице F , хранящей распределения температуры для последовательных моментов времени.

Видео:Метод Фурье для неоднородного уравнения теплопроводностиСкачать

11.5. Компьютерное моделирование периодических тепловых воздействий

Этап подготовки данных для расчета показан на REF _Ref510633371 h Рис. 11 . 7 . Предполагается исследовать температурные колебания в латунном стержне длиной 39 мм, если на правом торце задана температура, пульсирующая около среднего значения 800˚ C c амплитудой 320˚ C , а на левом торце – постоянная нулевая температура. Эти граничные условия имитируются заданием соответствующих температур жидкости и таких значений сеточного числа Био, которые соответствуют очень большим величинам коэффициентов теплоотдачи на правом и левом торцах объекта.

Подготовительные вычисления на REF _Ref510633371 h Рис. 11 . 7 понятны без комментариев. В конце этого фрагмента c с помощью функции matrix формируется вектор начальных температур T 0 , задающий равномерное начальное распределение 100˚ C , а также вектор значений относительного коэффициента теплопроводности Λ.

В данном примере величина Λ принимается постоянной. Если теплопроводность зависит от координат и/или температуры, то понадобится дополнительная процедура или вставка в программу Coef . Лучшим способом осреднения коэффициента теплопроводности будет вычисление среднегармонических значений по формуле REF _Ref510624678 h ( 11 . 4 ) . Мы оставляем эти усовершенствования программы как упражнение для читателей.

Рис. 11 . 7 . Пример расчета: исходные данные

Рис. 11 . 8 . Результаты расчета: температурные волны

Результаты вычислений показаны на REF _Ref510637802 h Рис. 11 . 8 . Правый график представляет серию распределений температуры вдоль оси x для различных моментов времени. Более наглядным является трехмерное представление на левой диаграмме. Плоскость в основании графика построена на продольной координате x (отметки 1..40) и временной координате (отметки номера временного шага в пределах 1..200). По вертикальной оси отложена температура. Отчетливо видны вынужденные колебания температуры на одном из торцов ( отметка 40 по координатной оси) и уменьшение амплитуды колебаний по мере приближения к противоположному торцу.

Эффектное представление нестационарного температурного поля можно получить с помощью анимации, как это описано в главе 9 .

Видео:6-1. Уравнение теплопроводностиСкачать

11.6. Заключение

Мы привели только один пример применения разработанной в этой главе программы для решения одномерного нестационарного уравнения теплопроводности. Численный эксперимент можно продолжить в следующих направлениях:

1. Исследовать влияние коэффициента теплоотдачи на температуру горячего торца ( на ее максимальное и среднее значение)
2. Исследовать влияние свойств материала на максимальную температуру поверхности и теплоотвод через стержень
3. Исследовать глубину проникания температурных волн, взяв длинный стержень и задав адиабатические условия на левом торце
4. Исследовать теплопередачу через стенку, если коэффициент теплоотдачи пульсирует во времени (для этого понадобится несколько модифицировать программу Coef , обеспечив вариации числа Био так же, как сейчас это сделано с температурой жидкости)
5. Исследовать температурные режимы стенок цилиндра двигателя внутреннего сгорания с воздушным охлаждением
6. Исследовать температурные режимы стен здания при погодных и сезонных изменениях температуры
7. И т.д.

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

Последние три главы книги посвящены реализации в Mathcad ’е численных методов для дифференциальных уравнений в частных производных. Достоинство разработанных на Mathcad ’е простых программ в том, что они демонстрируют основные алгоритмы численного анализа, будучи доступными и открытыми для модификации и экспериментов. Более широкое представление о численном анализе как фундаменте современных вычислительных программ дают включенные в список литературы книги: [ Ошибка! Источник ссылки не найден. , Ошибка! Источник ссылки не найден. , Ошибка! Источник ссылки не найден. , Ошибка! Источник ссылки не найден. , Ошибка! Источник ссылки не найден. , Ошибка! Источник ссылки не найден. , Ошибка! Источник ссылки не найден. , Ошибка! Источник ссылки не найден. , Ошибка! Источник ссылки не найден. , Ошибка! Источник ссылки не найден. , Ошибка! Источник ссылки не найден. , Ошибка! Источник ссылки не найден. , Ошибка! Источник ссылки не найден. , Ошибка! Источник ссылки не найден. , Ошибка! Источник ссылки не найден. ]

Далее мы приводим краткие сведения из справочной системы Mathcad по встроенным функциям для численного интегрирования дифференциальных уравнений. Мы покажем , как решается одна и та же система различными методами. В качестве примера используется уравнение Фолкнера –Скэн.

odesolve ( x , b , step )

Интегрирование одного дифференциального уравнения (первого, второго и т.д. порядка) при решении начальной или краевой задачи.

f = odesolve (x,b,step)

Решение возвращается как функция. Можно строить график функции f ( x ) или вычислять ее значения, например, f (6) (см. рисунок).

• x – переменная интегрирования
• b – конечная точка отрезка интегрирования
• step — необязательный параметр, число шагов интегрирования.

Параметры x , b , step должны быть действительными (не комплексными). Рекомендуется поэскпериментировать со значениями параметра step , чтобы оценить точность решения. Можно выбрать применяемый при работе odesolve метод интегрирования: с фиксированным шагом ( fixed step ) или с адаптирующимся шагом (а daptive step size ) , дважды щелкнув мышкой на слове odesolve и отмечая Fixed или Adaptive в контекстном меню.

Применение odesolve для решения дифференциального уравнения третьего порядка

rkfixed (IC, x1, x2, npoints, D)

Решение с фиксированным шагом интегрирования.

S = rkfixed (IC, x1, x2, npoints, D)

S — матрица, содержащая независимую переменную x в первом столбце, а в остальных столбцах — искомую функцию и ее первые n -1 производные.

• n –порядок системы дифференциальных уравнений
• IC — n –вектор зависимой переменной F
• D — определенная пользователем вектор-функция D ( x , F ), содержащая правые части дифференциальных уравнений
• x 1, x 2 — конечные точки отрезка интегрирования. Начальные значения в IC определены при x 1
• npoints – число шагов интегрирования; задает число строк (1 + npoints ) в матрице S .

Применение rkfixed для решения системы дифференциальных уравнений, эквивалентной уравнению Фолкнера-Скэн

rkadapt ( IC , x 1, x 2, acc , D , kmax , s )

Решение с адаптирующимся шагом интегрирования.

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

• IC – вектор начальных значений
• x 1, x 2 начальная и конечная точки интервала интегрирования
• acc – параметр, контролирующий точность решения. Если вычислительная погрешность на данном шаге превосходит асс, то шаг уменьшается
• D —n -элементная вектор-функция, задающая производные от неизвестных функций (т.е. правые части системы дифференциальных уравнений)
• kmax – максимальное число промежуточных точек, в которых будет определено решение. Значение kmax задает верхнюю границу числа строк матрицы S , возвращаемой функцией rkfixed .
• s — наименьшее допустимое расстояние между точками, в которых определяется решение. Эта величина определяет нижнюю границу различия между двумя значениями в первом столбце возвращаемой матрицы S .

Применение rkadapt для решения системы уравнений, эквивалентной уравнению Фолкнер-Скан

Rkadapt ( y , x 1, x 2, npoints , D )

Аргументы Rkadapt такие же, как у функции rkfixed

Rkadapt возвращает решение в равноотстоящих точках (всего в 1+ npoints точках). На каждом таком интервале Rkadapt вызывает функцию rkadapt , для которой параметр точности acc равен значению переменной среды Mathcad TOL .

Применение Rkadapt для решения системы дифференциальных уравнений, эквивалентной уравнению Фолкнер-Скан

bulstoer ( y , x 1, x 2, acc , D , kmax , s )

Bulstoer(y, x1, x2, npoints, D)

Решение гладких систем.

Аргументы и формат возвращаемого значения для этих функций такие же, как для rkadapt и Rkadapt соответственно.

В справочной системе Mathcad функция bulstoer рекомендована для гладких ( smooth ) систем, а rkadapt — для медленно меняющихся систем ( slowly varying systems ). Обычный пользователь вряд ли заметит разницу в применении этих функций. Сопоставление можно найти в книге [ Ошибка! Источник ссылки не найден. ]. Принципиальные различия методов интегрирования rkfixed , rkadapt , stiffr обсуждаются на вычислительных примерах в главе 4.

Stiffr (Y, x1, x2, npoints, D, J)

Stiffb(y, x1, x2, npoints, D, J)

stiffr (y,x1,x2,acc,D,J,kmax,s)

stiffb(y,x1,x2,acc,D,J,kmax,s)

Эти интеграторы специально разработаны для жестких систем.

Список аргументов аналогичен тому, какой имеют функции rkadapt ( y , x 1, x 2, acc , D , kmax , s ) или Rkadapt ( y , x 1, x 2, npoints , D ). Однако появляется дополнительный параметр J .

J есть функция, возвращающая n ·( n +1)-матрицу, первый столбец которой содержит вектор производных от правой части системы по независимой переменной, а остальные столбцы образуют матрицу Якоби для системы дифференциальных уравнений.

Stiffb использует метод Булирша — Штоера для жестких систем, Stiffr — метод Розенброка [ REF _ Ref 510683830 w h Ошибка! Источник ссылки не найден. ].

Решение жестких систем подробно рассматривается в главе 4. Там же подробно обсуждаются принципиальные различия методов интегрирования rkfixed , rkadapt , stiffr .

sbval (v, x1, x2, D, load, score)

Решение двухточечной краевой задачи.

sbval сводит двухточечную краевую задачу к задаче с начальными данными. Возвращает недостающие начальные условия на левой границе x 1.

• v ( ξ в примере) есть вектор начальных приближений для недостающих начальных условий при x 1
• x 1, x 2 — начальная и конечная точки интервала интегрирования
• D —n -элементная вектор-функция содержащая правые части дифференциальных уравнений
• Load (SetInit в примере) — n -вектор начальных значений при x 1. Некоторые значения будут константами, определенными пользователем при задании известных начальных условий (первые два элемента вектора SetInit в примере). Другие будут неизвестными вначале, но будут найдены функцией sbval . Если значение неизвестно, то используется соответствующее начальное приближение из v
• score ( discrepancy в примере) – вектор-функция, имеющая столько же элементов, что и v . Каждый элемент есть граничное условие при x 2, записанное в форме выражения, которое должно быть нулем ( F ₁-1 вместо F ₁=1, как в примере).

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

Функция sbval решает двухточечную краевую задачу методом пристрелки (или стрельбы), который в деталях обсуждается в главе 6.

Решение двухточечной краевой задачи для уравнения Фолкнера-Скэн методом sbval

bvalfit (v1, v2, x1, x2, xin, D, load1, load2, score)

Использует промежуточную точку для восстановления начальных условий. Возвращает вектор, содержащий те начальные условия, которые недоопределены в начальной точке x 1.

• v 1 – вектор начальных приближений для начальных значений, не заданных при x 1. v 2 — вектор начальных приближений для начальных значений, оставшихся не заданными при x 2.
• x 1, x 2 начальная и конечная точки интервала интегрирования
• xin – точка между x 1 и x 2, в которой траектории решений, начинающихся в x 1 и в x 2, должны совпасть.
• D – правые части уравнений, входящих в систему порядка n
• load 1 (в примере SetInitL ) – вектор-функция, n элементов которой являются значениями n искомых функций при x 1. Некоторые из этих значений предопределены заданными начальными условиями. Если же значение в начальной точке неизвестно (не задано), то используется начальное приближение из v 1.
• load 2(в примере SetInitR ) – аналогично load 1, но сказанное относится к точке x 2.
• score ( xin , y ) — n -элементная вектор-функция, определяющая условия согласования решений y (F в примере) в промежуточной точке xin . Обычно записывают score ( xin , y ) := y , чтобы приравнять значения всех искомых функций системы y по обеим сторонам границы xin .

Решение двухточечной XE «Двухточечная краевая задача» краевой задачи для уравнения Фолкнера XE «Фолкнера-Скэн уравнение» -Скан методом bvalfit XE » bvalfit «

Применение функции bvalfit XE » bvalfit » подробно рассмотрено в главе 10.

📹 Видео

1 Одно уравнениеСкачать

Решение неоднородного уравнения теплопроводностиСкачать

15. Решение уравнения теплопроводности в кругеСкачать

Решение задачи теплопроводности (Явная разностная схема)Скачать

Решение уравнения теплопроводности в одномерной постановке в Excel с применением неявной схемыСкачать

Решение первой краевой задачи для неоднородного уравнения теплопроводности.Скачать

Вычисление простых выражений в MathCAD 14 (4/34) Часть 1Скачать

Одномерное уравнение теплопроводности. Виды краевых задачСкачать