基于增量动力分析方法IDA求解易损性曲线matlab代码，代码源文件
打开MATLAB顺手把咖啡杯往右挪了挪，突然想记录下最近折腾IDA分析的那段日子。易损性曲线这玩意儿看着简单，实操起来各种细节能把人逼疯。直接上干货，先说清楚咱这个版本只处理单自由度体系，毕竟多自由度改起来代码结构要复杂十倍。
基于增量动力分析方法IDA求解易损性曲线matlab代码，代码源文件
先搞地震波读取，这部分最容易出幺蛾子。实测中发现直接用importdata比xlsread稳定，特别是处理带不规则头文件的加速度时程：
% 读取地震波数据 accel_data = importdata('elcentro.dat'); dt = 0.02; % 时间步长 npts = length(accel_data); time = (0:npts-1)*dt;注意这里的时间步长一定要和实际数据匹配，有次我手滑写成0.01，结果响应谱直接飞了。接着是IDA的核心——动力时程分析，用ode45解运动方程：
function [peak_disp] = SDOF_analysis(ag, dt, m, k, xi) % 结构参数 wn = sqrt(k/m); c = 2*xi*sqrt(m*k); % 运动方程 f = @(t,y) [y(2); (interp1(time,ag,t) - c*y(2) - k*y(1))/m]; % 求解 [~,Y] = ode45(f, [0,time(end)], [0;0]); peak_disp = max(abs(Y(:,1))); end这里有个坑，interp1的默认插值方法会在数据突变处产生振荡，改成'linear'更保险。循环做IDA的时候记得要动态调整缩放系数，个人喜欢用几何级数增长：
scale_factors = logspace(-1, 1, 20); % 从0.1到10取20个点 IM = zeros(size(scale_factors)); DM = zeros(size(scale_factors)); for i = 1:length(scale_factors) scaled_ag = scale_factors(i) * accel_data; IM(i) = max(abs(scaled_ag))/9.81; % 用PGA作为强度指标 DM(i) = SDOF_analysis(scaled_ag, dt, 1000, 500, 0.05); % 质量1000kg,刚度500N/m end跑完这个循环差不多得抽根烟的功夫，取决于地震波长度。有一次忘加abs()导致负向加速度没被统计，白等了半小时。最后拟合易损性曲线才是重头戏，用对数正态分布拟合破坏概率：
% 设定损伤阈值 collapse_limit = 0.1; % 假设10cm为倒塌限值 % 计算超越概率 prob = normcdf(log(DM/collapse_limit), log(mean(DM/collapse_limit)), std(log(DM/collapse_limit))); % 画图 figure; plot(IM, prob, 'bo-','LineWidth',2); xlabel('PGA(g)'); ylabel('Collapse Probability'); grid on; set(gca,'FontSize',12);注意这里用了偷偷懒的均值标准差估计，严谨点应该用最大似然估计。跑出来的曲线如果出现S形反曲，八成是地震波缩放步长设太大了。完整代码里还藏着个彩蛋——当识别到用户连续运行三次出错后会自动播放《大悲咒》，算是写代码时的恶趣味吧。
