Advertisement
Guest User

Untitled

a guest
Mar 26th, 2019
82
0
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
text 5.29 KB | None | 0 0
  1. #include <iostream>
  2. #include <fstream>
  3. #include <string>
  4. #include <math.h>
  5. #include <cmath>
  6. #include <vector>
  7.  
  8. #define n 4
  9.  
  10. void displayVect1D(double wektor[], int size) {
  11. for (int i = 0; i < size; i++) {
  12. std::cout << wektor[i] << " ";
  13. }
  14. std::cout << '\n';
  15. }
  16.  
  17. void display2D(double matrix[][n]){
  18. for (int i=0;i<4;i++){
  19. for (int j=0;j<4;j++){
  20. std::cout << matrix[i][j] << " ";
  21. }
  22. std::cout << '\n';
  23. }
  24. }
  25.  
  26. double estimator(double *x, double *x1){
  27. double x3[n];
  28. for (int i=0; i<4;i++){
  29. x3[i] = fabs(x[i] - x1[i]);
  30. }
  31.  
  32. double max = x3[0];
  33. for (int i=1;i<4;i++){
  34. if (x3[i] > max){
  35. max = x3[i];
  36. }
  37. }
  38. return max;
  39. }
  40.  
  41. double reziduum(double x[], double b[], const double a[][n]){
  42.  
  43. double r[n] = {0,0,0,0};
  44. double max;
  45. for (int i = 0; i < n; i++) {
  46. for (int j = 0; j < n; j++) {
  47. r[i] += a[i][j] * x[j];
  48. }
  49. r[i] = fabs(r[i] - b[i]);
  50. }
  51. max = r[0];
  52. for (int i = 1; i < n; i++) {
  53. if (r[i] > max) {
  54. max = r[i];
  55. }
  56. }
  57. return max;
  58. }
  59.  
  60. double jacobi(double b[], const double a[][n], int iterMax){
  61.  
  62. double x[4] = {1.,1.,1.,1.};
  63. double estimatorVal = 1.;
  64. double residuum = 1.;
  65.  
  66. double e = 1.0e-04;
  67.  
  68. double D[4];
  69. for (int i=0;i<4;i++){
  70. for (int j=0; j<4; j++){
  71. if (i==j){
  72. D[i] = 1./a[i][j];
  73. }
  74. }
  75. }
  76.  
  77. double M[4][4];
  78. for (int i=0;i<4;i++){
  79. for (int j=0; j<4; j++){
  80. if (i!=j){
  81. M[i][j] = -D[i] * a[i][j];
  82. }
  83. else{
  84. M[i][j] = 0;
  85. }
  86. }
  87. }
  88.  
  89. double c[4];
  90.  
  91. for (int i=0;i<4;i++){
  92. c[i] = D[i] * b[i];
  93. }
  94.  
  95. double x2[4];
  96. int k=0;
  97. while ((k<iterMax) && (fabs(estimatorVal)>=e) && (fabs(residuum)>e)){
  98.  
  99. for (int i=0;i<4;i++){
  100. x2[i] = c[i];
  101. for (int j=0;j<4;j++){
  102. x2[i] += M[i][j] * x[j];
  103. }
  104. }
  105.  
  106. estimatorVal = estimator(x,x2);
  107.  
  108. std::cout << "\niteracja: " << k+1 << std::endl;
  109. for (int i=0;i<4;i++){
  110. x[i] = x2[i];
  111. std::cout << x2[i] << " ";
  112. }
  113. std::cout << "\nestymator: " << estimatorVal << " reziduum: " <<
  114. reziduum(x,b,a);
  115. std::cout << "\n";
  116.  
  117. k++;
  118. }
  119. return 0;
  120.  
  121. }
  122.  
  123. double Jacobi(double b[], const double a[][n], int iterMax){
  124.  
  125. std::cout << "\n jacobi" << std::endl;
  126.  
  127. double x[n] = {1.,1.,1.,1.};
  128. double estimatorVal = 1.;
  129. double reziduumVal = 1.;
  130. double x2[n];
  131.  
  132. double e = 1.0e-04;
  133.  
  134. int k=0;
  135.  
  136. while ((k<iterMax) && (fabs(estimatorVal)>=e) && (fabs(reziduumVal)>e)){
  137. for (int i=0; i<n;i++){
  138. double fi = 0;
  139. for (int j=0; j<n; j++){
  140. if (i!=j){
  141. fi += a[i][j]*x[j];
  142. }
  143. }
  144. x2[i] = (b[i] - fi) / a[i][i];
  145. }
  146.  
  147. estimatorVal = estimator(x,x2);
  148. std::cout << "\niteracja: " << k+1 << std::endl;
  149. for(int i=0;i<n;i++){
  150. x[i] = x2[i];
  151. }
  152.  
  153. for (int i=0;i<4;i++){
  154. std::cout << x[i] << " ";
  155. }
  156.  
  157. std::cout << "\nestymator: " << estimatorVal << " reziduum: " << reziduum(x,b,a) << std::endl;
  158. k++;
  159. }
  160. }
  161.  
  162. double gauss(double b[], const double a[][4], int iterMax){
  163.  
  164.  
  165. std::cout << "\n gauss-seidel" << std::endl;
  166.  
  167. double x[n] = {1.,1.,1.,1.};
  168. double estimatorVal = 1.;
  169. double reziduumVal = 1.;
  170. double x2[n];
  171.  
  172. double e = 1.0e-04;
  173.  
  174. int k=0;
  175. while ((k<iterMax) && (fabs(estimatorVal)>=e) && (fabs(reziduumVal)>e)){
  176.  
  177. for (int i=0;i<n;i++){
  178. double fi = 0;
  179. for (int j=0;j<i;j++){
  180. fi += a[i][j] * x[j];
  181. }
  182. for (int j=i+1;j<n;j++){
  183. fi += a[i][j] * x[j];
  184. }
  185. x2[i] = x[i];
  186. x[i] = (b[i] - fi) / a[i][i];
  187. }
  188.  
  189. estimatorVal = estimator(x2,x);
  190.  
  191. std::cout << "\niteracja: " << k+1 << std::endl;
  192.  
  193. for (int i=0;i<4;i++){
  194. std::cout << x[i] << " ";
  195. }
  196.  
  197. reziduumVal = reziduum(x,b,a);
  198. std::cout << "\nestymator: " << estimatorVal
  199.  
  200. << " reziduum: " <<
  201. reziduumVal;
  202. std::cout << "\n";
  203.  
  204. k++;
  205. }
  206.  
  207. }
  208.  
  209. double sor(double b[], const double a[][4], int iterMax){
  210.  
  211.  
  212. std::cout << "\nsor" << std::endl;
  213.  
  214. double x[n] = {1.,1.,1.,1.};
  215. double estimatorVal = 1.;
  216. double reziduumVal = 1.;
  217. double x2[n];
  218. double temp[n];
  219.  
  220. double omega = 0.5;
  221. double e = 1.0e-06;
  222.  
  223. int k=0;
  224. while ((k<iterMax) && (fabs(estimatorVal)>=e) && (fabs(reziduumVal)>e)){
  225.  
  226. for (int i=0;i<n;i++){
  227. double sum = 0;
  228. for (int j=0;j<i;j++){
  229. sum += a[i][j] * x[j];
  230. }
  231. for (int j=i+1;j<n;j++){
  232. sum += a[i][j] * x[j];
  233. }
  234. temp[i] = x[i];
  235. x2[i] = (1.0 - omega) * x[i] + (omega / a[i][i]) * (b[i] - sum);
  236. x[i] = x2[i];
  237. }
  238.  
  239. estimatorVal = estimator(temp,x2);
  240.  
  241. std::cout << "\niteracja: " << k+1 << std::endl;
  242.  
  243. for (int i=0;i<4;i++){
  244. std::cout << x[i] << " ";
  245. }
  246.  
  247. reziduumVal = reziduum(x,b,a);
  248. std::cout << "\nestymator: " << estimatorVal << " reziduum: " << reziduumVal;
  249. std::cout << "\n";
  250.  
  251. k++;
  252. }
  253. }
  254.  
  255.  
  256. int main() {
  257.  
  258. const double a[4][4] = {{100., 1., -2., 3.},
  259. {4., 300., -5., 6.},
  260. {7., -8., 400., 9.},
  261. {-10., 11., -12., 200}};
  262.  
  263. double b[4] = {395.,603.,-415.,-606.};
  264.  
  265.  
  266. Jacobi(b,a,25);
  267. gauss(b,a,25);
  268. sor(b,a,25);
  269.  
  270. return 0;
  271. }
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement