Решить дифференциальное уравнение методом Рунге-Кутта 4 порядка (ode45)

@Moro1 · Регистрация: 31.12.2013

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

Добрый день!
Необходимо решить дифференциальное уравнение методом Рунге-Кутта 4 порядка.
Уравнение подготовил, а вот как реализовать его через функцию ode45, с построением графика не знаю.
Подскажите пожалуйста, какой код необходимо дописать?
Код программы:

Matlab M
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
e=15;
V1=10;
t=(1:0.1:30);
V2=1;
f=10000;
omega=2*pi*f;
k=1.36*10-23;
Lk=0.000001;
R=1600;
T=222;
B=222;
Vc=0.92837;
R=1;
Vn=sqrt(4*k*R*T*B);
function q(2)=solverDE(t,y);
f=@(t,q)([q(2); V1+V2*sin*(omega*t)-Lk*q(2)-R*q(1)-Vn-Vc*sin((2*pi*q(1))/2*e)/Lk]);

Зосима · 20.02.2014, 13:21

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

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

Обратите внимание, когда 0 умножается на синус (специально обвел линией), ну никак не может последний остаться.

я обратил внимание, а также и на то, что именно строится, поэтому и написал, то что написал

а разница порой бывает мало различима, особенно, если масштаб по осям берется автоматом

@Moro1 · 22.02.2014, 02:40 **[ТС]**

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

Matlab M
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
function dq = zod(t,q)
global L R Va Vc T
 
for i = 1:length(T)
    DQ(i,:)=zod(T(i),Q(i,:));
    U(i)=V(T(i),Q(i,1));
end
 
dq=zeros(2,1);
dq(1,1)=q(2);
dq(2,1)=1/L*((Va*t+Vc*sin(pi*q(1))/1.6*10^-19)-R*q(2));
end
 
function y=V(t,q);
global Vc Va
y=(Va*t+Vc*sin(pi*q)/1.6*10^-19);
end

программа расчета

Matlab M
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
clear, clc
global L R Va Vc
 
L=1*10^-7;
R=1*10^6;
Va=1*10^-10;
Vc=1*10^-4;
 
t=0,5e-13;
q0=[0 0];
[T, Q] = ode45('zod',t,q0);
 
subplot(4,1,1); plot(T,Q(:,1))
subplot(4,1,2); plot(T,DQ(:,1))
subplot(4,1,3); plot(T,DQ(:,2))
subplot(4,1,4); plot(U,DQ(:,1))
legend('dqdt(V(t))')

Добавлено через 4 часа 39 минут
Три графика адаптировать удалось, не получается получить зависимость I(V). Не знаю как адаптировать строчку U(i)=V(T(i),Q(i,1)); Ниже программа.
Функция:

Matlab M
1
2
3
4
5
6
function dq = zod(t,q)
global L R Va Vc 
dq=zeros(2,1);
dq(1)=q(2);
dq(2)=1/L*((Va*t+Vc*sin(pi*q(1))/1.6*10^-19))-R*q(2);
end

Рассчет:

Matlab M
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
function zaq
q0=[0 0];
global L R Va Vc
 
L=1*10^-7;
R=1*10^6;
Va=1*10^-10;
Vc=1*10^-4;
[T, Q] = ode45('zod',[0,5e-13],q0);
for i = 1:length(T)
    DQ(i,:)=zod(T(i),Q(i,:));
end
subplot(4,1,1); plot(T,Q(:,1))
subplot(4,1,2); plot(T,DQ(:,1))
subplot(4,1,3); plot(T,DQ(:,2))
subplot(4,1,4); plot(U,DQ(:,1))
legend('dqdt(V(t))')
end

Добавлено через 20 часов 54 минуты
Зосима, Спасибо, разобрался. Такой вопрос. Можно ли решить ДУ минуя вариант:

Matlab M
1
2
3
function dd=func1(x,y) % 
dd=-2*x*y; % 
>> [T,Y]=ode45(@func1,[0,2],1);

libo cherez

Matlab M
1
f=@(t,q)([q(2); V1+V2*sin*(omega*t)-Lk*q(2)-R*q(1)-Vn-Vc*sin((2*pi*q(1))/2*e)/Lk]);

В обоих случаях из решения ОДУ 2 порядка появляется синус, возможно ли как то исключить его?

Добавлено через 2 часа 10 минут
всеже с физическои точки зрения колебания не возможны, обьясните пожалуиста их природу. Откуда они возникают в уравнении второго порядка?

Добавлено через 8 часов 41 минуту
Я не понимаю, что я не так делаю, пытаюсь ввести интерполяцию на выходе
Error using chckxy (line 42)
There should be at least two data points.

Error in spline (line 54)
[x,y,sizey,endslopes] = chckxy(x,y);

Matlab M
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
function tempzx
[T, Q] = ode45(@zod, [0:1e-8:1e-7],[0 0]);
 
for j=1:length(T)
qi=0:1e-9:1e-7; 
yy=spline(T(j),Q(:,j),qi)         
end
 
for i = 1:length(T)
    DQ(i,:) = zod(T(i),Q(i,:));
    U(i) = V(T(i), Q(i,1));
end
 
subplot(4,1,1); plot(T,Q(:,1),T,YY)
subplot(4,1,2); plot(T,DQ(:,1))
subplot(4,1,3); plot(T,DQ(:,2))
subplot(4,1,4); plot(U,DQ(:,1))
legend('dqdt(V(t))')
end

