diff --git a/Final/io.cpp b/Final/io.cpp
index ae754b7..186b4a8 100644
--- a/Final/io.cpp
+++ b/Final/io.cpp
@@ -1,433 +1,523 @@
 #include <iostream>
 #include <cassert>
 #include <cmath>
 #include <fstream>
 
 #include "io.hpp"
 #include "LinSys.hpp"
 
 using namespace std;
 
 Lin_sys_struct::Lin_sys_struct(int msize) {
     A.SetSize(msize,msize);
     b.SetSize(msize);
     size=msize;
 }
 
 Lin_sys_struct::~Lin_sys_struct() {}
 
 void out_2_file(Matrix A,string filename,int n_sig){
     ofstream write_output(filename);
     assert((write_output.is_open())&& "Cant write because file is open. Choose another filename or close file.");
     write_output.precision(n_sig);
     write_output<<A.NumberOfRows()<<endl;
     for (int i=0;i<A.NumberOfRows();i++){
         for (int j=0;j<A.NumberOfColumns();j++){
             write_output<<A(i+1,j+1)<<"  ";
         }
         write_output<<endl;
     }
     cout<<"Successfully wrote file "<<filename<<endl;
     write_output.close();
 }
 
 void out_2_file(Vector b,string filename,int n_sig){
     ofstream write_output(filename);
     assert((write_output.is_open())&& "Cant write because file is open. Choose another filename or close file.");
     write_output.precision(n_sig);
     write_output<<b.GetSize()<<endl;
     for (int i=0;i<b.GetSize();i++){
             write_output<<b[i]<<endl;
         }
     cout<<"Successfully wrote file "<<filename<<endl;
     write_output.close();
 }
 
 Matrix read_from_file(string filename,int flag){
     ifstream read_file(filename);
     assert(read_file.is_open());
 
     int size;
     read_file>>size;
 
     Matrix A(size);
     for (int i=0;i<A.NumberOfRows();i++){
         for (int j=0;j<A.NumberOfColumns();j++){
             read_file>>A(i+1,j+1);
         }
     }
     cout<<"Successfully read file "<<filename<<endl;
     read_file.close();
 return A;
 }
 
 Vector read_from_file(string filename){
     ifstream read_file(filename);
     assert(read_file.is_open());
 
     int size;
     read_file>>size;
 
     Vector b(size);
     for (int i=0;i<b.GetSize();i++){
             read_file>>b(i+1);
         }
     cout<<"Successfully read file "<<filename<<endl;
     read_file.close();
     return b;
 }
 
 Lin_sys_struct Initialize(){
     // Initialization
     cout<<"\nThis program allows solving linear systems of the type Ax=b.\nStart by giving the size of the linear system : ";
     int size;
     cin>>size;
     Matrix A=Matrix(size);
     Vector b=Vector(size);
     Lin_sys_struct instance(size);
 
 
     // Matrix A
     cout<<"\nLets fill matrix A."<<endl;
     for (int i=0;i<size;i++){
         for (int j=0;j<size;j++){
             cout<<"Give value of A("<<i+1<<","<<j+1<<") : ";
             cin>>A(i+1,j+1);
         }
     }
     TerminalPrint(A,"Matrix A :");
 
     //Vector b
     cout<<"\nLets fill vector B."<<endl;
     for (int i=0;i<size;i++){
         cout<<"Give value of b("<<i+1<<") : ";
         cin>>b(i+1);
     }
     TerminalPrint(b,"Vector b :");
     instance.A=A;
     instance.b=b;
     return instance;
 }
 
 string method(){
     // Description of methods available
     cout<<"Description of available methods:\n----------------------\n------------------------\n";
     cout<<"Exact : (suitable for relatively small systems)\n---------------------------------------------\n";
     cout<<"Cholesky : Fast but for semipositive definite systems only."<<endl;
     cout<<"LU : Works for all definite systems."<<endl<<endl;
     cout<<"Approximative : (suitable for large systems)\n----------------------------------------\n";
     cout<<"Conjugate_gradient : very fast"<<endl;
-    cout<<"Gauss_Seidel : Works for all definite systems."<<endl;
-    cout<<"Jacobi : Fast but for semipositive definite systems only."<<endl;
+    cout<<"Gauss_Seidel : Works best for diagonally dominant systems."<<endl;
+    cout<<"Jacobi : Fast but for semipositive definite systems and diagonally dominant systems."<<endl;
     cout<<"Richardson : Works for all definite systems."<<endl<<endl;
 
     // Choosing of method
     string my_method;
     cout<<"\nChoose a method by typing a string corresponding to the name of the method, or Exit to abort program\n";
     cin>>my_method;
     return my_method;
 }
 
