改进以下代码a=0;b=1;h=(b-a);T1=0.920735492403948;e=1.5*10^-7;i=0; syms x f(x)=sin(x)./x; while 1 s=0;x=a+h/2; s=s+f(x);x=x+h; if x<b s=s+f(x);x=x+h; else T2=vpa(T1/2+h/2*s); end if abs(T2-T1)<e print(T2) break; else h=h/2; T1=T2; return end end使其为复化梯形公式算法
时间: 2023-09-18 12:10:33 浏览: 87
下面是改进后的代码,使用复化梯形公式算法:
a=0;
b=1;
n=1;
T1=0.920735492403948;
e=1.5e-7;
syms x;
f(x)=sin(x)./x;
while 1
h=(b-a)/n;
s=0;
for i=1:n
x=a+i*h;
s=s+f(x);
end
T2=(T1+h/2*s);
if abs(T2-T1)<e
disp(T2)
break;
else
n=n*2;
T1=T2;
end
end
对比原来的代码,主要改动如下:
1.引入变量n表示分割数目,即将区间[a,b]分割为n个小区间。
2.在for循环中,循环n次,对每个小区间进行计算,得到中点法的近似值。
3.将中点法的近似值T1替换为复化梯形公式的近似值T2。
4.将步长h除以2,即将n乘以2,进行迭代计算。
相关问题
hbin = 10 for i in tqdm(data['WindNumber'].unique()): col = 'WindSpeed' cond = (data.WindNumber==i) & (data.label == 0) temp_df = data[cond] h_bins = pd.cut(temp_df.Power, np.arange(-1000, 3000, hbin)) temp_df['hbins'] = h_bins groups = [] for index,temp in temp_df.groupby("hbins"): if temp.shape[0]==0: continue iqr = temp[col].quantile(0.75) - temp[col].quantile(0.25) t1 = temp[col].quantile(0.25) - 1.5 * iqr t2 = temp[col].quantile(0.75) + 1.5 * iqr temp = temp[((temp[col]<t1) | (temp[col]>t2))] groups.append(temp) groups = pd.concat(groups).reset_index(drop = True) cond = (data.WindNumber==i) & (data.Time.isin(groups.Time)) data.loc[cond,'label'] = 3 cond = (data.WindNumber==1) vis_result_2D(data[cond],data[cond]['label'])详细解释
这段代码是对数据集中的每个风机进行处理,首先根据风机编号筛选数据,然后根据功率值将数据分为若干个区间,每个区间内的数据按照风速进行离群值检测,如果检测出离群值,则将对应的数据标记为异常(标签为3),最后将标记为异常的数据可视化展示。
具体解释如下:
1. `hbin = 10`:将功率值分为10个区间,每个区间的宽度为1000W。
2. `for i in tqdm(data['WindNumber'].unique()):`:循环遍历数据集中的所有风机。
3. `cond = (data.WindNumber==i) & (data.label == 0)`:根据风机编号和标签筛选数据,其中标签为0表示数据未被标记为异常。
4. `temp_df = data[cond]`:将筛选后的数据赋值给临时变量temp_df。
5. `h_bins = pd.cut(temp_df.Power, np.arange(-1000, 3000, hbin))`:根据功率值将数据分为若干个区间,并将区间标签赋值给 h_bins 变量。
6. `temp_df['hbins'] = h_bins`:将区间标签添加到数据集中。
7. `groups = []`:创建一个空列表,用于存储检测出的离群值数据。
8. `for index,temp in temp_df.groupby("hbins"):`:根据区间标签将数据分组,循环遍历每个区间。
9. `if temp.shape[0]==0:`:如果区间内没有数据,则跳过本次循环。
10. `iqr = temp[col].quantile(0.75) - temp[col].quantile(0.25)`:计算数据的四分位距。
11. `t1 = temp[col].quantile(0.25) - 1.5 * iqr`:计算下界阈值。
12. `t2 = temp[col].quantile(0.75) + 1.5 * iqr`:计算上界阈值。
13. `temp = temp[((temp[col]<t1) | (temp[col]>t2))]`:筛选出离群值数据,即风速小于下界阈值或大于上界阈值的数据。
14. `groups.append(temp)`:将检测出的离群值数据添加到 groups 列表中。
15. `groups = pd.concat(groups).reset_index(drop = True)`:将 groups 列表中的数据合并为一个数据集。
16. `cond = (data.WindNumber==i) & (data.Time.isin(groups.Time))`:根据风机编号和时间标记异常数据。
17. `data.loc[cond,'label'] = 3`:将标记为异常的数据标签设置为3。
18. `cond = (data.WindNumber==1)`:根据风机编号筛选数据。
19. `vis_result_2D(data[cond],data[cond]['label'])`:将标记后的数据可视化展示,其中 vis_result_2D 是一个自定义的函数。
将以下代码改为matlab的语言:syms a b; c=[a b]'; A=[79627.43,88131.12,95186.11,106762.43,112557.43,112350.00,113801.00, 104818.00,110600.00,120477.00]; B=cumsum(A); n=length(A); for i=1:(n-1) C(i)=(B(i)+B(i+1))/2; end D=A;D(1)=[]; D=D'; E=[-C;ones(1,n-1)]; c=inv(E*E')*E*D; c=c'; a=c(1);b=c(2); F=[];F(1)=A(1); for i=2:(n+5) F(i)=(A(1)-b/a)/exp(a*(i-1))+b/a; end G=[];G(1)=A(1); for i=2:(n+5) G(i)=F(i)-F(i-1); end t1=1:10; t2=1:15; G h=plot(t1,A,'o',t2,G,'-'); set(h,'LineWidth',3);
syms a b;
c=[a b]';
A=[79627.43,88131.12,95186.11,106762.43,112557.43,112350.00,113801.00, 104818.00,110600.00,120477.00];
B=cumsum(A);
n=length(A);
for i=1:(n-1)
C(i)=(B(i)+B(i+1))/2;
end
D=A;
D(1)=[];
D=D';
E=[-C;ones(1,n-1)];
c=inv(E*E')*E*D;
c=c';
a=c(1);b=c(2);
F=[];
F(1)=A(1);
for i=2:(n+5)
F(i)=(A(1)-b/a)/exp(a*(i-1))+b/a;
end
G=[];
G(1)=A(1);
for i=2:(n+5)
G(i)=F(i)-F(i-1);
end
t1=1:10;
t2=1:15;
h=plot(t1,A,'o',t2,G,'-');
set(h,'LineWidth',3);
阅读全文