u结构可靠度分析与设计的编程实践.doc
【精品文档】如有侵权，请联系网站删除，仅供学习与交流u结构可靠度分析与设计的编程实践.精品文档.结构可靠度分析与设计的编程实践摘要：本文根据结构可靠度分析与设计的基本原理，选用了MATLAB语言，编写计算机程序，对课程中所给出的例题进行编程实现，并给出结果展示。1. 引言已建结构的可靠性理论产生于20世纪70年代发达国家维修改造业迅速发展时期，研究内容主要集中于已建结构的损伤评估、模式识别和可靠性评价。1994年美国的姚治平和德国的Natke对已建结构的损伤检测及可靠性评价作了更深入的研究。现有的国际标准和一些国家的标准都可用可靠性指标来度量结构的可靠性。它是结构可靠性分析中的又一重要概念。求解结构的可靠性指标是一个非常复杂的问题，这是由于影响结构可靠性的因素繁杂，而且各因素之间相互影响。文中引入了验算点法和JC法，主要是介绍用MATLAB来实现这两种算法的程序。2. 结构可靠度分析的基本原理 2.1结构可靠度分析建立的结构可靠与不可靠的界限，称为极限状态。我国将极限状态分为承载能力极限状态(包括条件极限状态)和正常使用极限状态两类。承载能力极限状态对应于结构或结构构件达到最大承载能力或出现不适于继续承载的变形；正常使用极限状态对应于结构或结构构件达到正常使用或耐久性能的某项规定限值；条件极限状态也称破坏安全极限状态，对应于已局部出现破坏的结构的最大承载能力。结构的极限状态可用下列极限状态方程描述：式中：，是指结构上各种作用或作用效应、材料性能、几何参数等. 其中，结构的功能函数或功效函数为： 对于承载努力极限状态，若令为结构抗力，为作用综合效应，则有式中：，结构处于可靠状态；，结构处于失效状态；，结构处于极限状态。 根据结构的极限状态和功能函数可得结构的可靠度(即可靠概率)和失效概率 用概率干涉法计算的公式为： 由于功能函数中有许多个若用概率干涉法，直接计算则十分复杂，而且不一定能计算出。若令，则可用近似计算，其中为结构可靠度指标，本文中基于验算点法和JC法用MATLAB进行编程计算。2.2 验算点法2.2.1 功能函数为线性函数 功能函数随机变量是一个正态随机变量，其概率密度函数和U的密度曲线如图1示 图1 一个随机变量时的可靠指标（左图为正态随机变量，右图为标准正态随机变量）假定存在n个相互独立的随机变量，其均值为，标准差为结构功能函数为： （1）其中为常数将随机变量变换为标准正态随机变量则由（1）表示的功能函数表示成从而功能函数的平均值和标准差表示为按照严格的可靠度指标定义 （3）可靠度指标和结构失效概率存在精确的对应关系对极限状态方程 两端同时除以得到： （4）与公式（3）比较，有 （5）令公式（5）可以写成： (6）公式（6）表示的是一法线式的直线方程，为法线与坐标轴夹角余弦 图2 可靠度指标的几何意义及验算点验算点在空间（标准正态空间）表示为：在空间表示为：两者之间的关系为：根据几何关系有：在空间，验算点坐标值：通常表示为：2.1.3功能函数为线性函数假定随机变量服从正态分布，但结构功能函数不再是线性函数，显然，这时精确 求解的平均值和标准差是非常困然的。同结构功能函数为非线性的情形一样，如果将可靠指标定义为标准正态坐标系中坐标原点到极限状态曲面的距离，垂足为验算点，则不管结构极限状态方程的数学表达形式如何，只要具有相同的力学或物理含义，在标准正态坐标系中，所表示的都将是同一个曲面，曲面上与坐标原点距离最近的点也只有一个。因而，所得到的可靠指标是唯一的，不像中心点法那样，随结构极限状态方程数学表达式的形式而变。 图3验算点取法如果验算点已知可以在该点一次项展开：其均值和标准差为：所以可靠度指标：实际上验算点不可知，需要补充条件：对比表达式得到：2.3 JC法验算算点法的特点是能够考虑非正态的随机变量，在计算工作量增加不多的条件下，可对可靠指标进行精度较高的计算。对于极限状态方程中包含非正态分布的随机变量的情形，在进行其可靠度分析时，一般要把非正态随机变量当量化为正态随机变量。当量正态化方法即为JC法。它的基本思想就是：在设计验算点处，当量正态随机变量 (其均值,标准差为)的分布函数值与原随机变量(其均值，标准差为)的分布函数值相等；在设计验算点处，当量正态随机变量 (其均值,标准差为)的概率密度函数值与原随机变量(其均值，标准差为)的概率密度函数值相等。 下面详细介绍一下验算点法的具体步骤。引入标准正态随机变量，令： i=1,2,n 极限状态方程可表示为 定义方向余弦： 其中，表示标准正态空间内极限状态曲面的切平面的法线垂足，又可称为设计验算点。表示功能函数对基本随机变量的偏导数在设计验算点处的值。根据方向余弦的定义可得： 因此根据，两式可将原随机变量表示为：，i=1,2,n 联合，三式可求解和原始设计验算点。上述计算过程通常需要迭代直到前后两次计算所得可靠度相差不大，且原始设计验算点满足极限状态方程式为止。因此， 和的迭代计算过程为：(1)设置原始随机向量初始迭代点，如取均值，即，i=1,2,n；(2)由式计算标准正态空间的设计验算点；(3)计算功能函数梯度和方向余弦在设计验算点处的值；(4)将方向余弦值和随机变量的均值、标准差代入，两式求出 值和原始设计验算点；(5)判断前后两次计算所得的值之间的误差是否满足精度要求，以及设计验算点是否满足极限状态方程。如果不满足要求则将(3)步中计算所得的代入步骤(2)，重复步骤(3)至步骤(5)，直到结果满足要求为止。3. 结构可靠性设计的基本原理基于可靠度的概率设计水准分为四个等级：半概率设计，全概率设计，最优概率设计，近似概率设计。直接基于可靠度的概率设计方法是：l 给定：荷载与抗力的概率分布类型和特征值l 寻找：设计向量l 校准：基本原理：l 给定：荷载与抗力的概率分布类型和特征值l 寻找：抗力均值l 满足：4. 例题的程序运行结果展示 4.1 Project2.1的结果展示 (1)112页例题5.4迭代次数113.182510200000000.00082-4.6783510.66946-3.2E+07-0.0003834.24159.836384115555890.00014843.9556510.06272116583830.00021653.7146510.0768796331330.00029663.3911510.0811267408980.00041173.208510.0650348036630.00057683.1815510.0474944509570.00068393.1805510.0442443970810.000709103.1805510.0437843754260.000713113.1805510.0435843698710.000714(2)142页作业5.3迭代次数13.15374375121.228437.523.1537236.8678121.228708.38 4.2 Project2.2的结果展示(1)123页例题5.9迭代次数14.425220010023.903123.2315123.231533.7274144.9522144.952243.7033156.0206156.020653.7015159.5028159.502863.7014160.4051160.405173.7014160.6263160.6263(2)127页作业5.11迭代次数1100402000299.9993423.423112342.285399.9995829.77532977.516499.9996732.388363238.825599.9996832.713523271.342699.9996932.732663273.256(3)142页作业5.4迭代次数12.382312.52422.24689.41848328.2554532.213910.5011131.5033442.213210.7060932.1182652.213210.7256232.17686 4.3 Project2.3的结果展示 令the number of simulation（N）=1e7 计算结果得：, 以下是各个变量的概率图：4.4 Project2.4的结果展示 (1)第六章例题6.3迭代次数1159.68692186.10723193.78654194.66215194.73346194.73927194.73978194.73979194.7397(2)第六章例题6.7迭代次数1-0.64020.76822-0.79640.60483-0.80.64-0.80.64.5 Project2.5的结果展示(1)第六章例题7.1时，迭代次数13.220141.2310.6724.025519.1702510.84418.32615333.838223.3355511.1385812.1969743.809925.5688210.9803114.5885153.808226.2806210.929315.3513163.808226.4194710.9167815.5026973.808226.4431410.9144715.52867时，迭代次数13.319530.92255.3723.776214.278515.3818678.89664433.592919.013715.44442613.5692843.583121.114335.39187815.7224553.58321.436435.38054816.0558863.58321.459095.37912616.079975. 结论与讨论由上面的分析和计算过程，我们可以知道，计算结构可靠度指标是一件非常繁琐的工作。而恰当地利用计算机进行编程计算能大大降低工作量。验算点法较为简单，适用于计算随机变量分布为正态分布的功能函数的结构可靠度指标，JC法较为复杂，擅长于计算随机变量为非正态分布的功能函数的结构可靠度指标。6. 参考文献1. 张明著. 结构可靠度分析 方法与程序. 北京：科学出版社，2009年.2. 李国强等编著. 工程结构荷载与可靠度设计原理（第三版）. 北京：中国建筑工业出版社，2005.3. 贡金鑫，魏巍巍. 工程结构可靠性设计原理. 北京：机械工业出版社，2007.4. Nowak, A. S., and K. R. Collins (2000). Reliability of Structures. McGraw-Hill, New York, NY. 张川导读，重庆：重庆大学出版社，2005年3月7. 附件 附件1：改进一次二阶矩方法的程序和一次可靠度方法（R-F方法）的程序clearclcdisplay('请选择计算模式')judge1=input('验算点法请输入1，jc法请输入2:');if judge1=1 display('计算任意功能函数的可靠指标请输入1'); display('计算112页例题5.4请输入2，计算142页作业5.3请输入3'); judge11=input(':'); %输入数据 if judge11=2 g=sym('x1/360-0.0069*x2*x14/(x3*x4)'); %mean1是均值矩阵 mean1=5 10 2*107 8*10(-4); %mean1=38 54; %std是标准差矩阵 std=0 0.4 0.5*10(7) 1.5*10(-4); %std=3.8 5.4; distribution=1 1 1 1; elseif judge11=3 g=sym('x1*x2-x3'); mean1=4375 121.2 28437.5; cov=0.3 0 0.13; std=mean1.*cov; distribution=1 1 1; elseif judge11=1 g=input('请将随机变量依次换成x1,x2,x3.后输入功能函数(形式如:x1+x2+x3):','s'); g=sym(g); mean1=input('请依次对应每一个变量输入均值矩阵(形式如1 2 3 4):'); judge111=input('选择输入变异系数矩阵请输入1，选择输入标准差矩阵请输入2:'); if judge111=1 cov=input('请依次对应每一个变量输入变异系数矩阵(形式如0.05 0.12 0.20 0.15):'); std=mean1.*cov; else std=input('请依次对应每一个变量输入标准差矩阵(形式如1 2 3 4):'); end distribution(1:length(mean1)=1; else endelse display('计算任意功能函数的可靠指标请输入1，计算123页例题5.9请输入2') display('计算127页作业5.11请输入3,计算142页作业5.4请输入4,计算第七章例7.1请输入5') judge22=input(':'); if judge22=2; g=sym('x1-x2'); mean1=200 100; std=20 12; distribution=1 3; elseif judge22=3; g=sym('x1*x2-x3'); mean1=100 40 2000; cov=0.0004 0.1 0.1; std=mean1.*cov; distribution=1 2 3; elseif judge22=5; g=sym('x1-x2-x3'); lamada=1.33 1.06 0.70; cov=0.17 0.07 0.29; rou=input('请输入：'); x3_k=10; x2_k=x3_k/rou; k=1.55; x1_k=k*(x2_k+x3_k); x_k=x1_k x2_k x3_k; mean1=lamada.*x_k; std=cov.*mean1; distribution=2 1 3; elseif judge22=4; g=sym('3*x1-x2'); mean1=12.5 24; cov=0.125 0.15; std=mean1.*cov; distribution=2 3; else g=input('请将随机变量依次换成x1,x2,x3.后输入功能函数(形式如:x1+x2+x3):','s'); g=sym(g); display('请依次对应每一个变量输入分布函数类型矩阵') distribution=input('对应每个随机变量，正态分布取1，对数正态分布取2，极值1型取3(形式如1 2 3 3):'); mean1=input('请依次对应每一个变量输入均值矩阵(形式如1 2 3 4):'); judge111=input('选择输入变异系数矩阵请输入1，选择输入标准差矩阵请输入2:'); if judge111=1 cov=input('请依次对应每一个变量输入变异系数矩阵(形式如0.05 0.12 0.20 0.15):'); std=mean1.*cov; else std=input('请依次对应每一个变量输入标准差矩阵(形式如1 2 3):'); end endendvar=symvar(g); %找出g中的随机变量%length(var) length_var=length(var); %计算随机变量的个数for i=1:length_var g1(i)=diff(g,var(i);end %分别对每个随机变量求偏导g=char(g);for i=1:length_var gg(i)=char(g1(i);endmean=zeros(1,10);sd=zeros(1,10);mean(1:length_var)=mean1;sd(1:length_var)=std;x=mean;beta1=1;beta=0;count=0;cos_theta=zeros(1,10);mean_=mean;sd_=sd;while(abs(beta-beta1)>1e-5) beta1=beta;%当量正态化 for i=1:length_var if distribution(i)=2 a1=sqrt(log(1+(sd_(i)./mean_(i).2); a2=log(mean_(i)-0.5*a1.2; F_x=normcdf(log(x(i)-a2)./a1,0,1); f_x=(1./sqrt(2*pi)./a1./x(i).*exp(-0.5*(log(x(i)-a2)./a1).2); sd(i)=normpdf(norminv(F_x,0,1),0,1)./f_x; mean(i)=x(i)-norminv(F_x,0,1).*sd(i); elseif distribution(i)=3 alpha=pi./sd_(i)./sqrt(6); u=mean_(i)-0.5772/alpha; a1=exp(-alpha*(x(i)-u); f_x=alpha.*a1.*exp(-a1); F_x=exp(-a1); sd(i)=normpdf(norminv(F_x,0,1),0,1)./f_x; mean(i)=x(i)-norminv(F_x,0,1).*sd(i); end endfor i=1:length_var canshu2=mean-x; x1=x(1);x2=x(2);x3=x(3);x4=x(4);x5=x(5);x6=x(6);x7=x(7);x8=x(8);x9=x(9);x10=x(10); canshu(i)=eval(ggi).*canshu2(i); canshu1(i)=(eval(ggi).*sd(i).2; canshu11(i)=-eval(ggi).*sd(i);end count=count+1;beta=(eval(g)+sum(canshu)./(sqrt(sum(canshu1);beta_m(count)=beta;for i=1:length_var cos_theta(i)=canshu11(i)./(sqrt(sum(canshu1); cos_theta_m(count,i)=cos_theta(i); xx(count,i)=x(i);end x=mean+beta.*sd.*cos_theta;endbeta beta_m=beta_m'cos_theta_mxx 附件2：结构构件可靠度校准分析的程序 6.5:clcclearg=sym('x1-x2-x3');mean1=0 50 70;std=0 3.5 20;distribution=2 1 3;var=symvar(g); %length(var) length_var=length(var); for i=1:length_var g1(i)=diff(g,var(i);end g=char(g);for i=1:length_var gg(i)=char(g1(i);endmean=zeros(1,10);sd=zeros(1,10);mean(1:length_var)=mean1;sd(1:length_var)=std;x=mean;x(1)=120;x11=1;count=0;cos_theta=zeros(1,10);mean_=mean;sd_=sd;beta=3.2;while(abs(x11-x(1)>1e-5) x11=x(1); for i=1:length_var if distribution(i)=2 sd(i)=x(i).*sqrt(log(1+0.172); elseif distribution(i)=3 alpha=pi./sd_(i)./sqrt(6); u=mean_(i)-0.5772/alpha; a1=exp(-alpha*(x(i)-u); f_x=alpha.*a1.*exp(-a1); F_x=exp(-a1); sd(i)=normpdf(norminv(F_x,0,1),0,1)./f_x; mean(i)=x(i)-norminv(F_x,0,1).*sd(i); end endfor i=1:length_var x1=x(1);x2=x(2);x3=x(3);x4=x(4);x5=x(5);x6=x(6);x7=x(7);x8=x(8);x9=x(9);x10=x(10); canshu1(i)=(eval(ggi).*sd(i).2; canshu11(i)=-eval(ggi).*sd(i);end count=count+1;%beta=(eval(g)+sum(canshu)./(sqrt(sum(canshu1);for i=1:length_var cos_theta(i)=canshu11(i)./(sqrt(sum(canshu1); cos_theta_m(count,i)=cos_theta(i);end x(2)=mean(2)+beta.*sd(2).*cos_theta(2); x(3)=mean(3)+beta.*sd(3).*cos_theta(3); x(1)=x(2)+x(3); xx(count)=x(1);endmean(1)=x(1)-sd(1)*beta*cos_theta(1)miu_r=sqrt(1+0.172)*exp(mean(1)/x(1)-1+log(x(1)Rk=miu_r/1.33agk=Rk/(0.01/0.31/0.02+1)/0.31 6.7:clearclc%g=sym('x1-x2');%g1=diff(g,x1);%g2=diff(g,x2);miu_q=1;miu_r=1;x1=1;x2=x1;x11=0;count=0;beta=3;while (abs(x1-x11)>1e-5) count=count+1; x11=x1; bili=x2/x1;gr=0.1*bili*miu_r;gq=-0.12*miu_q;ar=-gr/(sqrt(gr2+gq2);a111(count)=ar;%arr(count)=ar;aq=-gq/(sqrt(gq2+gr2);a222(count)=aq;%aqq(count)=aq;x1=1+ar*beta*0.1*1;x2=1+aq*beta*0.12*1;%miu_rr(count)=miu_r;endx11/x1x2 附件3：分项系数校准分析的程序 见附件1内部资料，请勿外传！9JWKffwvG#tYM*Jg&6a*CZ7H$dq8KqqfHVZFedswSyXTy#&QA9wkxFyeQ!djs#XuyUP2kNXpRWXmA&UE9aQGn8xp$R#͑GxGjqv$UE9wEwZ#QcUE%&qYpEh5pDx2zVkum&gTXRm6X4NGpP$vSTT#&ksv*3tnGK8!z89AmYWpazadNu#KN&MuWFA5uxY7JnD6YWRrWwcvR9CpbK!zn%Mz849GxGjqv$UE9wEwZ#QcUE%&qYpEh5pDx2zVkum&gTXRm6X4NGpP$vSTT#&ksv*3tnGK8!z89AmYWpazadNu#KN&MuWFA5uxGjqv$UE9wEwZ#QcUE%&qYpEh5pDx2zVkum&gTXRm6X4NGpP$vSTT#&ksv*3tnGK8!z89AmYWpazadNu#KN&MuWFA5uxY7JnD6YWRrWwcvR9CpbK!zn%Mz849GxGjqv$UE9wEwZ#QcUE%&qYpEh5pDx2zVkum&gTXRm6X4NGpP$vSTT#&ksv*3tnGK8!z89AmUE9aQGn8xp$R#͑GxGjqv$UE9wEwZ#QcUE%&qYpEh5pDx2zVkum&gTXRm6X4NGpP$vSTT#&ksv*3tnGK8!z89AmYWpazadNu#KN&MuWFA5uxY7JnD6YWRrWwcvR9CpbK!zn%Mz849GxGjqv$UE9wEwZ#QcUE%&qYpEh5pDx2zVkum&gTXRm6X4NGpP$vSTT#&ksv*3tnGK8!z89AmYWpazadNu#KN&MuWFA5uxGjqv$UE9wEwZ#QcUE%&qYpEh5pDx2zVkum&gTXRm6X4NGpP$vSTT#&ksv*3tnGK8!z89AmYWpazadNu#KN&MuWFA5uxY7JnD6YWRrWwcvR9CpbK!zn%Mz849GxGjqv$UE9wEwZ#QcUE%&qYpEh5pDx2zVkum&gTXRm6X4NGpP$vSTT#&ksv*3tnGK8!z8vG#tYM*Jg&6a*CZ7H$dq8KqqfHVZFedswSyXTy#&QA9wkxFyeQ!djs#XuyUP2kNXpRWXmA&UE9aQGn8xp$R#͑GxGjqv$UE9wEwZ#QcUE%&qYpEh5pDx2zVkum&gTXRm6X4NGpP$vSTT#&ksv*3tnGK8!z89AmYWpazadNu#KN&MuWFA5uxY7JnD6YWRrWwcvR9CpbK!zn%Mz849GxG89AmUE9aQGn8xp$R#͑GxGjqv$UE9wEwZ#QcUE%&qYpEh5pDx2zVkum&gTXRm6X4NGpP$vSTT#&ksv*3tnGK8!z89AmYWpazadNu#KN&MuWFA5uxY7JnD6YWRrWwcvR9CpbK!zn%Mz849GxGjqv$UE9wEwZ#QcUE%&qYpEh5pDx2zVkum&gTXRm6X4NGpP$vSTT#&ksv*3tnGK8!z89AmYWpazadNu#KN&MuWFA5uxGjqv$UE9wEwZ#QcUE%&qYpEh5pDx2zVkum&gTXRm6X4NGpP$vSTT#&ksv*3tnGK8!z89AmYWpazadNu#KN&MuWFA5uxY7JnD6YWRrWwcvR9CpbK!zn%Mz849GxGjqv$UE9wEwZ#QcUE%&qYpEh5pDx2zVkum&gTXRm6X4NGpP$vSTT#&ksv*3tnGK8!z8vG#tYM*Jg&6a*CZ7H$dq8KqqfHVZFedswSyXTy#&QA9wkxFyeQ!djs#XuyUP2kNXpRWXmA&UE9aQGn8xp$R#͑GxGjqv$UE9wEwZ#QcUE%&qYpEh5pDx2zVkum&gTXRm6X4NGpP$vSTT#&ksv*3tnGK8!z89AmYWpazadNu#KN&MuWFA5uxY7JnD6YWRrWwcvR9CpbK!zn%Mz849GxGjqv$UE9wEwZ#QcUE%&qYpEh5pDx2zVkum&gTXRm6X4NGpP$vSTT#&ksv*3tnGK8!z89AmYWpazadNu#KN&MuWFA5uxGjqv$UE9wEwZ#QcUE%&qYpEh5pDx2zVkum&gTXRm6X4NGpP$vSTT#&ksv*3tnGK8!z89AmYWpazadNu#KN&MuWFA5uxY7JnD6YWRrWwcvR9CpbK!zn%Mz849GxGjqv$UE9wEwZ#QcUE%&qYpEh5pDx2zVkum&gTXRm6X4NGpP$vSTT#&ksv*3tnGK8!z89AmUE9aQGn8xp$R#͑GxGjqv$UE9wEwZ#QcUE%&qYpEh5pDx2zVkum&gTXRm6X4NGpP$vSTT#&ksv*3tnGK8!z89AmYWpazadNu#KN&MuWFA5uxY7JnD6YWRrWwcvR9CpbK!zn%Mz849GxGjqv$UE9wEwZ#QcUE%&qYpEh5pDx2zVkum&gTXRm6X4NGpP$vSTT#&ksv*3tnGK8!z89AmYWpazadNu#KN&MuWFA5uxGjqv$UE9wEwZ#QcUE%&qYpEh5pDx2zVkum&gTXRm6X4NGpP$vSTT#&ksv*3tnGK8!z89AmYWpazadNu#KN&MuWFA5uxY7JnD6YWRrWwcvR9CpbK!zn%Mz849GxGjqv$UE9wEwZ#QcUE%&qYpEh5pDx2zVkum&gTXRm6X4NGpP$vSTT#&ksv*3tnGK8!z89AmYWv*3tnGK8!z89AmYWpazadNu#KN&MuWFA5uxY7JnD6YWRrWwcvR9CpbK!zn%Mz849GxGjqv$UE9wEwZ#QcUE%&qYpEh5pDx2zVkum&gTXRm6X4NGpP$vSTT#&ksv*3tnGK8!z89AmYWpazadNuGK8!z89AmYWpazadNu#KN&MuWFA5uxY7JnD6YWRrWwcvR9CpbK!zn%Mz849GxGjqv$UE9wEwZ#QcUE%&qYpEh5pDx2zVkum&gTXRm6X4NGpP$vSTT#&ksv*3tnGK8!z89AmYWpazadNu#KN&MuWFA5uxGjqv$UE9wEwZ#QcUE%&qYpEh5pDx2zVkum&gTXRm6X4NGpP$vSTT#&ksv*3tnGK8!z89AmYWpazadNu#KN&MuWFA5uxY7JnD6YWRrWwcvR9CpbK!zn%Mz849GxGjqv$UE9wEwZ#QcUE%&qYpEh5pDx2zVkum&gTXRm6X4NGpP$vSTT#&ksv*3tnGK8!z89AmYWv*3tnGK8!z89AmYWpazad
