Advertisement
Georgiy031

Untitled

Dec 14th, 2020
82
0
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
C++ 3.82 KB | None | 0 0
  1. #define _CRT_SECURE_NO_WARNINGS
  2. #include <bits/stdc++.h>
  3. using namespace std;
  4. typedef long long ll;
  5. typedef unsigned long long ull;
  6. #define all(x) x.begin(), x.end()
  7. #define rall(x) x.rbegin(), x.rend()
  8. //#define endl '\n'
  9. #define boostIO() ios_base::sync_with_stdio(0);cin.tie(0);cout.tie(0);
  10. ll gcd(ll a, ll b) { return (b == 0 ? a : gcd(b, a % b)); }
  11.  
  12. typedef long double ld;
  13. typedef vector<vector<ld>> matrix;
  14.  
  15. matrix mul(const matrix& a, const matrix& b) {
  16.     matrix ans(a.size(), vector<ld>(b[0].size()));
  17.     for (int i = 0; i < a.size(); ++i)
  18.         for (int j = 0; j < b[0].size(); ++j)
  19.             for (int k = 0; k < b.size(); ++k)
  20.                 ans[i][j] += a[i][k] * b[k][j];
  21.  
  22.     return ans;
  23. }
  24. void print(string name, const matrix& a) {
  25.     cout << name << " = \n";
  26.     for (auto row : a) {
  27.         for (auto x : row) {
  28.             cout << fixed << setw(15) << x;
  29.         }
  30.         cout << endl;
  31.     }
  32. }
  33.  
  34. const ld EPS = 1e-6;
  35. const ld INF = LONG_MAX;
  36.  
  37. int gauss(matrix& a, vector<ld>& ans) {
  38.     int n = (int)a.size();
  39.     int m = (int)a[0].size() - 1;
  40.  
  41.     vector<int> where(m, -1);
  42.     for (int col = 0, row = 0; col < m && row < n; ++col) {
  43.         int sel = row;
  44.         for (int i = row; i < n; ++i)
  45.             if (abs(a[i][col]) > abs(a[sel][col]))
  46.                 sel = i;
  47.         if (abs(a[sel][col]) < EPS)
  48.             continue;
  49.         for (int i = col; i <= m; ++i)
  50.             swap(a[sel][i], a[row][i]);
  51.         where[col] = row;
  52.  
  53.         for (int i = 0; i < n; ++i)
  54.             if (i != row) {
  55.                 double c = a[i][col] / a[row][col];
  56.                 for (int j = col; j <= m; ++j)
  57.                     a[i][j] -= a[row][j] * c;
  58.             }
  59.         ++row;
  60.     }
  61.     //print("a", a);
  62.     ans.assign(m, 0);
  63.  
  64.     for (int i = 0; i < m; ++i)
  65.         if (where[i] != -1)
  66.             ans[i] = a[where[i]][m] / a[where[i]][i];
  67.     for (int i = 0; i < n; ++i) {
  68.         double sum = 0;
  69.         for (int j = 0; j < m; ++j)
  70.             sum += ans[j] * a[i][j];
  71.         if (abs(sum - a[i][m]) > EPS)
  72.             return 0;
  73.     }
  74.  
  75.     for (int i = 0; i < m; ++i)
  76.         if (where[i] == -1)
  77.             return INF;
  78.     return 1;
  79. }
  80.  
  81. ld a = -2;
  82. ld b = -1;
  83. ld y_a = 0;
  84. ld y_b = 0;
  85.  
  86. ld p(ld x) {
  87.     return 0;
  88. }
  89. ld q(ld x) {
  90.     return -2.0 / (x * x);
  91. }
  92. ld f(ld x) {
  93.     return 3 * log(-x);
  94. }
  95.  
  96. vector<ld> get_y(int n) {
  97.     ld h = (b - a) / n;
  98.     vector<ld> x(n + 1), y(n + 1);
  99.     x[0] = a;
  100.     for (int i = 1; i <= n; ++i) {
  101.         x[i] = a + i * h;
  102.     }
  103.     vector<vector<ld>> mtrx(n - 1, vector<ld>(n));
  104.     vector<ld> ans(n - 1);
  105.  
  106.     mtrx[0][0] = h * h * q(x[1]) - 2;
  107.     mtrx[0][1] = 1 + h / 2 * p(x[1]);
  108.     mtrx[0].back() = h * h * f(x[1]) - (1 - h / 2 * p(x[1])) * y_a;
  109.  
  110.     for (int i = 1; i < n - 2; ++i) {
  111.         mtrx[i][i - 1] = (1 - h / 2 * p(x[i + 1]));
  112.         mtrx[i][i] = h * h * q(x[i + 1]) - 2;
  113.         mtrx[i][i + 1] = 1 + h / 2 * p(x[i + 1]);
  114.  
  115.         mtrx[i].back() = h * h * f(x[i]);
  116.     }
  117.  
  118.     mtrx[n - 2][n - 3] = (1 - h / 2 * p(x[n - 1]));
  119.     mtrx[n - 2][n - 2] = h * h * q(x[n - 1]) - 2;
  120.     mtrx[n - 2].back() = h * h * f(x[n - 1]) - (1 + h / 2 * p(x[n - 1])) * y_b;
  121.  
  122.     gauss(mtrx, ans);
  123.     vector<ld> result(n + 1);
  124.     result[0] = y_a;
  125.     result[n] = y_b;
  126.     for (int i = 0; i < n - 1; ++i)
  127.         result[i + 1] = ans[i];
  128.     return result;
  129. }
  130.  
  131. signed main() {
  132.     for (int n = 5;; ++n) {
  133.         auto y1 = get_y(n); //size == n + 1
  134.         auto y2 = get_y(2 * n); //size == 2 * n + 1
  135.         ld max_eps = 0;
  136.         for (int i = 0; i <= n; ++i) {
  137.             ld mx = max(y1[i], y2[i * 2]);
  138.             ld mn = min(y1[i], y2[i * 2]);
  139.             max_eps = max(max_eps, mx / mn - 1);
  140.         }
  141.         if (max_eps < 0.01) {
  142.             cout << "n = " << n << "\n";
  143.             cout << "ANS: \n";
  144.             for (auto y : y1) cout << y << "\n";
  145.             break;
  146.         }
  147.     }
  148. }
  149.  
  150. /*
  151. 0 //
  152. 1 =====
  153. 2
  154. 3
  155.  
  156.  
  157. n - 3
  158. n - 2
  159. n - 1 ====
  160. n //
  161.  
  162. */
  163.  
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement