Решить систему двух нелинейных ОДУ

@Dimarik21 · Регистрация: 18.02.2015

Студворк — интернет-сервис помощи студентам

Господа не нашел похожего в ваших темах. Может кто знает, как решается данная система методом Рунге-Кутта. Посмотрите мой код, может увидите, что там неправильно. Здесь коэффициенты A,B,P(большое), Лямбда, D и Альфа являются некими константами.
$https://www.cyberforum.ru/cgi-bin/latex.cgi?\left\{\begin{matrix}\frac{dQ}{dl}=\frac{-A+\sqrt{A^2+4*B*(P^2-p^2)}}{2B}\\ \frac{dp}{dl}=\(\lambda \frac{Q}{D}+2\frac{dQ}{dl})*\frac{\alpha *Q*p}{2*\alpha *Q^2-p}\\ \end{matrix}\right.$

Matlab M
1
2
3
4
5
6
7
8
9
10
11
12
13
function dydx=Func1(t, y)
%FUNC1 Summary of this function goes here
%   Detailed explanation goes here
ag = 0.013;
bg = 6.17*10^(-8);
P = 20.7;
L = 2.4;
D = 0.177;
Al = 1;
dydx = zeros(2,1);
dydx(1) = (-ag+(ag^2+4*bg*(P^2-y(2)^2))^0.5)/(2*bg);
dydx(2) = (L*y(1)/D+2*dydx(1))*Al*y(1)*y(2)/(2*Al*y(1)^2-y(2));
end

Matlab M
1
2
>> [T H] = ode45(@Func1,[0 100],[0 20]);
>> plot(T,H);

Начальные условия при L=0 Q=0 P=20, нужна графическая зависимость на L от 0 до 100.

Зосима · 18.02.2015, 20:28

Dimarik21, вроде все верно, но система неустойчива

Может в исходной формуле ошибка?

@Dimarik21 · 18.02.2015, 20:33 **[ТС]**

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

Зосима · 19.02.2015, 14:32

Dimarik21, странно... возможно дело в константах, единицы измерений перевели в СИ ?

@Dimarik21 · 19.02.2015, 15:24 **[ТС]**

Зосима, коэффициенты да, в системе СИ за исключением давления P. Тогда:

Matlab M
1
2
3
4
5
6
7
8
9
10
11
12
13
function dydx=Func1(t, y)
%FUNC1 Summary of this function goes here
%   Detailed explanation goes here
ag = 0.013;
bg = 6.17*10^(-8);
P = 20.7*10^6;
L = 2.4;
D = 0.177;
Al = 1;
dydx = zeros(2,1);
dydx(1) = (-ag+(ag^2+4*bg*(P^2-y(2)^2))^0.5)/(2*bg);
dydx(2) = (L*y(1)/D+2*dydx(1))*Al*y(1)*y(2)/(2*Al*y(1)^2-y(2));
end

вводя следующее:

Matlab M
1
ode45(@Func1,[0 100],[0 20]);

полученыый график улетает куда-то в космос. ВОбщем я не знаю что делать. СИжу ,, думаю который день уже

Зосима · 19.02.2015, 16:29

Dimarik21, странно...

если в функции давление P=20.7*10^6, то и нач.условия д.б. [0, 20e6] (правда это не особо помогает)
Я попробовал посмотреть начальные значения решения, взял первые 130 отсчетов:

Matlab M
1
2
3
[L U] = ode45(@func1,[0 100],[0 20e6]);
n = 130; %кол-во отсчетов
plot(L(1:n),U(1:n,2))

получается так:

т.е. давление падает, а потом начинает расти, причем матлаб ругается на комплексные данные (откуда они там??? комплексные числа могут получится только если берется корень отрицательного числа, но давление y(2) не превосходит P, поэтому ОТКУДА???

)

@Dimarik21 · 19.02.2015, 16:51 **[ТС]**

Зосима, да правильно давление должно падать о длине, а дебит наоборот расти. У меня тоже кстати ругается на комплексные числа. На самом деле достаточно точек 5 на длину в 100 метров.

@Dimarik21 · 19.02.2015, 17:03 **[ТС]**

Зосима, пробую сейчас в МатКаде написать программу еще. С МатКадом не знакомы случайно, ошибка вылезла не могу понять как исправить.

Зосима · 19.02.2015, 17:25

*попробуй в Odesolve написать y1(x) y2(x)

@Dimarik21 · 19.02.2015, 17:28 **[ТС]**

Не помогает, все также

Зосима · 19.02.2015, 17:45

увы, я с маткадом не дружу

@Dimarik21 · 19.02.2015, 18:17 **[ТС]**

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

