文档章节

python里求解物理学上的双弹簧质能系统

wangxuwei
 wangxuwei
发布于 07/21 22:50
字数 724
阅读 7
收藏 1

物理的模型如下:

在这个系统里有两个物体,它们的质量分别是m1和m2,被两个弹簧连接在一起,伸缩系统为k1和k2,左端固定。假定没有外力时,两个弹簧的长度为L1和L2。
由于两物体有重力,那么在平面上形成摩擦力,那么摩擦系数分别为b1和b2。所以可以把微分方程写成这样:

这是一个二阶的微分方程,为了使用python来求解,需要把它转换为一阶微分方程。所以引入下面两个变量:

这两个相当于运动的速度。通过运算可以改为这样:

这时可以线性方程改为向量数组的方式,就可以使用python定义了,代码如下:

# Use ODEINT to solve the differential equations defined by the vector field
from scipy.integrate import odeint
 
def vectorfield(w, t, p):
    """
    Defines the differential equations for the coupled spring-mass system.
    Arguments:
        w :  vector of the state variables:
                  w = [x1,y1,x2,y2]
        t :  time
        p :  vector of the parameters:
                  p = [m1,m2,k1,k2,L1,L2,b1,b2]
    """
    x1, y1, x2, y2 = w
    m1, m2, k1, k2, L1, L2, b1, b2 = p
 
    # Create f = (x1',y1',x2',y2'):
    f = [y1,
         (-b1 * y1 - k1 * (x1 - L1) + k2 * (x2 - x1 - L2)) / m1,
         y2,
         (-b2 * y2 - k2 * (x2 - x1 - L2)) / m2]
    return f
 
# Parameter values
# Masses:
m1 = 1.0
m2 = 1.5
# Spring constants
k1 = 8.0
k2 = 40.0
# Natural lengths
L1 = 0.5
L2 = 1.0
# Friction coefficients
b1 = 0.8
b2 = 0.5
 
# Initial conditions
# x1 and x2 are the initial displacements; y1 and y2 are the initial velocities
x1 = 0.5
y1 = 0.0
x2 = 2.25
y2 = 0.0
 
# ODE solver parameters
abserr = 1.0e-8
relerr = 1.0e-6
stoptime = 10.0
numpoints = 250
 
# Create the time samples for the output of the ODE solver.
# I use a large number of points, only because I want to make
# a plot of the solution that looks nice.
t = [stoptime * float(i) / (numpoints - 1) for i in range(numpoints)]
 
# Pack up the parameters and initial conditions:
p = [m1, m2, k1, k2, L1, L2, b1, b2]
w0 = [x1, y1, x2, y2]
 
# Call the ODE solver.
wsol = odeint(vectorfield, w0, t, args=(p,),
              atol=abserr, rtol=relerr)
 
with open('two_springs.dat', 'w') as f:
    # Print & save the solution.
    for t1, w1 in zip(t, wsol):        
        out = '{0} {1} {2} {3} {4}\n'.format(t1, w1[0], w1[1], w1[2], w1[3]);
        print(out)
        f.write(out);


在这里把结果输出到文件two_springs.dat,接着写一个程序来把数据显示成图片,就可以发表论文了,代码如下:

# Plot the solution that was generated
 
from numpy import loadtxt
from pylab import figure, plot, xlabel, grid, hold, legend, title, savefig
from matplotlib.font_manager import FontProperties
 
t, x1, xy, x2, y2 = loadtxt('two_springs.dat', unpack=True)
 
figure(1, figsize=(6, 4.5))
 
xlabel('t')
grid(True)
lw = 1
 
plot(t, x1, 'b', linewidth=lw)
plot(t, x2, 'g', linewidth=lw)
 
legend((r'$x_1$', r'$x_2$'), prop=FontProperties(size=16))
title('Mass Displacements for the\nCoupled Spring-Mass System')
savefig('two_springs.png', dpi=100)

最后来查看一下输出的png图片如下:

本文转载自:https://blog.csdn.net/caimouse/article/details/78070550

共有 人打赏支持
wangxuwei
粉丝 24
博文 332
码字总数 115454
作品 0
杭州
其他
私信 提问
5 分钟快速入门 Python

您是否正在寻找 一种容易学习的编程语言, 来帮助完成科学工作? 不必舍近求远, Python 就能办到。 Python 非常适合用来 快速实现和测试新旧算法, 并将多个物理学代码编组到一起。 它易于使...

IBMdW
04/25
0
0
【ZZ】Python下划线与命名规范

Python下划线与命名规范 以下分四种情况说明下划线的作用,python对成员域没有严格控制,大部份只是作为命名规范存在,以下英文部份摘自python官方网站 singleleadingunderscore: weak "int...

学徒1986
2013/09/03
0
0
人人都能学会的python编程教程1:第一行代码

前言 众所周知,现在IT行业很火,行业薪酬也很高,国家在2017年也发布了人工智能教育的推广计划,人人会编程的时代将要到来。不会编程、不懂编程可能有些跟不上时代的节奏,普通人不懂技术也...

编程老司机
05/10
0
0
Python基础教程系列:一、Python起步

毕竟不是写论文,也不是照本宣科,尽管写的内容主要还是来自于书本或者网络,但也尽量以自己理解的方式写出来吧,呵呵! 首先软件的安装我就不说了,建议尽量在Linux系统下玩,而且一般Linux...

东方神剑
2015/01/20
0
2
人人都能学会的python编程教程2:数据类型和变量

数据类型 了解一门编程语言最开始就是了解它的数据类型了,python基本的数据类型分为如下几类: 整数 Python可以处理任意大小的整数,当然包括负整数,在程序中的表示方法和数学上的写法一模...

编程老司机
05/10
0
0

没有更多内容

加载失败,请刷新页面

加载更多

MyBaties分页插件PageHelper的简单使用

抛出问题: 如果想要将现有的select语句改为支持分页功能的查询语句该怎么做呢? 最简单的一种做法就是将所有的select语句都加上limit来实现分页,这种做法有什么问题呢? 有没有一种简便方法...

嘴角轻扬30
21分钟前
1
0
创业公司如何快速构建高效的监控系统?

12 月 7 日,在 2018 ArchSummit 全球架构师峰会·运维与监控专场,七牛云资深运维开发工程师贺强带来了主题为《如何快速构建高效的监控系统》的内容分享。 本文是对演讲内容的实录整理。 大...

七牛云
21分钟前
0
0
Docker与自动化测试及其测试实践

Docker 与自动化测试 对于重复枯燥的手动测试任务,可以考虑将其进行自动化改造。自动化的成本在于自动化程序的编写和维护,而收益在于节省了手动执行用例的时间。简而言之,如果收益大于成本...

微笑向暖wx
21分钟前
1
0
Spring Cloud 分布式链路跟踪 Sleuth + Zipkin + Elasticsearch

随着业务越来越复杂,系统也随之进行各种拆分,特别是随着微服务架构的兴起,看似一个简单的应用,后台可能很多服务在支撑;一个请求可能需要多个服务的调用;当请求迟缓或不可用时,无法得知...

编程SHA
24分钟前
1
0
Swift-清除缓存

func removeCache (){ // 取出cache文件夹路径.如果清除其他位子的可以将cachesDirectory换成对应的文件夹 let cachePath = NSSearchPathForDirectoriesInDomains(FileMan...

west_zll
25分钟前
1
0

没有更多内容

加载失败,请刷新页面

加载更多

返回顶部
顶部