地震newmark法

function [deltaudeltadiffudeltaddiffudeltaPdeltadiffPdiffuddiffuu]=newmark（MCKu2diffu2timedtP2） 
clear all
M=input（‘输入质量M ：‘）;
C=input（‘输入阻尼C：‘）;
K=input（‘输入刚度K：‘）;
P=input（‘输入荷载P：‘）;
u2=input（‘输入初始位移u2： ‘）;
diffu2=input（‘输入初始速度diffu2： ‘）;
time=input（‘输入模拟时间time：‘）;
dt=input（‘输入时间步长dt ：‘）;
n=time/dt;
ddiffu2=（P（:1）-C*diffu2-K*u2）/M;
Ke=K+（3*C）/dt+（6*M）/（dt^2）;
a=（6*M）/dt+3*C;
b=3*M+（dt*C）/2;
for i=2:1:n+2;
    if i==2;
        u（i）=u2;