Форум программистов, компьютерный форум, киберфорум
Наши страницы
Численные методы
Войти
Регистрация
Восстановить пароль
 
 
Рейтинг 5.00/7: Рейтинг темы: голосов - 7, средняя оценка - 5.00
супер тупой
6 / 6 / 3
Регистрация: 29.08.2014
Сообщений: 89
Завершенные тесты: 1
#1

Численное решение задачи Дирихле для уравнения Лапласа на полусфере

26.03.2016, 01:14. Просмотров 1221. Ответов 23
Метки нет (Все метки)

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

Постановка задачи:
http://www.cyberforum.ru/cgi-bin/latex.cgi?\left\{\begin{matrix}& \\ \Delta u(r, \varphi, \psi ) = 0& \\ {u\mid}_{\psi  = 0} = {u}_{0}(\varphi)& \\ {u\mid}_{\psi  = \frac{\pi}{4}} = {u}_{1}(\varphi)& \\ {u\mid}_{\psi  = \frac{\pi}{2}} = {u}_{2}(\varphi)\end{matrix}\right.
Где http://www.cyberforum.ru/cgi-bin/latex.cgi?u = u(r, \varphi, \psi ) - искомая функция, http://www.cyberforum.ru/cgi-bin/latex.cgi?{u}_{0} \,,\, {u}_{1} \,,\, {u}_{2} - краевые условия на сфере, http://www.cyberforum.ru/cgi-bin/latex.cgi?\mathbf{R} - радиус сферы.
Т.е. заданы значения на поверхности сферы при http://www.cyberforum.ru/cgi-bin/latex.cgi?\psi = 0\,,\,\psi = \frac{\pi }{4}\,,\,\psi = \frac{\pi }{2}.

Нужно найти значения функции http://www.cyberforum.ru/cgi-bin/latex.cgi?\mathbf{u} на верхней полусфере. Если http://www.cyberforum.ru/cgi-bin/latex.cgi?u = u(x, y, z), то пользуемся такой заменой:
http://www.cyberforum.ru/cgi-bin/latex.cgi?\left\{\begin{matrix}& \\ x = Rcos(\varphi )cos(\psi )& \\ y = Rsin(\varphi )cos(\psi )& \\ z = Rsin(\psi )\end{matrix}\right.,0 \leq \varphi \leq 2\pi  ,  0 \leq \psi \leq \frac{\pi }{2}

Решение задачи:
Кликните здесь для просмотра всего текста

Пусть
http://www.cyberforum.ru/cgi-bin/latex.cgi?\begin{matrix}\\ \varphi = i{h}_{\varphi }\\ \psi = j{h}_{\psi }\end{matrix}

Распишем оператор Лапласа в наших координатах:
http://www.cyberforum.ru/cgi-bin/latex.cgi?\Delta u = \frac{1}{{R}^{2}{\cos(\psi)}^{2}}{u}_{\varphi \varphi} + \frac{1}{{R}^{2}}{u}_{\psi \psi} - \frac{tg(\psi)}{{R}^{2}}{u}_{\psi } = 0

Умножим на http://www.cyberforum.ru/cgi-bin/latex.cgi?\mathbf{R} и положим http://www.cyberforum.ru/cgi-bin/latex.cgi?tg(\psi ) = v(\psi ), тогда
http://www.cyberforum.ru/cgi-bin/latex.cgi?v'(\psi ){u}_{\varphi \varphi} + {u}_{\psi \psi} - v(\psi){u}_{\psi } = 0

Пусть
http://www.cyberforum.ru/cgi-bin/latex.cgi?\begin{matrix}\\ {u}_{\varphi} = \frac{{u}_{i+1;j} - {u}_{i-1;j}}{2{h}_{\varphi}},{u}_{\psi} = \frac{{u}_{i;j+1} - {u}_{i;j-1}}{2{h}_{\psi}}\\ {u}_{\varphi \varphi} = \frac{{u}_{i+1;j} - 2{u}_{i;j} + {u}_{i-1;j}}{{{h}^{2}}_{\varphi}},{u}_{\psi \psi} = \frac{{u}_{i;j+1} - 2{u}_{i;j} + {u}_{i;j-1}}{{{h}^{2}}_{\psi}}\end{matrix}

После подстановки в уравнение Лапласа получаем:
http://www.cyberforum.ru/cgi-bin/latex.cgi?{u}_{i;j} = \alpha {u}_{i+1;j} + \alpha {u}_{i-1;j} + \beta {u}_{i;j+1} + \gamma {u}_{i;j-1}
Где
http://www.cyberforum.ru/cgi-bin/latex.cgi?\alpha = \frac{d}{a}\,,\, \beta = \frac{b}{a}\,,\, \gamma = \frac{c}{a}
И
http://www.cyberforum.ru/cgi-bin/latex.cgi?a = 2\frac{v'(\psi){{h}^{2}}_{\psi} + {{h}^{2}}_{\varphi}}{{{h}^{2}}_{\varphi}{{h}^{2}}_{\psi}}\, ,\, b = \frac{2 - v(\psi){h}_{\psi}}{2{{h}^{2}}_{\psi}}\, ,\, c = \frac{2 + v(\psi){h}_{\psi}}{2{{h}^{2}}_{\psi}}\, , \, d = \frac{v'(\psi)}{{{h}^{2}}_{\varphi}}

Реализация на C++:
Кликните здесь для просмотра всего текста

main.cpp
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
#include "head.h"
 
