Advertisement
Guest User

Untitled

a guest
Oct 25th, 2014
118
0
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
C++ 1.87 KB | None | 0 0
  1. float max_nondiagonal_elem(int &l, int &m, matrix a)
  2. {
  3.     float max = abs(a[0][1]);
  4.     int i = 0, j = 1; l = i; m = j;
  5.     for (int i = 0; i < n - 1; i++)
  6.         for (int j = n - 1; j > i; j--)
  7.             if (abs(a[i][j]) > max)
  8.             {max = abs(a[i][j]); l = i; m = j;}
  9.     return a[l][m];
  10. }
  11.  
  12. float rot_condition(matrix a)
  13. {
  14.     float s = 0.0;
  15.     for (int i = 0; i < n; i++)
  16.         for (int j = 0; j < n; j++)
  17.             if (i != j) s += a[i][j] * a[i][j];
  18.     return sqrt(s);
  19. }
  20.  
  21. void iterative_method_of_rotations()
  22. {
  23.     int l = 0, m = 0; matrix h, ht, t, a1, u;
  24.     copy_matr(a, a1);
  25.     fprintf(fw, "Iterative method of rotations\n\n");
  26.     fprintf(fw, "Count of iterations = %d\n\n", 0);
  27.     float c = rot_condition(a1);
  28.     int k = 0;
  29.     printmatr(a1); e_matr(u);
  30.     fprintf(fw, "T(A)^(1/2) = %.8f\n\n", c);
  31.     while (c >= eps)
  32.     {
  33.         ++k;
  34.         float max = max_nondiagonal_elem(l, m, a1);
  35.         fprintf(fw, "max = %.3f\n\n", max);
  36.         float alpha = atan(2 * a1[l][m] / (a1[l][l] - a1[m][m])) / 2;
  37.         fprintf(fw, "alpha = %.8f\n\n", alpha);
  38.         for (int i = 0; i < n; i++)
  39.             for (int j = 0; j < n; j++)
  40.                 if (i == j) h[i][j] = ht[i][j] = 1.0;
  41.                 else h[i][j] = ht[i][j] = 0.0;
  42.         h[l][l] = ht[l][l] = h[m][m] = ht[m][m] = cos(alpha);
  43.         h[m][l] = ht[l][m] = sin(alpha); h[l][m] = ht[m][l] = -sin(alpha);
  44.         fprintf(fw, "Rotation matrix\n");
  45.         printmatr(h);
  46.         mult_matr(ht, a1, t);
  47.         mult_matr(t, h, a1);
  48.         mult_matr(u, h, t);
  49.         copy_matr(t, u);
  50.         fprintf(fw, "Count of iterations = %d\n\n", k);
  51.         fprintf(fw, "A matrix\n");
  52.         printmatr(a1);
  53.         c = rot_condition(a1);
  54.         fprintf(fw, "T(A)^(1/2) = %.8f\n\n", c);
  55.     }
  56.     float y[n];
  57.     for (int j = 0; j < n; j++) {
  58.         fprintf(fw, "l%d = %.8f\n\n", j + 1, a1[j][j]);
  59.         fprintf(fw, "vector x%d : ", j + 1);
  60.         for (int i = 0; i < n; i++)
  61.         {fprintf(fw, "%.8f ", u[i][j]); y[i] = u[i][j];}
  62.         fprintf(fw,"\n\n");
  63.         fprintf(fw, "norma(gamma) = %.8f\n\n", desp_max(a1[j][j], y, a));  
  64.     }
  65. }
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement