声振论坛

 找回密码
 我要加入

QQ登录

只需一步,快速开始

查看: 1406|回复: 0

[综合讨论] 用fmincon求相邻点距离不变问题超限

[复制链接]
发表于 2014-4-10 18:25 | 显示全部楼层 |阅读模式

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

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

x
本帖最后由 牛小贱 于 2014-4-10 18:39 编辑

问题描述:
在一个长(x向)0.2mm,宽(y向)0.4mm,高(z向)0.1mm的方体内有N个有序已知点(顺序相连可成一条曲线,可设N点坐标为一个N*3的矩阵X),要求除第一个点不动外,对其余N-1个点进行微调((N-1)*3的矩阵xd),使调整后的相邻点间距离为一定值lx,当然调整后的点不能出界。
fmincon进行最优化,但第一个点也被调整了,且不知是提取的数学问题不对还是程序有问题,得到的微调值很大,超限很多,请高手帮忙指点一下,先谢了。
数学问题:
目标函数:Σx(j)^2   x(j)为微调xd的各分量,j=1, 2……3N-3
约束条件:|X(i+1)+xd(i+1)-X(i)- xd(i)|=lx,    i=1, 2…..N
0<X(i,1)+xd(i,1)<0.2, 0<X(i,2)+xd(i,2)<0.4,0<X(i,3)+xd(i,3)<0.1
假设N=20,已知这20个点构成的坐标矩阵X, lx=1e-5;
程序:
  1. format long e   %因所给量很小,为避免舍入误差太大,用long e
  2. X=[1.000000000000000e-005                         0    5.000000000000000e-005
  3.     1.084212146258562e-004    1.374317617706318e-004    3.353201365358312e-005
  4.    1.743335211270057e-004   3.211374518483659e-004   3.183321977317479e-005
  5.    1.981769753590540e-004  -9.021471158199482e-004   4.409472569732121e-005
  6.    3.607778176281379e-005  -6.450571057147658e-004   5.362661879261433e-005
  7.    1.992600048898033e-005   1.115269281749919e-004   5.289233224855776e-005
  8.    1.830891687691136e-005   3.502329268561499e-004   5.041792172515980e-005
  9.    1.724894446411056e-005   2.082276257844818e-004   4.989077981552538e-005
  10.    1.442138466276086e-005   2.856263517830034e-004   4.968606079655937e-005
  11.    1.599026055319464e-005   2.196338795975163e-004   4.906672591045520e-005
  12.    2.649812284492222e-005   3.016661141552573e-004   4.797951158729175e-005
  13.    3.501905903500460e-005   1.947676718678416e-004    4.656622210738030e-005
  14.    3.715026059364602e-005   2.017765869090804e-004   4.685650795284104e-005
  15.    5.924772262711515e-005   1.078990962840916e-004   4.935954495744724e-005
  16.    8.305504484629322e-005   1.320501850534990e-004   5.104872456086212e-005
  17.    7.104610355191896e-005   1.015464768267407e-004   5.116284821760612e-005
  18.    9.016579851488848e-005   3.569370145761626e-004   5.109457596677840e-005
  19.    1.786477300524624e-004   1.542941709813852e-004   5.048413376542842e-005
  20.    2.247149253348480e-004   1.823588732474208e-004   4.997876808826507e-005
  21. lx=1e-5;
  22. m1=size(X,1);m2=size(X,2);        %m1,m2分别为各节点坐标矩阵的行数、列数
  23. %-------------数据传送到平台
  24. assignin('base','Xz',X);
  25. assignin('base','m1z',m1);
  26. assignin('base','m2z',m2);
  27. assignin('base','lxz',lx);
  28. %-------------程序主体
  29. x0=zeros(1,m1*m2);                     %x进行波动的初始值
  30. %%%求线性上、下限
  31. %%x,y,z边界-X(i,1)<=x(i,1)<=2e-4-rx-X(i,1);rx-X(i,2)<=x(i,2)<=4e-4-rx-X(i,2); rx-X(i,3)<=x(i,3)<=1e-4-rx-X(3,1)
  32. lb=zeros(1,m1*m2);ub=zeros(1,m1*m2);
  33. for i=2:m1
  34.    mm=3*(i-1);
  35.    lb(mm+1)=rx-X(i,1); ub(mm+1)=2e-4-rx-X(i,1);
  36.    lb(mm+2)=rx-X(i,2); ub(mm+2)=4e-4-rx-X(i,2);
  37.    lb(mm+3)=rx-X(i,3); ub(mm+3)=1e-4-rx-X(i,3);
  38. end
  39.    
  40. %options=optimset('LargeScale','off','display','off');
  41. options=optimset('LargeScale','off','display','off','Algorithm','active-set');
  42. [x,fval]=fmincon(@myfun,x0,[],[],[],[],lb,ub,@mycon,options);
  43. for i=1:m1
  44.    for j=1:m2
  45.        xd(i,j)=x((i-1)*m2+j);
  46.    end
  47. end
  48. %%%%%%%%%%防止调整的微小位移过大的限制程序
  49. for i=1:m1
  50.    ifsum(xd(i,:).^2)>1e-12   %既要求微小位移不能大于1e-6mm
  51.      disp('调整的微小距离大于0.001mm,不对')
  52.    end
  53. end
  54. xd,   %显示的xd值很怪,超限很多
  55. Jdxyz=X+xd;                     %满足各项条件的节点坐标
  56. return
  57. function f=myfun(x)
  58. f=sum(x.^2);
  59. return
  60. function [c,ceq]=mycon(x)
  61. %----------从平台引入数据
  62. X=evalin('base','Xz');
  63. m1=evalin('base','m1z');
  64. m2=evalin('base','m2z');
  65. lx=evalin('base','lxz');
  66. %----------节点限制条件的数学表示
  67. for i=1:m1-1
  68.    if m2==3
  69.        ceq(i)=((X(i+1,1)+x(i*m2+1))-(X(i,1)+x((i-1)*m2+1)))^2+...
  70.            ((X(i+1,2)+x(i*m2+2))-(X(i,2)+x((i-1)*m2+2)))^2+...
  71.            ((X(i+1,3)+x(i*m2+3))-(X(i,3)+x((i-1)*m2+3)))^2-lx^2;
  72.      
  73.    end
  74. end
  75. %防止调整的微小位移过大,限定矢量xd的模小于jx,但这步好像没起作用
  76. jx=5e-6;
  77. for i=1:m1
  78.    if m2==3
  79.        c(i)=x((i-1)*m2+1)^2+x((i-1)*m2+2)^2+x((i-1)*m2+3)^2-jx^2;
  80.    end
  81. end
  82. return
复制代码


回复
分享到:

使用道具 举报

您需要登录后才可以回帖 登录 | 我要加入

本版积分规则

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

GMT+8, 2024-11-13 10:59 , Processed in 0.074565 second(s), 17 queries , Gzip On.

Powered by Discuz! X3.4

Copyright © 2001-2021, Tencent Cloud.

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