int main() {
 
    Data data;//класс данных, описан в .h файле
 
    size_t iter     = 0;//номер итерации
    bool   go       = true;//индикатор продолжения счета
    double max_diff = .0;//максимальная разница на данной итерации
    double h_phi2   = data.h_phi * data.h_phi, h_psi2 = data.h_psi * data.h_psi;//квадраты некоторых величин
 
    double a, b, c, d, alpha, beta, gamma, phi, psi, v, dv, new_val;//нужные переменные
 
    do {
        iter++;//начинаем с первой итерации
        max_diff = .0;//каждый раз обнуляем разницу
 
        for (int i = 0; i < data.num_of_units - 1; i++) {//цикл по phi
            for (int j = 1; j < data.num_of_units - 1; j++) {//цикл по psi
                if (j == (data.num_of_units - 1) / 2) continue;//не пересчитываем значения для psi / 4 
 
                phi = i * data.h_phi;
                psi = j * data.h_psi;
 
                v  = tan(psi);
                dv = 1.0 / (cos(psi) * cos(psi));
                
                a = (2.0 * (dv * h_psi2 + h_phi2)) / (h_phi2 * h_psi2);
                b = (2.0 - v * data.h_psi) / (2.0 * h_psi2);
                c = (2.0 + v * data.h_psi) / (2.0 * h_psi2);
                d = dv / h_phi2;
 
                alpha = d / a;
                beta  = b / a;
                gamma = c / a;
 
                //считаем точку phi = 2pi за точку phi = 0
                new_val = alpha * (data.u[((i + 1) == (data.num_of_units - 1)) ? 0 : (i + 1)][j] + data.u[((i - 1) > data.num_of_units) ? (data.num_of_units - 2) : (i - 1)][j]) + beta * data.u[i][j + 1] + gamma * data.u[i][j - 1];
 
                if (fabs(new_val - data.u[i][j]) > max_diff && iter != 1)//если не первая итерация и данная разница больше, чем максимальная на предыдущем шаге
                    max_diff = fabs(new_val - data.u[i][j]);//заменяем ее
                data.u[i][j] = new_val;//записываем новое значение
            }
        }
 
        if (max_diff < data.EPS && iter != 1) go = false;//если максимальная разница стала меньше EPS, то завершаем счет
 
    } while(iter <= 150000 && go);//если итераций больше, чем заданное число, то завершаем счет
 
    data.write_debug();//записываем данные в файл
 
    system("pause");
    return 0;
}
head.h
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
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
#pragma once
#include <iostream>
#include <fstream>
 
inline double orig(double r, double phi, double psi) {//оригинальная функция
    return r * cos(phi) * cos(psi) + r * sin(phi) * cos(psi) + r * sin(psi);
}
 
class Data {
public:
    Data() {//конструктор, вызывающийся при создании объекта класса
        std::cout << "part = "; std::cin >> part;//спрашиваем величину разбиения
        std::cout << "R = ";    std::cin >> R;//и радиус сферы
 
        if (part % 2 != 0) part++;//если величина разбиения нечетна, то делаем ее четной
        //нужно для того, чтобы узел psi = pi / 4 был ровно в середине массива
 
        h_r   = R      / part;//считаем шаги
        h_phi = PI_m_2 / part;
        h_psi = PI_d_2 / part;
 
        num_of_units = part + 1;//количество узлов
        //нужно сохранять на 1 значение больше, т.к. не будет включаться граничная точка
        //например part = 2 и R = 1
        //получаем h_r = 0.5
        //должны быть узлы 0, 0.5, 1
        //но т.к. part = 2, то массив будет только из 2х элементов => увеличиваем на 1
 
        u = new double*[num_of_units - 1];//элементы для phi, т.к. 0 == 2PI, то делаем на один элемент меньше и считаем нулевой за 0 и 2PI
        for (size_t i = 0; i < num_of_units - 1; i++)
            u[i] = new double[num_of_units];//элементы для psi
 
        //j == 0                      <=> psi = 0
        //j == (num_of_units - 1) / 2 <=> psi = pi / 4
        //j == num_of_units - 1       <=> psi = pi / 2
        for (size_t i = 0; i < num_of_units - 1; i++) {
            for (size_t j = 0; j < num_of_units; j++) {
                if (j == 0 || j == (num_of_units - 1) / 2 || j == num_of_units - 1) u[i][j] = orig(R, i * h_phi, j * h_psi);//на границах заполняем массив известными данными
                else u[i][j] = .0;//внутренние точки заполняем нулями
                //else u[i][j] = orig(R, i * h_phi, j * h_psi);
            }
        }
    }
    ~Data() {//деструктор
        for (size_t i = 0; i < num_of_units - 1; i++)//удаляем массив при удалении объекта
            delete[] u[i];
        delete[] u;
    }
 
    void write_debug() {//записываем данные в файл
        std::ofstream dfout;
        dfout.open("data.txt");
 
        size_t num = 0;//количество погрешностей(для средней)
        double mid_diff = 0, max_diff = 0;//средняя погрешность и максимальная погрешность
 
        dfout << "part = "        <<      part     << " <- величина разбиения\n";
        dfout << "num_of_units = " << num_of_units << " <- количество узлов\n";
        dfout << "R = "            <<       R      << " <- радиус\n\n";
 
        dfout << "h_r   = " << h_r   << " <- мелкость разбиения по переменной r\n";
        dfout << "h_phi = " << h_phi << " <- мелкость разбиения по переменной phi\n";
        dfout << "h_psi = " << h_psi << " <- мелкость разбиения по переменной psi\n\n\n";
 
        for (size_t i = 0; i < num_of_units - 1; i++) {
            dfout << "\n";
            for (size_t j = 0; j < num_of_units; j++) {
                double diff = fabs(u[i][j] - orig(R, i * h_phi, j * h_psi));//считаем погрешность в точке
                if (diff > max_diff) max_diff = diff;//находим максимальную
                mid_diff += diff;//складываем все погрешности
                num++;//увеличиваем количество погрешностей
                dfout << "u[" << i << ", " << j << "] = " << u[i][j] << " || u[" << i * h_phi << ", " << j * h_psi << "] = " << u[i][j] << " || orig(this) = " << orig(R, i * h_phi, j * h_psi) << " || diff = " << diff << "\n";
            }
        }
 
        mid_diff /= num;//вычисляем среднюю погрешность
        dfout << "\n\nmax_diff = " << max_diff << "\nmid_diff = " << mid_diff;
 
        dfout.close();
    }
 
