matlab过电压仿真

clear

fid=fopen（‘n.txt‘‘r‘）;
[Ncount]=fscanf（fid‘%d‘inf）;
fclose（fid）;%读取节点数

fid=fopen（‘t.txt‘‘r‘）;
[TTcount]=fscanf（fid‘%f‘inf）;
fclose（fid）;%读取时间常数
tt=TT（1）;tmin=TT（2）;tmax=TT（3）;

fid=fopen（‘et.txt‘‘r‘）;
[ETcount]=fscanf（fid‘%f‘inf）;
fclose（fid）;%读取入侵电压

fid=fopen（‘ruqin.txt‘‘r‘）;
[Z0count]=fscanf（fid‘%d‘inf）;
fclose（fid）;%读取入侵波节点及波阻抗
fid=fopen（‘jd_c.txt‘‘r‘）;

[Ccount]=fscanf（fid‘%d %f‘[2inf]）;%C第一行表示节点，第二行表示该节点接地电容的大小（pf）
[mcnc]=size（C）;
fclose（fid）;%读取节点对地电容

fid=fopen（‘zhl_z_l.txt‘‘r‘）;
[Zcount]=fscanf（fid‘%d %d %f %f‘[4inf]）;%Z第一二行表示支路的节点，第三行阻抗（欧姆），第四行长度（m）
[mznz]=size（Z）;
fclose（fid）;%读取支路波阻抗

fid=fopen（‘blq.txt‘‘r‘）;
if fid==1
    [BLQcount]=fscanf（fid‘%d %f %f %f‘[4inf]）;%第一行：节点；第二行：放电电压；第三行：C；第四行：a
    [mbnb]=size（BLQ）;
    fclose（fid）;%读取避雷器
else
    nb=0;
end

T=Z（4:）/300;
Z=[Z;T];%Z第五行为每条支路波每一次传播的时间us

Y=zeros（NN）;
Y（Z0（1）Z0（1））=Y（Z0（1）Z0（1））+1/Z0（2）;%加入入侵支路阻抗
for ii=1:nz
    Y（Z（1ii）Z（1ii））=Y（Z（1ii）Z（1ii））+1/Z（3ii）;
    Y（Z（2ii）Z（2ii））=Y（Z（2ii）Z（2ii））+1/Z（3ii）;
end%加入支路阻抗
for ii=1:nc
    Y（C（1ii）C（1ii））=Y（C（1ii）C（1ii））+2*C（2ii）*（1e-6）/tt;
end%加入接地电容阻抗
Y;
y=inv（Y）;
%…………………………………………………………以上求节点导纳矩阵…………………………………………………………%

k0=-1/tt*tmin+1;kmax=1/tt*（tmax-tmin）+1;
%设定t最大小值，为使数组下标为正整数，在t和k间建立映射关系，k=1/tt*（t-tmin）+1即t=（k-1）*tt+tmin，求t=0时k0及tmax时kmax

Iz=zeros（NNkmax）;Ic=zeros（Nkm