1 / 1 / 0
Регистрация: 28.01.2017
Сообщений: 32
1

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

28.01.2017, 11:20. Показов 1554. Ответов 33
Метки нет (Все метки)

Здравствуйте, мне нужно сделать программу, которая будет решать различные системы дифуров, заданные мною же.
У меня есть, мною же написанная программа, решающая дифуры по отдельности, методом рунге кута 4го порядка, программа выводит точки удовлетворяющие решению и строит по ним график функции.
Вот пример для y'=y :
Добавить функционал в готовую программу решения систем дифферинциальных уравнений

Решением является функция C*e^x , как видно на скрине, точки удовлетворяют этой функции, ну в общем, это и есть решение.

Как можно этот метод прикрутить для решения сразу системы?

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

В общем, может я что-то не догоняю или вообще мыслю не правильно, помогите пожалуйста, направьте на путь истинный!

Если что, вот функция решающая, возвращает массив точек, которые и выводятся в таблице:
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
double *RungeKutt(int n,double x0,double y0,double h,double (*f)(double x,double y))
{
  register int i;
  double *res,s;
  double k1,k2,k3,k4;
  res=new double [n+1];
  res[0]=y0;
  s=x0;
  for(i=1;i<=n;i++)
  {
    k1=h*f(s,res[i-1]);
    k2=h*f(s+h/2,res[i-1]+k1/2);
    k3=h*f(s+h/2,res[i-1]+k2/2);
    k4=h*f(s+h,res[i-1]+k3);
 
    res[i]=res[i-1]+(k1+2*k2+2*k3+k4)/6;
    s+=h;
  }
  return res;
}
__________________
Помощь в написании контрольных, курсовых и дипломных работ, диссертаций здесь
0
Лучшие ответы (1)
Programming
Эксперт
94731 / 64177 / 26122
Регистрация: 12.04.2006
Сообщений: 116,782
28.01.2017, 11:20
Ответы с готовыми решениями:

Помогите найти программу решения систем алгебраических уравнений методом Зейделя
Ребята, помогите найти программу решения систем алгебраических уравнений методом Зейделя. Код нужен...

Готовые библиотеки для решения систем уравнений
Методом Гаусса. Приведением к диагональному виду. Разрядность 32 бита, при обращении элементов...

Библиотеки для решения (недоопределённых) систем линейных уравнений
Знает кто-нибудь сабж?. Если система недоопределена, то нужно в некотором формате отдавать её общее...

Разработать функцию для решения систем линейных уравнений методом Зейделя
Ребят, срочно, выручайте, задание на зачет... Прототип функции int slau(int n, double e, const...

33
1269 / 1026 / 470
Регистрация: 25.12.2016
Сообщений: 3,333
28.01.2017, 11:33 2
Цитата Сообщение от andrew95qq Посмотреть сообщение
Как можно этот метод прикрутить для решения сразу системы?
Формулы будут те же самые, только теперь h1, h2, h3, h4, y, f будут не скалярами, а векторами. Длина векторов равна числу уравнений. Подробности есть в википедии.
2
1 / 1 / 0
Регистрация: 28.01.2017
Сообщений: 32
28.01.2017, 11:36  [ТС] 3
Я так понимаю, не h1,h2,h3,h4, а k1,2,3,4 ? Или что вы имели ввиду?
0
1269 / 1026 / 470
Регистрация: 25.12.2016
Сообщений: 3,333
28.01.2017, 11:41 4
Цитата Сообщение от andrew95qq Посмотреть сообщение
не h1,h2,h3,h4, а k1,2,3,4
Да, разумеется.
1
1 / 1 / 0
Регистрация: 28.01.2017
Сообщений: 32
28.01.2017, 12:13  [ТС] 5
Более менее разобрался вот по такому примеру:
Добавить функционал в готовую программу решения систем дифферинциальных уравнений

Но, не могу понять, в моем примере, где решалось просто одно ДУ, функция возвращала 1 массив чисел,
ну и в итоге была таблица из 2х колонок, Х, где значения были от Х0 и далее +h, ну и Y, там у нас был тот самый массив чисел, полученных при решении, получались точки (x;y) и мы по ним строили график функции решения.
А тут получается, если система состоит из 2х ДУ, то функция будет возвращать 2 массива чисел, из примера, Z и Y, соответственно, если будет 3 ДУ, то будет 3 массива чисел? Как тогда строить график? и что выводить?
Или я чего то не догоняю?)
0
1269 / 1026 / 470
Регистрация: 25.12.2016
Сообщений: 3,333
28.01.2017, 12:21 6
Вот, нашёл свой старый исходник. Здесь использован собственный тип VectorD с перегрузкой операции умножения на скаляр. F имеет тип VectorFn - вектор функций с перегруженной операцией вызова функции.
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
void ProcessRungeKuthMethod3( VectorFn& F, VectorD& Y,
                              double x, double xmax, double h )
{
    VectorD K1, K2, K3, K4; double x0 = x;
    while (x < xmax)
    {
        K1 = h * F( x, Y );
        K2 = h * F( x + h/2, Y + 0.5*K1 );
        K3 = h * F( x + h/2, Y + 0.5*K2 );
        K4 = h * F( x + h, Y + K3 );
        Y += 1./6*( K1 + 2*K2 + 2*K3 + K4 );
        x += h;
        printf("%6.2f    ", x);
        PrintVectorD( Y );
        printf("\n");
    }
}
2
1 / 1 / 0
Регистрация: 28.01.2017
Сообщений: 32
28.01.2017, 14:39  [ТС] 7
Извиняюсь, можно по подробнее про VectorD и VectorFn?