    const double PI_d_2 = 1.57079632679,// pi / 2
                 PI     = 3.14159265359,// pi
                 PI_m_2 = 6.28318530718,// pi * 2
                 EPS    = 0.00000000001;
 
    size_t  part;//величина разбиения
    double  R;//радиус сферы
 
    double  h_r, h_phi, h_psi;//шаг
    size_t  num_of_units;//количество узлов
 
    double** u;//массив данных
};
0
Надоела реклама? Зарегистрируйтесь и она исчезнет полностью.
Similar
Эксперт
41792 / 34177 / 6122
Регистрация: 12.04.2006
Сообщений: 57,940
26.03.2016, 01:14
Ответы с готовыми решениями:

Решение задачи Дирихле для уравнения Пуассона
Если есть у кого решение данной задачи, выложите срочно!! не получается в своей...

Численное решение уравнения Лапласа
Очень нужна ваша помощь. Никак не могу разобраться с решением уравнения...

Численное решение краевой задачи для ОДУ
доброго времени суток есть задача по теме &quot;Численное решение краевой задачи...

Найти численное решение задачи Коши для системы дифференциальных уравнений
Добрые день.Дана вот такая задача: Электронная схема во временном интервале...

Найти численное решение задачи Коши для ДУ-1 (по методу Эйлера). Построить график интегральной кривой
Найти численное решение задачи Коши для ДУ-1 (по методу Эйлера) на отрезке с...

23
bobah16
370 / 340 / 42
Регистрация: 14.07.2015
Сообщений: 2,882
30.03.2016, 22:15 #2
Цитата Сообщение от супер тупой Посмотреть сообщение
Буду рад советам, хотя бы указывающим, куда капать.
Причина скорее всего проста - мелкая ошибка в коде. Ищите.
0
супер тупой
6 / 6 / 3
Регистрация: 29.08.2014
Сообщений: 89
Завершенные тесты: 1
30.03.2016, 23:11  [ТС] #3
bobah16, код переписывался, проверялись отдельные части, если бы все было так просто
Конечно не отрицаю такой возможности, но не понимаю тогда где собака зарыта.
0
cmath
Модератор
2490 / 1714 / 145
Регистрация: 11.08.2012
Сообщений: 3,293
Завершенные тесты: 6
31.03.2016, 07:24 #4
Для теста рекомендую принять R=1 и считать, что функция не зависит от r, i.e. расстояния до центра координат. Это раз.
Два - не вижу, где применяются граничные условия.
Три - у вас схема-то неявная - тем не менее, вы все равно используете формулу с альфами-бетами непосредственно в вычислениях.
***
main.cpp и head.h можно отправлять на помойку, в них ничего ценного нет, алгоритм некорректен. Слово "мелкость" рекомендую заменить словом "шаг", ибо общеупотребительная терминология.
***
Алгоритм решения состоит в применении связки уравнение + граничные условия непосредственно для построения СЛАУ размерности n*m, где n - кол-во узлов по "фи", m - по "пси". Учитывая, что http://www.cyberforum.ru/cgi-bin/latex.cgi?\psi = \frac{\pi}{4} \in [0, \frac{\pi}{2}], то шаг по этой переменной нужно задавать, отталкиваясь от величины http://www.cyberforum.ru/cgi-bin/latex.cgi?\frac{\pi}{4} а не http://www.cyberforum.ru/cgi-bin/latex.cgi?\frac{\pi}{2}, чтобы часть наших узлов обязательно имела координату http://www.cyberforum.ru/cgi-bin/latex.cgi?\psi = \frac{\pi}{4}. Нужно также задать соответствие для индексов [J] вектора-решения СЛАУ и индексами [i, j] матрицы, в которой будем хранить значения функции u.
На вашем месте я бы также вынес в отдельные функции вычисление коэффициентов a,b,c,d.
Заполнять матрицу коэффициентов СЛАУ нужно так (считаем что все её элементы при инициализации = 0): пробегаем по узлам [i,j], используем их координаты для вычисления коэффициента при u[i, j] в уравнении, прибавляем полученное значение к элементу матрицы СЛАУ с индексом [J, J]. Если координата "пси" узла соответствует граничному условию - используем граничное условие вместо уравнения, не забывая записывать в элемент вектора свободных членов с индексом [J] значение граничного условия в узле [i, j] (если используем уравнение, то элемент вектора равен нулю).
По окончанию процедуры полученную СЛАУ нужно решить, например, методом Гаусса (можно попробовать прогонку, если звезды сойдутся :-) ). Затем вектор решения записываете в матрицу u, используя связь между индексом J и индексами i,j.
***
Профит

P.S. Почитайте еще про god object, и почему его следует избегать. И еще - в конструкторе класса не стоит писать код для ввода значений, это нужно делать отдельно в какой-нибудь процедуре.
0
супер тупой
6 / 6 / 3
Регистрация: 29.08.2014
Сообщений: 89
Завершенные тесты: 1
01.04.2016, 01:09  [ТС] #5
cmath, спасибо большое за ответ!

Цитата Сообщение от cmath Посмотреть сообщение
Для теста рекомендую принять R=1 и считать, что функция не зависит от r, i.e. расстояния до центра координат. Это раз.
Так я же и решаю задачу для фиксированного R или я вас не понял?

Цитата Сообщение от cmath Посмотреть сообщение
Два - не вижу, где применяются граничные условия.
Граничные условия я применяю для заполнения значений точек на границе области(38я строка в .h).

Цитата Сообщение от cmath Посмотреть сообщение
Три - у вас схема-то неявная - тем не менее, вы все равно используете формулу с альфами-бетами непосредственно в вычислениях.
Не понял о чем речь(видимо не хватает каких-то знаний, буду признателен, если подскажите, где почитать).