Зосима · 19.02.2015, 20:27

твою функцию func1(последний вариант) я не менял, а чтобы получить такой график:

Matlab M
1
2
3
[L U] = ode45(@func1,[0 100],[0 20e6]);
n = 130; %кол-во отсчетов
plot(L(1:n),U(1:n,2))

Будь другом, напиши размерности ВСЕХ величин и констант!

@Dimarik21 · 19.02.2015, 23:46 **[ТС]**

Размерности коэффициентов:
ag МПа²*сут/тыс.м³
bg МПа²*сут²/(тыс.м³)²
P МПа
D метр
L безразмерный
Al безразмерный, его можно взять 10^8

Зосима · 21.02.2015, 11:11

Сообщение от Dimarik21

ag МПа2*сут/тыс.м3
bg МПа2*сут2/(тыс.м3)2
P МПа
D метр
L безразмерный
Al безразмерный

О! замечательно, теперь можно проверить формулы

Только я для проверки переведу единицы в СИ:

ag $https://www.cyberforum.ru/cgi-bin/latex.cgi?\left[ \frac{Pa^2s}{m^3}\right]$

bg $https://www.cyberforum.ru/cgi-bin/latex.cgi?\left[ \frac{Pa^2s^2}{m^6}\right]$

P $https://www.cyberforum.ru/cgi-bin/latex.cgi?\left[ Pa \right]$

D, l $https://www.cyberforum.ru/cgi-bin/latex.cgi?\left[ m \right]$

Q $https://www.cyberforum.ru/cgi-bin/latex.cgi?\left[ \frac{m^3}{s} \right]$

Тогда для первого уравнения:

$https://www.cyberforum.ru/cgi-bin/latex.cgi?\large \frac{\frac{m^3}{s}}{m} = \frac{\frac{Pa^2s}{m^3} + \sqrt{\frac{Pa^4s^2}{m^6}+\frac{Pa^2s^2}{m^6}Pa^2} }{\frac{Pa^2s^2}{m^6}}; \\\frac{m^2}{s} = \frac{\frac{Pa^2s}{m^3}}{\frac{Pa^2s^2}{m^6}};\\\frac{m^2}{s} \neq \frac{m^3}{s};$

Видишь, правая и левая части не сходятся! (Правую часть нужно еще поделить на метр)

И для второго:

$https://www.cyberforum.ru/cgi-bin/latex.cgi?\large \frac{Pa}{m} = \left( \frac{\frac{m^3}{s}}{m} + \frac{\frac{m^3}{s}}{m}\right) \cdot \frac{\frac{Pa^2s}{m^3}\cdot\frac{m^3}{s}\cdot Pa}{\frac{Pa^2s}{m^3}\cdot\frac{m^6}{s^2}-Pa}; \\\frac{Pa}{m} \neq \left(\frac{m^2}{s}\right) \cdot \frac{Pa^3}{\frac{Pa^2m^3}{s}-Pa};$

Тут в знаменателе единицы совсем не пляшут

одинокое p все портит.

В качестве подгониана я бы предложил такую запись: $https://www.cyberforum.ru/cgi-bin/latex.cgi?2\alpha Q^2 - Qp^2$
тогда единицы справа и слева сойдутся, а дробь примет вид:

$https://www.cyberforum.ru/cgi-bin/latex.cgi?\frac{dp}{dl} = \left( \lambda \frac{Q}{D} + 2\frac{dQ}{dl}\right) \cdot \frac{\alpha p}{2\alpha Q - p^2}$

но не уверен, что такая запись корректна

@Dimarik21 · 21.02.2015, 12:27 **[ТС]**

Зосима, сейчас перепроверю все размерности коэффициентов по литературам.

Зосима · 21.02.2015, 14:01

Проверка лишней не будет, но лично мое заключе: ошибка в формуле (причем не "просто знак перепутали")

@Dimarik21 · 21.02.2015, 14:54 **[ТС]**

В разной литературе, разные формулы. Нашел правильную ГОСТовскую запись этих уравнений.
$https://www.cyberforum.ru/cgi-bin/latex.cgi?\left\{\begin{matrix}\frac{dQ}{dl}=\frac{-A+\sqrt{A^2+4*B*(P^2-p^2)}}{2*B}\\ \frac{dP}{dl}=(\lambda \frac{Q}{D}+2\frac{dQ}{dl})*\frac{\alpha *Q*p}{2*\alpha *Q^2-p^2}\end{matrix}\right.$
Здесь \alpha =20 безразмерный параметр
A=8 безразмерный
В=0,5 1/м^2
D=0.177 м
\lambda =0,003 безразмерный
P=20 МПа
l метры
Q м^3/сут

