R中的非线性摆：获取角度随时间的精确演变的问题

问题描述

我正在研究一个简单的摆式问题，以便在知道数据集后的时间间隔内获得角度θ的变化。为了得出解决方案，我以Meléndez et alt的工作为基础，得出θ（t）的时间演化的以下精确表达式。

我在R中的位置是这样的：

PositionAng <- function(t,AngIn,FreAng)
{
  2 * acos(sin(AngIn/2) * sn(K.fun((sin(AngIn/2)^2)-FrecAng * t),(sin(AngIn/2)^2))
}

函数θ（t）的值的计算变得越来越复杂，当我运行程序时，基于梅伦德斯的工作，得到的图非常不同。

我的图没有椭圆曲率，是线性的。

我想知道我的代码是否有错误，或者我可以用其他方式解决该问题。

从以下角度开始：角度pi / 4（45度），绳索长度为0.8m，重力为9.8m / s ^ 2。

我已经使用库（椭圆形）来找到这样的结果，但是我不知道图形会发生什么。

我建立了以下代码：

library(elliptic)

L <- 0.8
g <- 9.8
AngIn <- pi/4

#Angular frequency
FrequencyAng <- function(g,L)
{
  2*pi*(sqrt(g/L))
}

#Obtain the value of the angle in a certain second t.
PositionAng <- function(t,(sin(AngIn/2)^2))
}

#Function that applies the prevIoUs formula using the different intervals (t).
Test <- function(L,g,AngIn)
{
  FrecAng <- 0
  t = seq(from= 0,to= 10,by= 0.1)
  PosAng <- PositionAng(t,FrecAng)
  return (PosAng)
}

#Evolution of the angle as a function of time
EvolAng <- Test(L,AngIn)

plot(EvolAng,type = "l",main = "Evolution of the angle in time",ylab = "Angle",xlab = "Time
)

sn 是Jacobi的椭圆函数。

我想要获得如下所示的内容

that

。

解决方法

正如奥利弗（Oliver）指出的那样，您的问题是您正在FrecAng中将Test设置为0。由于PositionAng仅使用t一次，并乘以FrecAng（即0），因此该项始终为0，因此此函数会产生一个常数。我认为应该是：

Test <- function(L,g,AngIn)
{
  FrecAng <- FrequencyAng(g,L)
  t = seq(from = 0,to = 10,by = 0.1)
  PosAng <- PositionAng(t,AngIn,FrecAng)
  return (PosAng)
}

哪个会产生：

EvolAng <- Test(L,AngIn)

plot(EvolAng,type = "l",main = "Evolution of the angle in time",ylab = "Angle",xlab = "Time")

这似乎给出了明智的结果-对于更长的摆，我们看起来像简单的谐波运动：

plot(Test(20,9.8,pi/4),xlab = "Time")

但是对于较短的长度，我们会看到非线性影响：

plot(Test(0.4,xlab = "Time")

pendulum plot r r