function out=Fspring(u,Force)

% The inputs for this funtion are the (50 by 1) displacement vector 'u' and
% the scalar 'Force', which represents the force applied at the right-most
% node.  The output is the (50 by 1) vector that represents the system of
% force balance equations evaluated at the current displacements u.

% To set up the function one must consider various cases that depend on
% whether we are balancing forces around the first node or
% last node and on whether the force in the spring included in the force
% balance equation is governed by F_i = k_i(du_i)^0.9 or F_i=k(du_i)^0.5.


% First set up array of spring constants.  Note that what I am doing here
% is generally not considered good programming practice.  Since the array
% of spring constants is the same for all iterations, it would be best to
% set it up only one time outside of this function and pass it in.
% However, to keep the functions self-contained for the purpose of this
% assignment, the spring constants will be set up each time the function is
% called.

% Constants for sprigns 1-25
for i=1:25
    k(i) = 1.0+0.01*i;
end
% Spring constant for springs 26-35
k(26:35)=2;
% Constants for springs 36-50
for i=35:50
   k(i) = 1.0+0.01*i; 
end

% Force balance about first node
out(1,1) = k(2)*sign(u(2)-u(1))*(abs(u(2)-u(1)))^(9/11) - k(1)*sign(u(1))*(abs(u(1)))^(9/11);
% Force balance about nodes 2-24
for i=2:24
   out(i,1)= k(i+1)*sign(u(i+1)-u(i))*(abs(u(i+1)-u(i)))^(9/11) - k(i)*sign(u(i)-u(i-1))*(abs(u(i)-u(i-1)))^(9/11); 
end
% Force balance about node 25
out(25,1) = k(26)*sign(u(26)-u(25))*(abs(u(26)-u(25))) - k(25)*sign(u(25)-u(24))*(abs(u(25)-u(24)))^(9/11);
% Force balace about nodes 26-34
for i=26:34
   out(i,1)=k(i+1)*sign(u(i+1)-u(i))*(abs(u(i+1)-u(i))) - k(i)*sign(u(i)-u(i-1))*(abs(u(i)-u(i-1)));
end
% Force balance about node 35
out(35,1) = k(36)*sign(u(36)-u(35))*(abs(u(36)-u(35)))^(9/11) - k(35)*sign(u(35)-u(34))*(abs(u(35)-u(34)));
%Force balance about nodes 36-49
for i=36:49
   out(i,1) = k(i+1)*sign(u(i+1)-u(i))*(abs(u(i+1)-u(i)))^(9/11) - k(i)*sign(u(i)-u(i-1))*(abs(u(i)-u(i-1)))^(9/11); 
end
% Force balance about 50th node
out(50,1)=Force-k(50)*sign(u(50)-u(49))*(abs(u(50)-u(49)))^(9/11);