Добавлено через 1 час 43 минуты
Измучился уже пробовать, ничего не получается

Добавлено через 5 минут
Сделал вот так:
C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
double **RungeKutt(int n,double x0,double y0,double h,double (*f1)(double x,double y),double (*f2)(double x,double y))
{
  register int i;
  double **res,s;
  double k11,k21,k31,k41;
  double k12,k22,k32,k42;
  res=new double* [n+1];
  for(i=0;i<n+1;i++)
  res[i]=new double[2];
  res[0][0]=x0;
  res[0][1]=y0;
  s=x0;
  for(i=1;i<=n;i++)
  {
    k11=h*f1(s,res[i-1][0]);
    k21=h*f1(s+h/2,res[i-1][0]+k11/2);
    k31=h*f1(s+h/2,res[i-1][0]+k21/2);
    k41=h*f1(s+h,res[i-1][0]+k31);
 
    k12=h*f2(s,res[i-1][1]);
    k22=h*f2(s+h/2,res[i-1][1]+k12/2);
    k32=h*f2(s+h/2,res[i-1][1]+k22/2);
    k42=h*f2(s+h,res[i-1][1]+k32);
 
    res[i][0]=res[i-1][0]+(k11+2*k21+2*k31+k41)/6; // точки по Х
    res[i][1]=res[i-1][1]+(k12+2*k22+2*k32+k42)/6; // точки по У
    s+=h;
  }
  return res;
}
считает всё, получаются точки, но не те, ибо вбиваю систему, решение которой знаю и график знаю, по полученным точкам, получается совершенно не тот график
0
1269 / 1026 / 470
Регистрация: 25.12.2016
Сообщений: 3,333
28.01.2017, 14:57 8
Цитата Сообщение от andrew95qq Посмотреть сообщение
Сделал вот так:
На первый взгляд ошибок не видно.
Цитата Сообщение от andrew95qq Посмотреть сообщение
считает всё, получаются точки, но не те, ибо вбиваю систему, решение которой знаю и график знаю, по полученным точкам, получается совершенно не тот график
Возможно, дело не в программе, а в известном решении. Может оно от другой системы или для других начальных условий. Я бы проверил работу программы с помощью какого-нибудь хорошо известного простого уравнения, решение которого точно известно. Например, уравнения свободных колебаний.
1
1 / 1 / 0
Регистрация: 28.01.2017
Сообщений: 32
28.01.2017, 15:16  [ТС] 9
Беру просто ДУ:
Название: 1.PNG
Просмотров: 21

Размер: 1.1 Кб
Его решение:
Название: 2.PNG
Просмотров: 21

Размер: 1.1 Кб
График у:
Добавить функционал в готовую программу решения систем дифферинциальных уравнений

