分段线性滞回的MATLAB语言实现

guokaiyuan

[img]file:///C:/DOCUME~1/凯@豪/LOCALS~1/Temp/Z)[{DR5TX)G6U%X{5AZ_EAH.jpg[/img] [img]file:///C:/DOCUME~1/凯@豪/LOCALS~1/Temp/UY4(]4I`C274GG]$$T_1K@R.jpg[/img]
上图是一个分段线性滞回的力学模型，右面是它的数学表达式，我将它嵌入到Newmark-beta法当中，得到相应的响应，然后画出上面的图形出来，可是发现得到的图形是错误的，我的编程如下：
PS：问题之一：本来应该是执行if语句，而且也感觉是符合if的判断语句，结果偏偏执行第一个elseif语句。这是怎么回事呢？

clear
clc
M=1;       %输入质量矩阵
K=1*10^4;  %输入刚度矩阵
C=20;    %输入阻尼矩阵
%%%%%%%%%%%%%%%%%%%%%%%%以上计算结构的物理特性参数
[row,column]=size(M);  
h=0.001;         %输入时间步长  
t=4;             %输入最终时间
N=t/h;           %求出计算点数
x=0:h:(N-1)*h;
P = 100*sin(25*x);  %简谐激励
M_ni=inv(M);                 %对质量矩阵求逆
S(:,1)=zeros(row,1);         %取初始位移为零
V(:,1)=zeros(row,1);         %取初始速度为零
F_SMA(:,1) = 0;              %令初始力为0
%% 接下来求线性部分的位移、速度和加速度作为非线性迭代的初始值
A(:,1)=M_ni*(P(:,1)-C*V(:,1)-K*S(:,1)-F_SMA(:,1));  
equK=K+3*C/h+6*M/(h^2);      %计算等效刚度
equK_ni=inv(equK);      
for i=1:N-1;
    equP=P(:,i+1)-F_SMA(:,i)+M*(6*S(:,i)/(h^2)+6*V(:,i)/h+2*A(:,i))+C*(3*S(:,i)/h+2*V(:,i)+h*A(:,i)/2);    %计算等效荷载
    S(:,i+1)=equK_ni*equP;            %计算位移迭代初始值
    V(:,i+1)=3/h*(S(:,i+1)-S(:,i))-2*V(:,i)-h*A(:,i)/2;     %计算速度迭代初始值
   
%%%%%%%%%%%%%%%%%%%%%%%%%%求解非线性恢复力
K1=5000; K2=0.5*K1; X_b=0.004; X_a=0.5*X_b;

    X_d = 0.0071;
    X_c = X_d-0.002;
if (S(:,i+1)-S(:,i)>=0)&&(0<=S(:,i+1)<=X_b)
    F_SMA(:,i+1) = K1*S(:,i+1);                 %oab
elseif (S(:,i+1)-S(:,i)>0)&&(X_b<S(:,i+1)<=X_d)
    F_SMA(:,i+1) = (K1-K2)*X_b+K2*S(:,i+1);       %bd
elseif (S(:,i+1)-S(:,i)<=0)&&(X_c<S(:,i+1)<X_d)
    F_SMA(:,i+1) = (K1-K2)*(X_b-X_d)+K1*S(:,i+1);             %dc
elseif (S(:,i+1)-S(:,i)<=0)&&(X_a<=S(:,i+1)<=X_c)
    F_SMA(:,i+1) =  K2*S(:,i+1)+(K1-K2)*X_a;      %ca
elseif (S(:,i+1)-S(:,i)<=0)&&(-X_b<=S(:,i+1)<X_a)
    F_SMA(:,i+1) = K1*S(:,i+1);                 %b'ao
elseif (S(:,i+1)-S(:,i)<=0)&&(-X_b>S(:,i+1)>=-X_d)
    F_SMA(:,i+1) = -(K1-K2)*X_b+K2*S(:,i+1);       %b'd'
elseif (S(:,i+1)-S(:,i)>=0)&&(-X_c>S(:,i+1)>=-X_d)
    F_SMA(:,i+1) =  -(K1-K2)*(X_b-X_d)+K1*S(:,i+1);            %d'c'
elseif (S(:,i+1)-S(:,i)>=0)&&(-X_c<=S(:,i+1)<=-X_a)
    F_SMA(:,i+1) = K2*S(:,i+1)-(K1-K2)*X_a;       %c'a'
else
    F_SMA(:,i+1) = K1*S(:,i+1);                 %a'o
end
A(:,i+1)=M_ni*(P(:,i+1)-F_SMA(:,i+1)-C*V(:,i+1)-K*S(:,i+1));         %计算加速度迭代初始值
end

%% figure
figure(1);
plot(S,F_SMA);
grid on;

yufeng

条件判断 不能联写 要一个个写 然后用and

guokaiyuan

回复 2 # yufeng 的帖子

谢谢您，不能连写，要一个个写，然后用&&具体怎么操作呀？可能没完全理解您的意思。
以if (S(:,i+1)-S(:,i)>=0)&&(0<=S(:,i+1)<=X_b)为例，该怎么改写呢？十分感谢

guokaiyuan

回复 2 # yufeng 的帖子

yufeng，我非常感谢你给我的建议，真的是太感谢你了，再次感谢，哈哈