Ayuda con Codigo de Atractor y uso de if
Publicado por Alejandra (7 intervenciones) el 16/05/2015 23:14:27
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
#include <cmath>
#include <cstdlib>
#include <fstream>
#include <iostream>
#include <math.h>
using namespace std;
#include "vector.hpp"using namespace cpl;
double sigma = 10; // Lorenz model constants in textbook
double b = 8.0 / 3.0;
double r = 28;Vector f(Vector txyz) {
double t = txyz[0];
double x = txyz[1];
double y = txyz[2];
double z = txyz[3];
Vector f(4);
f[0] = 1;
f[1] = - sigma * x + sigma * y;
f[2] = - x * z + r * x - y;
f[3] = x * y - b * z;
return f;
}void RK4Step( // 4th order Runge-Kutta
Vector& y, // extended solution vector
double h) // step size
{Vector k1 = h * f(y);
Vector k2 = h * f(y + 0.5 * k1);
Vector k3 = h * f(y + 0.5 * k2);
Vector k4 = h * f(y + k3);
y += (k1 + 2 * k2 + 2 * k3 + k4) / 6.0;
}Vector txyz(4); // global variable to hold t,x,y,z
void initialize() { // initial conditions in textbook
txyz[0] = 0.0;
txyz[1] = 0.0;
txyz[2] = 1.0;
txyz[3] = 0.0;
}double dt = 0.001; // time step for integration
void findNextCrossing() {
Vector txyzOld = txyz; // save the old point
while (true) {
RK4Step(txyz, dt); // take a step
// check whether y changes sign, i.e., crosses y = 0
if (txyz[2] * txyzOld[2] < 0)
break;
else
txyzOld = txyz; }// use linear interpolation like in projectile.cpp to find intersection
double r = txyzOld[2] / txyz[2];
for (int i = 0; i < 4; i++)
txyz[i] = (txyzOld[i] + r * txyz[i]) / (r + 1);
}int main(int argc, char *argv[]) {
cout << " Trajectory and Poincare Section for the Lorenz Attractor\n"
<< " using 4th order Runge-Kutta with time step dt = " << dt << "\n"
<< " sigma = " << sigma << ", b = " << b << ", r = " << r << endl;
initialize();
cout << " initial conditions: x = " << txyz[1] << "\t"
<< ", y = " << txyz[2] << "\t" << ", z = " << txyz[3] << endl;
// transient trajectory
double t = 50;string fileName = "transient.data";
cout << " Integrating to time t = " << t << "\n"
<< " trajectory in file " << fileName << endl;
ofstream dataFile(fileName.c_str());
dataFile << txyz[0] << "\t" << txyz[1] << "\t"
<< txyz[2] << "\t" << txyz[3] << "\n";
int step = 0; int skip = 5;while (txyz[0] < t) {
RK4Step(txyz, dt);
if (++step % skip != 0) // record every skip steps
continue;
dataFile << txyz[0] << "\t" << txyz[1] << "\t"
<< txyz[2] << "\t" << txyz[3] << "\n";
}dataFile.close();
// Poincare section
fileName = "section.data";
int points = 1000;cout << " Finding " << points << " Poincare section points at y = 0\n"
<< " section data in file " << fileName << endl;
dataFile.open(fileName.c_str());
for (int point = 0; point < points; point++) {
findNextCrossing();
dataFile << txyz[0] << "\t" << txyz[1] << "\t"
<< txyz[2] << "\t" << txyz[3] << "\n";
} t = 1000;for (int point = 0; point < points; point++) {
RK4Step(txyz, dt);
if ((txyz[2] = sqrt(b * (r - 1))) && (txyz[1] >= sqrt(b *(r - 1))))
{ cout << " z = " << txyz[3] << "," << "1" << endl;
}if ((txyz[2] = (-1) * sqrt(b * (r - 1))) && (txyz[1] <= (-1) * (sqrt(b * (r - 1)))))
{ cout << " z = " << txyz[3] << "," << "0" << endl; }
}}Buenas ese es mi código pero no se que ocurre al final con los ciclos for e if tengo un problema no me entra a los ciclos, entonces no se que estoy colocando mal.
2do) En el segundo if, si yo coloco else if me arroja un resultado distinto, alguien podria ayudarme?
Gracias de antemano
Valora esta pregunta
0