Билиберда которая получается у меня:
Добавить функционал в готовую программу решения систем дифферинциальных уравнений
0
1 / 1 / 0
Регистрация: 28.01.2017
Сообщений: 32
28.01.2017, 15:19  [ТС] 10
Вот еще график частного решения, тоже по У:
Добавить функционал в готовую программу решения систем дифферинциальных уравнений

Но это тоже не то.
0
1269 / 1026 / 470
Регистрация: 25.12.2016
Сообщений: 3,333
28.01.2017, 15:28 11
А программа какой график строит? x(t), y(t) или y(x)?
0
1 / 1 / 0
Регистрация: 28.01.2017
Сообщений: 32
28.01.2017, 15:35  [ТС] 12
по полученным точкам в таблице
0
1269 / 1026 / 470
Регистрация: 25.12.2016
Сообщений: 3,333
28.01.2017, 15:39 13
Цитата Сообщение от andrew95qq Посмотреть сообщение
по полученным точкам в таблице
То есть y(x). А надо строить два графика: x(t) и y(t), где t = 0, h, 2h, 3h, ...
0
1 / 1 / 0
Регистрация: 28.01.2017
Сообщений: 32
28.01.2017, 16:01  [ТС] 14
пф, а где я возьму эти функции? суть решения ведь в том, что мы получаем точки и по ним строим график

Добавлено через 6 минут
В 1м самом посте, где я описывал решение просто ду, там я тоже строю просто по точкам, которые получаю из таблицы и всё верно получается
0
1269 / 1026 / 470
Регистрация: 25.12.2016
Сообщений: 3,333
28.01.2017, 16:11 15
Цитата Сообщение от andrew95qq Посмотреть сообщение
В 1м самом посте, где я описывал решение просто ду, там я тоже строю просто по точкам, которые получаю из таблицы и всё верно получается
Вот в той таблице первый столбец - это и есть значение времени с заданным шагом, а второй столбец - найденное решение в данных точках.

Иными словами, задача сводится к нахождению значений функции y = y(t) в точках t0, t0+h, t0+2h, ...
В общем случае y является вектором-строкой, то есть для каждого момента времени мы вычисляем значения нескольких функций. И для каждой найденной функции будет свой график yi = yi(t).
0
1 / 1 / 0
Регистрация: 28.01.2017
Сообщений: 32
28.01.2017, 16:23  [ТС] 16
C++
1
2
3
4
5
6
 for(i=0;i<=n;i++)
  {
    Series2->AddXY(t,Yres[i][0],"",Series2->SeriesColor); // Yres[i][0] тут наши решения для х
   Series3->AddXY(t,Yres[i][1],"",Series3->SeriesColor); // Yres[i][1] тут для у
    t+=h;
  }
Так? Фигня какая то получается.
Добавить функционал в готовую программу решения систем дифферинциальных уравнений
0
1 / 1 / 0
Регистрация: 28.01.2017
Сообщений: 32
28.01.2017, 16:50  [ТС] 17
Да и опять же, это получается, тоже самое, что я бы решал по отдельности два этих ДУ, у меня просто 2 массива с числами, а мне то нужно решение для системы
0
1269 / 1026 / 470
Регистрация: 25.12.2016
Сообщений: 3,333
28.01.2017, 17:06 18
Лучший ответ Сообщение было отмечено andrew95qq как решение

Решение

C++
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
#include <iostream>
#include <iomanip>
using namespace std;
 
typedef double (*TFunc) (double t, double x, double y);
 
double** RungeKutt(int n, double t0, double x0, double y0, double h, TFunc f1, TFunc f2)
{
    register int i;
    double **res, s;
    double k11, k21, k31, k41;
    double k12, k22, k32, k42;
    res=new double* [n+1];
    for(i=0; i<n+1; i++)
        res[i]=new double[2];
    res[0][0] = x0;
    res[0][1] = y0;
    s = 0;
    for(i=1; i<=n; i++)
    {
        double x = res[i-1][0];
        double y = res[i-1][1];
 
        k11=h*f1(s, x, y);
        k12=h*f2(s, x, y);
 
        k21=h*f1(s+h/2, x+k11/2, y+k12/2);
        k22=h*f2(s+h/2, x+k11/2, y+k12/2);
 
        k31=h*f1(s+h/2, x+k21/2, y+k22/2);
        k32=h*f2(s+h/2, x+k21/2, y+k22/2);
 
        k41=h*f1(s+h,   x+k31,   y+k32);
        k42=h*f2(s+h,   x+k31,   y+k32);
 
        res[i][0] = x + (k11+2*k21+2*k31+k41)/6;
        res[i][1] = y + (k12+2*k22+2*k32+k42)/6;
 
        s+=h;
    }
 
    return res;
}
 