Добавлено через 2 минуты
Зосима, вроде правильная, но все может быть. Уравнение написано так, чтобы не обязательно было в систему СИ все переводить. МПа надо оставить также, Q тоже, вобщем все коэффициенты так как написаны должны быть без всяких домножений на 10^какой-то степени

Зосима · 21.02.2015, 15:22

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

- в первом под корнем сумма безразмерной величины и Па^2/м^2 (а в результате еще должны откуда то получится сутки!)
- во втором опять же разность в знаменателе (м^3/сут)^2 и Па^2

так что "воз и ныне там"

Не по теме:

я не вредный, но стремлюсь найти истину :)

@Dimarik21 · 21.02.2015, 15:27 **[ТС]**

Зосима, вот эти безразмерные параметры они как раз и учитывают различные размерности внутри себя.

Новые блоги и статьи Все статьи Все блоги /
PhpStorm 2025.3: WSL Terminal всегда стартует в ~ and_y87 14.12.2025 PhpStorm 2025. 3: WSL Terminal всегда стартует в ~ (home), игнорируя директорию проекта Симптом: После обновления до PhpStorm 2025. 3 встроенный терминал WSL открывается в домашней директории. . .	Access VikBal 11.12.2025 Помогите пожалуйста !! Как объединить 2 одинаковые БД Access с разными данными.	Новый ноутбук volvo 07.12.2025 Всем привет. По скидке в "черную пятницу" взял себе новый ноутбук Lenovo ThinkBook 16 G7 на Амазоне: Ryzen 5 7533HS 64 Gb DDR5 1Tb NVMe 16" Full HD Display Win11 Pro	Музыка, написанная Искусственным Интеллектом volvo 04.12.2025 Всем привет. Некоторое время назад меня заинтересовало, что уже умеет ИИ в плане написания музыки для песен, и, собственно, исполнения этих самых песен. Стихов у нас много, уже вышли 4 книги, еще 3. . .	От async/await к виртуальным потокам в Python IndentationError 23.11.2025 Армин Ронахер поставил под сомнение async/ await. Создатель Flask заявляет: цветные функции - провал, виртуальные потоки - решение. Не threading-динозавры, а новое поколение лёгких потоков. Откат?. . .
Поиск "дружественных имён" СОМ портов Argus19 22.11.2025 Поиск "дружественных имён" СОМ портов На странице: https:/ / norseev. ru/ 2018/ 01/ 04/ comportlist_windows/ нашёл схожую тему. Там приведён код на С++, который показывает только имена СОМ портов, типа,. . .	Сколько Государство потратило денег на меня, обеспечивая инсулином. Programma_Boinc 20.11.2025 Сколько Государство потратило денег на меня, обеспечивая инсулином. Вот решила сделать интересный приблизительный подсчет, сколько государство потратило на меня денег на покупку инсулинов. . . .	Ломающие изменения в C#.NStar Alpha Etyuhibosecyu 20.11.2025 Уже можно не только тестировать, но и пользоваться C#. NStar - писать оконные приложения, содержащие надписи, кнопки, текстовые поля и даже изображения, например, моя игра "Три в ряд" написана на этом. . .	Мысли в слух kumehtar 18.11.2025 Кстати, совсем недавно имел разговор на тему медитаций с людьми. И обнаружил, что они вообще не понимают что такое медитация и зачем она нужна. Самые базовые вещи. Для них это - когда просто люди. . .	Создание Single Page Application на фреймах krapotkin 16.11.2025 Статья исключительно для начинающих. Подходы оригинальностью не блещут. В век Веб все очень привыкли к дизайну Single-Page-Application . Быстренько разберем подход "на фреймах". Мы делаем одну. . .

Зосима 5245 / 3573 / 379 Регистрация: 02.04.2012 Сообщений: 6,477 Записей в блоге: 18
	18.02.2015, 20:28
	Dimarik21, вроде все верно, но система неустойчива Может в исходной формуле ошибка? 0

@Dimarik21 0 / 0 / 0 Регистрация: 18.02.2015 Сообщений: 18
	18.02.2015, 20:33 [ТС]
	Система верна, сомнений нет. Это даже упрощенная версия системы, в настоящих формулах те постоянные коэффициенты на самом деле не постоянные, а также изменяются, но для упрощения они приняты пока постоянными. На графике выходит ерунда полнейшая. Уравнения представляют собой распределение давления и Дебита газовой скважины по стволу горизонтального ствола. 0

