function retval = rhs(el,M,dt,n,U,u,gmma,bta,pnodes,qnodes)
% usage: b = rhs(el,M,dt,n,U,u,gmma,bta,pnodes,qnodes)
% description: returns the time dependent right hand side
% vector for this problem, at time t = n*dt. Input the 
% prior values of U, u at (n-1)*dt.

% local variables: dx, xnodes, Fnodes, qnodes,

retval = zeros(M,1);
dx = el/M;
xnodes = linspace(dx,el,M)';
Fnodes = F(xnodes,n*dt);
retval = dx*dx*(Fnodes + (gmma/exp(gmma*dt)*qnodes).*(U + dt/2*u));
% fix the first entry
retval(1,1) = retval(1,1) + (1+pnodes(1)*dx/2)*mu(n*dt);
% fix the last entry by centered differences, etc
retval(M,1) = retval(M,1) + (1-pnodes(M)*dx/2)*2*dx/bta*alpha(n*dt);
% or do last row by a one-sided derivative formula (more accurate)
% retval(M,1) = 2*dx*alpha(n*dt);