Модераторы: Daevaorn
  

Поиск:

Ответ в темуСоздание новой темы Создание опроса
> Алгоритм Рунге-Кутта 
:(
    Опции темы
larev01
Дата 25.1.2011, 14:44 (ссылка) | (нет голосов) Загрузка ... Загрузка ... Быстрая цитата Цитата


Новичок



Профиль
Группа: Участник
Сообщений: 2
Регистрация: 11.7.2010

Репутация: нет
Всего: нет



Добрый день. Столкнулся с проблемой. Необходимо решить уравнение методом Рунге-Кутта четвертого порядка с точностью 0.0001 (для достижения точности использую метод двойного пересчета). Написал программу, но интервалы получаются какими-то подозрительно маленькими. Помогите пожалуйста найти ошибку (если она есть).
Уравнение:
user posted image

Исходник:
Код

#include <iostream>
#include <cmath>
#include <conio.h>
 
inline double f(double x, double y)
{
        return ((6-y*y)*std::cos(x)+2*y);
}
 
void runge_kutta(double h, double* x, double* y, int n)
{
        double k1, k2, k3, k4;
 
        for(int i=1; i<=n; i++)
                x[i] = x[i-1]+h;
 
        for(int i=0; i<n; i++) {
                k1 = f(x[i], y[i]);
                k2 = f(x[i]+h/2, y[i]+(h/2)*k1);
                k3 = f(x[i]+h/2, y[i]+(h/2)*k2);
                k4 = f(x[i]+h, y[i]+h*k3);
                y[i+1] = y[i]+(h/6)*(k1+2*k2+2*k3+k4);
        }
}
 
int main()
{
        double h = 0.5;
        const double precision = 0.0001;
 
        int n;
        double* x;
        double* y;
        while(true) {
                double* x_temp;
                double* y_temp;
                n = (static_cast<int>(1/h)+1)*2;
                x_temp = new double[n+1];
                y_temp = new double[n+1];
                x = new double[n+1];
                y = new double[n+1];
                x_temp[0] = 0;
                y_temp[0] = 0.3;
                x[0] = 0;
                y[0] = 0.3;
                runge_kutta(h, x_temp, y_temp, n);
                runge_kutta(h/2, x, y, n);
 
                bool need_break = false;
                for(int i=1; i<=n; i++)
                        if(std::fabs(y_temp[i]-y[i])<precision*15)
                                need_break = true;
                if(need_break) {
                        delete[] x_temp;
                        delete[] y_temp;
                        break;
                }
 
                delete[] x_temp;
                delete[] y_temp;
                delete[] x;
                delete[] y;
                h/=2;
        }
 
        std::cout.precision(4);
        for(int i=0; i<=n; i++)
                std::cout << x[i] << '\t' << y[i] << '\n';
 
        delete[] x;
        delete[] y;
 
        _getch();
        return 0;
}


PM MAIL   Вверх
Peter
Дата 7.3.2011, 19:29 (ссылка) | (нет голосов) Загрузка ... Загрузка ... Быстрая цитата Цитата


Опытный
**


Профиль
Группа: Участник
Сообщений: 771
Регистрация: 28.7.2003
Где: Ставрополь

Репутация: -1
Всего: 1



Мне кажется, ошибка здесь:
Код
                runge_kutta(h, x_temp, y_temp, n);
                runge_kutta(h/2, x, y, n);

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


--------------------
всё, что делаете, делайте от души, как для Господа (Послание апостола Павла колоссянам, 3:23).
PM MAIL WWW   Вверх
  
Ответ в темуСоздание новой темы Создание опроса
Правила форума "С++:Общие вопросы"
Earnest Daevaorn

Добро пожаловать!

  • Черновик стандарта C++ (за октябрь 2005) можно скачать с этого сайта. Прямая ссылка на файл черновика(4.4мб).
  • Черновик стандарта C (за сентябрь 2005) можно скачать с этого сайта. Прямая ссылка на файл черновика (3.4мб).
  • Прежде чем задать вопрос, прочтите это и/или это!
  • Здесь хранится весь мировой запас ссылок на документы, связанные с C++ :)
  • Не брезгуйте пользоваться тегами [code=cpp][/code].
  • Пожалуйста, не просите написать за вас программы в этом разделе - для этого существует "Центр Помощи".
  • C++ FAQ

Если Вам понравилась атмосфера форума, заходите к нам чаще! С уважением, Earnest Daevaorn

 
0 Пользователей читают эту тему (0 Гостей и 0 Скрытых Пользователей)
0 Пользователей:
« Предыдущая тема | C/C++: Общие вопросы | Следующая тема »


 




[ Время генерации скрипта: 0.2048 ]   [ Использовано запросов: 22 ]   [ GZIP включён ]


Реклама на сайте     Информационное спонсорство

 
По вопросам размещения рекламы пишите на vladimir(sobaka)vingrad.ru
Отказ от ответственности     Powered by Invision Power Board(R) 1.3 © 2003  IPS, Inc.