Зосима 5245 / 3573 / 379 Регистрация: 02.04.2012 Сообщений: 6,477 Записей в блоге: 18
	19.02.2015, 14:32
	Dimarik21, странно... возможно дело в константах, единицы измерений перевели в СИ ? 0

@Dimarik21 0 / 0 / 0 Регистрация: 18.02.2015 Сообщений: 18
	19.02.2015, 16:51 [ТС]
	Зосима, да правильно давление должно падать о длине, а дебит наоборот расти. У меня тоже кстати ругается на комплексные числа. На самом деле достаточно точек 5 на длину в 100 метров. 0

@Dimarik21 0 / 0 / 0 Регистрация: 18.02.2015 Сообщений: 18
	19.02.2015, 17:03 [ТС]
	Зосима, пробую сейчас в МатКаде написать программу еще. С МатКадом не знакомы случайно, ошибка вылезла не могу понять как исправить. Миниатюры 0

Зосима 5245 / 3573 / 379 Регистрация: 02.04.2012 Сообщений: 6,477 Записей в блоге: 18
	19.02.2015, 17:25
	попробуй в Odesolve написать y1(x) y2(x)* 0

@Dimarik21 0 / 0 / 0 Регистрация: 18.02.2015 Сообщений: 18
	19.02.2015, 17:28 [ТС]
	Не помогает, все также 0

Зосима 5245 / 3573 / 379 Регистрация: 02.04.2012 Сообщений: 6,477 Записей в блоге: 18
	19.02.2015, 17:45
	увы, я с маткадом не дружу 0

@Dimarik21 0 / 0 / 0 Регистрация: 18.02.2015 Сообщений: 18
	19.02.2015, 18:17 [ТС]
	Зосима, а не могли бы написать исходный код, работая с которым получили такие графики 0

@Dimarik21 0 / 0 / 0 Регистрация: 18.02.2015 Сообщений: 18
	19.02.2015, 23:46 [ТС]
	Размерности коэффициентов: ag МПа²сут/тыс.м³ bg МПа²сут²/(тыс.м³)² P МПа D метр L безразмерный Al безразмерный, его можно взять 10^8 0

@Dimarik21 0 / 0 / 0 Регистрация: 18.02.2015 Сообщений: 18
	21.02.2015, 12:27 [ТС]
	Зосима, сейчас перепроверю все размерности коэффициентов по литературам. 0

Зосима 5245 / 3573 / 379 Регистрация: 02.04.2012 Сообщений: 6,477 Записей в блоге: 18
	21.02.2015, 14:01
	Проверка лишней не будет, но лично мое заключе: ошибка в формуле (причем не "просто знак перепутали") 0

@Dimarik21 0 / 0 / 0 Регистрация: 18.02.2015 Сообщений: 18
	21.02.2015, 14:54 [ТС]
	В разной литературе, разные формулы. Нашел правильную ГОСТовскую запись этих уравнений. $https://www.cyberforum.ru/cgi-bin/latex.cgi?\left\{\begin{matrix}\frac{dQ}{dl}=\frac{-A+\sqrt{A^2+4B(P^2-p^2)}}{2B}\\ \frac{dP}{dl}=(\lambda \frac{Q}{D}+2\frac{dQ}{dl})\frac{\alpha Qp}{2\alpha Q^2-p^2}\end{matrix}\right.$ Здесь \alpha =20 безразмерный параметр A=8 безразмерный В=0,5 1/м^2 D=0.177 м \lambda =0,003 безразмерный P=20 МПа l метры Q м^3/сут Добавлено через 2 минуты Зосима, вроде правильная, но все может быть. Уравнение написано так, чтобы не обязательно было в систему СИ все переводить. МПа надо оставить также, Q тоже, вобщем все коэффициенты так как написаны должны быть без всяких домножений на 10^какой-то степени 0

Зосима 5245 / 3573 / 379 Регистрация: 02.04.2012 Сообщений: 6,477 Записей в блоге: 18
	21.02.2015, 15:22
	да, можно оставить МПа и сутки, но речь сейчас не об этом. При таком раскладе единиицы тоже не пляшут, нельзя складывать килограммы и секунды! - в первом под корнем сумма безразмерной величины и Па^2/м^2 (а в результате еще должны откуда то получится сутки!) - во втором опять же разность в знаменателе (м^3/сут)^2 и Па^2 так что "воз и ныне там" Не по теме: я не вредный, но стремлюсь найти истину :) 0

@Dimarik21 0 / 0 / 0 Регистрация: 18.02.2015 Сообщений: 18
	21.02.2015, 15:27 [ТС]
	Зосима, вот эти безразмерные параметры они как раз и учитывают различные размерности внутри себя. 0

Решить систему двух нелинейных ОДУ

Решение