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

Метод монте карло для решения СЛАУ - C++

Восстановить пароль Регистрация
 
Рейтинг: Рейтинг темы: голосов - 10, средняя оценка - 4.90
olyaolua
0 / 0 / 0
Регистрация: 08.02.2012
Сообщений: 4
19.11.2012, 20:23     Метод монте карло для решения СЛАУ #1
не понимаю почему программа идеально работающая для системы 2χ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
#include "stdafx.h"
#include "stdafx.h"
#include <iostream>
#include <stdlib.h>
#include <time.h>
 
using namespace std; 
 
void main()
{
    int n; //Размерность системы
    float x=0,y=0; //Решение системы
    float **a; //Исходная матрица
    float *f; //Правая часть системы
    float *h;
    float *pi; //Вектор нач. вероятностей цепи Маркова
    float **p; //Матрица переходных состояний цепи Маркова
    int N; //Длина цепи Маркова
    int *i; //Цепь Маркова
    float *Q; //Веса состояний цепи Маркова
    float *ksi; //СВ
    int m; //Количество реализаций цепи Маркова
    float alpha; //блуждающаяСВ
 
    n=2;
 
    a=new float*[n];
 
    for(int k=0;k<n;k++)
    {
    a[k]=new float[n];
    }
    
    
    f=new float[n];
    h=new float[n];
    pi=new float[n];
    p=new float*[n];
 
    for(int k=0;k<n;k++)
    {
    p[k]=new float[n];
    }
    
 
    N=1000;
    i=new int[N+1];
    Q=new float[N+1];
    m=10000;
    
    ksi=new float[m];
    
    for(int j=0;j<m;j++)
        
            ksi[j]=0;
 
        a[0][0]=0.7; a[0][1]=-0.1;
        a[1][0]=-0.2; a[1][1]=0.4;
        
        f[0]=0.4;
        f[1]=0.5;
        
        h[0]=1;
        h[1]=0;
    
        pi[1]=0.5;
        pi[0]=0.5;
        
        p[0][0]=0.5; p[0][1]=0.5;
        p[1][0]=0.5; p[1][1]=0.5;
        //p[1][0]=0.0; p[1][1]=1.0;
 
    
 
    //Моделируем m цепей Маркова длины N
    for(int j=0;j<m;j++)
    {
        alpha=rand()/float(RAND_MAX);//случайное число от 0 до 1
        
        if(alpha<pi[0]) 
        
        i[0]=0; //реализуется 1-е состояние
        
        else 
        
        i[0]=1; //реализуется 2-е состояние
    
        
        for(int k=1;k<=N;k++)
        {
        
        alpha=rand()/float(RAND_MAX);
        
        if(alpha<0.5)
                
        i[k]=0;
            
        else 
        
        i[k]=1;
            
        }
    
 
        //Вычисляем веса цепи Маркова
        if(pi[i[0]]>0)
        
            Q[0]=h[i[0]]/pi[i[0]];
        
        else 
        
            Q[0]=0;
        
        
        for(int k=1;k<=N;k++)
        {
            if(p[i[k-1]][i[k]]>0)
            
            Q[k]=Q[k-1]*a[i[k-1]][i[k]]/p[i[k-1]][i[k]];
            
            else
            
            Q[k]=0;
            
        }
        
 
        for(int k=0;k<=N;k++)
        
            ksi[j]+=Q[k]*f[i[k]];
        
        
        
    }
        for(int k=0;k<m;k++)
        
        x=x+ksi[k];
        x=x/m;
        cout<<"x=";
        cout<<x;
        
    }
Добавлено через 8 минут
#include "stdafx.h"
#include "stdafx.h"
#include <iostream>
#include <stdlib.h>
#include <time.h>

using namespace std;

void main()
{
int n; //Размерность системы
float x=0,y=0; //Решение системы
float **a; //Исходная матрица
float *f; //Правая часть системы
float *h;
float *pi; //Вектор нач. вероятностей цепи Маркова
float **p; //Матрица переходных состояний цепи Маркова
int N; //Длина цепи Маркова
int *i; //Цепь Маркова
float *Q; //Веса состояний цепи Маркова
float *ksi; //СВ
int m; //Количество реализаций цепи Маркова
float alpha; //блуждающаяСВ

n=3;

a=new float*[n];

for(int k=0;k<n;k++)
{
a[k]=new float[n];
}


f=new float[n];
h=new float[n];
pi=new float[n];
p=new float*[n+1];

for(int k=0;k<n+1;k++)
{
p[k]=new float[n+1];
}


N=1;
i=new int[N+1];
Q=new float[N+1];
m=100;

ksi=new float[m];

for(int j=0;j<m;j++)

ksi[j]=0;

a[0][0]=0.6; a[0][1]=0.2;a[0][2]=0.1;
a[1][0]=-0.1; a[1][1]=-0.1;a[1][2]=0.3;
a[2][0]=0.5; a[2][1]=-0.1;a[2][2]=-0.3;

f[0]=0.9;
f[1]=0.7;
f[2]=0.4;

h[0]=0;
h[1]=0;
h[2]=1;

pi[0]=0.3333;
pi[1]=0.3333;
pi[2]=0.3333;

p[0][0]=0.3333; p[0][1]=0.3333;p[0][2]=0.3333;
p[1][0]=0.3333; p[1][1]=0.3333;p[1][2]=0.3333;
p[2][0]=0.3333; p[2][1]=0.3333;p[2][2]=0.3333;
p[3][0]=0.0; p[3][1]=0.0; p[3][2]=1.0;



//Моделируем m цепей Маркова длины N
for(int j=0;j<m;j++)
{
alpha=rand()/float(RAND_MAX);//случайное число от 0 до 1

if(alpha<pi[0])
{
i[0]=0; //реализуется 1-е состояние
}
else
{
i[0]=1; //реализуется 2-е состояние
}

if((pi[0]<alpha) && (alpha<pi[1]))

i[1]=0; //реализуется 1-е состояние

else

i[1]=2; //реализуется 2-е состояние

if((pi[1]< alpha) && (alpha<pi[2]))

i[2]=0; //реализуется 1-е состояние

else

i[2]=3; //реализуется 2-е состояние


for(int k=1;k<=N;k++)
{

alpha=rand()/float(RAND_MAX);

if(alpha<0.3333)

i[k]=0;

else

i[k]=1;

}


//Вычисляем веса цепи Маркова
if(pi[i[0]]>0)

Q[0]=h[i[0]]/pi[i[0]];

else

Q[0]=0;


for(int k=1;k<=N;k++)
{
if(p[i[k-1]][i[k]]>0)

Q[k]=Q[k-1]*a[i[k-1]][i[k]]/p[i[k-1]][i[k]];

else

Q[k]=0;

}


for(int k=0;k<=N;k++)

ksi[j]+=Q[k]*f[i[k]];



}
for(int k=0;k<m;k++)

x=x+ksi[k];
x=x/m;
cout<<"x=";
cout<<x;

}
Similar
Эксперт
41792 / 34177 / 6122
Регистрация: 12.04.2006
Сообщений: 57,940
19.11.2012, 20:23     Метод монте карло для решения СЛАУ
Посмотрите здесь:

метод Монте-Карло C++
C++ метод Монте-Карло
C++ Метод Монте-Карло
Решения кратного интеграла методом Монте Карло на С++ C++
C++ Метод монте Карло
C++ вроде метод монте карло
Метод Монте-Карло. Объем сферы C++
C++ Метод Монте-Карло

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

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

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