Цитата Сообщение от cmath Посмотреть сообщение
Учитывая, что http://www.cyberforum.ru/cgi-bin/latex.cgi?\psi =\frac{\pi }{4}\in \left[0,\frac{\pi }{2} \right], то шаг по этой переменной нужно задавать, отталкиваясь от величины http://www.cyberforum.ru/cgi-bin/latex.cgi?\frac{\pi }{4} а не http://www.cyberforum.ru/cgi-bin/latex.cgi?\frac{\pi }{2}, чтобы часть наших узлов обязательно имела координату http://www.cyberforum.ru/cgi-bin/latex.cgi?\psi = \frac{\pi }{4}.
Я просто требую четное количество узлов по http://www.cyberforum.ru/cgi-bin/latex.cgi?\psi, этого будет достаточно?

Цитата Сообщение от cmath Посмотреть сообщение
алгоритм некорректен
А что именно некорректно, как я вычисляю значения u[i;j]?

И последнее, не совсем понял как задавать СЛАУ и про какие индексы [J] вы говорите.

P.S.
Цитата Сообщение от cmath Посмотреть сообщение
Почитайте еще про god object, и почему его следует избегать. И еще - в конструкторе класса не стоит писать код для ввода значений, это нужно делать отдельно в какой-нибудь процедуре.
Этот код был написан чисто для того, чтобы выложить его на форум, старался сделать максимально просто и кратко
Но про god object почитал, спасибо!
0
bobah16
370 / 340 / 42
Регистрация: 14.07.2015
Сообщений: 2,882
01.04.2016, 01:51 #6
Цитата Сообщение от супер тупой Посмотреть сообщение
Так я же и решаю задачу для фиксированного R или я вас не понял?
Короче, вы составили разностную схему, написали код, теперь нужно проверить правильность написания схемы и кода, для этого подбираете решение, которое будет вычислено с нулевой погрешностью. Для этого свое решение подставляете в выражение для погрешности аппроксимации и смотрите на погрешность...делаете выводы. В общем метод пробных функций надо бы вам освоить.
Ладно, еще короче, у меня есть код для решения уравнения Пуассона на фортране, могу скинуть, посмотрите че к чему, если захотите.
Ну а так, пользуйтесь общеизвестной литературой по численным методам, там все есть. Кое-где даже код для вас вроде написан...правда не на срр.
0
cmath
Модератор
2490 / 1714 / 145
Регистрация: 11.08.2012
Сообщений: 3,293
Завершенные тесты: 6
01.04.2016, 02:54 #7
Цитата Сообщение от супер тупой Посмотреть сообщение
Так я же и решаю задачу для фиксированного R или я вас не понял?
У вас при решении R ушло в знаменатель (вы еще домножили уравнение на R^2, чтобы избавиться от знаменателя). Если R=1, то умножать на R^2 не нужно и в последствии его можно не учитывать при окончательной сборке решения и проверке погрешности. Это для теста, чтобы малой кровью убедиться, что все работает. Затем можно сделать R произвольным.
Цитата Сообщение от супер тупой Посмотреть сообщение
Я просто требую четное количество узлов
Нужно наоборот, нечетное.
Цитата Сообщение от супер тупой Посмотреть сообщение
Не понял о чем речь
У вас значение функции в текущем узле зависит не только от значений в предыдущих узлах, но и от последующих тоже, и ладно бы, если бы таких узлов было 1 шт, так их 2 шт, т.е. вычислить значение в [i, j] с помощью этой формулы непосредственно нельзя.
Цитата Сообщение от супер тупой Посмотреть сообщение
Граничные условия я применяю для заполнения значений точек на границе области
Вижу. Это бы прокатило бы при явной схеме. Но увы. Да и одних граничных условий при построении решения не достаточно. Нужно в добавок к этому использовать условие непрерывности. Вы чертеж строили?
Цитата Сообщение от супер тупой Посмотреть сообщение
И последнее, не совсем понял как задавать СЛАУ и про какие индексы [J] вы говорите.
Когда составляете СЛАУ, там нужно, как бы это попроще - нужно уметь переписать значения матрицы в вектор и обратно без разрушения структуры, по этому нужно задать соответствие между элементами матрицы и элементами вектора.
Например:
http://www.cyberforum.ru/cgi-bin/latex.cgi?<br />
\begin{pmatrix}2 & 3\\ -3 & 10\end{pmatrix} \Rightarrow \{2,3,-3,10\}\\<br />
[0,0]\rightarrow 0, [0,1]\rightarrow 1, [1,0]\rightarrow 2, [1,1]\rightarrow 3
Т.е. J=i*m+j где m - кол-во столбцов матрицы

Добавлено через 7 минут
Цитата Сообщение от супер тупой Посмотреть сообщение
видимо не хватает каких-то знаний, буду признателен, если подскажите, где что почитать
Самарского. У него пара толстеньких книжек есть, посвященных именно численному решению ОДУ и УЧП

Добавлено через 30 минут

Не по теме:

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

0
супер тупой
6 / 6 / 3
Регистрация: 29.08.2014
Сообщений: 89
Завершенные тесты: 1
02.04.2016, 01:26  [ТС] #8
bobah16, так я беру известное u, по нему задаю краевые условия и решаю задачу. А потом сравниваю мое решение с точными значениями функции. Вы об этом?
Цитата Сообщение от bobah16 Посмотреть сообщение
у меня есть код для решения уравнения Пуассона на фортране, могу скинуть
Хуже от этого явно не будет, если вам не сложно, буду благодарен!

cmath,
Цитата Сообщение от cmath Посмотреть сообщение
Если R=1, то умножать на R^2 не нужно и в последствии его можно не учитывать при окончательной сборке решения и проверке погрешности.
Я может не правильно делал, но у меня уравнение Лапласа и если http://www.cyberforum.ru/cgi-bin/latex.cgi?r = const \neq 0, то я от него могу спокойно избавиться.

Цитата Сообщение от cmath Посмотреть сообщение
Нужно наоборот, нечетное.
Опять не понял, у меня отрезок http://www.cyberforum.ru/cgi-bin/latex.cgi?\left[0;\frac{\pi }{2} \right]. Беру четное количество узлов http://www.cyberforum.ru/cgi-bin/latex.cgi?n = 2k,\, k = 1, 2, ..., тогда шаг будет равен http://www.cyberforum.ru/cgi-bin/latex.cgi?\frac{\left( \frac{\pi }{2} - 0\right)}{n} = \frac{\pi }{2n} = \frac{\pi }{4k} = \left( \frac{\pi }{4}\right)\left( \frac{1}{k}\right), т.е. у меня всегда будет узел при http://www.cyberforum.ru/cgi-bin/latex.cgi?\psi = \frac{\pi }{4}.

Цитата Сообщение от cmath Посмотреть сообщение
Нужно в добавок к этому использовать условие непрерывности. Вы чертеж строили?
Да, строил. Вы говорите про равенство значений функции в точках http://www.cyberforum.ru/cgi-bin/latex.cgi?\varphi = 0 и http://www.cyberforum.ru/cgi-bin/latex.cgi?\varphi = 2\pi? Т.е. http://www.cyberforum.ru/cgi-bin/latex.cgi?u\left(r,\varphi ,\psi  \right) = u\left(r,\varphi + 2\pi  ,\psi  \right). Если да, то я считаю точки http://www.cyberforum.ru/cgi-bin/latex.cgi?\varphi = 0 = 2\pi внутренними и, грубо говоря, считаю точку 0 за точку http://www.cyberforum.ru/cgi-bin/latex.cgi?2\pi.

Про индексы [J] понял, а вот как заполнить элемент [i;j] не понял. Но есть подозрение, что нужно в решение мне использовать другие разностные схемы, так? Просто не вижу смысл заполнения [i;j] элемента СЛАУ значением коэф. при u[i;j]...
0
cmath
Модератор
2490 / 1714 / 145
Регистрация: 11.08.2012
Сообщений: 3,293
Завершенные тесты: 6
02.04.2016, 03:08 #9
Цитата Сообщение от супер тупой Посмотреть сообщение
тогда шаг будет равен
Нет, не будет. Вы не правильно считаете. Конфигурация такая:

(x[0]=0)----------(x[1]=pi/4)----------(x[2]=pi/2), шаг = pi/4, 3 узла.

(x[0]=0)----------(x[1]=pi/8)----------(x[2]=pi/4)----------(x[3]=3pi/8)----------(x[4]=pi/2), шаг=pi/8 5 узлов

и т.д. Координата узла x[i] = i * h, счет идет с нуля, при этом h=pi/2(n-1), где n-количество узлов
Цитата Сообщение от супер тупой Посмотреть сообщение
т.е. у меня всегда будет узел при
Как видите, не будет.
Цитата Сообщение от супер тупой Посмотреть сообщение
то я от него могу спокойно избавиться
Нет, не можете. Если у вас радиус фигурирует в функции каким-либо боком, то его нужно будет учитывать при сборке решения после вычислений.
Цитата Сообщение от супер тупой Посмотреть сообщение
грубо говоря
Только грубо. Это понадобиться при проходе по http://www.cyberforum.ru/cgi-bin/latex.cgi?\varphi=2\pi и http://www.cyberforum.ru/cgi-bin/latex.cgi?\varphi=0
Цитата Сообщение от супер тупой Посмотреть сообщение
как заполнить элемент [i;j]
i = J mod m, j = J - i*m
Цитата Сообщение от супер тупой Посмотреть сообщение
Просто не вижу смысл заполнения [i;j] элемента СЛАУ значением коэф. при u[i;j]...
Вы краевые задачи для обычных ОДУ порядка выше первого хоть раз решали? Далеко не все краевые задачи для ОДУ, а уж тем более для УЧП решаются численно с помощью явных вычислений, когда последующие значения вычисляются через предыдущие (это вообще говоря, возможно только для задач Коши или метода пристрелки, когда краевая задача сводится к задаче Коши, но это не ваш случай) значения. Простенький пример:
http://www.cyberforum.ru/cgi-bin/latex.cgi?y''=-\sin x, y(0)=y(\pi)=0 решение этой задачи известно. y=sin x. А вы постройте-ка разностную схему и решите численно. Посмотрим, что у вас получится. Задача простецкая и её можно решить с неплохой точностью даже на листе бумаги.
0
супер тупой
6 / 6 / 3
Регистрация: 29.08.2014
Сообщений: 89
Завершенные тесты: 1
03.04.2016, 02:17  [ТС] #10
cmath,
Цитата Сообщение от cmath Посмотреть сообщение
Нет, не будет. Вы не правильно считаете.
Считаю правильно, а вот выражаюсь не правильно
Я требую, чтобы четным числом было не количество узлов, а количество отрезков. Узлов соответственно будет нечетное число, конечно.

Цитата Сообщение от cmath Посмотреть сообщение
Нет, не можете. Если у вас радиус фигурирует в функции каким-либо боком, то его нужно будет учитывать при сборке решения после вычислений.
В какой именно функции? И про какую сборку вы говорите?

Я сейчас полистал несколько книжек и во всех именно такой способ, которым я решал задачу Дирихле в прямоугольнике. По аналогии решал в круге. Теперь по аналогии хочу решить на сфере, а вы мне говорите, что это совершенно не верно
Где я повернул не туда?
И я кажется начинаю понимать, какой способ решения вы пытаетесь до меня донести. Такое я тоже видел в одной книге, но после этого метода была глава с итерационным методом(который делает по сути тоже самое), которым я собственно и пользовался.

