Форум программистов, компьютерный форум CyberForum.ru

Решение СЛАУ методом сопряженных градиентов - C++

Восстановить пароль Регистрация
 
Рейтинг: Рейтинг темы: голосов - 26, средняя оценка - 4.62
ftarr
0 / 0 / 0
Регистрация: 30.11.2012
Сообщений: 5
23.12.2012, 18:25     Решение СЛАУ методом сопряженных градиентов #1
Всем доброго времени суток!

Суть вопроса. есть код (ниже) по решению СЛАУ методом сопряженных градиентов. делал по алгоритму, описанному тут

http://www.hpcc.unn.ru/?dir=847

пример для проверки

3x0 -x1 = 3,
-x0 +3x1 = 7.

На первой итерации градиент g1=(-3,-7), вектор направления d1=(3, 7), значение величины смещения s1=0,439. очередное приближение к точному решению системы x1=(1,318, 3,076).

На второй итерации градиент g2=(-2,121, 0,909), вектор направления d2=(2,397, -0,266), а величина смещения – s2=0,284. Очередное приближение совпадает с точным решением системы x2=(2, 3).

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

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
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
#include <iostream>
#include <math.h>
#include <fstream>
 
using namespace std;
 
typedef double arr[20];
 
int main()
{
    bool f;
    int n,i,j,k;
    double a[20][20];
    arr b,xk,xkk,gk,gkk,dk,dkk,c;
    ifstream inp("input.txt"); // создаём объект класса ofstream для записи и связываем его с файлом cppstudio.txt
    inp >> n;
    double t,e,p1,p2,p,sk,xkp,xkkp,z;
    for(i=0;i<n;i++)
        for(j=0;j<n;j++)
        {
            inp >> t;
            a[i][j]=t;
            //cout << t << " ";
        }
    //cout << endl;
    for(i=0;i<n;i++)
    {
        inp >> t;
        //cout << t << " ";
        b[i]=t;
        gk[i]=0;
        gkk[i]=0;
        dk[i]=0;
        dkk[i]=0;
        xk[i]=0;
    }
    cout << endl;
    inp.close();
    cout << "Введите точность: ";
    cin >> e;
    k=0;
    f = false;
    do
    {
        //-------------------шаг 1 gk= A*xkk - b, k - текущая итерация, kk - предыдущая
        if (k==1) f=true;
        //cout << "шаг 1 " << k;
        for(i=0;i<n;i++)
        {
            c[i]=0;
            for(j=0;j<n;j++)
            c[i]=c[i]+a[i][j]*xk[i];
        }
        for(i=0;i<n;i++)
            gk[i]=c[i]-b[i];
        //cout << "вектор gk ";
        //for(i=0;i<n;i++)
        //  cout << gk[i] << " ";
        //cout << endl;
 
        //------------------шаг 2 dk=...
        //cout << "шаг 2 ";
        p1=0; p2=0; p=0;
        for(i=0;i<n;i++)
        {
            p1=p1+gk[i]*gk[i];
            p2=p2+gkk[i]*gkk[i];
        //  cout << "p1 " << p1 << " p2 " << p2 << endl;
        }
        if (p2 !=0)
        p=p1/p2;
        for(i=0;i<n;i++)
            dk[i]= -gk[i] + p*dkk[i];
 
        //cout << "вектор dk ";
        //for(i=0;i<n;i++)
         //cout << dk[i] << " ";
        //cout << endl;
 
        //-----------------шаг 3 sk=...
        //cout << "шаг 3 ";
        p1=0;p2=0;
        for(i=0;i<n;i++)
            p1=p1+dk[i]*gk[i];
 
        for(i=0;i<n;i++)
        {
            c[i]=0;
            for(j=0;j<n;j++)
            {
                //cout << endl;
                //cout << " считаем С в SK ";
                c[i]=c[i]+dk[j]*a[i][j];
                //cout << " i,j " << i << "," << j << " dkj " << dk[j] << " aij " << a[i][j] << " ci " << c[i] << endl;
            }
        }
 
        for(i=0;i<n;i++)
            p2=p2+c[i]*dk[i];
 
        //cout << "п1 и п2" << p1 << " " << p2;
        sk=p1/p2;
 
        //cout << "sk " << sk << endl;
        //------------------шаг 4 xk=
        //cout << "шаг 4 ";
        for(i=0;i<n;i++)
            xk[i]=xkk[i]+sk*dk[i];
 
        //cout << "вектор xk ";
            //for(i=0;i<n;i++)
            //  cout << xk[i] << " ";
        //cout << endl;
 
        //-----------------постшаг
        for(i=0;i<n;i++)
        {
            dkk[i]=dk[i];
            gkk[i]=gk[i];
            xkk[i]=xk[i];
        }
 
        // --------------счет для точности
        p1=0;p2=0;
        for(i=0;i<n;i++)
        {
            p1=p1+xk[i]*xk[i];
            p2=p2+xkk[i]*xkk[i];
        }
        xkp =sqrt(p1);
        xkkp =sqrt(p2);
        if (xkp>xkkp)
            z=xkp-xkkp;
        else
            z=xkkp-xkp;
        cout << z <<" ";
        k=k+1;
        if (k<2) z = e*2;
    }
    while (z>=e);
 
    for(i=0;i<n;i++)
        cout << xk[i] << " ";
    return 0;
}
ЗЫ: если есть советы по оптимизации или исключению ненужных костылей - буду премного блгодарен.
Similar
Эксперт
41792 / 34177 / 6122
Регистрация: 12.04.2006
Сообщений: 57,940
23.12.2012, 18:25     Решение СЛАУ методом сопряженных градиентов
Посмотрите здесь:

C++ Решение СЛАУ методом Крамера
Решение СЛАУ методом Гаусса C++
C++ Решение СЛАУ методом Якоби
C++ Решение СЛАУ большой размерности методом сопряженных градиентов
C++ Решение СЛАУ методом Гаусса
Решение СЛАУ методом прогонки C++
Решение СЛАУ методом вращения C++
C++ Решение слау методом релаксации

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

Или воспользуйтесь поиском по форуму:
После регистрации реклама в сообщениях будет скрыта и будут доступны все возможности форума.
Ответ Создать тему
Опции темы

Текущее время: 00:22. Часовой пояс GMT +3.
КиберФорум - форум программистов, компьютерный форум, программирование
Powered by vBulletin® Version 3.8.9
Copyright ©2000 - 2016, vBulletin Solutions, Inc.
Рейтинг@Mail.ru