| 查看: 399 | 回复: 0 | |||
[交流]
请大家帮我看下我的MATLAB卷积分问题出在哪里？谢谢！
|
|
用MATLAB自带卷积分函数计算弹性半空间在一三角脉冲荷载下自由表面的竖向位移，与文献正确结果相差10的5次方倍，但自己编写卷积分代码，反而得到文献结果。请大家帮我看下问题出在哪？谢谢！ 另外，为何MATLAB自带卷积分函数结果矩阵的维数会是被卷积两矩阵维数之和减1？ 1. 问题描述，如下图1所示：  2. 格林函数，如下图2所示：  3. 参考文献正确解，如下图3所示：  4. 本人直接用自编卷积分MATLAB代码算的出解（如下图4）及相应代码：  clear; clc; %loading history dt=0.01; ti=0.0008; te=12; t=ti:dt:te; pt=(10*t./1.5).*(t>=0 & t<=1.5)+(20-10*t./1.5).*(t>1.5 & t<=3)+0.*(t>3); m=length(pt); %load distribution in space rp=0.1; ri=0; rc=0.1; dr=0.001; r=ri:dr:rc; pr=1*(r<=rp)+0*(r>rp); n=length(pr); %load function with respect to t and r p=pr.'*pt; %green's function G=1; cs=1; for i=1:1:m for j=1:1:n u(j,i)=heaviside(cs*t(i)-r(j))/(pi*G*sqrt(t(i)^2-(r(j)/cs)^2)); end end %convolution and response of displacement for i=1:1:m for j=1:1:n v(j,i)=0; for k=1:1:i for g=1:1:j v(j,i)=v(j,i)+p(g,k)*u(j-g+1,i-k+1)*dr*dt; end; end; end; end; %plot the response history plot(t,v(n,  );xlabel('t/s'); ylabel('v/m'); grid on; title('response of point A or C'); 5. 直接用MATLAB自带卷积分算的出解（如下图5）及相应代码：  clear; clc; %loading history dt=0.01; ti=0.0008; te=12; t=ti:dt:te; pt=(10*t./1.5).*(t>=0 & t<=1.5)+(20-10*t./1.5).*(t>1.5 & t<=3)+0.*(t>3); m=length(pt); %load distribution in space rp=0.1; ri=0; rc=0.1; dr=0.001; r=ri:dr:rc; pr=1*(r<=rp)+0*(r>rp); n=length(pr); %load function with respect to t and r p=pr.'*pt; %green function G=1; cs=1; for i=1:1:m for j=1:1:n u(j,i)=heaviside(cs*t(i)-r(j))/(pi*G*sqrt(t(i)^2-(r(j)/cs)^2)); end end %convolution and response of displacement v=conv2(u,p); %plot the response history rr=r(n); [tpu,rpu]=meshgrid(ti:dt:2*te-dt,ri:dr:2*rc); [X,Y,Z]=meshgrid(linspace(min(tpu( ),max(tpu()),linspace(min(rpu(),max(rpu()),linspace(min(v(),max(v()));V=Y; h=contourslice(X,Y,Z,V,tpu,rpu,v,[0 0]+rr); set(h,'edgecolor','k'); contourslice(X,Y,Z,V,tpu,rpu,v,[0 0]+rr); xlabel('t/s'); ylabel('r/m'); zlabel('v/m'); axis([0 12 0 0.1 0 12e4]); view(0,0); grid on; title('displacement response history of point A or C'); |
回复此楼» 猜你喜欢
售SCI文章，我:8O.5.5.1O.54,科目全,可十急
已经有4人回复
-
售SCI一区文章，我:8.O.55.1.O.54,科目齐全,可伽急
已经有5人回复
-
我就是申请一个面上项目而已，这评审意见是按照杰青的条件评的吧？
已经有4人回复
-
售SCI一区文章，我:8.O.55.1.O.54,科目齐全,可伽急
已经有4人回复
-
售SCI-T0P文章，我:8O.5.5.1.O.54,科目齐全,可+急
已经有4人回复
-
有多少人是今天查系统知道结果的？
已经有18人回复
-
国自然评审意见
已经有15人回复
-
为什么资助数各大高校都创新高，自己申请怎么就这么难
已经有13人回复
-
国自然面上复盘~欢迎讨论
已经有14人回复
-
系统查不到
已经有11人回复
30