Цитата Сообщение от cmath Посмотреть сообщение
Вы краевые задачи для обычных ОДУ порядка выше первого хоть раз решали?
Решал, второго например, разбивая его на систему из двух ОДУ первого порядка

Цитата Сообщение от cmath Посмотреть сообщение
А вы постройте-ка разностную схему и решите численно.
Построил и решил.
0
cmath
Модератор
2490 / 1714 / 145
Регистрация: 11.08.2012
Сообщений: 3,293
Завершенные тесты: 6
03.04.2016, 04:21 #11
Цитата Сообщение от супер тупой Посмотреть сообщение
В какой именно функции? И про какую сборку вы говорите?
Я сейчас полистал несколько книжек и во всех именно такой способ, которым я решал задачу Дирихле в прямоугольнике. По аналогии решал в круге. Теперь по аналогии хочу решить на сфере, а вы мне говорите, что это совершенно не верно
У вас функция u, в общем случае, зависит от ТРЕХ переменных, среди которых есть r - расстояние до центра координат.
Цитата Сообщение от супер тупой Посмотреть сообщение
Решал, второго например, разбивая его на систему из двух ОДУ первого порядка
Что решали? Краевую задачу? В которой нет начальных условий для производных, как в той, что я вам дал? Вы НЕ можете её таким способом решить, разве что методом стрельбы, но он не идеален и годен только для линейных систем и линейных условий, если закроется хоть какая-либо нелинейщина - все, пиши пропало. Да и для УЧП метод тоже не пригоден, если применять его "в лоб".
***
Что за итерационный метод? Вы хоть дайте ссылку на свой источник. Описывать здесь не надо, у вас это плохо получается, из кода ничего не понятно - там я вижу применение в лоб разностной схемы, которую применять нельзя.
Цитата Сообщение от супер тупой Посмотреть сообщение
Построил и решил.
И, каким методом? Что делали? СЛАУ составляли? до сих пор непонятно, о какой СЛАУ речь?
0
супер тупой
6 / 6 / 3
Регистрация: 29.08.2014
Сообщений: 89
Завершенные тесты: 1
04.04.2016, 01:27  [ТС] #12
cmath,
Цитата Сообщение от cmath Посмотреть сообщение
Вы хоть дайте ссылку на свой источник.
Н.В. Копченова, И.А. Марон - Вычислительная математика в примерах и задачах(стр. 266)

Цитата Сообщение от cmath Посмотреть сообщение
И, каким методом? Что делали?
Как-то так:
Кликните здесь для просмотра всего текста

http://www.cyberforum.ru/cgi-bin/latex.cgi?\begin{matrix}\\ {y}^{''} = -\sin \left(x \right), \:  {y}^{''} = \frac{{y}_{i+1} - 2{y}_{i} + {y}_{i-1}}{{h}^{2}}\\=> \, \frac{{y}_{i+1} - 2{y}_{i} + {y}_{i-1}}{{h}^{2}} = -\sin \left({x}_{i} \right)\\\\...\\\\{y}_{i} = \frac{1}{2}\,\left({y}_{i+1}\,+\,{y}_{i-1}\,-\,{h}^{2}\sin \left({x}_{i} \right) \right)\: ,\\\\{x}_{i} = a + ih\, ,\, y(a) = {y}_{0}\, ,\, y(b) = {y}_{1}\end{matrix}
Получилось
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
#include <iostream>
#include <vector>
#include <cmath>
 
double orig(double x) {//искомая функция
    return sin(x);
}
double g(double x) {//y'' = g(x)
    return -sin(x);
}
 
int main() {
    const size_t N = 101;//количество узлов
 
    const double a   = .0,//левый край
                 b   = 3.14159265359,//правый край
                 y0  = .0,//y(a)
                 y1  = .0,//y(b)
                 h   = (b - a) / (N - 1),//шаг
                 h2  = h * h,//h^2
                 EPS = 1e-10;
 
 
    std::vector<double> y(N, .0);//вектор приближений
    y[0]     = y0;//заполняем граничные элементы
    y[N - 1] = y1;
 
    double max_diff, diff, new_val;//максимальная разница на итерации, разница, новое значение
    do {
        max_diff = 0;
        for(size_t i = 1; i < N - 1; i++) {//проходим по всем узлам
            new_val = 0.5 * (y[i + 1] + y[i - 1] - h2 * g(a + i * h));//вычисляем новое приближение
 
            diff = fabs(new_val - y[i]);//вычисляем разницу нового и старого приближения
            if (diff > max_diff) max_diff = diff;//если она максимальна, запоминаем ее
 
            y[i] = new_val;//записываем новое приближение
        }
    } while(max_diff > EPS);//считаем, пока максимальная разница не станет меньше заданной точности
 
    for (size_t i = 0; i < N; i++)//для проверки
        std::cout << "y[" << i << "] = " << y[i] << "; diff = " << fabs(orig(a + i * h) - y[i]) << ";\n";
 
    return 0;
}


Цитата Сообщение от cmath Посмотреть сообщение
до сих пор непонятно, о какой СЛАУ речь?
Как видите, обошелся без СЛАУ.
0
cmath
Модератор
2490 / 1714 / 145
Регистрация: 11.08.2012
Сообщений: 3,293
Завершенные тесты: 6
04.04.2016, 03:03 #13
Цитата Сообщение от супер тупой Посмотреть сообщение
Как видите, обошелся без СЛАУ.
И зря. В вашей книге, про ваш же метод, в частности, написано, что перед применением процесса усреднения Либмана нужно сначала решить задачу обычным сеточным методом на какой-нибудь крупной сетке либо с помощью интерполяции по граничным условиям. Ни того, ни другого ваш код не делает. При этом интерполировать с нулевыми граничными условиями смысла нет. По этой причине имеется условие при http://www.cyberforum.ru/cgi-bin/latex.cgi?\psi=\frac{\pi}{4} - как опорная точка для интерполяции многочленами степени 1 или 2, но в общем случае неизвестно, что это будут за условия, поэтому следует вычислять, на первой итерации, с помощью обычного сеточного метода. Впрочем, в конкретно вашей задаче этого можно и не делать, тут я погорячился, признаю. Но интерполяцию сделать нужно обязательно. Хотя бы линейную, но лучше квадратичную, тогда ваш процесс быстрее сойдется. Код из поста #1 советую все-таки переписать начисто, а не переделывать.
1
bobah16
370 / 340 / 42
Регистрация: 14.07.2015
Сообщений: 2,882
04.04.2016, 18:53 #14
супер тупой
Код внутри отчета.
0
Вложения
Тип файла: docx отчёт(2 задание).docx (42.7 Кб, 13 просмотров)
VSI
04.04.2016, 20:17
  #15

