diff --git a/Jacobi_iteration.m b/Jacobi_iteration.m
new file mode 100644
index 0000000000000000000000000000000000000000..fd175d20693fd1ebd9482cd2ec7ff15d8cf8528b
--- /dev/null
+++ b/Jacobi_iteration.m
@@ -0,0 +1,77 @@
+% Jacobi iteration
+
+% INITIALIZE MATLAB
+
+close all;
+clc;
+clear all;
+
+%%%%%%%%%%%%%%%%%%%%%%
+%% DashBoard
+%%%%%%%%%%%%%%%%%%%%%%
+
+%Defrine [a] and [b]
+
+A = [ 8 2 3 1 ; ...
+      0 6 4 0 ; ...
+      2 3 9 3 ; ...
+      1 2 3 7 ];
+
+B = [25;24;47;42];
+
+%Tolerance
+
+ tol = 1e-6;
+
+ 
+
+%%%%%%%%%%%%%%%%%%%%%%
+%% Solve Ax-B using Jacobi iteration
+%%%%%%%%%%%%%%%%%%%%%%
+
+% First we should check for square matrix
+
+% Determine Matrix size
+[M,N] = size(A);
+if M~=N
+    error('Not a square matrix');
+end
+
+%Check that diagonally dominant
+
+for m=1 : M
+    row = abs(A(m,:));
+    d = sum(row)-row(m);
+    if row(m)<=d
+        error('Not diagonally dominant');
+    end
+end
+
+%Initial guess
+
+x=ones(M,1);
+
+% Calculate diagonal matrix D
+D=diag(diag(A));
+
+%Main Loop
+iter = 0;
+err=inf;
+while err > tol
+
+    %update [x]
+    dx = D\ (B-A*x);
+    x=x+dx;
+    iter = iter+1;
+
+    %error
+    err= max(abs(dx./x));
+
+end
+
+%Display answer
+disp(['Number of iterations:' num2str(iter)]);
+x
+
+
+