Advertisement
Not a member of Pastebin yet?
Sign Up,
it unlocks many cool features!
- #include <iostream>
- #include <fstream>
- #include <string>
- #include <math.h>
- #include <cmath>
- #include <vector>
- #define n 4
- void displayVect1D(double wektor[], int size) {
- for (int i = 0; i < size; i++) {
- std::cout << wektor[i] << " ";
- }
- std::cout << '\n';
- }
- void display2D(double matrix[][n]){
- for (int i=0;i<4;i++){
- for (int j=0;j<4;j++){
- std::cout << matrix[i][j] << " ";
- }
- std::cout << '\n';
- }
- }
- double estimator(double *x, double *x1){
- double x3[n];
- for (int i=0; i<4;i++){
- x3[i] = fabs(x[i] - x1[i]);
- }
- double max = x3[0];
- for (int i=1;i<4;i++){
- if (x3[i] > max){
- max = x3[i];
- }
- }
- return max;
- }
- double reziduum(double x[], double b[], const double a[][n]){
- double r[n] = {0,0,0,0};
- double max;
- for (int i = 0; i < n; i++) {
- for (int j = 0; j < n; j++) {
- r[i] += a[i][j] * x[j];
- }
- r[i] = fabs(r[i] - b[i]);
- }
- max = r[0];
- for (int i = 1; i < n; i++) {
- if (r[i] > max) {
- max = r[i];
- }
- }
- return max;
- }
- double jacobi(double b[], const double a[][n], int iterMax){
- double x[4] = {1.,1.,1.,1.};
- double estimatorVal = 1.;
- double residuum = 1.;
- double e = 1.0e-04;
- double D[4];
- for (int i=0;i<4;i++){
- for (int j=0; j<4; j++){
- if (i==j){
- D[i] = 1./a[i][j];
- }
- }
- }
- double M[4][4];
- for (int i=0;i<4;i++){
- for (int j=0; j<4; j++){
- if (i!=j){
- M[i][j] = -D[i] * a[i][j];
- }
- else{
- M[i][j] = 0;
- }
- }
- }
- double c[4];
- for (int i=0;i<4;i++){
- c[i] = D[i] * b[i];
- }
- double x2[4];
- int k=0;
- while ((k<iterMax) && (fabs(estimatorVal)>=e) && (fabs(residuum)>e)){
- for (int i=0;i<4;i++){
- x2[i] = c[i];
- for (int j=0;j<4;j++){
- x2[i] += M[i][j] * x[j];
- }
- }
- estimatorVal = estimator(x,x2);
- std::cout << "\niteracja: " << k+1 << std::endl;
- for (int i=0;i<4;i++){
- x[i] = x2[i];
- std::cout << x2[i] << " ";
- }
- std::cout << "\nestymator: " << estimatorVal << " reziduum: " <<
- reziduum(x,b,a);
- std::cout << "\n";
- k++;
- }
- return 0;
- }
- double Jacobi(double b[], const double a[][n], int iterMax){
- std::cout << "\n jacobi" << std::endl;
- double x[n] = {1.,1.,1.,1.};
- double estimatorVal = 1.;
- double reziduumVal = 1.;
- double x2[n];
- double e = 1.0e-04;
- int k=0;
- while ((k<iterMax) && (fabs(estimatorVal)>=e) && (fabs(reziduumVal)>e)){
- for (int i=0; i<n;i++){
- double fi = 0;
- for (int j=0; j<n; j++){
- if (i!=j){
- fi += a[i][j]*x[j];
- }
- }
- x2[i] = (b[i] - fi) / a[i][i];
- }
- estimatorVal = estimator(x,x2);
- std::cout << "\niteracja: " << k+1 << std::endl;
- for(int i=0;i<n;i++){
- x[i] = x2[i];
- }
- for (int i=0;i<4;i++){
- std::cout << x[i] << " ";
- }
- std::cout << "\nestymator: " << estimatorVal << " reziduum: " << reziduum(x,b,a) << std::endl;
- k++;
- }
- }
- double gauss(double b[], const double a[][4], int iterMax){
- std::cout << "\n gauss-seidel" << std::endl;
- double x[n] = {1.,1.,1.,1.};
- double estimatorVal = 1.;
- double reziduumVal = 1.;
- double x2[n];
- double e = 1.0e-04;
- int k=0;
- while ((k<iterMax) && (fabs(estimatorVal)>=e) && (fabs(reziduumVal)>e)){
- for (int i=0;i<n;i++){
- double fi = 0;
- for (int j=0;j<i;j++){
- fi += a[i][j] * x[j];
- }
- for (int j=i+1;j<n;j++){
- fi += a[i][j] * x[j];
- }
- x2[i] = x[i];
- x[i] = (b[i] - fi) / a[i][i];
- }
- estimatorVal = estimator(x2,x);
- std::cout << "\niteracja: " << k+1 << std::endl;
- for (int i=0;i<4;i++){
- std::cout << x[i] << " ";
- }
- reziduumVal = reziduum(x,b,a);
- std::cout << "\nestymator: " << estimatorVal
- << " reziduum: " <<
- reziduumVal;
- std::cout << "\n";
- k++;
- }
- }
- double sor(double b[], const double a[][4], int iterMax){
- std::cout << "\nsor" << std::endl;
- double x[n] = {1.,1.,1.,1.};
- double estimatorVal = 1.;
- double reziduumVal = 1.;
- double x2[n];
- double temp[n];
- double omega = 0.5;
- double e = 1.0e-06;
- int k=0;
- while ((k<iterMax) && (fabs(estimatorVal)>=e) && (fabs(reziduumVal)>e)){
- for (int i=0;i<n;i++){
- double sum = 0;
- for (int j=0;j<i;j++){
- sum += a[i][j] * x[j];
- }
- for (int j=i+1;j<n;j++){
- sum += a[i][j] * x[j];
- }
- temp[i] = x[i];
- x2[i] = (1.0 - omega) * x[i] + (omega / a[i][i]) * (b[i] - sum);
- x[i] = x2[i];
- }
- estimatorVal = estimator(temp,x2);
- std::cout << "\niteracja: " << k+1 << std::endl;
- for (int i=0;i<4;i++){
- std::cout << x[i] << " ";
- }
- reziduumVal = reziduum(x,b,a);
- std::cout << "\nestymator: " << estimatorVal << " reziduum: " << reziduumVal;
- std::cout << "\n";
- k++;
- }
- }
- int main() {
- const double a[4][4] = {{100., 1., -2., 3.},
- {4., 300., -5., 6.},
- {7., -8., 400., 9.},
- {-10., 11., -12., 200}};
- double b[4] = {395.,603.,-415.,-606.};
- Jacobi(b,a,25);
- gauss(b,a,25);
- sor(b,a,25);
- return 0;
- }
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement