[程序分享] 用deSolve作Lorenz混沌系统的仿真 [推广有奖]

0关注
0粉丝

硕士生

68%

还不是VIP/贵宾

威望: 0 级
论坛币: 451 个
通用积分: 8.4091
学术水平: 2 点
热心指数: 1 点
信用等级: 0 点
经验: 1554 点
帖子: 95
精华: 0
在线时间: 215 小时
注册时间: 2014-11-15
最后登录: 2023-6-3

楼主

onczhao 发表于 2015-9-26 16:11:41 |只看作者 |坛友微信交流群|倒序 |AI写论文

相似文件

换一批

是否 +2 论坛币

k人参与回答

经管之家送您一份

应届毕业生专属福利!

求职就业群

赵安豆老师微信：zhaoandou666

经管之家联合CDA

送您一个全额奖学金名额~ !

立即领取

感谢您参与论坛问题回答

经管之家送您两个论坛币！

+2 论坛币

library(deSolve);
parameters <- c(a = -8/3,
            b = -10,
            c = 28);
state <- c(X = 1,
         Y = 1,
         Z = 1);
Lorenz<-function(t, state, parameters)
{
  with(as.list(c(state, parameters)),
   {
      # rate of change
      dX <- a*X + Y*Z
      dY <- b * (Y-Z)
      dZ <- -X*Y + c*Y - Z

      # return the rate of change
      list(c(dX, dY, dZ))
   }) # end with(as.list ...
};
Span=0.006
SimTime=100;
IterTimes=seq(1,SimTime/Span,by=1);
out=matrix(data = NA, nrow = length(IterTimes), ncol = 4)
for (i in IterTimes){
  times <- seq(0, Span, by = Span);
  out1 <- ode(y = state, times = times, func = Lorenz, parms = parameters);
  nrosOut1=dim(out1)[1];
  state <- c(out1[nrosOut1,2],out1[nrosOut1,3], out1[nrosOut1,4]);
  out[i,]=c(Span*(i-1),state);
}

head(out);
tail(out);
par(oma = c(0, 0, 3, 0),lty="dotdash",col=rgb(0, 0, 1))#Set or Query Graphical Parameters
plot(out[,1], out[,2],xlab = "time", ylab = "x(t)",type="l",col=rgb(0, 0, 1))
mtext(outer = TRUE, side = 3, "Lorenz model", cex = 1.5)
plot(out[,1], out[,3],xlab = "time", ylab = "y(t)",type="l",col=rgb(0, 1, 0))
mtext(outer = TRUE, side = 3, "Lorenz model", cex = 1.5)
plot(out[,1], out[,4],xlab = "time", ylab = "y(t)",type="l",col=rgb(1, 0, 0))
#plot(out[, "3"], out[, "4"], pch = ".")
mtext(outer = TRUE, side = 3, "Lorenz model", cex = 1.5)