Не по теме:

bobah16, при попытке открыть Ваш файл возникает "ошибка при открытии файла". Может быть это только у меня так, но лучше Вам выложить код отдельно от doc'овского документа, тем более, что если Вы внимательно прочитаете пункт 5.18. Правил форума... :rtfm:

0
bobah16
370 / 340 / 42
Регистрация: 14.07.2015
Сообщений: 2,882
04.04.2016, 20:41 #16

Не по теме:

Цитата Сообщение от VSI Посмотреть сообщение
ошибка при открытии файла
У меня все ок.
Цитата Сообщение от VSI Посмотреть сообщение
если Вы внимательно прочитаете пункт 5.18. Правил форума...
А вы предлагаете мне перепечатать весь отчет с кучей формул сюда? :)
Это ж ведь не решение и не задание, а всего лишь пример с объяснениями+небольшое исследование.


В док формате должен открыться у всех.
Fortran
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
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
      program poiss_eq
      implicit none
      include 'link_fnl_static.h'
      
      parameter nx=20,ny=20   !число узлов по х и у 
      integer i,j,k
      parameter lda=nx+1-2,ldfact=nx+1-2,ncoda=nx-2
      real*8 k1,k2,g1,g2,g3,g4,f,kappa3,kappa4,u
      real*8 fact(nx+1-2,(nx-2)*ny),rcond,eps
      real*8 hx,hy,x(nx),y(nx),ab(nx-2+1,(nx-2)*ny),v((nx-2)*ny),
     &  b((nx-2)*ny),err((nx-2)*ny)
                               !ab матрица коэфф интегро-интерп метода
                               !v вектор искомых значений ф-ции
                               !b вектор правой части
      !задание интервалов интегрирования и узлов
      open(1,file='out.dat')
      kappa3=1.0d0
      kappa4=1.0d0
!      kappa3=-12.0d0
!      kappa4=24.0d0
      x(1)=1.0d0
      x(nx)=4.0d0
      y(1)=4.0d0
      y(nx)=8.0d0
      hx=(x(nx)-x(1))/(nx-1)
      hy=(y(nx)-y(1))/(ny-1)
!      write(1,*) '    x(i)    ','       y(j)'
!      write(1,*) x(1),y(1)
      do i=2,nx
       x(i)=x(i-1)+hx
       y(i)=y(i-1)+hy
!       write(1,*) x(i),y(i)
      end do
      write(1,*)
      
      !задание коэфф во внутр точках с учётом гр усл 1 рода по х 
      do j=2,ny-1
       do i=2,nx-1
        k=(j-1)*(nx-2)+(i-1) 
        ab(nx-2+1,k)=hy/hx*(k1(x(i)+0.5d0*hx,y(j))+k1(x(i)-0.5d0*
     &  hx,y(j)))+hx/hy*(k2(x(i),y(j)+0.5d0*hy)+k2(x(i),y(j)-0.5d0*hy))
        if(i.ne.(nx-1)) then
        ab(nx-2,k+1)=-hy/hx*k1(x(i)+0.5d0*hx,y(j))
        end if
        ab(1,k+nx-2)=-hx/hy*k2(x(i),y(j)+0.5d0*hy)
        if((3-(nx-2)).le.0) then
        b(k)=hx*hy*f(x(i),y(j))
        end if
        if(i.eq.2) then
        b(k)=hx*hy*f(x(i),y(j))+hy/hx*k1(x(i)-0.5d0*hx,y(j))*g1(y(j))
        end if
        if(i.eq.(nx-1)) then
        b(k)=hx*hy*f(x(i),y(j))+hy/hx*k1(x(i)+0.5d0*hx,y(j))*g2(y(j))
        end if
       end do
      end do
          
      !гр усл 3 рода по у c учётом гран условий 1 рода по х
      do i=2,nx-1
       j=1
       k=(j-1)*(nx-2)+(i-1) 
       ab(nx-2+1,k)=0.5d0*hy/hx*(k1(x(i)+0.5d0*hx,y(j))+k1(x(i)-0.5d0*
     & hx,y(j)))+hx/hy*k2(x(i),y(j)+0.5d0*hy)+hx*kappa3
       if(i.ne.(nx-1)) then
       ab(nx-2,k+1)=-0.5d0*hy/hx*k1(x(i)+0.5*hx,y(j))                   
       end if          
       ab(1,k+nx-2)=-hx/hy*k2(x(i),y(j)+0.5d0*hy) 
       if((3-(nx-2)).le.0) then
       b(k)=0.5d0*hx*hy*f(x(i),y(j))+hx*g3(x(i))
       end if
       if(i.eq.2) then
       b(k)=0.5d0*hx*hy*f(x(i),y(j))+hx*g3(x(i))+0.5d0*hy/hx*k1(x(i)-
     & 0.5d0*hx,y(j))*g1(y(j))
       end if
       if(i.eq.(nx-1)) then
       b(k)=0.5d0*hx*hy*f(x(i),y(j))+hx*g3(x(i))+0.5d0*hy/hx*k1(x(i)+
     & 0.5d0*hx,y(j))*g2(y(j))
       end if
       j=ny 
       k=(j-1)*(nx-2)+(i-1)           
       ab(nx-2+1,k)=0.5d0*hy/hx*(k1(x(i)+0.5d0*hx,y(j))+k1(x(i)-0.5d0*
     & hx,y(j)))+hx/hy*k2(x(i),y(j)-0.5d0*hy)+hx*kappa4
       if(i.ne.(nx-1)) then
       ab(nx-2,k+1)=-0.5d0*hy/hx*k1(x(i)+0.5d0*hx,y(j))
       end if
       if((3-(nx-2)).le.0) then
       b(k)=0.5d0*hx*hy*f(x(i),y(j))+hx*g4(x(i))
       end if
       if(i.eq.2) then
       b(k)=0.5d0*hx*hy*f(x(i),y(j))+hx*g4(x(i))+
     & 0.5d0*hy/hx*k1(x(i)-0.5d0*hx,y(j))*g1(y(j))
       end if
       if(i.eq.(nx-1)) then
       b(k)=0.5d0*hx*hy*f(x(i),y(j))+hx*g4(x(i))+
     & 0.5d0*hy/hx*k1(x(i)+0.5d0*hx,y(j))*g2(y(j))
       end if      
      end do  
      
      call dlfcqs((nx-2)*ny,ab,lda,ncoda,fact,ldfact,rcond)  
      call dlfsqs((nx-2)*ny,fact,ldfact,ncoda,b,v)
     
      do j=1,ny
       do i=2,nx-2
        k=(j-1)*(nx-2)+i-1
        err(k)=dabs(v(k)-u(x(i),y(j)))
       end do
      end do
      eps=0.0d0
      do i=1,(nx-2)*ny
       if(err(i).gt.eps) then
       eps=err(i)
       end if
      end do
 
      write(1,*)
      write(1,*)'rcond=',1.0d0/rcond
      write(1,*) eps
      end
       
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!      
      double precision function k1(x,y)                                 
      real*8 x,y
      k1=1.0d0
