Изменение имён переменных

Для более точного соответствия отчёту

include/Header.hpp

@@ -2,7 +2,7 @@
 #include <vector>
 #include <cmath>
 #include <fstream>
-#include <iomanip>
+#include <chrono>
 
 #define M_PI 3.14159265358979323846
 

src/main.cpp

@@ -9,6 +9,13 @@ double t_cur = 0;
 using namespace std;
 using namespace Eigen;
 
+VectorXd Phi = VectorXd::Zero(DIMENTION);
+VectorXd delta_Phi(DIMENTION);
+VectorXd prev_Phi = VectorXd::Zero(DIMENTION);
+VectorXd prev_prev_Phi = VectorXd::Zero(DIMENTION);
+MatrixXd Y = MatrixXd::Zero(DIMENTION, DIMENTION);
+VectorXd I(DIMENTION);
+
 double E1(double t) {
     const double amplitude = 10.0;
     const double period = 1e-4;
@@ -18,8 +25,8 @@ double E1(double t) {
 }
 
 // Заполнение матрицы
-void fill_Jacobian(MatrixXd& Y, const VectorXd& cur_X) {
+void fill_Jacobian() {
-    double d_Id = It / MFt * exp(cur_X[U_Cb] / MFt);
+    double d_Id = It / MFt * exp(Phi[U_Cb] / MFt);
 
     for (int i = 0; i < 3; i++) {
         Y(i, i) = 1.0;
@@ -65,58 +72,54 @@ void fill_Jacobian(MatrixXd& Y, const VectorXd& cur_X) {
     Y(I_E2, phi_4) = -1.;
 }
 // Заполнение вектора невязок
-void fill_Nevyzka(VectorXd& b, const VectorXd& cur_X, const VectorXd& prev_X) {
+void fill_Nevyzka() {
-    double Id = It * exp(cur_X[U_Cb] / MFt);
+    double Id = It * exp(Phi[U_Cb] / MFt);
 
     // Первые три уравнения (производные)
-    b[dU_C1] = cur_X[dU_C1] - (cur_X[U_C1] - prev_X[U_C1]) / dt;
+    I[dU_C1] = Phi[dU_C1] - (Phi[U_C1] - prev_Phi[U_C1]) / dt;
-    b[dU_Cb] = cur_X[dU_Cb] - (cur_X[U_Cb] - prev_X[U_Cb]) / dt;
+    I[dU_Cb] = Phi[dU_Cb] - (Phi[U_Cb] - prev_Phi[U_Cb]) / dt;
-    b[dI_L] = cur_X[dI_L] - (cur_X[I_L] - prev_X[I_L]) / dt;
+    I[dI_L] = Phi[dI_L] - (Phi[I_L] - prev_Phi[I_L]) / dt;
 
     // Уравнения для переменных состояния
-    b[U_C1] = cur_X[U_C1] - cur_X[phi_5];
+    I[U_C1] = Phi[U_C1] - Phi[phi_5];
-    b[U_Cb] = cur_X[U_Cb] - (cur_X[phi_3] - cur_X[phi_4]);
+    I[U_Cb] = Phi[U_Cb] - (Phi[phi_3] - Phi[phi_4]);
-    b[I_L] = L1 * cur_X[dI_L] - (cur_X[phi_2] - cur_X[phi_5]);
+    I[I_L] = L1 * Phi[dI_L] - (Phi[phi_2] - Phi[phi_5]);
 
     // Уравнения Кирхгофа для узлов
-    b[phi_1] = cur_X[I_E1] + (cur_X[phi_1] - cur_X[phi_2]) / R1;
+    I[phi_1] = Phi[I_E1] + (Phi[phi_1] - Phi[phi_2]) / R1;
-    b[phi_2] = -(cur_X[phi_1] - cur_X[phi_2]) / R1 + cur_X[I_L] + (cur_X[phi_2] - cur_X[phi_3]) / Rb;
+    I[phi_2] = -(Phi[phi_1] - Phi[phi_2]) / R1 + Phi[I_L] + (Phi[phi_2] - Phi[phi_3]) / Rb;
-    b[phi_3] = -(cur_X[phi_2] - cur_X[phi_3]) / Rb + Id + (cur_X[phi_3] - cur_X[phi_4]) / Ru + Cb * cur_X[dU_Cb];
+    I[phi_3] = -(Phi[phi_2] - Phi[phi_3]) / Rb + Id + (Phi[phi_3] - Phi[phi_4]) / Ru + Cb * Phi[dU_Cb];
-    b[phi_4] = -Id - (cur_X[phi_3] - cur_X[phi_4]) / Ru - Cb * cur_X[dU_Cb] + cur_X[I_E2];
+    I[phi_4] = -Id - (Phi[phi_3] - Phi[phi_4]) / Ru - Cb * Phi[dU_Cb] + Phi[I_E2];
-    b[phi_5] = -cur_X[I_L] + C1 * cur_X[dU_C1] + cur_X[phi_5] / R2;
+    I[phi_5] = -Phi[I_L] + C1 * Phi[dU_C1] + Phi[phi_5] / R2;
 
     // Уравнения для источников
-    b[I_E1] = E1(t_cur) - cur_X[phi_1];
+    I[I_E1] = E1(t_cur) - Phi[phi_1];
-    b[I_E2] = E2 - cur_X[phi_4];
+    I[I_E2] = E2 - Phi[phi_4];
 
-    b = -1 * b;
+    I = -1 * I;
 }
 
 int main() {
     ofstream result_file(RESULT);
-    VectorXd delta_X(DIMENTION);
-    VectorXd cur_X = VectorXd::Zero(DIMENTION);
-    VectorXd prev_X = VectorXd::Zero(DIMENTION);
-    VectorXd prev_prev_X = VectorXd::Zero(DIMENTION);
-    MatrixXd Y = MatrixXd::Zero(DIMENTION, DIMENTION);
-    VectorXd b(DIMENTION);
     double step_count = 0;              // Счетчик числа шагов
     double it_count = 0;                // Счетчик числа итераций метода Ньютона
     int counter = 0;                    // Счетчик числа шагов
+    
+    auto start = std::chrono::steady_clock::now();
     while (t_cur < T) {
         bool is_end = false;            // Индикатор завершения итераций метода Ньютона
         size_t n = 0;                   // Счетчик числа итераций метода Ньютона
         while (!is_end) {
             // Заполнение матрицы узловых проводимостей и вектора невязок значениями
-            fill_Jacobian(Y, cur_X);
+            fill_Jacobian();
-            fill_Nevyzka(b, cur_X, prev_X);
+            fill_Nevyzka();
 
             // Решение СЛАУ
-            delta_X = Y.partialPivLu().solve(b);
+            delta_Phi = Y.partialPivLu().solve(I);
-            cur_X += delta_X;
+            Phi += delta_Phi;
 
             // Проверка точности
-            is_end = delta_X.cwiseAbs().maxCoeff() <= ACR;
+            is_end = delta_Phi.cwiseAbs().maxCoeff() <= ACR;
             n++;
             // Если требуемая точность не достигнута
             if (!is_end) {
@@ -124,7 +127,7 @@ int main() {
                 if (n > N_iter) {
                     n = 0;
                     dt *= 0.5;
-                    cur_X = prev_X;
+                    Phi = prev_Phi;
                     if (dt < SMN) {
                         throw domain_error("Решение не сходится, dt < SMN");
                     }
@@ -135,42 +138,46 @@ int main() {
         // Расчет текущей ошибки метода Ньютона
         for (int i = dU_C1; i < DIMENTION; i++) {
             //if (i == dU_Cb) continue;
-            double tmp = 0.5 * dt * dt * abs(((cur_X(i) - prev_X(i)) / dt - (prev_X(i) - prev_prev_X(i)) / prev_dt));
+            double tmp = 0.5 * dt * dt * abs(((Phi(i) - prev_Phi(i)) / dt - (prev_Phi(i) - prev_prev_Phi(i)) / prev_dt));
             cur_delta = (tmp > cur_delta) ? tmp : cur_delta;
         }
         // Если интегрирование неудовлетворительно по точности
         if (cur_delta > EPS_MAX && dt > SMN) {
             // Уменьшаем шаг по времени и отбрасываем результаты текущего шага
             dt *= 0.5;
-            cur_X = prev_X;
+            Phi = prev_Phi;
         }
 
         else {
             // Сохранение значений с предыдущего шага
-            prev_prev_X = prev_X;
+            prev_prev_Phi = prev_Phi;
-            prev_X = cur_X;
+            prev_Phi = Phi;
             prev_dt = dt;
             // Вывод значения переменной на текущем временном шаге
             result_file
                 << t_cur << ' '
-                << cur_X(phi_1) << ' '
+                << Phi(phi_1) << ' '
-                << cur_X(phi_2) << ' '
+                << Phi(phi_2) << ' '
-                << cur_X(phi_3) << ' '
+                << Phi(phi_3) << ' '
-                << cur_X(phi_4) << ' '
+                << Phi(phi_4) << ' '
-                << cur_X(phi_5) << std::endl;
+                << Phi(phi_5) << std::endl;
             // Шаг по времени
             t_cur += dt;
             counter++;
             step_count += dt;
             it_count += n;
-            cout << counter << " n_it=" << n << " t=" << t_cur << endl;
+            //cout << counter << " n_it=" << n << " t=" << t_cur << endl;
             if (cur_delta < EPS_MIN) {
                 dt *= 2;
             }
         }
     }
-    cout << "Число шагов по времени: " << counter << endl;
+    auto end = std::chrono::steady_clock::now();
-    cout << "Средний шаг по времени: " << step_count / counter << endl;
+    auto sec = std::chrono::duration<double>(end - start);
-    cout << "Среднее число итераций: " << it_count / counter << endl;
+
+    std::cout << "Time: " << sec.count() << " sec\n";
+    cout << "Time steps: " << counter << endl;
+    cout << "Average step by time: " << step_count / counter << endl;
+    cout << "Average iteration count: " << it_count / counter << endl;
     return 0;
 }