大家好，给大家分享一下python数学建模基础教程，很多人还不知道这一点。下面详细解释一下。现在让我们来看看！
1️⃣python 语言快速入门
- 1.7 输入一个整数，判断它是否是水仙花数：水仙花数指的是各位数字的立方和等于该数本身的的3位整数
def func(n): a, b, c = map(int, n) # 对字符串序列进行映射 return (a**3 + b**3 + c**3) == int(n) x = input("请输入三位整数") if func(x): print("这是水仙花数") else: print("这不是水仙花数")
- 1.8随机产生一个3位整数，将它的十位数字变为0.（假设生成的三位整数为738，则输出为708）
import random x = random.randint(100, 1000) # 随机生成的三位数 def handle(x): X = str(x) return int(X[0] + '0' + X[-1]) print("x=%d" % x, "handle x=%d" % handle(x))
- 1.9 输入整数 x , y x,y x,y 和 z z z,若 x 2 + y 2 + z 2 x^2+y^2+z^2 x2+y2+z2大于1000，则输出 x 2 + y 2 + z 2 x^2+y^2+z^2 x2+y2+z2千位以上的数字，否则输出三个数之和
def func(*args): res = sum(map(lambda x:x*x,args)) if res>1000: return int(res/1000) else: return res print("result:",func(10,12,20))#644 print("another:",func(20,20,30))#1
- 1.10 某运输公司在计算运费时，按运输距离 s ~s~ s 给运费给定一定的折扣率 d ~d d,折扣率 d ~d~ d 的标准输入如表1.10所示。输入基本运费 p ~p~ p ,货物的重量 w ~w~ w ,距离 s ~s~ s ,计算总运费 f ~f~ f ,其中，总运费的计算公式为 f = p w s ( 1 − d ) f=pws(1-d) f=pws(1−d)
运输距离 s ~s~ s 折扣率$ d/$% s < 250 ~s~<250 s <250 0 250 ≤ s < 500 250\leq~s~<500 250≤ s <500 2.5 500 ≤ s < 1000 500\leq~s~<1000 500≤ s <1000 4.5 1000 ≤ s < 2000 1000\leq~s~<2000 1000≤ s <2000 7.5 2000 ≤ s < 2500 2000\leq~s~<2500 2000≤ s <2500 9.0 2500 ≤ s < 3000 2500\leq~s~<3000 2500≤ s <3000 12.0 3000 ≤ s 3000\leq~s~ 3000≤ s 15.0 f=lambda p,w,s,d:p*w*s*(1-d) lis=[i/100for i in [2.5,4.5,7.5,9.0,12.0,15.0]]#百分数转换为小数 def much(p,w,s): if s<250: return 0 elif s<500: return f(p,w,s,lis[0]) elif s<1000: return f(p,w,s,lis[1]) elif s<2000: return f(p,w,s,lis[2]) elif s<2500: return f(p,w,s,lis[3]) elif s<3000: return f(p,w,s,lis[4]) else: return f(p,w,s,lis[5])#考察程序流程控制的if 分支语句
- 1.11 编写一个Python 程序，将日期作为输入并打印该日期是一周当中的周几，用户输入有3个： m ( 月 ) ， d ( 日 ) ， y ( 年 ) m(月)，d(日)，y(年) m(月)，d(日)，y(年).对于 m ~m~ m ,用1表示一月…以此类推python用turtle画简单树形图。对于输出，0表示周日，1表示周一，2表示周二，以此类推 计算阳历日期对应的周几可以用一下的公式:
y 0 = y − ( 14 − m ) / / 12 y_0=y-(14-m)//12 y0=y−(14−m)//12
x = y 0 + y 0 / / 4 − y 0 / / 100 + y 0 / / 400 x=y_0+y_0//4-y_0//100+y_0//400 x=y0+y0//4−y0//100+y0//400
m 0 = m + 12 × ( ( 14 − m ) / / 12 ) − 2 m_0=m+12\times((14-m)//12)-2 m0=m+12×((14−m)//12)−2
d 0 = ( d + x + ( 31 × m 0 ) / / 12 ) d_0=(d+x+(31\times m_0)//12)%7 d0=(d+x+(31×m0)//12)
def time(day,month,year): y=year-(14-month)//12 x=y+y//4-y//100+y//400 m=month+12*((14-month)//12)-2 d=(day+x+(31*m)//12)%7 return d print("今天是周=",time(9,7,2022))
- 1.12编写一个python程序，在给定年限 N ~N~ N 和年利率 r ~r~ r 的情况下，计算当贷款金额为 P ~P~ P 时，每月需要还贷的金额，每月的还贷公式为** P r ′ ( 1 + r ′ ) N ′ ( 1 + r ′ ) N ′ − 1 \frac{Pr'(1+r')^{N'}}{(1+r')^{N'}-1} (1+r′)N′−1Pr′(1+r′)N′** ,其中 N ′ = 12 N , r ′ = r / 12 N'=12N,r'=r/12 N′=12N,r′=r/12为月利息
f=lambda N,r,P:(P*r/12*(1+r/12)**(12*N))/((1+r/12)**(12*N)-1)#lambda 表达式即可
- 1.13 设地球表面A点的经度和纬度分别为 x 1 和 y 1 x1和y1 x1和y1，B点为 x 2 和 y 2 x2和y2 x2和y2，编写程序打印地球上A,B两点的大圆弧距离 d ~d~ d .大圆弧计算公式如下：
d = R arccos [ cos ( x 1 − x 2 ) cos ( y 1 ) cos ( y 2 ) + sin ( y 1 ) sin ( y 2 ) ] d=R \arccos[\cos(x_1-x_2)\cos(y_1)\cos(y_2)+\sin(y_1)\sin(y_2)] d=Rarccos[cos(x1−x2)cos(y1)cos(y2)+sin(y1)sin(y2)]
from math import * def distance(x1,x2,y1,y2): R=6370 lis=[radians(i) for i in [x1,x2,y1,y2]] x1,x2,y1,y2=lis return R*acos(cos(x1-x2)*cos(y1)*cos(y2)+sin(y1)*sin(y2))
- 1.4 如果一个数等于它的因子之和（不包括数本身），则称该数为完数.如6的因子为1,2,3 且6=1+2+3
找出1000以内的所有完数
#找出所有的因子 def all(x): a=1;s=0 while a<x: if x%a==0 : s+=a a+=1 return s #打印出所有的完数 for i in range(1,1001): if all(i)==i: print(i)# 6，28,496
- 1.16 编写函数，接受一个字符串，分别统计大写字母，小写字母，数字，其他字符的个数，并以元祖的形式返回结果
#利用字符串的判断方法，返回的bool值和int存在隐式转换最后进行求和统计 def judge(x:str): a= [ i.isupper() for i in x ] b= [i.islower()for i in x] c=[i.isdigit() for i in x] d=[ not (i.isdigit()|i.isalpha() )for i in x ] return map(sum,[a,b,c,d]) res=judge("aBNDE*&cs12")#测试 print(tuple(res))#(4, 3, 2, 2)
- 1.17 编写一个Python程序，将用户输入的一个1~999的整数转换成其对应的英文表示，例如，729将转换成 seven hundred and twenty nine.(程序中要尽可能使用函数封装的一些常用转换)
import num2word a=num2word.word(1290) print(a) #直接利用第三方库
2️⃣数据处理与可视化
- 2.1 考虑下面的 3 × 4 ~3 \times 4~ 3×4 矩阵
[ 1 2 3 4 5 6 7 8 9 10 11 12 ] \begin{bmatrix} 1 & 2 & 3&4 \\ 5&6&7&8 \\9&10&11&12\end{bmatrix} \quad ⎣ ⎡159261037114812⎦ ⎤
- 使用 a r r a y ~array~ array 函数在 P y t h o n ~Python~ Python 中构键该矩阵.
np.arange(1,13).reshape(3,4)
- 使用 a r a n g e ~arange~ arange 函数构造该矩阵.
a = list(range(1, 13)) m=np.array([a[0:4],a[4:8],a[8:12]]) print(m)
- 表达式 A[2,:] 的结果是什么？类似的表达式 A[2:]的结果是什么？
a = list(range(1, 13)) m=np.array([a[0:4],a[4:8],a[8:12]]) print(m[2:].shape) print(m[2:]) print('-'*30) print(m[2,:]) print(m[2,:].shape) #(1, 4) #[[ 9 10 11 12]] ------------------------------ #[ 9 10 11 12] #(4,)
- 2.2 构造范德蒙德矩阵 （ V a n d e r m o n d e m a t r i x ） ~（Vandermonde~matrix） （Vandermonde matrix）
V = [ 1 1 1 1 1 2 4 8 1 3 9 27 1 4 16 64 ] V=\begin{bmatrix} 1 & 1 & 1&1 \\ 1&2&4&8 \\1&3&9&27\\1&4&16&64\end{bmatrix} \quad V=⎣ ⎡1111123414916182764⎦ ⎤np.vander([1,2,3,4],increasing=True)#通过numpy命令直接构造
- 2.3 令 v = [ 1 , − 1 , 1 ] T ~v= [1,-1,1]^T~ v=[1,−1,1]T ,构造如下投影矩阵：
P = v v T v T v P =\frac{vv^T}{v^Tv}~ P=vTvvvT 和 Q = I − P ~Q=I-P Q=I−P
v=np.mat([1,-1,1]).T#创建矩阵（二维） P=(v*v.T)/(v.T*v) print(P) Q=np.eye(3)-P#三阶单位阵eye（3） print(Q)结果： [[ 0.33333333 -0.33333333 0.33333333] [-0.33333333 0.33333333 -0.33333333] [ 0.33333333 -0.33333333 0.33333333]] [[ 0.66666667 0.33333333 -0.33333333] [ 0.33333333 0.66666667 0.33333333] [-0.33333333 0.33333333 0.66666667]]
- 2.4 编写程序，生成1000个0~100内的随机整数，并统计每个元素出现的次数
#创建Series数据对象，进行值统计 num=pd.Series(np.random.randint(0,100,1000)) print(num.value_counts())
- 2.5 编写程序，生成包含20个随机数的数组，然后将前10个元素升序排列，后10个元素降序排列，并输出结果
num=np.random.randint(0,100,20)#生成0-100 以内的元素共20个的数组 print(num) index=np.argsort(-num[10:])+10#对数组的后十个元素进行反向排序返回索引 index1=np.argsort(num[0:10])#对前十个元素进行正向排序 num=np.concatenate((num[index1],num[index]),axis=None)#拼接两个数组 print(num) #result： [47 86 55 51 36 35 69 46 86 79 86 37 97 46 65 63 68 15 29 15] [35 36 46 47 51 55 69 79 86 86 97 86 68 65 63 46 37 29 15 15]
- 2.6 求矩阵 A = [ 9 80 205 40 90 − 60 96 1 210 − 3 101 80 ] A=\begin{bmatrix} 9 & 80 & 205&40 \\ 90&-60&96&1 \\210&-3&101&80\end{bmatrix} \quad A=⎣ ⎡99021080−60−32059610140180⎦ ⎤ 的鞍点，即该位置上元素是该行上的最大值，是该列上的最小值. （矩阵可能存在多个鞍点或没有鞍点）
A=np.array([[9,80,205,40],[90,-60,96,1],[210,-3,101,89]]) a=np.max(A,axis=1) b=np.min(A,axis=0) x=np.intersect1d(a,b)#intersect1d返回两个数组中的重复的元素 print(x) # 96
- 2.7 假设有一个英文文本文件，编写程序读取其内容，并将其中的大写字母变为小写字母，小写字母变为大写字母
f=open('new.txt','w') with open('文本', 'r') as file: for line in file: f.write(line.swapcase())#swapcase 将字符串中的大小写字母反转 file.close()
- 2.8 在一个图形界面中分别画出6条曲线
y = k x 2 sin ( x ) + 2 k + cos ( x 3 ) , k = 1 , 2 , … , 6 y=kx^2\sin(x)+2k+\cos(x^3), k=1,2,\dots ,6 y=kx2sin(x)+2k+cos(x3),k=1,2,…,6
x = np.linspace(-10, 10, 100) col = ['tomato', 'orange', 'red', 'c', 'lime', 'm'] for i in range(6): k = i + 1 plt.plot(x, (k *x**2*np.sin(x) + 2 * k + np.cos(x ** 3)), color=col[i], label='k={0}'.format(k)) plt.legend(ncol=2) plt.show()
- 9 把屏幕开成 2行 3 列 6 个子窗口，每个子窗口画一条曲线，画出曲线：
y = k x 2 sin ( x ) + 2 k + cos ( x 3 ) , k = 1 , 2 , … , 6 y=kx^2\sin(x)+2k+\cos(x^3), k=1,2,\dots ,6 y=kx2sin(x)+2k+cos(x3),k=1,2,…,6
# 利用循环绘制 fig,pic=plt.subplots(2,3,sharex=True,sharey=True) x = np.linspace(-10, 10, 1000) k=1 for i in range(2): for j in range(3): picture=pic[i,j] picture.plot(x,(k*x**2*np.sin(x)+2*k+np.cos(x**3)),'b',label='k={0}'.format(k)) picture.legend(loc='upper left') k+=1 plt.show()
- 2.10 分别画出下列二次曲面
- 单叶双曲面 x 2 8 + y 2 10 − z 2 6 = 1 ~\frac{x^2}{8}+\frac{y^2}{10}-\frac{z^2}{6}=1 8x2+10y2−6z2=1 ;
参数方程： x 2 a 2 + y 2 b 2 − z 2 c 2 = 1 \frac{x^2}{a^2}+\frac{y^2}{b^2}-\frac{z^2}{c^2}=1 a2x2+b2y2−c2z2=1
{ x = a cosh θ cos ϕ y = b cosh θ sin ϕ z = c sinh θ \left\{ \begin{aligned} x = & a\cosh \theta\cos\phi\\ y = & b\cosh \theta\sin\phi \\ z = & c \sinh\theta \end{aligned} \right. ⎩ ⎨ ⎧x=y=z=acoshθcosϕbcoshθsinϕcsinhθ
import numpy as np import matplotlib.pyplot as plt from mpl_toolkits import mplot3d plt.rcParams['font.sans-serif'] = ['SimHei'] plt.rcParams['axes.unicode_minus'] = False theta1 = np.linspace(-10, 10, 100) phi1 = np.linspace(-np.pi, np.pi, 100) theta, phi = np.meshgrid(theta1, phi1) x = np.sqrt(8) * np.cosh(theta) * np.cos(phi) y = np.sqrt(10) * np.cosh(theta) * np.sin(phi) z = np.sqrt(6) * np.sinh(theta) aex = plt.subplot(1, 1, 1, projection='3d') aex.plot_surface(x, y, z, cmap='autumn') plt.title('单叶双曲面') plt.show()
- 双叶双曲面 x 2 8 − y 2 12 − z 2 8 = 1 ~\frac{x^2}{8}-\frac{y^2}{12}-\frac{z^2}{8}=1~ 8x2−12y2−8z2=1 ;
参数方程： x 2 a 2 + y 2 b 2 − z 2 c 2 = − 1 \frac{x^2}{a^2}+\frac{y^2}{b^2}-\frac{z^2}{c^2}=-1 a2x2+b2y2−c2z2=−1
{ x = a sinh θ cos ϕ y = b sinh θ sin ϕ z = c cosh θ \left\{ \begin{aligned} x = & a\sinh \theta\cos\phi\\ y = & b\sinh \theta\sin\phi \\ z = & c \cosh\theta \end{aligned} \right. ⎩ ⎨ ⎧x=y=z=asinhθcosϕbsinhθsinϕccoshθ
#也可以尝试用参数方程来绘制 import matplotlib.pyplot as plt from mpl_toolkits import mplot3d plt.rcParams['font.sans-serif']=['SimHei'] plt.rcParams['axes.unicode_minus']=False z=np.linspace(-6,6,30) y=np.linspace(-6,6,30) Z,Y=np.meshgrid(z,y) X=np.sqrt(2/3*Y**2+Z**2+8) ax=plt.subplot(1,1,1,projection='3d') ax.plot_surface(-X, Y, Z, cmap='summer')# X取负值也可以-对称 plt.title("双叶双曲面") plt.show()
- 椭圆双曲面 x 2 10 + y 2 6 = z ~\frac{x^2}{10}+\frac{y^2}{6}=z~ 10x2+6y2=z .
import matplotlib.pyplot as plt from mpl_toolkits import mplot3d plt.rcParams['font.sans-serif']=['SimHei'] plt.rcParams['axes.unicode_minus']=False x=np.linspace(-6,6,30) y=np.linspace(-6,6,30) X,Y=np.meshgrid(x,y) Z= X**2/10+Y**2/6 ax=plt.subplot(1,1,1,projection='3d') ax.plot_surface(X,Y,Z,cmap='rainbow') plt.title("椭圆抛物面") plt.show()
- 2.11 默比乌斯带是一种拓扑学结构，他只有一个面和一个边界，是1858年由德国数学家，天文学假默比乌斯和约翰.李斯丁独立发现的。其参数方程为
{ x = ( 2 + s 2 cos x 2 ) cos t y = ( 2 + s 2 cos x 2 ) sin t z = s 2 sin t 2 \left\{ \begin{aligned} x = &(2+\frac{s}{2}\cos\frac{x}{2})\cos t \\ y = & (2+\frac{s}{2}\cos\frac{x}{2})\sin t \\ z = & \frac{s}{2}\sin\frac{t}{2} \end{aligned} \right. ⎩ ⎨ ⎧x=y=z=(2+2scos2x)cost(2+2scos2x)sint2ssin2t
其中， 0 ≤ t ≤ 2 π , − 1 ≤ s ≤ 1 0\le t\le 2\pi,~-1\le s \le 1 0≤t≤2π, −1≤s≤1.绘制默比乌斯带import numpy as np import matplotlib.pyplot as plt from mpl_toolkits import mplot3d plt.rcParams['font.sans-serif'] = ['SimHei'] plt.rcParams['axes.unicode_minus'] = False t = np.linspace(0, 2 * np.pi, 60) s = np.linspace(-1, 1, 60) S, T = np.meshgrid(s, t) x = (2 + (S / 2) * np.cos(T / 2)) * np.cos(T) y = (2 + (S / 2) * np.cos(T / 2)) * np.sin(T) z = (s / 2) * np.sin(T / 2) aex=plt.subplot(1,1,1,projection='3d') aex.plot_surface(x,y,z,cmap='winter') plt.title('默比乌斯带') plt.show()
- 2.11 画出如下函数的等高线，并进行标注：
- z = x e − x 2 − y 2 , 2 ≤ x ≤ 2 , − 2 ≤ y ≤ 3 z=xe^{-x^2-y^2},~2\le x \le 2,-2 \le y \le3 z=xe−x2−y2, 2≤x≤2,−2≤y≤3
import numpy as np import matplotlib.pyplot as plt plt.rcParams['font.sans-serif'] = ['SimHei'] plt.rcParams['axes.unicode_minus'] = False x = np.arange(-2, 2.001, 0.01) y = np.arange(-2, 3.001, 0.01) # 将x,y扩充为等维度的二维数组 X = np.array([x.tolist()] * y.shape[0]) print(y) print('-' * 30) Y = np.array([y.tolist()] * x.shape[0]) Y = Y.T # 根据函数构造z Z = X * np.e ** (-(X ** 2 + Y ** 2)) # 绘制等高线 con = plt.contour(X, Y, Z) plt.clabel(con) # 标注 plt.xlabel('$x$') plt.ylabel('$y$') plt.ylim([-2,2]) plt.xlim([-1.8,1.8]) plt.show()
- z = ( 1 − x 2 − y 2 ) e − y 2 / 3 , − 1.5 ≤ x ≤ 2 , − 1.5 ≤ y ≤ 2 z=(1-x^2-y^2)e^{-y^2/3},~-1.5\le x \le 2,-1.5\le y \le2 z=(1−x2−y2)e−y2/3, −1.5≤x≤2,−1.5≤y≤2
#将上述的X,Y,Z进行修改 Z =(1-X**2-Y**2)*np.e**(-y**3/3)
- 2.13 画出附件1：区域高程数据.xlsx 画出该区域的三维网格图和等高线图，
在 A ( 30 , 0 ) 和 B （ 43 ， 30 ） A(30,0)和B（43，30） A(30,0)和B（43，30）单位（km）点出建立了两个基地，在等高线上标注出 这两点
#绘制高程图 plt.rcParams['font.sans-serif'] = ['SimHei'] plt.rcParams['axes.unicode_minus'] = False data=pd.read_excel(r'D:\《python数学建模》程序数据\程序及数据\程序及数据\02第2章 数据处理与可视化\附件1：区域高程数据.xlsx') shape=data.shape i,j=shape[0],shape[1] x=np.linspace(0,43.65,j) y=np.linspace(0,58.2,i) z=np.array(data.values.tolist()) #绘制等高线图 plt.subplot(1,1,1) con=plt.contour(x,y,z) plt.clabel(con) plt.text(30,0,'A',fontsize=10,color='r' ) plt.text(43,30,'B',fontsize=10, color='b') plt.title('1.区域高程图') plt.show()
#绘制三维曲面图(利用上述处理好的数据) aex=plt.subplot(1,1,1,projection='3d') X,Y=np.meshgrid(x,y) aex.plot_surface(X,Y,z,cmap='viridis') plt.title("2.三维高程图") plt.show()
9460
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
