声振论坛

 找回密码
 我要加入

QQ登录

只需一步,快速开始

查看: 1395|回复: 11

[非线性振动] 打耙法——活动耙点!!!!!!

[复制链接]
发表于 2013-4-26 16:07 | 显示全部楼层 |阅读模式

马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。

您需要 登录 才可以下载或查看,没有账号?我要加入

x
对于活动耙点,并且具有两个初始条件(ICs)的打耙法,有大神们做过吗?若采用牛顿迭代,其雅克比矩阵中的未知量该采用哪个初始条件计算的结果呢?

回复
分享到:

使用道具 举报

 楼主| 发表于 2013-4-26 16:07 | 显示全部楼层
急求........
 楼主| 发表于 2013-4-26 17:09 | 显示全部楼层

以下是三自由度两点边值问题的拟牛顿迭代打耙法,这个算法是传统的,因为只具有一个初始条件;对于两个初始条件的两点边值问题,大家遇到过吗?


%%%%%%%%%%%%%%%%%%%%主程序
clear all                  %采用Fx(差值)作为精度依据
clc
tspan=0:0.1:2;
x0=[-0.538 0.288 0.49883]';  %定义初始耙点
x=x0;
norm(x);
y0=[1.076 0 0 x']';          %定义初始条件
A0=eye(3);

e1=0.0000001;
[t,y]=ode45('shooting_li2',tspan,y0);
Fx0=[y(end,1)-0;y(end,2)-0.576;y(end,3)-0.997661];
%s=-inv(A)*Fx;              % 由此可以看出,非点乘,点乘则要求矩阵的行列是一样的
c0=norm(Fx0,1);
c=c0;
%Y=zeros();
    A=A0;                    %定义初始迭代矩阵(单位矩阵)
    Fx=Fx0;                  %定义初始精度条件
while c>e1
    A;
    s=-inv(A)*Fx;
    x=x+s;                   %耙点迭代
    %Y=[Y' c]';
    y0=[1.076 0 0 x']';
    [t,y]=ode45('shooting_li2',tspan,y0);                   %解初始问题
    Fx2=[y(end,1)-0;y(end,2)-0.576;y(end,3)-0.997661];      %精度条件
    y(:,1);
    c=norm(Fx2,1);
    y=Fx2-Fx;
    A=A+(y-A*s)*s'/(s'*s);   %矩阵迭代
    Fx=Fx2;                  
    %C=[t y(end,1) y(end,2) y(end,3)];
while c<=e1
       [t,y]=ode45('shooting_li2',tspan,y0);               %解初始问题
        y1=y(:,1);
        %y2=y(2)
        %y3=y(3)
        C=[t y(:,1) y(:,2) y(:,3)]                         %输入结果
        break                                              %防止不断循环
end
end
%%%%%%%%%%%%%%%%%%%%%%%三自由度、6维函数
function hs=shooting_li2(t,y)
hs(1)=y(4); %y1'=y3
hs(2)=y(5); %y2'=y4
hs(3)=y(6);
hs(4)=-y(1)./(y(4).^2+y(5).^2+y(6).^2).^(3/2);  %y3'=
hs(5)=-y(2)./(y(4).^2+y(5).^2+y(6).^2).^(3/2);  %y4'=
hs(6)=-y(3)./(y(4).^2+y(5).^2+y(6).^2).^(3/2);
hs=hs(:);
end
 楼主| 发表于 2013-4-27 10:41 | 显示全部楼层
clear %%%%%%%%%%%%%%%%%%%%主程序
clc
x0=0:1:10;
怎么没有回复呢?好吧,再给大家分享一下secant method of shooting 。目前想解决多ICs的打耙法,或耙点活动的打耙法,,,亲,有支持的不

y01=8;
y0=[0 y01];
[t,y]=ode45('shooting_li2_secantmethod',x0',y0);
y10=y(end,1);

y00=0;
%y0=[0 7];
%[t,y]=ode45('shooting_li2_secantmethod',x0',y0);
%y20=y(end,1);


c1=y01;
c0=y00;

y1=y10;
y2=0;
e1=0.000001;
yy=y10;

while abs(yy)>e1
c2=c1-y1*(c1-c0)/(y1-y2)
y0=[0 c2];
[t,y]=ode45('shooting_li2_secantmethod',x0',y0);
y30=y(end,1)
c0=c1;
c1=c2;
y2=y1;
y1=y30;
yy=y30;

while abs(yy)<=e1
   
break
end
end


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%一自由度二维函数
function hs=shooting_li2_secantmethod(t,y)
hs(1)=y(2); %y1'=y3
hs(2)=8-1/4*y(1); %y2'=y4
hs=hs(:);
end
 楼主| 发表于 2013-4-27 15:41 | 显示全部楼层
双初初始条件,尘埃耙点的打耙法,讨论讨论呗......
 楼主| 发表于 2013-4-27 15:43 | 显示全部楼层
following is two DOFs secant method of shooting ......讨论讨论双初始条件,且具有活动耙点的打耙法呗.......

clear %%%%%%%%%%%%%%%%%%%%主程序
clc
x0=0:1:10;
y01=[8 1];
y0=[0 0 y01];
[t,y]=ode45('shooting_li2_secantmethod_two_dof',x0',y0);
y10=[y(end,1) y(end,2)]


y00=[1 0]
%y0=[0 7];
%[t,y]=ode45('shooting_li2_secantmethod',x0',y0);
%y20=y(end,1);


c1=y01
c0=y00

y1=y10
y0=[0 0];
e1=0.000001;
yy=y10;
c=norm(yy,1)
while c>e1
c2=c1-y1*((c1-c0)/(y1-y0))
y0=[0 0 c2];
[t,y]=ode45('shooting_li2_secantmethod_two_dof',x0',y0);
y30=[y(end,1) y(end,2)]
c0=c1;
c1=c2;
y0=y1;
y1=y30;
yy=y30;
c=norm(yy,1)

while abs(yy)<=e1
   
break
end
end


%%%%%%%%%%%%%%%%%%%%%%%
function hs=shooting_li2_secantmethod_two_dof(t,y)
hs(1)=y(3);                 %y1'=y3  
hs(2)=y(4);
hs(3)=-y(3)-(y(1)-y(2)).^3; %y2'=y4
hs(4)=-(y(2)-y(1)).^3;
hs=hs(:);
end

 楼主| 发表于 2013-4-27 15:53 | 显示全部楼层
将这种双初始条件(ICs),且具有活动耙点的例子,拿出来大家讨论讨论,我也不会啊。.....见图1、图2


图1

图1

图2

图2
 楼主| 发表于 2013-4-28 09:28 | 显示全部楼层
期待......
 楼主| 发表于 2013-4-28 10:26 | 显示全部楼层
打耙法中为什么使用两个初始条件(我的回答):传统打耙法只是一个求根工具,根据不断调整初始速度,以求射中耙点,最终可求得一个根(初始速度);但从理论情况来看,满足耙点的初始速度有很多,那么第二个初始条件就是决定了众多初始速度中的某一个,那么符合两个初始条件的特定解就解出来了。这是单从打耙法的思想解释为什么采用两个初始条件。

盟友们给点见解......
 楼主| 发表于 2013-4-28 10:27 | 显示全部楼层
若这种想法正确。 将这种想法如何用程序的语言进行描述呢?........
 楼主| 发表于 2013-4-28 22:12 | 显示全部楼层
终于找到了这种方法的专有名词,称为"multiple shooting method"。。。各位老师,请指教....
 楼主| 发表于 2013-4-29 11:32 | 显示全部楼层
但是这个“multiple shooting method" 范围还是太大,对我的问题针对性效果较差。目前可解决我的问题的方法的名称,好像是“拟合点打耙法”。。。。痛哭 也无泪啊,,,,怎么这么难呢
您需要登录后才可以回帖 登录 | 我要加入

本版积分规则

QQ|小黑屋|Archiver|手机版|联系我们|声振论坛

GMT+8, 2024-5-5 02:38 , Processed in 0.103667 second(s), 21 queries , Gzip On.

Powered by Discuz! X3.4

Copyright © 2001-2021, Tencent Cloud.

快速回复 返回顶部 返回列表