Ytukyg

Using a code like this, I am having a hard time applying a Robin Boundary condition for a instead of a dirichlet for the following problem:

IMG OF QUESTION FOR POISSON ROBIN BC

Nx = 10; % Number of sub-segments in x
Ny = 10; % Number of sub-segments in y
a = 0; % Location of boundary 'a' for 'x'
b = 1; % Location of boundary 'b' for 'x'
c = 0; % Location of boundary 'c' for 'y'
d = 1; % Location of boundary 'd' for 'y'
f = @(x,y) x + y; % Defining RH-side function 'f(x,y)'

function [gxy] = g(x,y,a,b,c,d) % Defining boundary condition 'g(x,y)'

gxy = 0; % Default value

if (x==a) % g(a,y), Left BC
gxy = 0;
end

if (x==b) % g(b,y), Right BC
gxy = 1;
end

if (y==c) % g(x,y), Bottom BC
gxy = 0;
end

if (y==d) % g(x,y), Top BC
gxy = 1;
end

end

%%% Setting up System of Eq's Au=B: %%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
x = linspace(a,b,Nx+1); % Defining locations 'x'
y = linspace(c,d,Ny+1); % Defining locations 'x'
h = (b-a)/Nx; % Sub-interval size, assumed hx=hy=h

% Defining "D" Block
D(1) = 4;
for ii = 2:(Nx-1)
D(ii,ii-1) = -1;
D(ii-1,ii) = -1;
D(ii,ii) = 4;
end

% Defining "I" Block
I = -eye([Nx-1,Nx-1]);

A(1:(Nx-1),1:(Nx-1)) = D;
for ii = 2:(Ny-1)
A((ii-1)(Nx-1)+1:(ii-1)(Nx-1)+Nx-1,...
(ii-1)(Nx-1)+1:(ii-1)(Nx-1)+Nx-1) = D;
A((ii-2)(Nx-1)+1:(ii-2)(Nx-1)+Nx-1,...
(ii-1)(Nx-1)+1:(ii-1)(Nx-1)+Nx-1) = I;
A((ii-1)(Nx-1)+1:(ii-1)(Nx-1)+Nx-1,...
(ii-2)(Nx-1)+1:(ii-2)(Nx-1)+Nx-1) = I;
end

% Defining RHS
B = zeros((Ny-1)*(Nx-1),1);
for ii = 1:(Ny-1)
for jj = 1:(Nx-1)

B((ii-1)*(Nx-1)+jj) = -h^2*f(x(jj+1),y(ii+1));

% Adding BC to 'B'
if (ii==1)
B((ii-1)(Nx-1)+jj) = B((ii-1)(Nx-1)+jj) + g(x(jj+1),y(ii),a,b,c,d);
end

% Adding BC to 'B'
if (ii==(Ny-1))
B((ii-1)(Nx-1)+jj) = B((ii-1)(Nx-1)+jj) + g(x(jj+1),y(ii+2),a,b,c,d);
end

% Adding BC to 'B'
if (jj==1)
B((ii-1)(Nx-1)+jj) = B((ii-1)(Nx-1)+jj) + g(x(jj),y(ii+1),a,b,c,d);
end

% Adding BC to 'B'
if (jj==(Nx-1))
B((ii-1)(Nx-1)+jj) = B((ii-1)(Nx-1)+jj) + g(x(jj+2),y(ii+1),a,b,c,d);
end

end
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

u = AB; % Solving Au=b for unknowns'u'

% Placing solution 'u' inside larger array with BCs
uu = zeros(Nx+1,Ny+1);%
for ii = 1:(Nx+1)
for jj = 1:(Ny+1)
uu(ii,jj) = g(x(jj),y(ii),a,b,c,d);
end
end
uu(2:(end-1),2:(end-1)) = reshape(u,(Nx-1),(Ny-1))';

% Plotting Solution
contourf(x,y,uu);
title('BVP Solution u(x,y)'); % Create plot title
xlabel('x'); % Label x-axis
ylabel('y'); % Label y-axis
colorbar;
set(findall(gcf,'-property','fontsize'),'fontsize',20); % Set font size