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
|
#include <cmath>
#include <iostream>
#include <ranges>
#include "include/praktable.hpp"
#include "include/prakphys.hpp"
using table = prak::table<f64>;
using f64p = prak::pvalue<f64>;
using f64v = std::vector<f64>;
using vecarg = const std::vector<f64> &;
using argvec = const std::vector<f64> &;
f64p add_nu(table& t) {
t.add_column("nu").apply([](vecarg a) { return a[0] * a[1] / (prak::R<f64> * a[2]); }, {"P", "V", "T"}, "nu");
return {t.col_avg("nu"), t.col_stddev("nu")};
}
void ex1(table &t1, const char *name) {
// t = const
f64 DA = prak::discrete_integral_trapezoid(t1.cbegin("V"), t1.cend("V"), t1.cbegin("P"), t1.cend("P"));
f64 DS = DA / t1.col_avg("T");
f64p nu = add_nu(t1);
t1 .add_columns({"DAi", "DSi", "lnViV0"})
.apply([&t1](vecarg a) -> f64 {
return prak::discrete_integral_trapezoid(
t1.begin("V"), std::find(t1.begin("V"), t1.end("V"), a[1])+1,
t1.begin("P"), std::find(t1.begin("P"), t1.end("P"), a[0])+1
);
}, {"P", "V"}, "DAi")
.apply([&t1](vecarg a){return a[0] / a[1];}, {"DAi", "T"}, "DSi")
.apply([&t1](vecarg a){ return std::log(a[0] / t1["V", 0]); }, "V", "lnViV0")
.apply([](vecarg a) {return std::abs(a[0]); }, "lnViV0", "lnViV0");
f64p DS_teor = nu * prak::R<f64> * std::abs(std::log(*(t1.end("V")-1) / *t1.begin("V")));
f64p A = std::abs(t1.least_squares_prop("lnViV0", "DSi", std::nullopt, 0.1));
f64p R = A / nu;
std::cout << "Упражнение 1\n" var(DA) var(DS) var(nu) var(DS_teor) << "доп данные (упражнение 4):\n"
var(A) var(R) << std::endl;
t1 .write_plot(std::string("PV") + name + ".plot", "V", "P", std::nullopt)
.write_plot(std::string("PVV") + name + ".plot", "V", "x", std::nullopt)
.write_plot(std::string("DSilnViV0_Tconst") + name + ".plot", "lnViV0", "DSi", std::nullopt);
std::ofstream f(std::string("ex1_") + name + ".print");
f << "Данные упражнения 1 (T = const)\n";
t1.print(f);
f.close();
}
void ex2(table &t) {
// p = const
std::cout << "Упражнение 2\n";
std::ofstream f("ex2.print");
f << "Данные упражнения 2 (P = const)\n";
f64p nu = add_nu(t);
t .add_columns({"lnTiT0", "DSi"})
.apply([&t](vecarg a) -> f64 {return std::log(a[0] / t["T", 0]); }, "T", "lnTiT0")
.apply([&t, &nu](vecarg a) -> f64 {
return nu.val * (prak::R<f64> * std::log(a[1] / t["V", 0]) + 5.0/2 * prak::R<f64> * a[2]);
}, {"T", "V", "lnTiT0"}, "DSi")
.write_plot("DSilnTiT0_Pconst.plot", "lnTiT0", "DSi", std::nullopt)
.write_plot("VTV.plot", "V", "x", std::nullopt)
.print(f);
f64p A = t.least_squares_prop("lnTiT0", "DSi", std::nullopt, 0.1);
f64p Cp = A / nu;
std::cout var(nu) var(A) var(Cp);
}
void ex3(table &t) {
// v = const
std::cout << "Упражнение 3\n";
std::ofstream f("ex3.print");
f << "Данные упражнения 3 (V = const)\n";
f64p nu = add_nu(t);
t .add_columns({"lnTiT0", "DSi"})
.apply([&t](vecarg a) -> f64 {return std::log(a[0] / t["T", 0]); }, "T", "lnTiT0")
.apply([&t, &nu](vecarg a) -> f64 {
return nu.val * (-prak::R<f64> * std::log(a[0] / t["P", 0]) + 7.0/2 * prak::R<f64> * a[1]);
}, {"P", "lnTiT0"}, "DSi")
.write_plot("DSilnTiT0_Vconst.plot", "lnTiT0", "DSi", std::nullopt)
.write_plot("PTP.plot", "P", "x", std::nullopt)
.print(f);
auto [A, B] = t.least_squares_linear("lnTiT0", "DSi", std::nullopt, 0.1);
f64p Cv = A / nu;
std::cout var(nu) var(A) var(Cv);
}
void ex(void) {
table tables[] = {table("data1"), table("data2"), table("data3"), table("data4")};
// approximate how the volume should have changed in ex3 using parabola
std::vector<f64> mrxdata; mrxdata.reserve(9);
std::vector<f64> mrydata; mrydata.reserve(3);
for (size_t i = 0; i < tables[2].rows; ++i) {
if (f64 v; !std::isnan(v = tables[2]["V", i])) {
f64 n = tables[2]["n", i];
mrxdata.push_back(n*n);
mrxdata.push_back(n);
mrxdata.push_back(1);
mrydata.push_back(v);
}
}
prak::matrix<f64> xs(3, 3, std::move(mrxdata)), ys(3, 1, std::move(mrydata));
prak::matrix<f64> poly = xs.inv().value() * ys;
for (size_t i = 0; i < tables[2].rows; ++i)
tables[2]["V", i] = (poly.tr() * prak::matrix<f64>(3, 1, {(double)(i+1)*(i+1), (double)(i+1), 1}))[0, 0];
std::function<f64(vecarg)> fs[] = {
[](vecarg a) { return a[0] * a[1]; },
[](vecarg a) { return a[0] / a[1]; },
[](vecarg a) { return a[0] / a[1]; },
[](vecarg a) { return a[0] * a[1]; },
};
std::vector<std::string> vs[] = {
{"P", "V"}, {"V", "T"}, {"P", "T"}, {"P", "V"},
};
for (auto elem : std::views::zip(tables, fs, vs))
std::get<0>(elem).apply(std::get<1>(elem), std::get<2>(elem), "x");
ex1(tables[0], "1");
ex1(tables[3], "4");
ex2(tables[1]);
ex3(tables[2]);
}
int main() {
ex();
return 0;
}
|