!      k1=y*x**2.0d0
      end
      
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! 
      double precision function k2(x,y)
      real*8 x,y
      k2=1.0d0
!      k2=y**2.0d0*x
      end 
      
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!      
      double precision function g1(y)
      real*8 y
      g1=y
!      g1=y**3.0d0
      end
      
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!      
      double precision function g2(y) 
      real*8 y
      g2=4.0d0*y
!      g2=64.0d0*y**3.0d0
      end
      
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
      double precision function g3(x)
      real*8 x
      g3=3.0d0*x
!      g3=-768.0d0*x**3.0d0*(1.0d0+x)
      end
      
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!  
      double precision function g4(x) 
      real*8 x
      g4=9.0d0*x
!      g4=12288.0d0*x**3.0d0*(1.0d0+x)
      end
      
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!  
      double precision function f(x,y) 
      real*8 x,y
      f=0.0d0
!      f=-12.0d0*x**3.0d0*y**3.0d0*(x+y)
      end
      
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!  
      double precision function u(x,y) 
      real*8 x,y
      u=x*y
!      u=x**3.0d0*y**3.0d0 
      end
      
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1
Вложения
Тип файла: rar отчёт(2 задание).rar (128.2 Кб, 2 просмотров)
супер тупой
6 / 6 / 3
Регистрация: 29.08.2014
Сообщений: 89
Завершенные тесты: 1
05.04.2016, 00:13  [ТС] #17
bobah16, спасибо большое!

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

Цитата Сообщение от cmath Посмотреть сообщение
По этой причине имеется условие при - как опорная точка для интерполяции многочленами степени 1 или 2
Эта точка имеется по некоторым другим причинам.

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

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

Не по теме:

Цитата Сообщение от cmath Посмотреть сообщение
Код из поста #1 советую все-таки переписать начисто, а не переделывать.
Да что же этот код вас так напугал :jokingly:
Он нигде не используется и я не просто не буду его переделывать, а вовсе удалю, как только разберусь с проблемой.

0
супер тупой
6 / 6 / 3
Регистрация: 29.08.2014
Сообщений: 89
Завершенные тесты: 1
07.04.2016, 18:08  [ТС] #18
К слову, вопрос остается открытым.
0
bobah16
370 / 340 / 42
Регистрация: 14.07.2015
Сообщений: 2,882
07.04.2016, 22:11 #19
супер тупой, ну так чем вам помочь? Перепилите свой код, и проверяйте на тестовых задачах. Для начала воспроизведите хотя бы константу с нулевой погрешностью. Если код ее не воспроизводит значит ошибка в схеме. Читайте Самарский А.А. Теория разностных схем. Мне эта книга больше всего нравится.
0
супер тупой
6 / 6 / 3
Регистрация: 29.08.2014
Сообщений: 89
Завершенные тесты: 1
07.04.2016, 23:04  [ТС] #20
bobah16,
Вот этим и помочь! Константу воспроизводит отлично. Самарского пошел читать.
0
07.04.2016, 23:04
MoreAnswers
Эксперт
37091 / 29110 / 5898
Регистрация: 17.06.2006
Сообщений: 43,301
07.04.2016, 23:04

Численное решение трансцендентного уравнения
Помогите решить e^(-x)-(x-1)=0 различными Методами: хорд, касательной,...

Численное решение уравнения Шредингера
По работе столкнулся со следующей задачей. Необходимо численно решить...

Численное решение задачи Коши dy/dx = xy ; y(0)=0
Здравствуйте, постановка задачи: dy/dx=xy y(0)=0 Допустим шаг h = 0.1,...


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

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

КиберФорум - форум программистов, компьютерный форум, программирование
Powered by vBulletin® Version 3.8.9
Copyright ©2000 - 2018, vBulletin Solutions, Inc.
Рейтинг@Mail.ru