Зосима · 03.03.2014, 16:40

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

for j=1:length(T)
qi=0:1e-9:1e-7;
yy=spline(T(j),Q(:,j),qi)
end

не, не! тут цикл не нужен:

Matlab M
qi=0:1e-9:1e-7; 
yy=spline(T,Q(:,1),qi);

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

Откуда они возникают в уравнении второго порядка?

Если характеристическое уравнение L^2+p*L+q=0 имеет сопряженные комплексные корни L_1,2 = a+-i*b, (дискриминант D<0), то общее решение однородного уравнения принимает вид:
y(x) = exp(a*x)*( C1*cos(b*x)+C2*sin(b*x) ), где C1, C2 – константы.
Вот так!

А если ты напряжение считаешь сразу внутри системы: dq(2,1)=1/L*((Va*t+Vc*sin(pi*q(1))/1.6*10^-19)-R*q(2)); то никак его потом не получить, кроме как опять переписать:

Matlab M
[T, Q] = ode45('zod',t,q0);
U = (Va*T+Vc*sin(pi*Q(:,1))/1.6*10^-19);

@Moro1 · 25.03.2014, 17:37 **[ТС]**

Скажите а зависимость DV/DI(I) каким образом можно построить? Я уге по всякому пробовал: и в цикле пытался записать

Matlab M
1
2
3
4
5
for i = 1:length(T)
    DQ(i,:)=zod(T(i),Q(i,:));
    U(i)=V(T(i),Q(i,1));
    DV(i)=(DQ(:,1),U(i));
end

и отделно формулой
U = ((Va*T+Vc*sin(pi*Q(:,1))/1.6*10^-19))/DQ(:,1);
Ничего не получается

(

@Moro1 · 26.03.2014, 18:53 **[ТС]**

Пытаюсь построить зависимость так, но она оказывается идентична зависимости I/V

Matlab M
1
2
3
4
5
6
7
8
9
10
11
for i = 1:length(T)
    DQ(i,:)=zod(T(i),Q(i,:)); 
end
 
U=Va*T+Vc*sin((pi*Q(:,1))/1.6*10^-19);
 
DV=U/DQ(i,1);
 
subplot(4,1,3); plot(U,DQ(:,1))
legend('U(T)');
subplot(4,1,4); plot(DV,DQ(:,1))

@Moro1 · 26.03.2014, 21:05 **[ТС]**

Вопрос почему зависимость U/T монотонно возрастающая при sin?

Новые блоги и статьи Все статьи Все блоги /
Программный домашний кинотеатр russiannick 27.09.2026 Сподобился на программный домашний кинотеатр. В качестве ЯВУ по традиции выбрал js. В помощники взял Яндекс-Алису. Было создано три зала на разные интересы. исторические и ретро сериал Хичкок. . .	Беседа с ИИ о программистах, недопускающих к созданию и правке кода генеративные ИИ и причины этого zorxor 21.09.2026 Раньше я радовался или получал некоторые эмоции, пусть небольшие, но всё же, от самого процесса написания кода, рекомпиляции и запуска, видя постепенное развитие программы и прочее. А теперь лень. . .	Мобильное приложение ColorStep pavlinmavlin 17.09.2026 Реализовал приложение Красный, Зеленый, Синий в Unity3d + c#. Название изменил на ColorStep. Приложение прошло модерацию и теперь доступно для скачивания. Делал его сам, шаг за шагом — и вот,. . .	Запрет дублирования строк в табличной части Maks 13.09.2026 Реализация из решения ниже выполнена на нетиповом справочнике "Нормы ТО" с табличной часть "Виды ТО", разработанного в КА2, со следующими реквизитами: - ВидТО (СправочникСсылка. ВидыТО); - ВидГСМ. . .
Скрипты Tampermonkey для CyberForum, ChatGPT, Claude и пр. Jin X 06.09.2026 Скрипты Tampermonkey для CyberForum, ChatGPT, Claude и пр. Работая с форумом и нейросетями в браузере часто хочется что-то подкорректировать или добавить какого-то функционала. Ниже прикреплён. . .	Программа опроса у.з. расходомера SLS-720F Argus19 02.09.2026 Программа опроса у. з. расходомера SLS-720F Программа опрашивает один раз в минуту три ультразвуковых расходомера SLS-720F через интерфейс RS-485 по протоколу Modbus RTU. Опрашиваются регистры. . .	Hyper-V: Компьютер должен поддерживать доверенный платформенный модуль 2.0. Maks 31.08.2026 При установке Windows 11 на виртуальную машину Hyper-V 2-го поколения вылезла такая ошибка: Решение: в параметрах виртуальной машины, в разделе "Безопасность" (Security) активировать флаг. . .	Архитектура биовида Стива в Майнкрафте: Зачем бонобо кубический каннибализм anaschu 30.08.2026 Кубический Вагинокапитализм в Minecraft: Математический инвариант ОДУ и рок Стивов-бонобо Главная задача разработанной «Модели Всего» — наглядно продемонстрировать наличие системной «судьбы». . .

@Moro1 0 / 0 / 0 Регистрация: 31.12.2013 Сообщений: 48
	26.03.2014, 21:05 [ТС]
	Вопрос почему зависимость U/T монотонно возрастающая при sin? Миниатюры 0