8×8二维载荷谱编制

Tiki ^

已于 2022-07-15 18:03:57 修改

阅读量823

点赞数 6

分类专栏：载荷数据处理文章标签： matlab 算法

于 2022-07-15 17:49:06 首次发布

本文链接：https://blog.csdn.net/Tikixxx/article/details/125809437

版权

载荷数据处理专栏收录该内容

2 篇文章 0 订阅

订阅专栏

一般均值服从正态分布，幅值服从威布尔分布，拟合求得正态分布均值muHat，标准差sigmaHat；三参数威布尔分布尺度参数A，形状参数B，位置参数C。然后根据极值概率 $P=1\times 10^-6$ ，求得均值、幅值的最大值Xmax、Ymax，最后根据联合概率密度函数计算形成8×8二维载荷谱。

$f(x,y)=\frac{1}{\sqrt{2\pi }\delta } e^{-\frac{(x-\mu)^{2}}{2\delta ^{2}}}\times \frac{\alpha }{\beta }(\frac{y-\varepsilon}{\beta })^{\alpha }e^{-(\frac{y-\varepsilon }{\beta })^{\alpha }}$

x为载荷均值，y为载荷幅值， $\alpha =B$ ， $\beta =A$ ， $\varepsilon =C$ ， $\mu =muHat$ , $\delta =sigmHat$ 。

N为8×8二维载荷谱：

N=nn(muHat,sigmaHat,A,B,C,Xmin,Xmax,Ymax);

主程序：

function  NN = nn(muHat,sigmaHat,A,B,C,Xmin,Xmax,Ymax)
% 8×8二维载荷谱计算
muHat=muHat; % 均值
sigmaHat=sigmaHat; % 标准差
A = A; % 尺度参数
B = B; % 形状参数
C = C; % 阈值参数
Xmin=Xmin;
Xmax=Xmax;
Ymax=Ymax;

% 均值均分为8级
for i = 2:8
   M(:,1) = Xmin;
   M(:,i) = M(:,i-1) + (Xmax-Xmin)/7;
end

% 幅值分级（不等间距）分为8级
pc = [0.125 0.275 0.425 0.575 0.725 0.85 0.95 1];
R = Ymax.*pc;

n=1e6;% 极值概率

% 联合概率密度函数
 fun = @(x,y) (1/(sqrt(2*pi)*sigmaHat)*exp(-(x-muHat).^2/(2*sigmaHat^2))).*(B/A*((y-C)/A).^(B-1).*exp(-((y-C)/A).^B));
%fun = @(x,y) 1/(sqrt(2*pi)*sigmaHat).*exp(-(x-muHat).^2/(2*(sigmaHat^2))) .* (B/A*(y/A).^(B-1).*exp(-(y/A).^B));

% 均值x的下、上限u_w,u_w+1; 幅值y的下、上限v_j,v_j+1;
% u_1=Xmin; u_9=Xmax; v_1=Ymin=0.125*Ymax; v_9=Ymax;

% 主程序： 
for i=1:8
    if i==1
        ymin = 0.125*Ymax;
        ymax = (R(1)+R(2))/2;
        for j=1:8
            if j==1
                xmin = Xmin;
                xmax = (M(1)+M(2))/2;
                NN(i,j)=n*integral2(fun,xmin,xmax,ymin,ymax);
            elseif j==8 
                xmin = (M(7)+M(8))/2;
                xmax = Xmax;
                NN(i,j)=n*integral2(fun,xmin,xmax,ymin,ymax);
            else
                xmin = (M(j-1)+M(j))/2;
                xmax = (M(j)+M(j+1))/2;
                NN(i,j)=n*integral2(fun,xmin,xmax,ymin,ymax);
            end
        end
    elseif i==8
        ymin = (R(7)+R(8))/2;
        ymax = Ymax;
        for j=1:8
            if j==1
                xmin = Xmin;
                xmax = (M(1)+M(2))/2;
                NN(i,j)=n*integral2(fun,xmin,xmax,ymin,ymax);
            elseif j==8 
                xmin = (M(7)+M(8))/2;
                xmax = Xmax;
                NN(i,j)=n*integral2(fun,xmin,xmax,ymin,ymax);
            else
                xmin = (M(j-1)+M(j))/2;
                xmax = (M(j)+M(j+1))/2;
                NN(i,j)=n*integral2(fun,xmin,xmax,ymin,ymax);
            end
        end
    else
        for j=1:8
            if j==1
                xmin = Xmin;
                xmax = (M(1)+M(2))/2;
                ymin = (R(i-1)+R(i))/2;
                ymax = (R(i)+R(i+1))/2;
                NN(i,j)=n*integral2(fun,xmin,xmax,ymin,ymax);
            elseif j==8
                xmin = (M(7)+M(8))/2;
                xmax = Xmax;
                ymin = (R(i-1)+R(i))/2;
                ymax = (R(i)+R(i+1))/2;
                NN(i,j)=n*integral2(fun,xmin,xmax,ymin,ymax);
            else
                xmin = (M(j-1)+M(j))/2;
                xmax = (M(j)+M(j+1))/2;
                ymin = (R(i-1)+R(i))/2;
                ymax = (R(i)+R(i+1))/2;
                NN(i,j)=n*integral2(fun,xmin,xmax,ymin,ymax);
            end
        end
    end
end

% 四舍五入
NN = round(NN);
end