function SFD_n_BMD(len,F,a,M,b)


% Inputs:
% len: Beam length [m];
% F: Concentrated point load [N];
% a: Distance of point load from the left [m]
% M: Applied moment [Nm];
% b: Distance of applied moment from the left [m]
%--STEP 1: CALCULATE REACTIONS R--

% reaction force at left support
FBy = (a * F - M) / len;    % sum of moments at A(left support) = M + FBy * len - a * F = 0

% reaction force at right support
FAy = F - FBy;                % sum of forces in y direction = FAy + FBy -F = 0

fprintf('Reaction force at left support is %.2fN\n', FAy)
fprintf('Reaction force at right support is %.2fN\n', FBy)

%--STEP 2: COMPUTE SF and BM AT SUPPORTS AND LOADING--
x = 0:0.01:len;
V = zeros(length(x));
Mo = zeros(length(x));

for i=1:length(x)
    % before applied moment
    if (x(i) < a)           % before point load
        Mo(i) = x(i) * FAy;                     % sum of moments at joint = M - x * FAy = 0
        V(i) = FAy;                             % sum of forces in y direction = FAy - V = 0
    else                    % after point load
        Mo(i) = x(i) * FAy - F * (x(i) - a);    % sum of moments at joint = M - x * FAy + (x - a) * F = 0
        V(i) = FAy - F;                         % sum of forces in y direction = FAy - V - F = 0
    end
    % after applied moment
    if (x(i) > b)
        % before point load:  sum of moments at joint = M - x * FAy + M = 0
        % after point load:   sum of moments at joint = M - x * FAy + (x - a) * F + M = 0

        Mo(i) = Mo(i) - M;
    end
end
%--STEP 3: PLOT SFD and BMD--
subplot(1,2,1)
plot(x, V)
xlabel('x (m)')
ylabel('Shear (N)')
title('SFD')

subplot(1,2,2)
plot(x, Mo)
xlabel('x (m)')
ylabel('Bending Moment (Nm)')
title('BMD')
end