function  [dis_new,vel_new,acc_new] = Newmark(M,C,K_inv,R_new,dis_old,vel_old,acc_old,a0,a1,a2,a3,a4,a5,a6,a7) 
%--(1)--OUTPUT VARIBALES
%-------dis_new---new displcement vector at time (t + delt t);
%-------vel_new---new velocity vector at time (t + delt t);
%-------acc_new---new acceloration vector at time (t + delt t);
%--(2)--INPUT VARIBALES
%-------M---------mass matrix;
%-------C---------damping matrix;
%-------K---------stiffness matrix;
%-------dis_old---old displcement vector at time t;
%-------vel_old---old velocity vector at time t;
%-------acc_old---old acceloration vector at time t;
%-------K_inv-----inverse of effective stiffness matrix;
%-------R_new-----new load vector at time (t + delt t);

%---(3)---Calculate effective loads at time (t + delt t)
R_eff = R_new + M*(a0*dis_old + a2*vel_old + a3*acc_old)+ C*(a1*dis_old + a2*vel_old + a5*acc_old );
%---(4)----Calculate new displacement at time (t + delt t)
dis_new = K_inv*R_eff;
%---(5)----Calculate new velocity at time (t + delt t)
vel_new = vel_old + a6*acc_old + a7*acc_new;
%---(6)----Calculate new acceleration at time (t + delt t)
acc_new = a0*(dis_new - dis_old )- a2*vel_old - a3 * acc_old ;