求简谐下滑曲线的方程

Ickiverar

从能量守恒进行验算：你的解的动能 $1/2(\dot x^2+\dot y^2)=2{cos^2 2θ}/cos^2θ$ ，重力势能 $y$ 中却含有$ln cos^2θ$，这两个东西加起来不可能成为常数。

Ickiverar

我把方程简化了一下，最终需要求解的方程是$\frac{\text{d}y}{\text{d}t}=\sqrt{\sin^2t-ay}$：

1. Clear[x,equ,Energy,EnergyConst,EnergyPhi];
   
2. x=A Sin[\[Omega] t];
   
3. (*能量，实际上是两倍能量=vx^2+vy^2+2gy*)
   
4. Energy=D[x,t]^2+D[y[t],t]^2+2g y[t]
   
5. (*将能量表达式自变量由t换为\[Phi]=\[Omega]t，相位\[Phi]和t是线性关系，可以看作无量纲时间。*)
   
6. EnergyPhi=Energy/.{y[t]->y[\[Phi]],y'[t]->y'[\[Phi]]\[Omega]}/.t->\[Phi]/\[Omega]
   
7. (*根据能量守恒，上式应为常数。考虑t\[Equal]0时，x\[Equal]0，小球位于水平方向简谐运动平衡点，x方向速度最大，位于y方向最低点y\[Equal]0，且y方向速度为0，即y'[\[Phi]]\[Equal]0。则该常数为*)
   
8. EnergyConst=EnergyPhi/.{y[\[Phi]]->0,y'[\[Phi]]->0}/.\[Phi]->0
   
9. (*现在已经得到微分方程*)
   
10. EnergyPhi-EnergyConst
    
11. (*令\[Eta]为 y/A，a\[Equal](2 g)/(A \[Omega]^2)，可以将上式无量纲化：*)
    
12. equ=FullSimplify[(EnergyPhi-EnergyConst)/(A \[Omega])^2/.{y[\[Phi]]->A \[Eta][\[Phi]],y'[\[Phi]]->A \[Eta]'[\[Phi]],g->a A \[Omega]^2/2}]
    
13. (*考虑 0<\[Phi]<\[Pi]/2的区间，此时小球在上升，取上式中 \[Eta]'[\[Phi]]>0 的解，则最终需解的方程为*)
    
14. DSolve[{\[Eta]'[\[Phi]]==Sqrt[Sin[\[Phi]]^2-a \[Eta][\[Phi]]],\[Eta][0]==0},\[Eta][\[Phi]],\[Phi]]

复制代码

hujunhua

@wayne 在8#和 @Ickiverar在10#、11#的质疑是对的，纠正了我一直以来的错误印象，可能@mathe有跟我一样的印象。
wayne在3#得出的方程是正确的, 不过应用边界条件能进一步得到 $H=1/4$

设滑行速度为 `v`, 曲线切线倾角为  `φ`, 切线斜率`y'=\tan φ`。
沿用3#的符号， `y',y''`表示对 `x` 求导, $dot x,ddot x, dot y, ddot y$表示对时间 `t` 求导。
则切向加速度为 `\dot v=-g \sin φ`, 水平速度为 $\dot x=vcosφ$.
所以水平加速度为\[
\ddot x=(v\cos φ)'_t=\dot v\cos φ-v\sin φ·\dot φ=-g\sin φ\cos φ-v\sin φ·\dot φ\tag{1}
\]8#的水平加速度与之相比，少了 `-v\sin φ·\dot φ`项。即楼上指出的轨道曲率的作用。
在端点 `v=0`, 所以水平方向最大加速度为 $-g sin φ cos φ$, 其大小不超 `g/2`, 这一点倒是没有变.
所以仍应设定端点加速度大小为 $g/2$, 即端点坡度为 45°. Ickiverar在7#也指出了端点斜率不能大太。
  $x=sin ωt,\dot x=ωcosωt,\ddot x=-ω^2x$
可见水平方向速度最大值平方和加速度最大值均为 $ω^2$, 故有 `ω^2=2gh=g/2`→ `h=1/4, ω=\sqrt{g/2}`
将$H=1/4$代入3#的 (3) 式得
\[
(1-x^2)y'^2+4y=x^2\tag{3}
\]试了一下，Mathematica14都解不出这个非线性方程。
将`y'=\dot y/\dot x`代入可得参数方程\[
\begin{cases}x=\sin ωt\\\dot y^2=ω^2(\sin^2 ωt-4y)\end{cases}\tag{4}
\]取 `g=2→ω=1`, 代入得\[
\begin{cases}x=\sin t\\\dot y^2=\sin^2 t-4y\end{cases}\tag{5}
\]第2个微分方程与楼上是一致的，用Mathematica14也解不出来。

wayne

关于8楼的加速度分解的物理原理,可以说不存在, 也可以说存在, 既然@mathe 也来凑热闹了, 那我认真补充一下,都说的没问题,只是表述不同, ^_^
1) 受力分析,滑块只受到两个力的作用,一个是重力$G$, 一个是滑道对滑块的斜向内的支持力$N$,反正就是跟重力的方向对抗着来就是了,是曲线的法线方向,向内. 根据牛顿第二运动定律 , $G$和$N$的合力$F =ma$,  方向是曲线的切线方向.
2) 我们引入虚力,$F'=-F$,方向相反, 那么$F',N,G$三个力构成的矢量三角形是等效的静力学平衡.  容易得到 $F=Gsin(\phi)$,也就是说,切向加速度是 $a=\dot v=gsin\phi$

所以,只根据1),我可以认为加速度分解的物理基础是不存在的. 但是依据2),我们又可以继续使用加速度分解(引入虚拟力,产生的等效的静力学平衡)

wayne

\[
(1-x^2) y'^2+4y=x^2
\]解得的级数解是

1. FullSimplify[url=home.php?mod=space&uid=6175]@[/url] Normal[AsymptoticDSolveValue[{y'[x]^2 + 1 == (1 - 4 y[x])/(1 - x^2), y[0] == 0}, y[x], {x, 0, 15}]]

复制代码

\[-\frac{5 \left(539096+398507 \sqrt{2}\right) x^{10}}{60082624}-\frac{\left(991+741 \sqrt{2}\right) x^8}{15778}-\frac{1}{392} \left(38+29 \sqrt{2}\right) x^6-\frac{1}{28} \left(5+4 \sqrt{2}\right)
   x^4-\frac{1}{2} \left(1+\sqrt{2}\right) x^2\]

\[
\begin{array}{l}
\{2\}\to \frac{1}{2} \left(\sqrt{2}-1\right) \\
\{4\}\to \frac{1}{28} \left(4 \sqrt{2}-5\right) \\
\{6\}\to \frac{1}{392} \left(29 \sqrt{2}-38\right) \\
\{8\}\to \frac{741 \sqrt{2}-991}{15778} \\
\{10\}\to \frac{5 \left(398507 \sqrt{2}-539096\right)}{60082624} \\
\{12\}\to \frac{720697 \sqrt{2}-982152}{28815136} \\
\{14\}\to \frac{5 \left(1556988441871 \sqrt{2}-2132772637848\right)}{394632047321344} \\
\{16\}\to \frac{3 \left(146153667670419 \sqrt{2}-200965448001404\right)}{27278940271087904} \\
\{18\}\to \frac{5577416513496157005 \sqrt{2}-7691649054995007296}{415512818209210953728} \\
\{20\}\to \frac{382274845963828979775 \sqrt{2}-528417116345422804928}{33448781865841481775104} \\
\{22\}\to \frac{8825496040763685318196437 \sqrt{2}-12222733902541405671775552}{892948680690504197468176384} \\
\{24\}\to \frac{15815453984392449125957646369 \sqrt{2}-21938161381073333601837761344}{1826749763522598961970521837568} \\
\{26\}\to \frac{526864747435226624256851925681610149 \sqrt{2}-731814521130374706391111765711823360}{68728873174872637944010733879574003712} \\
\{28\}\to \frac{2572133743225521495968838539079515036807 \sqrt{2}-3576809893268935045151999051418931038208}{375500198590916657407102644551052569280512} \\
\{30\}\to \frac{1301511075185300333676450611035735994747343571 \sqrt{2}-1811692688679492890043701881732847702383593984}{210974035577909342130865806635719783945117106176} \\
\{32\}\to \frac{11442145613755372645211283277249453402129263652007 \sqrt{2}-15941343912322709126027496733168265258106473208288}{2045439425498113739627335122228504869050648338743296} \\
\{34\}\to \frac{3 \left(423972869450180339192265214202670695402023814502197329 \sqrt{2}-591142647652215400577478832801243270222614081887686656\right)}{249249066633498147856028548654276689323035803965903077376}
   \\
\{36\}\to \frac{3 \left(241497139595400960755646913953543561605064168603825289 \sqrt{2}-336950670131267026730932132073807106065799019401875456\right)}{154810406438453001869055794060355244752542999349064368128}
   \\
\{38\}\to \frac{5 \left(2146090975517365160727709924590836799632418482056181035890600347
   \sqrt{2}-2996206299652137648892611263684715781412476601205754915738005504\right)}{2488437095863444102559947377260233969067374517738407168513427898368} \\
\{40\}\to \frac{155676033619008589548102244675766500616812671859230872185363458155
   \sqrt{2}-217464497704669351857310146858272216863264307836920345226375483392}{39015138753001855750850603522044382586449193331684312392049815977984} \\
\{42\}\to \frac{21699163138967256491129534029822458107820422496591195515953069859115591918041
   \sqrt{2}-30327047276537941519099465298535227467088457922360313459008343183779789553664}{5854629608312031730215436522489073181478384065195519534867591894383111570456576} \\

\end{array}
\]

1. a[1] = 0;
   
2. a[n_] := If[Mod[n, 2] == 1, 0, Subscript[a, n]]; Do[
   
3. fx = Sum[If[i == 2 m, Subscript[a, 2 m], a[i]] x^i, {i, 1, 2 m + 1}];
   
4.   sol = SolveValues[
   
5.    Coefficient[1 - 4 fx - (D[fx, x]^2 + 1) (1 - x^2), x^(2 m)] == 0,
   
6.    Subscript[a, 2 m]];
   
7. If[Length[sol] > 0, a[2 m] = FullSimplify[Last@sol]], {m, 20}];
   
8. fx

复制代码

hujunhua

设 $y=a_1x^2+a_2x^4+a_3x^6+...+a_rx^{2r}+...$
验算
`y(1)≈a_1+a_2+...+a_{10}≈0.247524,\\ y'(1)=2a_1+4a_2+...+20a_{10}≈0.656957`
收敛好慢，是不是算错了。

hujunhua

上述级数解，y(1)应该能收敛到1/4, 但 y'(1)没有希望收敛到1（太慢）, y''(1)更不可能收敛到 0. 所以没有希望从它得到曲线正确的形貌。
原因，应该是多项式形式不适用，因为函数有3个奇点，在x=-1, 0, 1处。
10#的方程虽然不正确，但11#的图像却反映了曲线的概貌。事实上，将 y 坐标 除以 4(1-ln2), 压缩的结果与上述微分方程的数值解极其接近。

从微分方程\[
\frac{1-4y}{1+y'^2}=1-x^2
\] 可以得到曲线的一些特征值，比如 `y(0)=0, -∞, y'(0)=0, ∞, y(1)=1/4, y''(1)=0, y''(0)=\sqrt2-1`。

曲线有两个尖点，也许可以用一条 6 次曲线来逼近它。方程如下：\[
(a_1x^4+b_1x^2+c_1)y^2+(a_2x^4+b_2x^2+c_2)y+a_3x^4+b_3x^2=0
\]待有空了来试试这些待定系数。
这个方程，即使作变换 `x^2->u`, 也只能降为4次曲线，不能直接运用比较熟悉的三次曲线方法。

Jack315

微信视频号“喵喵喵世界”中的视频：
简谐运动视频

获取视频号中的视频或其链接的方法：
视频号视频下载神器

hujunhua

mathe在8#的评论中说“加速度分解本质上就是力的分解”， 是符合牛顿第二定律的。
如图，滑块受重力G和曲面支撑力N的作用，合力F=ma, 在曲面曲率 κ 不为零时，F不在曲面切向。
因为滑块对曲面的压力除了重力的法向分力 G cos φ, 还有离心力 mκv², 所以曲面支撑力 N=G cosφ + mκv²。
N - G cosφ = mκv²， 即支撑力减去重力的法向分量，剩余部分提供了向心力。
滑块受到的作用力在水平方向的分量 = N sinφ = G cosφ sinφ + mκv² sinφ, 故水平加速度为
`\ddot x=-(g\cosφ\sin φ + κv² \sinφ)`
8#的错误就在于忽视了向心力的水平分量κv² sinφ。

hujunhua

简谐下滑曲线的数值解是这样的，它以x=±1为渐近线。
左图是局部，右图是下伸更深的图。
轨道解有3支，蓝色的像一段圆弧的容易理解，它也是实际可用于盘子母线的部分。
绿色和橙色的关于 y 轴对称，两者交点是曲线的一个奇点（二重点），曲线在此附近近似于蓝色虚线所示的两条抛物线。
这意味着滑块经过此二重点时呈自由抛落状态，轨道支撑力为零。继续下滑，若滑块还在轨道上方，将在离心力作用下甩脱。
为了保持在轨道上，向下经过二重点后，滑块进入轨道下方，被轨道压住以免甩脱。这是单程运动，永坠深渊，无路返回。（又断头了？）


如此一来，则17#的6次曲线就不对了，应该用一条7次曲线来近似：\[
(a_1x^4+2b_1x^2+c_1)y^3+3(a_2x^4+2b_2x^2+c_2)y^2+3(a_3x^4+2b_3x^2+c_3)y+a_4x^4+2b_4x^2=0
\]因为，它与竖直线（x=常数）一般有3个交点，与水平线（y=常数）有4个交点（在x轴以下有两个虚交点）。
可以取`a_1=1`.