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
| #include <iostream>
#include <fstream>
#include <iomanip>
#include <cmath>
#include <string>
#define W_E1 3.1415926535897932384626433832795
#define E2 5.0
#define Cb 2e-12
#define C1 1e-6
#define C2 1e-6
#define I 1e-12
#define R 1e+6
#define Rb 20.0
#define R1 1000.0
#define R2 1000.0
#define R3 1000.0
#define L 1e-5
using namespace std;
double euler(double, double, double, double, double, double, double, double&, double&, double&, double&); //передача параметров по ссылке
int main()
{
double E1i,E10,ti, tf, dt, tmax, Ucbf, Uc1f, Uc2f, Ilf, Ucbi, Uc1i, Uc2i, Ili;
const string method = "Simple Euler";
/* output: file and formats */
ofstream file;
file.open ("ode.txt");
file.precision(6);
file.setf(ios::fixed | ios::showpoint);
cout.precision(6);
cout.setf(ios::fixed | ios::showpoint);
/* initial information */
ti = 0.0; // initial value for variable
dt = 0.00001; // step size for integration
tmax = 0.001 ; // inegrate from ti till tmax
Ucbi = 0.0;
Uc1i = 0.0;
Uc2i = 0.0;
Ili = 0.0;
E10 = 10*sin(20000*W_E1*ti);
file << setw(24) << method << endl;
file << setw(16) << "t" << setw(16) << "ucb" << setw(16) << "uc1" << setw(16) << "uc2" << setw(16) << "il" << setw(16) << "E1" << endl;
file << setw(16) << ti << setw(16) << Ucbi << setw(16) << Uc1i << setw(16) << Uc2i << setw(16) << Ili << setw(16) << E10 << endl;
/* step 2: integration of ODE */
while (ti <= tmax)
{
tf = ti + dt;
E1i = 10*sin(20000*W_E1*(ti+dt));
euler(E1i,ti,tf,Ucbi,Uc1i,Uc2i,Ili,Ucbf,Uc1f,Uc2f,Ilf);
file << setw(16) << tf << setw(16) << Ucbf << setw(16) << Uc1f << setw(16) << Uc2f << setw(16) << Ilf << setw(16) << E1i << endl;
ti = tf;
Ucbi = Ucbf;
Uc1i = Uc1f;
Uc2i = Uc2f;
Ili = Ilf;
}
cout << "Hello World!" << endl;
return 0;
}
double euler(double E1i, double ti, double tf, double Ucbi, double Uc1i, double Uc2i, double Ili, double& Ucbf, double& Uc1f, double& Uc2f, double& Ilf)
{
//Ucf = Uci + Ili/C*(tf-ti);
//Ilf = Ili + (E-Uci)/L*(tf-ti);
Ucbf = Ucbi + ((1/Rb)*(E1i - E2 - Uc1i - Ucbi) - (1/R)*Ucbi - I*exp(Ucbi/0.026-1))/Cb*(tf-ti);
Uc1f = Uc1i + (Ili + (1/Rb)*(E1i - E2 - Uc1i - Ucbi) - (1/R1)*Uc1i)/C1*(tf-ti);
Uc2f = Uc2i + (Ili - (1/R3)*Uc2i)/C2*(tf-ti);
Ilf = Ili - (Uc1i + Uc2i - E1i + R2*Ili)/L*(tf-ti);
} |