function [dsp,vel,acc] = newmark(M,C,K,Ft,Gama,Beta,DeltT,Npoint,NNDOF)
%================================================================
%function [X,t] = newmark(M,C,K,f,dt,gamma,beta,Xi,Xdi)
% This function is intended to perform the numerical
% integration of a structural system subjected to an
% external dynamic excitation.
% The structural model is assumed to be a lumped mass shear
% model.The integration scheme utilized for this analysis
% is the newmark alpha-beta method. The newmark alpha-beta
% method is an implicit time steping scheme so stability of
% the system need not be considered.
%
% Input Variables:
% [M] = Mass Matrix (nxn)
% [C] = Damping Matrix (nxn)
% [K] = Stiffness Matrix (nxn)
% {Ft} = Excitation Vector (mx1)
% dt = Time Stepping Increment
% beta= Newmark Const (1/6 or 1/4 usually)
% gamma = Newmark Const (1/2)
% Xi = Initial Displacement Vector (nx1)
% Xdi = Initial Velocity Vector (nx1)
%
% Output Variables:
% [dsp] = Displacement Response Matrix (mxn)
% [vel] = Velocity Response Matrix (mxn)
% [acc] = Acceleration Response Matrix (mxn)
%================================================================
a0=1/(Beta*DeltT^2);
a1=Gama/(Beta*DeltT);
a2=1/(Beta*DeltT);
a3=1/(2*Beta)-1;
a4=Gama/Beta-1;
a5=(Gama/(2*Beta)-1)*DeltT;
KK=K+a0*M+a1*C; %等效刚度矩阵
dsp=zeros(NNDOF,Npoint);vel=zeros(NNDOF,Npoint); acc=zeros(NNDOF,Npoint); %赋初值
% UPP_M=chol(M);
% UPP_KK=chol(KK);
dsp(:,1)=0; vel(:,1)=0;
% temp1=UPP_M'\Ft(:,1);
%acc(:,1)=M\Ft(:,1); %初始值
acc(:,1)=0;
for i=1:Npoint-1
FFt(:,i+1)=Ft(:,i+1)+M*(a0*dsp(:,i)+a2*vel(:,i)+a3*acc(:,i))+C*(a1*dsp(:,i)+a4*vel(:,i)+a5*acc(:,i)); %t+DeltaT 时刻的等效荷载
%temp2=UPP_KK'\FFt(:,i+1);
dsp(:,i+1)=KK\FFt(:,i+1); %位移
%acc(:,i+1)=a0*(dsp(:,i+1)-dsp(:,i)-DeltT*vel(:,i)-(0.5-Beta)*DeltT*DeltT*acc(:,i)); %加速度
acc(:,i+1)=a0*(dsp(:,i+1)-dsp(:,i))-a2*vel(:,i)-a3*acc(:,i); %加速度
vel(:,i+1)=vel(:,i)+(1-Gama)*DeltT*acc(:,i)+Gama*DeltT*acc(:,i+1); %速度
end
return