物流公司网站建设,9377传奇手游盒子,wordpress 免费 博客主题,重庆网站推广效果问题三
最小螺距单目标优化模型的建立
问题二考虑了在螺距固定的条件下计算舞龙队盘入的终止时间#xff0c;问题三在第二问的基础提出了改变螺距的要求#xff0c;即求解在螺距最小为多少时#xff0c;龙头前把手能够沿着相应的螺线盘入到调头空间的边界。故可将其转换为…问题三
最小螺距单目标优化模型的建立
问题二考虑了在螺距固定的条件下计算舞龙队盘入的终止时间问题三在第二问的基础提出了改变螺距的要求即求解在螺距最小为多少时龙头前把手能够沿着相应的螺线盘入到调头空间的边界。故可将其转换为一个单目标优化问题目标函数为求解最小的螺距将板凳龙模型中具体物理要求与龙头前把手能够沿着相应的螺线盘入到调头空间的边界作为约束条件由于该目标函数非线性且较为复杂故本文采取了带罚函数的双重模拟退火算法与变步长搜索算法结合进行求解。
目标函数与约束条件
由于此时螺距的取值受到约束条件的限制故将螺距 a a a的表达修改为变量 A A A可写出目标函数的表达为 min A \min A minA 其中约束条件如下如下
1. 根据题意龙头能够进入调头空间的边界即此时结束位置的极径应该小于4.5m此时螺线方程为 r A 2 π ⋅ θ r\frac{A}{2\pi}\cdot \theta r2πA⋅θ故可列约束条件为 A 2 π ⋅ θ ≤ 4.5 ( m ) \frac{A}{2\pi}\cdot \theta \leq4.5(m) 2πA⋅θ≤4.5(m) 2.根据题意在龙头能够进入调头空间前板凳之间不能相撞即此时的决策变量应该满足的约束条件为 ∑ i 1 223 ∑ j 1 , 223 δ i j 0 \sum_{i1}^{223}\sum_{j1,}^{223}\delta_{ij}0 i1∑223j1,∑223δij0 3.板凳的宽度为15cm螺距的宽度能够支持两个板凳并肩移动即此时螺距应该满足的约束条件为 A ≥ 30 ( c m ) A \geq30(cm) A≥30(cm) 4.调头空间意味着舞龙队的总长度不能长于螺线的总弧长此时应该满足的约束条件为 A ∫ 0 32 π 1 θ 2 d θ ≥ 220 × 222 341 A\int_{0}^{32\pi}\sqrt{1\theta^2}d\theta \geq220\times222341 A∫032π1θ2 dθ≥220×222341 5.此时碰撞时间的求解模型基于问题一中对位置的迭代求解模型给出问题一中的迭代方程为 { ( x i 1 − x i ) 2 ( y i 1 − y i ) 2 16 5 2 x i 1 55 2 π ⋅ θ i 1 ⋅ cos ( θ i 1 ) y i 1 55 2 π ⋅ θ i 1 ⋅ sin ( θ i 1 ) \begin{cases} (x_{i1}-x_i)^2(y_{i1}-y_i)^2165^2\\ x_{i1}\frac {55}{2\pi}\cdot\theta_{i1}\cdot\cos(\theta_{i1})\\y_{i1}\frac {55}{2\pi}\cdot\theta_{i1}\cdot\sin(\theta_{i1})\\ \end{cases} ⎩ ⎨ ⎧(xi1−xi)2(yi1−yi)21652xi12π55⋅θi1⋅cos(θi1)yi12π55⋅θi1⋅sin(θi1) 此处 ( x i 1 , y i 1 ) (x_{i1},y_{i1}) (xi1,yi1)为下一龙身前孔的坐标 ( x i , y i ) (x_{i},y_{i}) (xi,yi)为上一龙身前孔的坐标 θ i 1 \theta_{i1} θi1为下一龙身前孔此时的弧度。
综上所述给出问题三最小螺距单目标优化模型如下 min A s.t. { A 2 π ⋅ θ ≤ 4.5 ( m ) ∑ i 1 223 ∑ j 1 , 223 δ i j 0 A ≥ 30 ( c m ) A ∫ 0 32 π 1 θ 2 d θ ≥ 220 × 222 341 ( x i 1 − x i ) 2 ( y i 1 − y i ) 2 16 5 2 x i 1 55 2 π ⋅ θ i 1 ⋅ cos ( θ i 1 ) y i 1 55 2 π ⋅ θ i 1 ⋅ sin ( θ i 1 ) \min A\\ \text{s.t.}\begin{cases} \frac{A}{2\pi}\cdot \theta \leq4.5(m)\\ \sum_{i1}^{223}\sum_{j1,}^{223}\delta_{ij}0\\ A \geq30(cm)\\ A\int_{0}^{32\pi}\sqrt{1\theta^2}d\theta \geq220\times222341\\ (x_{i1}-x_i)^2(y_{i1}-y_i)^2165^2\\ x_{i1}\frac {55}{2\pi}\cdot\theta_{i1}\cdot\cos(\theta_{i1})\\y_{i1}\frac {55}{2\pi}\cdot\theta_{i1}\cdot\sin(\theta_{i1})\\ \end{cases} minAs.t.⎩ ⎨ ⎧2πA⋅θ≤4.5(m)∑i1223∑j1,223δij0A≥30(cm)A∫032π1θ2 dθ≥220×222341(xi1−xi)2(yi1−yi)21652xi12π55⋅θi1⋅cos(θi1)yi12π55⋅θi1⋅sin(θi1)
模型的求解
本题的方程比较复杂且非线性但是最终结果为一维且搜索空间较小于是本文在求解的时候采取了带罚函数的模拟退火算法与定步长遍历算法结合的办法启发式算法可以很好的在解空间内探查到最优解而变步长遍历算法则可以在解空间较小的情况下遍历空间得到整个函数的分布情况进而查找到最优解。下面是两种算法的求解过程
带罚函数的双重模拟退火算法
双重模拟退火算法(Dual Annealing)是一种启发式算法适合用于求解复杂的非线性优化问题。在该问题中由于目标函数和约束条件都比较复杂使用双重模拟退火算法可以有效搜索到接近全局最优解的结果。为了处理约束条件本文引入罚函数将约束条件违背的程度引入到目标函数中确保求解过程中尽量满足约束条件。
Step1 设定初始解初始温度和罚函数
设定初始螺距 A 0 30 c m A_030cm A030cm初始温度 T 0 100 T_0100 T0100设定罚函数因子 λ 1 \lambda_1 λ1 λ 2 \lambda_2 λ2与 λ 3 \lambda_3 λ3其大小为某远大于目标值的数定义罚函数为 P ( A ) λ 1 ⋅ max ( 0 , A 2 π ⋅ θ − 4.5 ) λ 2 ⋅ ( ∑ i 1 223 ∑ j 1 , 223 δ i j ) λ 3 ⋅ max ( 0 , 220 × 222 341 − A ∫ 0 32 π 1 θ 2 d θ ) P(A) \lambda_1 \cdot \max(0, \frac{A}{2\pi}\cdot \theta - 4.5) \lambda_2 \cdot \left(\sum_{i1}^{223}\sum_{j1,}^{223}\delta_{ij}\right) \lambda_3 \cdot \max(0, 220\times 222 341 - A\int_{0}^{32\pi}\sqrt{1\theta^2}d\theta) P(A)λ1⋅max(0,2πA⋅θ−4.5)λ2⋅(i1∑223j1,∑223δij)λ3⋅max(0,220×222341−A∫032π1θ2 dθ) 由此实现对约束条件的处理
Step2 全局搜索
在当前解 x i x_i xi的邻域内生成新解 x i 1 x_{i1} xi1计算新解的函数值 f ( x i 1 ) f(x_{i1}) f(xi1)如果 f ( x i 1 ) f ( x i ) f(x_{i1})f(x_i) f(xi1)f(xi)则接受新解作为当前解。若不满足该条件则以概率 P exp ( f ( x i ) − f ( x i 1 ) T ) P\exp\left(\frac{f(x_i)-f(x_{i1})}T\right) Pexp(Tf(xi)−f(xi1))接受新解随着温度的降低接受劣解的概率逐渐减小。
Step3 局部搜索
在接受新解后使用局部优化算法如最速下降法或梯度下降法对当前解进行局部优化以快速收敛到邻域内的局部最优解。
Step4 温度下降
更新温度 T i 1 α ⋅ T i T_{i1}α⋅T_i Ti1α⋅Ti,逐步降低温度。
Step5 搜索结束输出结果
当最终结果收敛时将搜索得到的结果输出获得此时的螺距长度。
变步长遍历算法
变步长求解是一种结合了大步长的快速搜索和小步长的精细搜索的优化策略。在开始阶段使用较大的步长快速缩小解的范围后面逐步缩小步长通过类似二分查找的方式精确找到最优解。
Step1 设定初始解和初始步长
设定初始螺距 A 0 30 c m A_030cm A030cm最大螺距为在通过观察模拟退火算法求解过程基础上确定的边界解 A m a x 50 c m A_{max}50cm Amax50cm初始步长为 Δ A 2 c m \Delta A2cm ΔA2cm
Step2 大步长快速搜索
从初始螺距 A 0 A_0 A0开始搜索以 Δ A 2 c m \Delta A2cm ΔA2cm为步长开始遍历快速搜索解空间直达搜索到 A m a x A_{max} Amax计算每个搜索点的目标函数值若当前螺距下的目标函数值优于前面搜索到的最优值则更新最优值重复上述过程直到搜索完初始步长的范围。
Step3 二分精细搜索
经过大步长搜索找到最优螺距 A b e s t A_{best} Abest在其左右邻域 [ A b e s t − Δ A , A b e s t Δ A ] [A_{best}-\Delta A,A_{best}\Delta A] [Abest−ΔA,AbestΔA]内进行搜索此时设定步长为 Δ A 2 \frac{\Delta A}{2} 2ΔA。重复Step2与Step3的步骤。
Step4 导出搜索结果
当搜索步长缩小到 Δ A 1024 \frac{\Delta A}{1024} 1024ΔA时其值以达到误差容限内输出此时的目标函数值 A b e s t A_{best} Abest。
问题四
模型的建立
问题四设置盘入螺线的螺距为1.7m盘入螺线与盘出螺线为中心对称关系并且舞龙队需要在调头空间内实现调头其中调头位置在盘入螺线上调出位置在盘出螺线即现在确定了调头空间的调入点与调出点以及调入圆弧与调出圆弧之间的数量关系在未知圆弧进入点与圆弧出去点的情况下对圆弧曲线的长度进行优化求解出最短的调头曲线。为了求解这个问题首先需要求解调头螺线的数学方程进而通过几何关系建立起调头圆弧的数学模型根据圆弧的弧度和调头的位置确定整个调头曲线长度的表达式作为模型的目标函数然后通过数学模型以及实际板凳的物理模型给出目标函数的约束条件
调头弧线模型的建立
首先已知盘入螺线的螺距为 1.7 m,盘出螺线与盘入螺线关于螺线中心呈中心对称故可先求解出盘入螺线 C 1 C_1 C1与盘出螺线 C 2 C_2 C2的表达式分别为 { C 1 : r 170 2 π ⋅ θ C 2 : r − 170 2 π ⋅ θ \begin{cases} C1:r\frac{170}{2\pi}\cdot \theta\\ C2:r-\frac{170}{2\pi}\cdot \theta\\ \end{cases} {C1:r2π170⋅θC2:r−2π170⋅θ 由于目前只可推出进入调头空间的点与离开调头空间的点由于两点分别位于中心对称的两螺线与调头区域圆的交点上所以两点坐标存在中心对称关系,可设进入点为 A ( x A , y A ) A(x_A,y_A) A(xA,yA)离开点为 B ( x B , y B ) B(x_B,y_B) B(xB,yB)则两点满足的方程为 { x A x B 0 y A y B 0 x 2 y 2 45 0 2 \begin{cases} x_Ax_B0\\ y_Ay_B0\\ x^2y^2450^2 \end{cases} ⎩ ⎨ ⎧xAxB0yAyB0x2y24502 由于开始调头的点在盘入螺线和前段圆弧上且经过该点后龙头前把手运动方程从螺线方程变为圆方程结束调头的点则相反故可设开始调头的点为点 C C C坐标为 ( x C , y C ) (x_C,y_C) (xC,yC)结束调头的点为点 D D D坐标为 ( x D , y D ) (x_D,y_D) (xD,yD)。其示意图如下图所示
故可表示调头弧线的总长度 L L L可表示为 L A C ⌢ C D ⌢ D B ⌢ L\overset{\LARGE{\frown}}{AC}\overset{\LARGE{\frown}}{CD}\overset{\LARGE{\frown}}{DB} LAC⌢CD⌢DB⌢ 为确定两圆的圆心位置需计算出点 C C C处螺线的导数将其带入曲率半径的公式即可计算出前段圆弧的半径进而确定后段圆弧半径即可确定两圆心的位置以及调头弧线的长度。
由于点 C C C位于螺线上故可先求解螺线的 d x d θ \frac{dx}{d\theta} dθdx与 d y d θ \frac{dy}{d\theta} dθdy再通过链式法则即可求解得到 d y d x \frac{dy}{dx} dxdy计算过程如下 d x d θ 170 2 π ⋅ ( cos ( θ ) − θ ⋅ sin ( θ ) ) d y d θ 170 2 π ⋅ ( sin ( θ ) θ ⋅ cos ( θ ) ) d y d x d y d θ d x d θ sin ( θ ) θ ⋅ cos ( θ ) cos ( θ ) − θ ⋅ sin ( θ ) \frac{dx}{d\theta}\frac{170}{2\pi}\cdot(\cos(\theta)-\theta\cdot\sin(\theta))\\ \frac{dy}{d\theta}\frac{170}{2\pi}\cdot(\sin(\theta)\theta\cdot\cos(\theta))\\ \frac{dy}{dx}\frac{\frac{dy}{d\theta}}{\frac{dx}{d\theta}}\frac{\sin(\theta)\theta\cdot\cos(\theta)}{\cos(\theta)-\theta\cdot\sin(\theta)} dθdx2π170⋅(cos(θ)−θ⋅sin(θ))dθdy2π170⋅(sin(θ)θ⋅cos(θ))dxdydθdxdθdycos(θ)−θ⋅sin(θ)sin(θ)θ⋅cos(θ) 故点 C C C处的导数为 y C ′ ( x ) d y C d x C d y C d θ C d x C d θ C sin ( θ C ) θ C ⋅ cos ( θ C ) cos ( θ C ) − θ ⋅ sin ( θ C ) y_C(x)\frac{dy_C}{dx_C}\frac{\frac{dy_C}{d\theta_C}}{\frac{dx_C}{d\theta_C}}\frac{\sin(\theta_C)\theta_C\cdot\cos(\theta_C)}{\cos(\theta_C)-\theta\cdot\sin(\theta_C)} yC′(x)dxCdyCdθCdxCdθCdyCcos(θC)−θ⋅sin(θC)sin(θC)θC⋅cos(θC) 根据互相垂直的直线之间存在 k 1 ⋅ k 2 − 1 k_1\cdot k_2-1 k1⋅k2−1可知 C C C点切线的垂线斜率 k w 1 k_{w1} kw1为 k w 1 d y C d θ C d x C d θ C − cos ( θ C ) − θ ⋅ sin ( θ C ) ) sin ( θ C ) θ C ⋅ cos ( θ C k_{w1}\frac{\frac{dy_C}{d\theta_C}}{\frac{dx_C}{d\theta_C}}-\frac{\cos(\theta_C)-\theta\cdot\sin(\theta_C))}{\sin(\theta_C)\theta_C\cdot\cos(\theta_C} kw1dθCdxCdθCdyC−sin(θC)θC⋅cos(θCcos(θC)−θ⋅sin(θC)) 同理 D D D点切线的垂线斜率 k w 2 k_{w2} kw2为 k w 2 d y D d θ D d x D d θ D − cos ( θ D ) − θ ⋅ sin ( θ D ) ) sin ( θ D ) θ D ⋅ cos ( θ D ) k_{w2}\frac{\frac{dy_D}{d\theta_D}}{\frac{dx_D}{d\theta_D}}-\frac{\cos(\theta_D)-\theta\cdot\sin(\theta_D))}{\sin(\theta_D)\theta_D\cdot\cos(\theta_D)} kw2dθDdxDdθDdyD−sin(θD)θD⋅cos(θD)cos(θD)−θ⋅sin(θD)) 作出切线图后调头空间中存在如图的几何约束示意图如下图所示
其中 r 1 r_1 r1为前段圆的半径 r 2 r_2 r2为后段圆的半径设两圆圆心与坐标分别为 O 1 ( x O 1 , y O 1 ) O_1(x_{O1},y_{O1}) O1(xO1,yO1)与 O 2 ( x O 2 , y O 2 ) O_2(x_{O2},y_{O2}) O2(xO2,yO2)此时可依图列出方程组为 { x O 1 x C 1 k w 1 2 1 ⋅ r 1 y O 1 y C k w 1 k w 1 2 1 ⋅ r 1 x O 2 x D 1 k w 2 2 1 ⋅ r 2 y O 2 y D k w 2 k w 2 2 1 ⋅ r 2 r 1 2 r 2 ( x O 1 − x O 2 ) 2 ( y O 1 − y O 2 ) 2 r 1 r 2 \begin{cases} x_{O1}x_C\frac1{\sqrt{k_{w1}^21}}\cdot r_1\\ y_{O1}y_C\frac{k_{w1}}{\sqrt{k_{w1}^21}}\cdot r_1\\ x_{O2}x_D\frac1{\sqrt{k_{w2}^21}}\cdot r_2\\ y_{O2}y_D\frac{k_{w2}}{\sqrt{k_{w2}^21}}\cdot r_2\\ r_12r_2\\ \sqrt{(x_{O1}-x_{O2})^2(y_{O1}-y_{O2})^2}r_1r_2\\ \end{cases} ⎩ ⎨ ⎧xO1xCkw121 1⋅r1yO1yCkw121 kw1⋅r1xO2xDkw221 1⋅r2yO2yDkw221 kw2⋅r2r12r2(xO1−xO2)2(yO1−yO2)2 r1r2 求解上述方程组即可计算得到两圆圆心的坐标。
得到圆心坐标后则需要再计算出两圆的切点 E E E的坐标由于已知两圆的圆心坐标和圆心到 E E E点的直线的斜率故需要通过简单的比例计算就可得到 E E E点的坐标如下 { x E x O 1 r 1 r 1 r 2 ⋅ ( x O 2 − x O 1 ) y E y O 1 r 1 r 1 r 2 ⋅ ( y O 2 − y O 1 ) \begin{cases} x_Ex_{O1}\frac{r1}{r1r2}\cdot (x_{O2}-x_{O1})\\ y_Ey_{O1}\frac{r1}{r1r2}\cdot (y_{O2}-y_{O1}) \end{cases} {xExO1r1r2r1⋅(xO2−xO1)yEyO1r1r2r1⋅(yO2−yO1) 下面依据 E E E点的坐标注意到前段圆弧均为顺时针方向后段圆弧为逆时针方向在导出公式中将角度的顺序进行调整可导出两圆转过的弧度 ϕ 1 \phi_1 ϕ1与 ϕ 2 \phi_2 ϕ2为 { ϕ 1 arctan ( y C − y O 1 x C − x O 1 ) − arctan ( y E − y O 1 x E − x O 1 ) ϕ 2 arctan ( y D − y O 2 x D − x O 2 ) − arctan ( y E − y O 2 x E − x O 2 ) \begin{cases} \phi_1\arctan(\frac{y_C-y_{O1}}{x_C-x_{O1}})-\arctan (\frac{y_E-y_{O1}}{x_E-x_{O1}})\\ \phi_2\arctan(\frac{y_D-y_{O2}}{x_D-x_{O2}})-\arctan (\frac{y_E-y_{O2}}{x_E-x_{O2}})\\ \end{cases} {ϕ1arctan(xC−xO1yC−yO1)−arctan(xE−xO1yE−yO1)ϕ2arctan(xD−xO2yD−yO2)−arctan(xE−xO2yE−yO2) 即可求解出 L L L的表达式为 L 170 2 π ∫ θ C θ A 1 θ 2 d θ r 1 ⋅ ϕ 1 r 2 ⋅ ϕ 2 170 2 π ∫ θ D θ B 1 θ 2 d θ L\frac{170}{2\pi}\int_{\theta_C}^{\theta_A}\sqrt{1\theta^2}d\thetar_1\cdot\phi_1r_2\cdot\phi_2\frac{170}{2\pi}\int_{\theta_D}^{\theta_B}\sqrt{1\theta^2}d\theta L2π170∫θCθA1θ2 dθr1⋅ϕ1r2⋅ϕ22π170∫θDθB1θ2 dθ
建立最短弦长单目标优化模型
前文已经求解出 L L L的表达式根据题意需要通过调整圆弧保持各部分相切在条件下使得调头曲线变短故可基于此构建最短调头曲线的优化模型其目标函数为最短调头曲线弧长约束条件则基于调头区域内需要保持各部分相切以及圆弧之间的几何关系约束。
目标函数与约束条件
目标函数为最短调头曲线弧长由于前文的推导过程均基于 C C C点与 D D D点的不同坐标即 θ C 与 θ D \theta_C与\theta_D θC与θD故可视 L L L为 θ C \theta_C θC与 θ D \theta_D θD的函数故目标函数数学表达式为 min L ( θ C , θ D ) 170 2 π ∫ θ C θ A 1 θ 2 d θ r 1 ⋅ ϕ 1 r 2 ⋅ ϕ 2 170 2 π ∫ θ D θ B 1 θ 2 d θ \min L(\theta_C,\theta_D)\frac{170}{2\pi}\int_{\theta_C}^{\theta_A}\sqrt{1\theta^2}d\thetar_1\cdot\phi_1r_2\cdot\phi_2\frac{170}{2\pi}\int_{\theta_D}^{\theta_B}\sqrt{1\theta^2}d\theta minL(θC,θD)2π170∫θCθA1θ2 dθr1⋅ϕ1r2⋅ϕ22π170∫θDθB1θ2 dθ 其中 r 1 , r 2 , ϕ 1 , ϕ 2 r_1,r_2,\phi_1,\phi_2 r1,r2,ϕ1,ϕ2均可由 θ C \theta_C θC与 θ D \theta_D θD确定 θ A \theta_A θA与 θ B \theta_B θB为定值。
约束条件在前文调头弧线模型的建立已经进行了详细的推导得出了方程组(44)(45)(46)47与(48)方程组本问中需要的约束条件。
综上所述最短弦长单目标优化模型建立如下 min θ C , θ D L ( θ C , θ D ) 170 2 π ∫ θ C θ A 1 θ 2 d θ r 1 ⋅ ϕ 1 r 2 ⋅ ϕ 2 170 2 π ∫ θ D θ B 1 θ 2 d θ s.t. { k w 1 − cos ( θ C ) − θ ⋅ sin ( θ C ) ) sin ( θ C ) θ C ⋅ cos ( θ C k w 2 − cos ( θ D ) − θ ⋅ sin ( θ D ) ) sin ( θ D ) θ D ⋅ cos ( θ D ) x O 1 x C 1 k w 1 2 1 ⋅ r 1 y O 1 y C k w 1 k w 1 2 1 ⋅ r 1 x O 2 x D 1 k w 2 2 1 ⋅ r 2 y O 2 y D k w 2 k w 2 2 1 ⋅ r 2 r 1 2 r 2 ( x O 1 − x O 2 ) 2 ( y O 1 − y O 2 ) 2 r 1 r 2 x E x O 1 r 1 r 1 r 2 ⋅ ( x O 2 − x O 1 ) y E y O 1 r 1 r 1 r 2 ⋅ ( y O 2 − y O 1 ) ϕ 1 arctan ( y C − y O 1 x C − x O 1 ) − arctan ( y E − y O 1 x E − x O 1 ) ϕ 2 arctan ( y D − y O 2 x D − x O 2 ) − arctan ( y E − y O 2 x E − x O 2 ) \min_{\theta_C,\theta_D} L(\theta_C,\theta_D)\frac{170}{2\pi}\int_{\theta_C}^{\theta_A}\sqrt{1\theta^2}d\thetar_1\cdot\phi_1r_2\cdot\phi_2\frac{170}{2\pi}\int_{\theta_D}^{\theta_B}\sqrt{1\theta^2}d\theta \\ \text{s.t.} \begin{cases} k_{w1}-\frac{\cos(\theta_C)-\theta\cdot\sin(\theta_C))}{\sin(\theta_C)\theta_C\cdot\cos(\theta_C}\\ k_{w2}-\frac{\cos(\theta_D)-\theta\cdot\sin(\theta_D))}{\sin(\theta_D)\theta_D\cdot\cos(\theta_D)}\\ x_{O1}x_C\frac1{\sqrt{k_{w1}^21}}\cdot r_1\\ y_{O1}y_C\frac{k_{w1}}{\sqrt{k_{w1}^21}}\cdot r_1\\ x_{O2}x_D\frac1{\sqrt{k_{w2}^21}}\cdot r_2\\ y_{O2}y_D\frac{k_{w2}}{\sqrt{k_{w2}^21}}\cdot r_2\\ r_12r_2\\ \sqrt{(x_{O1}-x_{O2})^2(y_{O1}-y_{O2})^2}r_1r_2\\ x_Ex_{O1}\frac{r1}{r1r2}\cdot (x_{O2}-x_{O1})\\ y_Ey_{O1}\frac{r1}{r1r2}\cdot (y_{O2}-y_{O1})\\ \phi_1\arctan(\frac{y_C-y_{O1}}{x_C-x_{O1}})-\arctan (\frac{y_E-y_{O1}}{x_E-x_{O1}})\\ \phi_2\arctan(\frac{y_D-y_{O2}}{x_D-x_{O2}})-\arctan (\frac{y_E-y_{O2}}{x_E-x_{O2}})\\ \end{cases} θC,θDminL(θC,θD)2π170∫θCθA1θ2 dθr1⋅ϕ1r2⋅ϕ22π170∫θDθB1θ2 dθs.t.⎩ ⎨ ⎧kw1−sin(θC)θC⋅cos(θCcos(θC)−θ⋅sin(θC))kw2−sin(θD)θD⋅cos(θD)cos(θD)−θ⋅sin(θD))xO1xCkw121 1⋅r1yO1yCkw121 kw1⋅r1xO2xDkw221 1⋅r2yO2yDkw221 kw2⋅r2r12r2(xO1−xO2)2(yO1−yO2)2 r1r2xExO1r1r2r1⋅(xO2−xO1)yEyO1r1r2r1⋅(yO2−yO1)ϕ1arctan(xC−xO1yC−yO1)−arctan(xE−xO1yE−yO1)ϕ2arctan(xD−xO2yD−yO2)−arctan(xE−xO2yE−yO2)
位置速度计算模型基于调头弧线的优化
在问题一与问题二中求解板凳龙各个圆孔的位置与速度按照本文的思路需要计算螺线的积分但是第四问中路线修改为螺线与双圆弧进行拼接无法直接套用问题一中对于位置的求解模型于是在第一问的基础上进行了优化以 C C C点处的弧度 θ C \theta_C θC与 D D D点弧度 θ D \theta_D θD作为盘入圆弧与盘出圆弧的范围限制在从 [ θ D , θ C ] [\theta_D,\theta_C] [θD,θC]的弧度范围内设定映射关系 g ( θ ) g(\theta) g(θ)表述双圆弧段的板凳位置计算则整个板凳龙所处位置与弧度的关系可表示为函数 f ( x , y ) f(x,y) f(x,y)如下 f ( x , y ) { r ( θ ) , if θ ∈ [ θ C 32 π ] g ( θ ) , if θ ∈ [ θ D , θ C ] − r ( θ ) , if θ ∈ [ − 32 π θ D ] f(x,y)\begin{cases}r(\theta),\text{if}\theta\in[\theta_C32\pi]\\g(\theta),\text{if}\theta\in[\theta_D,\theta_C]\\-r(\theta),\text{if}\theta\in[-32\pi\theta_D]\end{cases} f(x,y)⎩ ⎨ ⎧r(θ),g(θ),−r(θ),ifθ∈[θC32π]ifθ∈[θD,θC]ifθ∈[−32πθD] 其中 r ( θ ) r(\theta) r(θ)与 − r ( θ ) -r(\theta) −r(θ)分别对应前段螺线和后段螺线的位置计算函数可沿用第一问模型对于 g ( θ ) g(\theta) g(θ)求解过程如下
由于驶入圆弧前与盘出圆弧后的位置 C ( x C , y C ) C(x_C,y_C) C(xC,yC)与 D x D , y D ) Dx_D,y_D) DxD,yD)可根据最短路径 L L L求解导出由此可导出对应的 θ C \theta_C θC与 θ D \theta_D θD进而求解出圆弧部分总长度 S ϕ 1 ⋅ r 1 ϕ 2 ⋅ r 2 S\phi_1\cdot r_1\phi_2\cdot r_2 Sϕ1⋅r1ϕ2⋅r2.
其中在整个圆弧路线中盘入点弧度 θ C \theta_C θC为其中弧度的最大值即 θ m a x \theta_{max} θmax同理盘出点 θ D \theta_D θD为其中弧度的最小值即 θ m i n \theta_{min} θmin故可定义映射系数 k k k表示螺旋线开始到当前角度 θ \theta θ 的路径长度相对于整个双圆弧路径长度的比例其公式为 k θ − θ m i n θ m a x − θ m i n k\frac{\theta-\theta_{min}}{\theta_{max}-\theta_{min}} kθmax−θminθ−θmin 基于比例系数与映射关系可计算出此时在圆弧部分前进路程与后段圆弧长度的差值 Δ x \Delta x Δx及其含义为 Δ x k ⋅ S − ϕ 2 ⋅ r 2 \Delta xk\cdot S-\phi_2\cdot r_2 Δxk⋅S−ϕ2⋅r2 其中
如果 Δ x 0 \Delta x 0 Δx0则当前点在前段圆弧上如果 Δ x ≤ 0 \Delta x \leq 0 Δx≤0则当前点在后段圆弧上。
若该点在前段圆弧上
首先需求解出该点在前段圆弧上的角度占比 η \eta η的表达式为 η k S − r 2 ϕ 2 2 π r 1 \eta\frac{kS-r_2\phi_2}{2\pi r_1} η2πr1kS−r2ϕ2 最终导出该点的极坐标下的极角 θ 0 \theta_0 θ0表达为 θ 0 arctan ( y E − y O 1 x E − x O 1 ) 2 π η \theta_0\arctan (\frac{y_E-y_{O1}}{x_E-x_{O1}})2\pi\eta θ0arctan(xE−xO1yE−yO1)2πη 最终将 θ 0 \theta_0 θ0与圆心 O 1 ( x O 1 y O 1 ) O_1(x_{O1}y_{O1}) O1(xO1yO1)的坐标结合即可推导出此时该点的位置坐标为 f ( x , y ) { x x O 1 r 1 cos θ 0 y y O 1 r 1 sin θ 0 f(x,y) \begin{cases} xx_{O1}r_1\cos\theta_0\\ yy_{O1}r_1\sin\theta_0 \end{cases} f(x,y){xxO1r1cosθ0yyO1r1sinθ0 若该点在后段圆弧上
推导过程与前段圆弧类似只是将其中相对的位置计算点改为 O 2 ( x 2 , y 2 ) O_2(x_2,y_2) O2(x2,y2)并且由于其在后段圆弧上投影范围从 [ θ D , θ C ] [\theta_D,\theta_C] [θD,θC]则其圆弧长度小于后段圆弧的总长故公式表达修改如下 η r 2 ϕ 2 − k S 2 π r 1 \eta\frac{r_2\phi_2-kS}{2\pi r_1} η2πr1r2ϕ2−kS θ 0 arctan ( y E − y O 2 x E − x O 2 ) 2 π η \theta_0\arctan (\frac{y_E-y_{O2}}{x_E-x_{O2}})2\pi\eta θ0arctan(xE−xO2yE−yO2)2πη f ( x , y ) { x x O 2 r 2 cos θ 0 y y O 2 r 2 sin θ 0 f(x,y) \begin{cases} xx_{O2}r_2\cos\theta_0\\ yy_{O2}r_2\sin\theta_0 \end{cases} f(x,y){xxO2r2cosθ0yyO2r2sinθ0
最终导出位置计算的全曲线表达式为 f ( x , y ) { r ( θ ) i f θ ∈ [ θ C , 32 π ] { x x O 1 r 1 cos ( θ 0 ) , y y O 1 r 1 sin ( θ 0 ) , i f θ ∈ [ θ D , θ C ] a n d Δ x 0 { x x O 2 r 2 cos ( θ 0 ) , y y O 2 r 2 sin ( θ 0 ) , i f θ ∈ [ θ D , θ C ] a n d Δ x ≤ 0 − r ( θ ) i f θ ∈ [ − 32 π , θ D ] \begin{aligned}f(x,y)\begin{cases}r(\theta)\mathrm{if~}\theta\in[\theta_C,32\pi]\\\begin{cases}xx_{O1}r_1\cos(\theta_0),\\yy_{O1}r_1\sin(\theta_0),\end{cases}\mathrm{if~}\theta\in[\theta_D,\theta_C]\mathrm{~and~}\Delta x0 \\\begin{cases}xx_{O2}r_2\cos(\theta_0),\\yy_{O2}r_2\sin(\theta_0),\end{cases}\mathrm{if~}\theta\in[\theta_D,\theta_C]\mathrm{~and~}\Delta x\leq0\\-r(\theta)\mathrm{if~}\theta\in[-32\pi,\theta_D]\end{cases}\end{aligned} f(x,y)⎩ ⎨ ⎧r(θ){xxO1r1cos(θ0),yyO1r1sin(θ0),{xxO2r2cos(θ0),yyO2r2sin(θ0),−r(θ)if θ∈[θC,32π]if θ∈[θD,θC] and Δx0if θ∈[θD,θC] and Δx≤0if θ∈[−32π,θD] 由于数学公式过于繁琐故在本问中速度的求解模型的建立只采取了问题一中以极小时间步长 Δ t \Delta t Δt进行数值模拟的方法具体过程将置于模型的求解中进行。
模型的求解
本题中需要求解的部分是单目标优化模型最优值的求解与速度计算小步长数值模拟的求解。对于本题最短弦长单目标优化模型的求解由于主要决策变量仅为 θ C \theta_C θC与 θ D \theta_D θD计算量较小故采取定步长搜索策略具体步骤如下
Step1 设定初始解向量和初始步长
设定初始解 [ θ C 0 , θ D 0 ] [\theta_{C0},\theta_{D0}] [θC0,θD0]为 [ π , π ] [\pi,\pi] [π,π],步长均为 2 π 360 \frac{2\pi}{360} 3602π即角度制中的 1 ° 1° 1°。
Step2 遍历解空间
采用双层循环的方式遍历解空间通过固定其中一个变量改变另一个变量并计算每一对 [ θ C , θ D ] [\theta_C,\theta_D] [θC,θD]的目标函数值若比之前的最优值更好则更新最优解和最优值。
Step3 得到最优解
将最终保存的最优解与最优值导出。
对于问题四速度的求解仍沿用问题一中提出的求解办法具体步骤参照问题一模型的求解部分。