double f1(double t, double x, double y)
{
    return -2*x + 4*y;
}
 
double f2(double t, double x, double y)
{
    return -x + 3*y;
}
 
 
int main()
{
    int n = 60;
    double **res = RungeKutt(n, 0, 0, 1, 0.1, f1, f2);
 
    for (int i=0; i<=n; i++)
    {
        double t = i*0.1;
        cout << setw(3) << t << "  " << setw(12) << res[i][0] << "  " << res[i][1] << endl;
    }
 
}
Добавлено через 9 минут
А вот как должны выглядеть графики x(t) и y(t):
https://www.wolframalpha.com/i... rom+0+to+2
https://www.wolframalpha.com/i... rom+0+to+2
1
1 / 1 / 0
Регистрация: 28.01.2017
Сообщений: 32
28.01.2017, 17:19  [ТС] 19
Капеееец, почему так? вроде тоже самое все написано, а цифры получаются совсем другие!!!
Круто, спасибо!!
Добавить функционал в готовую программу решения систем дифферинциальных уравнений

Решения правильные, но я так понимаю, опять же, это мы решаем по отдельности функции нашей системы? а общее решение как?)))))
Просто вот это вот типо график общего решения:
Добавить функционал в готовую программу решения систем дифферинциальных уравнений



з.ы.
ну и там в вашей функции, я так понимаю, t0 - лишняя?

зз.ы.
а нет, в вашем решении не тоже самое, увидел где, спасибо еще раз
0
1269 / 1026 / 470
Регистрация: 25.12.2016
Сообщений: 3,333
28.01.2017, 17:28 20
Цитата Сообщение от andrew95qq Посмотреть сообщение
вроде тоже самое все написано
Всё дело в порядке вычисления коэффициентов k1, k2, k3, k4.

Цитата Сообщение от andrew95qq Посмотреть сообщение
t0 - лишняя?
Нет, не лишняя. t0 - это начальный момент времени, x0 - это значение x при t=0, аналогично y0.
В данном случае функции f1 и f2 не зависят от t (в программе s), поэтому его в принципе можно не передавать и не вычислять. Но в общем случае такая зависимость может иметь место.

Добавлено через 3 минуты
Замечу, что я добавил ещё один параметр функциям f1 и f2. Если его убрать (он реально не используется), то вызовы функций должны выглядеть например так:
C++
1
k32=h*f2(x+k21/2, y+k22/2);
Добавлено через 2 минуты
Цитата Сообщение от andrew95qq Посмотреть сообщение
Решения правильные, но я так понимаю, опять же, это мы решаем по отдельности функции нашей системы? а общее решение как?
Не понял оба вопроса. Что значит по отдельности? При чём тут общее решение?
1
IT_Exp
Эксперт
87844 / 49110 / 22898
Регистрация: 17.06.2006
Сообщений: 92,604
28.01.2017, 17:28
Помогаю со студенческими работами здесь

Исправить код метод Ньютона для решения систем нелинейных уравнений под нужное условие
Данный код для решения системы ax+tg(xy)=0; (y^2-b^2)+lnx=0 Перепишите его,пожалуйста для...

Подскажите пожалуйста как добавить конструктор копирования в готовую программу
Вот код программы: # include &lt;iostream&gt; # include &lt;conio.h&gt; # include &lt;string&gt; # include...

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

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


Искать еще темы с ответами

Или воспользуйтесь поиском по форуму:
20
Ответ Создать тему
Опции темы

КиберФорум - форум программистов, компьютерный форум, программирование
Powered by vBulletin
Copyright ©2000 - 2022, CyberForum.ru