+//Evaluate if the Matrix A is diagonally dominant
+string diagonallyDominant(Matrix A){
+    string m_diagonallyDominant = "True";
+    int size = A.NumberOfColumns();
+    double sum = 0;
+    for (int i = 0; i < size; ++i) {
+        for (int j = 0; j < size; ++j) {
+            sum += A(j+1, i+1);
+        }
+        sum -= A(i+1, i+1);
+        if (A(i+1, i+1) < sum) {
+            m_diagonallyDominant = "False";
+        }
+    }
+    return m_diagonallyDominant;
+}
+
+
 void run_algo(Matrix A,Vector b,string my_method){
 
     if (my_method=="LU") {
         LU mySys=LU(A,b);
         string command;
 
         cout<<"\nYou chose the LU method. Type a command from the ones available below:\n";
         cout<<"U_matrix: outputs the upper triangular decomposition of A.\n";
         cout<<"L_matrix: outputs the lower triangular decomposition of A.\n";
         cout<<"P_matrix: outputs the permutation matrix of A.\n";
         cout<<"Solve: outputs the solution of the system Ax=b\n\n";
         cin >> command;
 
         while(command!="Exit") {
             if (command == "U_Matrix") { TerminalPrint(mySys.get_U(), "U decomposition of A."); }
             else if (command == "L_Matrix") { TerminalPrint(mySys.get_L(), "L decomposition of A."); }
             else if (command == "P_Matrix") { TerminalPrint(mySys.get_P(), "Permutation matrix of A."); }
             else if (command == "Solve" )
             {
                 Vector x=mySys.Solve(); int counter=0;
                 for (int i=0;i<b.GetSize();i++){
                     if (b[i]==0) counter++;
                 }
                 if (counter!=b.GetSize())  TerminalPrint(x, "Solution of system Ax=b.");
             }
             else {cout<<"Command not supported.";}
             cout<<"Type Exit to quit the algorithm or type a new command :\n\n";
             cin >> command;
         }
 
     } else if (my_method=="Cholesky") {
         if (A.Check_positive_definite() == false) {
             cout << "Matrix is not definite semi-positive. Try LU method instead" << endl;
         } else {
             Cholesky mySys = Cholesky(A, b);
             string command;
 
             cout << "\nYou chose the Cholesky method. Type a command from the ones available below:\n";
             cout << "L_matrix: outputs the lower triangular decomposition of A.\n";
             cout << "U_matrix: outputs the upper triangular decomposition of A, which is the transpose of L.\n";
             cout << "P_matrix: outputs the permutation matrix of A.\n";
             cout << "Solve: outputs the solution of the system Ax=b\n\n";
             cin >> command;
 
             while (command != "Exit") {
                 if (command == "U_Matrix") { TerminalPrint(mySys.get_U(), "U decomposition of A."); }
                 else if (command == "L_Matrix") { TerminalPrint(mySys.get_L(), "L decomposition of A."); }
                 else if (command == "P_Matrix") { TerminalPrint(mySys.get_P(), "Permutation matrix of A."); }
                 else if (command == "Solve") { TerminalPrint(mySys.Solve(), "Solution of system Ax=b."); }
                 else { cout << "Command not supported."; }
                 cout << "Type Exit to quit the algorithm or type a new command :\n\n";
                 cin >> command;
             }
         }
     }
      else if (my_method=="Conjugate_gradient"){
         Conjugate_Gradient mySys=Conjugate_Gradient(A,b);
         string command;
 
         cout<<"\nYou chose the Conjugate gradient method. Type a command from the ones available below:\n";
         cout<<"Set_tol: Sets the convergence criteria to the value given (1e-12 by default). \n";
         cout<<"Set_x0 : Sets the initial solution. \n";
         cout<<"Set_max_iter : Sets the max number of iterations possible to the value given (5000 by default). \n";
 
         cout<<"Solve: outputs the solution of the system Ax=b. Run this first to unlock following methods:\n\n";
         cout<<"iter: outputs the number of iterations needed for convergence. \n";
         cout<<"Convergence: outputs a vector with each index corresponding to the error at the iteration. \n";
 
         cin >> command;
 
         while(command!="Exit") {
             if (command == "Solve") { TerminalPrint(mySys.Solve(), "Solution of system Ax=b."); }
             else if (command == "iter") { TerminalPrint(mySys.get_n_iter(), "Number of iterations needed for convergence."); }
             else if (command == "Convergence") {
                 for (int i=0;i<mySys.get_n_iter();i++){
                     cout<<mySys.get_convergence()[i]<<endl;
                 }
                 cout<<endl;
             }
             else if (command == "Set_tol") {
                 double tol;
                 cout<<"Type tolerance value you wish to set:\n\n";
                 cin >> tol;
                 mySys.set_tol(tol); }
             else if (command == "Set_x0") {
                 Vector x0(b.GetSize());
                 cout<<"Type x0(1),x0(2),...,x0(n) values you wish to start algorithm with:\n\n";
                 for (int i=0;i<b.GetSize();i++){
                     cout<<"Give value of x0("<<i+1<<") : ";
                     cin>>x0(i+1);
                 }
                 TerminalPrint(x0,"Initial solution x0 :"); }
             else if (command == "Set_max_iter") {
                 int max_iter;
                 cout<<"Type the maximum number of iterations criteria:\n\n";
                 cin >> max_iter;
                 mySys.set_iter_max(max_iter); }
             else {cout<<"Command not supported.";}
             cout<<"Type Exit to quit the algorithm or type a new command :\n\n";
             cin >> command;
         }
     } else if (my_method=="Gauss_Seidel"){
         Gauss_Seidel mySys=Gauss_Seidel(A,b);
         string command;
 
-        cout<<"\nYou chose the Gauss_Seidel method. Type a command from the ones available below:\n";
-        cout<<"Set_tol: Sets the convergence criteria to the value given (1e-12 by default). \n";
-        cout<<"Set_x0 : Sets the initial solution. \n";
-        cout<<"Set_max_iter : Sets the max number of iterations possible to the value given (5000 by default). \n";
+        string isDD;
+        isDD = diagonallyDominant(A);
 
-        cout<<"Solve: outputs the solution of the system Ax=b. Run this first to unlock following methods:\n\n";
-        cout<<"iter: outputs the number of iterations needed for convergence. \n";
-        cout<<"Convergence: outputs a vector with each index corresponding to the error at the iteration. \n";
+        cout << "\nYou chose the Jacobi method.";
+
+        string incorrect = "True";
+        while (incorrect == "True") {
+            incorrect = "False";
+            if (isDD == "False") {
+                cout << "The matrix isn't diagonally dominant therefore the algorithm might not converge.\n";
+                cout << "You can either continue or change algorithm. Type a command from the ones available below:\n";
+                cout << "Continue: if you want to procceed \n";
+                cout << "Change: if you want to change method \n";
+
+                cin >> command;
+                if (command == "Continue") {
+                    cout << "Type a command from the ones available below:\n";
+                    cout<<"Set_tol: Sets the convergence criteria to the value given (1e-12 by default). \n";
+                    cout<<"Set_x0 : Sets the initial solution. \n";
+                    cout<<"Set_max_iter : Sets the max number of iterations possible to the value given (5000 by default). \n";
+                    cout<<"Solve: outputs the solution of the system Ax=b. Run this first to unlock following methods:\n\n";
+                    cout<<"iter: outputs the number of iterations needed for convergence. \n";
+                    cout<<"Convergence: outputs a vector with each index corresponding to the error at the iteration. \n";
+
+                    cin >> command;
+
+                } else if (command == "Change") {
+                    command = "Exit";
+                } else {
+                    cout << "incorrect command\n";
+                    incorrect = "True";
+
+                }
+            }else {
+                cout << " Type a command from the ones available below:\n";
+                cout<<"Set_tol: Sets the convergence criteria to the value given (1e-12 by default). \n";
+                cout<<"Set_x0 : Sets the initial solution. \n";
+                cout<<"Set_max_iter : Sets the max number of iterations possible to the value given (5000 by default). \n";
+                cout<<"Solve: outputs the solution of the system Ax=b. Run this first to unlock following methods:\n\n";
+                cout<<"iter: outputs the number of iterations needed for convergence. \n";
+                cout<<"Convergence: outputs a vector with each index corresponding to the error at the iteration. \n";
+
+                cin >> command;
+            }
+        }
 
-        cin >> command;
 
         while(command!="Exit") {
             if (command == "Solve") { TerminalPrint(mySys.Solve(), "Solution of system Ax=b."); }
             else if (command == "iter") { TerminalPrint(mySys.get_n_iter(), "Number of iterations needed for convergence."); }
             else if (command == "Convergence") {
                 for (int i=0;i<mySys.get_n_iter();i++){
                     cout<<mySys.get_convergence()[i]<<endl;
                 }
                 cout<<endl;
             }
             else if (command == "Set_tol") {
                 double tol;
                 cout<<"Type tolerance value you wish to set:\n\n";
                 cin >> tol;
                 mySys.set_tol(tol); }
             else if (command == "Set_x0") {
                 Vector x0(b.GetSize());
                 cout<<"Type x0(1),x0(2),...,x0(n) values you wish to start algorithm with:\n\n";
                 for (int i=0;i<b.GetSize();i++){
                     cout<<"Give value of x0("<<i+1<<") : ";
                     cin>>x0(i+1);
                 }
                 TerminalPrint(x0,"Initial solution x0 :"); }
             else if (command == "Set_max_iter") {
                 int max_iter;
                 cout<<"Type the maximum number of iterations criteria:\n\n";
                 cin >> max_iter;
                 mySys.set_iter_max(max_iter); }
             else {cout<<"Command not supported.";}
             cout<<"Type Exit to quit the algorithm or type a new command :\n\n";
             cin >> command;
         }
     } else if (my_method=="Jacobi"){
         Jacobi mySys=Jacobi(A,b);
         string command;
 
-        cout<<"\nYou chose the Jacobi method. Type a command from the ones available below:\n";
-        cout<<"Set_tol: Sets the convergence criteria to the value given (1e-12 by default). \n";
-        cout<<"Set_x0 : Sets the initial solution. \n";
-        cout<<"Set_max_iter : Sets the max number of iterations possible to the value given (5000 by default). \n";
+        string isDD;
+        isDD = diagonallyDominant(A);
 
-        cout<<"Solve: outputs the solution of the system Ax=b. Run this first to unlock following methods:\n\n";
-        cout<<"iter: outputs the number of iterations needed for convergence. \n";
-        cout<<"Convergence: outputs a vector with each index corresponding to the error at the iteration. \n";
+        cout << "\nYou chose the Jacobi method.";
+
+        string incorrect = "True";
+        while (incorrect == "True") {
+            incorrect = "False";
+            if (isDD == "False") {
+                cout << "The matrix isn't diagonally dominant therefore the algorithm might not converge.\n";
+                cout << "You can either continue or change algorithm. Type a command from the ones available below:\n";
+                cout << "Continue: if you want to procceed \n";
+                cout << "Change: if you want to change method \n";
+
+                cin >> command;
+                if (command == "Continue") {
+                    cout << "Type a command from the ones available below:\n";
+                    cout<<"Set_tol: Sets the convergence criteria to the value given (1e-12 by default). \n";
+                    cout<<"Set_x0 : Sets the initial solution. \n";
+                    cout<<"Set_max_iter : Sets the max number of iterations possible to the value given (5000 by default). \n";
+                    cout<<"Solve: outputs the solution of the system Ax=b. Run this first to unlock following methods:\n\n";
+                    cout<<"iter: outputs the number of iterations needed for convergence. \n";
+                    cout<<"Convergence: outputs a vector with each index corresponding to the error at the iteration. \n";
+
+                    cin >> command;
+
+                } else if (command == "Change") {
+                    command = "Exit";
+                } else {
+                    cout << "incorrect command\n";
+                    incorrect = "True";
+
+                }
+            }else {
+                cout << " Type a command from the ones available below:\n";
+                cout<<"Set_tol: Sets the convergence criteria to the value given (1e-12 by default). \n";
+                cout<<"Set_x0 : Sets the initial solution. \n";
+                cout<<"Set_max_iter : Sets the max number of iterations possible to the value given (5000 by default). \n";
+                cout<<"Solve: outputs the solution of the system Ax=b. Run this first to unlock following methods:\n\n";
+                cout<<"iter: outputs the number of iterations needed for convergence. \n";
+                cout<<"Convergence: outputs a vector with each index corresponding to the error at the iteration. \n";
+
+                cin >> command;
+            }
+        }
 
-        cin >> command;
 
         while(command!="Exit") {
             if (command == "Solve") { TerminalPrint(mySys.Solve(), "Solution of system Ax=b."); }
             else if (command == "iter") { TerminalPrint(mySys.get_n_iter(), "Number of iterations needed for convergence."); }
             else if (command == "Convergence") {
                 for (int i=0;i<mySys.get_n_iter();i++){
                     cout<<mySys.get_convergence()[i]<<endl;
                 }
                 cout<<endl;
             }
             else if (command == "Set_tol") {
                 double tol;
                 cout<<"Type tolerance value you wish to set:\n\n";
                 cin >> tol;
                 mySys.set_tol(tol); }
             else if (command == "Set_x0") {
                 Vector x0(b.GetSize());
                 cout<<"Type x0(1),x0(2),...,x0(n) values you wish to start algorithm with:\n\n";
                 for (int i=0;i<b.GetSize();i++){
                     cout<<"Give value of x0("<<i+1<<") : ";
                     cin>>x0(i+1);
                 }
                 TerminalPrint(x0,"Initial solution x0 :"); }
             else if (command == "Set_max_iter") {
                 int max_iter;
                 cout<<"Type the maximum number of iterations criteria:\n\n";
                 cin >> max_iter;
                 mySys.set_iter_max(max_iter); }
             else {cout<<"Command not supported.";}
             cout<<"Type Exit to quit the algorithm or type a new command :\n\n";
             cin >> command;
         }
     } else if (my_method=="Richardson"){
         Richardson mySys=Richardson(A,b);
         string command;
 
         cout<<"\nYou chose the Richardson method. Type a command from the ones available below:\n";
         cout<<"Set_tol: Sets the convergence criteria to the value given (1e-12 by default). \n";
         cout<<"Set_x0 : Sets the initial solution. \n";
         cout<<"Set_max_iter : Sets the max number of iterations possible to the value given (5000 by default). \n";
 
         cout<<"Solve: outputs the solution of the system Ax=b. Run this first to unlock following methods:\n\n";
         cout<<"iter: outputs the number of iterations needed for convergence. \n";
         cout<<"Convergence: outputs a vector with each index corresponding to the error at the iteration. \n";
 
         cin >> command;
 
         while(command!="Exit") {
             if (command == "Solve" )
             {
                 Vector x=mySys.Solve(); int counter=0;
                 for (int i=0;i<b.GetSize();i++){
                     if (x[i]==0) counter++;
                 }
                 if (counter!=b.GetSize())  TerminalPrint(x, "Solution of system Ax=b.");
             }
             else if (command == "iter") { TerminalPrint(mySys.get_n_iter(), "Number of iterations needed for convergence."); }
             else if (command == "Convergence") {
                 for (int i=0;i<mySys.get_n_iter();i++){
                     cout<<mySys.get_convergence()[i]<<endl;
                 }
                 cout<<endl;
             }
             else if (command == "Set_tol") {
                 double tol;
                 cout<<"Type tolerance value you wish to set:\n\n";
                 cin >> tol;
                 mySys.set_tol(tol); }
             else if (command == "Set_x0") {
                 Vector x0(b.GetSize());
                 cout<<"Type x0(1),x0(2),...,x0(n) values you wish to start algorithm with:\n\n";
                 for (int i=0;i<b.GetSize();i++){
                     cout<<"Give value of x0("<<i+1<<") : ";
                     cin>>x0(i+1);
                 }
                 TerminalPrint(x0,"Initial solution x0 :"); }
             else if (command == "Set_max_iter") {
                 int max_iter;
                 cout<<"Type the maximum number of iterations criteria:\n\n";
                 cin >> max_iter;
                 mySys.set_iter_max(max_iter); }
             else {cout<<"Command not supported.";}
             cout<<"Type Exit to quit this method and try a new one or type a new command :\n\n";
             cin >> command;
         }
     }else if (my_method=="Exit"){
         cout<<"Leaving program .. BYYE\n";
     } else {
         cout<<"\nMethod not recognized. Select another one.\n\n";
     }
      }
 
 void run_algo_2(Matrix A,Vector b,std::string my_method){
     if (my_method=="LU") {
         LU mySys=LU(A,b);
         out_2_file(mySys.get_U(),"U.dat",3);
         out_2_file(mySys.get_L(),"L.dat",3);
         out_2_file(mySys.get_P(),"P.dat",1);
         out_2_file(mySys.get_b_star(),"b*.dat",3);
         out_2_file(mySys.Solve(),"x.dat",3);
 
     } else if (my_method=="Cholesky"){
         if (A.Check_positive_definite() == false) {
             cout << "Matrix is not definite semi-positive. Try LU method instead" << endl;
         } else {
             Cholesky mySys = Cholesky(A, b);
             out_2_file(mySys.get_U(),"U.dat",3);
             out_2_file(mySys.get_L(),"L.dat",3);
             out_2_file(mySys.get_P(),"P.dat",1);
             out_2_file(mySys.get_b_star(),"b*.dat",3);
             out_2_file(mySys.Solve(),"x.dat",3);
         }
     } else if (my_method=="Conjugate_Gradient"){
         cout<<"System solved with default tol (1e-6), iter_max(5000) and initial solution set to 0."<<endl;
         Conjugate_Gradient mySys=Conjugate_Gradient(A,b);
         out_2_file(mySys.Solve(),"x.dat",3);
         TerminalPrint(mySys.get_n_iter(),"Number of iteration required for convergence.");
         out_2_file(mySys.get_convergence(),"Convergence.dat",3);
     } else if (my_method=="Gauss_Siedel"){
         cout<<"System solved with default tol (1e-6), iter_max(5000) and initial solution set to 0."<<endl;
         Gauss_Seidel mySys=Gauss_Seidel(A,b);
         out_2_file(mySys.Solve(),"x.dat",3);
         TerminalPrint(mySys.get_n_iter(),"Number of iteration required for convergence.");
         out_2_file(mySys.get_convergence(),"Convergence.dat",3);
     } else if (my_method=="Jacobi"){
         cout<<"System solved with default tol (1e-6), iter_max(5000) and initial solution set to 0."<<endl;
         Jacobi mySys=Jacobi(A,b);
         out_2_file(mySys.Solve(),"x.dat",3);
         TerminalPrint(mySys.get_n_iter(),"Number of iteration required for convergence.");
         out_2_file(mySys.get_convergence(),"Convergence.dat",3);
     } else if (my_method=="Richardson"){
         cout<<"System solved with default tol (1e-6), iter_max(5000) and initial solution set to 0."<<endl;
         Richardson mySys=Richardson(A,b);
         out_2_file(mySys.Solve(),"x.dat",3);
         TerminalPrint(mySys.get_n_iter(),"Number of iteration required for convergence.");
         out_2_file(mySys.get_convergence(),"Convergence.dat",3);
     }else if (my_method=="Exit"){
         cout<<"Leaving program .. BYYE\n";
     } else {
         cout<<"\nMethod not recognized. Rerun executable with arguments : string A_filename string b_filename string method.\n\n";
     }
 }