function [w, K, F, lambda, yA, MH] = EliptikKDD(a,b,c,d,N,M)

% u_xx + u_yy = f(x,y) , u(a,y)=fxa(x,y), u(b,y)=fxb(x,y), u(x,c)=fyc(x,y),
% u(x,d)=fyd(x,y)

f = @(x,y) -(cos(x+y) + cos(x-y));
fxa = @(x,y) cos(y);
fxb = @(x,y) -cos(y);
fyc = @(x,y) cos(x);
fyd = @(x,y) 0;

h = (b-a)/N;
k = (d-c)/M;

x = a:h:b;
y = c:k:d;

% Sınır koşulları

for j = 1:M+1
    w(1,j) = fxa(x(1),y(j));
    w(N+1,j) = fxb(x(N+1),y(j));
end

for i = 1:N+1
    w(i,1) = fyc(x(i),y(1));
    w(i,M+1) = fyd(x(i),y(M+1));
end

lambda = (h^2)/(k^2);
mu = 2*(1+lambda);


% Katsayılar matrisi

A = zeros(M-1);
B = zeros(M-1);

for i = 1:M-1
    A(i,i) = mu;
    B(i,i) = -1;
end

for i = 1:M-2
    A(i,i+1) = -lambda;
    A(i+1,i) = -lambda;
end

% K matrisi

K = zeros((N-1)*(M-1));

for i = 1:N-1
    
    K(1+(i-1)*(M-1):i*(M-1), 1+(i-1)*(M-1):i*(M-1)) = A;
    
end

for i = 1:N-2
    
    K(1+(i-1)*(M-1):i*(M-1), 1+i*(M-1):(i+1)*(M-1)) = B;
    K(1+i*(M-1):(i+1)*(M-1), 1+(i-1)*(M-1):i*(M-1)) = B;
    
end

% Sağ taraf vektörü F


F = zeros((N-1)*(M-1),1);

F(1,1) = -(h^2)*f(x(2),y(2)) + fxa(a,y(2)) + lambda*fyc(x(2),c); 

for i = 2:M-2
    
    F(i,1) = -(h^2)*f(x(2),y(i+1)) + fxa(a,y(i+1));
    
end

F(M-1,1) = -(h^2)*f(x(2),y(M)) + fxa(a,y(M)) + lambda*fyd(x(2),d);


for i = 3:N-1
    
    F((M-1)*(i-2)+1,1) = -(h^2)*f(x(i),y(2)) + lambda*fyc(x(i),c);
    F((M-1)*(i-1),1) = -(h^2)*f(x(i),y(M)) + lambda*fyd(x(i),d);
    
    for j = 3:M-1
        
        F((M-1)*(i-2) + j-1,1) = -(h^2)*f(x(i),y(j));
        
    end
end


F((N-2)*(M-1)+1,1) = -(h^2)*f(x(N),y(2)) + fxb(b,y(2)) + lambda*fyc(x(N),c); 

for i = 2:M-2
    
    F((N-2)*(M-1)+i,1) = -(h^2)*f(x(N),y(i+1)) + fxb(b,y(i+1));
    
end

F((N-1)*(M-1),1) = -(h^2)*f(x(N),y(M)) + fxb(b,y(M)) + lambda*fyd(x(N),d);

% Bilinmeyenler vektörü

ww = K\F;
    

for i = 2:N
    for j = 2:M
        
        w(i,j) = ww((j-1) + (M-1)*(i-2),1);
        
    end
end

% Analitik çözüm

[Y,X] = meshgrid(y,x);

yA = cos(X).*cos(Y);
    
MH = abs(yA-w);