pinyue / May 11 2019
Remix of Python by Nextjournal

数据溢出的解决过程

Remix this to get started with Python .

import sys; sys.version.split()[0]
0.6s
Python
Python 3
'3.6.8'

1. 数据溢出的问题

1.1. 查看运行环境

查看本文件的环境配置（在Bash环境中进行）

pip list
1.7s
Bash in Python
Python 3

1.2. 编写程序实现算法时出现数据溢出问题

要实现的算法

x=0.777rijωμρx = 0.777r_{i} \sqrt{\frac{j \omega \mu}{\rho}}
result=coth(x)=cosh(x)sinh(x)result = coth(x) = \frac{cosh(x)}{sinh(x)}
import numpy as np
def coth(r_i,f,mu,Rho):
    omega = 2*np.pi*f
    temp = omega*mu/Rho
    temp1 = 1j*temp
    m = np.sqrt(temp1)
    x = 0.777*m*r_i
    result = np.cosh(x)/np.sinh(x)
    return result
0.7s
Python
Python 3
f = 10**10
r_i = 5.9*0.001
mu = 4*np.pi*10**(-7)
Rho = 172.41*10**(-10)
result = coth(r_i,f,mu,Rho)
result
0.7s
Python
Python 3
(nan+nanj)

当f的取值大于10**9时，就会发生溢出

2. 解决数据溢出问题

安装mpmath解决数据溢出问题

pip install mpmath
33.5s
Bash in Python
Python 3

再次查看本文件的环境配置（在Bash环境中进行），说明安装mpmath已经成功安装在此运行环境中。

pip list
1.5s
Bash in Python
Python 3
import numpy as np
import mpmath as mp
def coth(r_i,f,mu,Rho):
    omega = 2*np.pi*f
    temp = omega*mu/Rho
    temp1 = 1j*temp
    # m = math.sqrt(temp1)
    m = np.sqrt(temp1)
    print(m)
    x = 0.777*m*r_i
    print(x)
    print(mp.cosh(x))
    print(mp.sinh(x))
    result = mp.cosh(x)/mp.sinh(x)
    return result
0.7s
Python
Python 3
r_i = 5.9*0.001
f = 10**20
mu = 4*np.pi*10**(-7)
Rho = 172.41*10**(-10)
result = coth(r_i,f,mu,Rho)
result
0.7s
Python
Python 3
mpc(real='1.0', imag='0.0')

以上说明数据溢出的问题已经解决。

3. numpy和mpmath的结合问题

当计算的数据需要以numpy的array的形式输入和输出时，会出现报错：TypeError:cannot create mpf from array

import numpy as np
import mpmath as mp
def coth(x):
    return (mp.exp(x)+mp.exp(-x))/(mp.exp(x)-mp.exp(-x))


def calc_Zc(f,r,mu,rho):
    m=np.sqrt(1j*2*np.pi*f*mu/rho)
    # n=np.shape(c_xy)[0]
    n = 14
    Zc=np.empty(n,np.complex128)
    Zc=rho*m/(2*np.pi*r)*coth(0.777*m*r)+0.356*rho/(np.pi*r**2)
    return Zc
0.5s
Python
Python 3
f=5000000
conductors_calc_radius=0.001*np.array([5.9,7.00,9.5,109.1,109.1,7.60,5.35,5.9,7.00,9.5,109.1,109.1,7.60,5.35])
r=conductors_calc_radius          # 单位 m
#mu=4*np.pi*10**-7     #
mu=np.array([4*np.pi*10**-7,4*np.pi*10**-7,4*np.pi*10**-7,8.575*4*np.pi*10**-7,8.575*4*np.pi*10**-7,4*np.pi*10**-7,4*np.pi*10**-7,4*np.pi*10**-7,4*np.pi*10**-7,4*np.pi*10**-7,8.575*4*np.pi*10**-7,8.575*4*np.pi*10**-7,4*np.pi*10**-7,4*np.pi*10**-7])

conductors_rho=10**-10*np.array([172.41,146.15,298.96,1666.7,1666.7,352.74,251.77,172.41,146.15,298.96,1666.7,1666.7,352.74,251.77])

rho=conductors_rho
Zc=calc_Zc(f,r,mu,rho)
print(Zc*10**3)
print(Zc.tolist())
0.4s
Python
Python 3

4. 解决问题--np.frompyfunc

pip install mpmath
4.5s
Bash in Python
Python 3
import mpmath as mp
import numpy as np
def coth(x):
    #coth = (mp.exp(x)+mp.exp(-x))/(mp.exp(x)-mp.exp(-x))
    coth = mp.cosh(x)/mp.sinh(x)
    return coth

coth_vec=np.frompyfunc(coth,1,1)

def calc_Zc(f,r,mu,rho):
    m=np.sqrt(1j*2*np.pi*f*mu/rho)
    # n=np.shape(c_xy)[0]
    n = 14
    Zc=np.empty(n,np.complex128)
    Zc=rho*m/(2*np.pi*r)*coth_vec(0.777*m*r)+0.356*rho/(np.pi*r**2)
    return Zc
0.4s
Python
Python 3
f=5000000
conductors_calc_radius=0.001*np.array([5.9,7.00,9.5,109.1,109.1,7.60,5.35,5.9,7.00,9.5,109.1,109.1,7.60,5.35])
r=conductors_calc_radius          # 单位 m
#mu=4*np.pi*10**-7     #
mu=np.array([4*np.pi*10**-7,4*np.pi*10**-7,4*np.pi*10**-7,8.575*4*np.pi*10**-7,8.575*4*np.pi*10**-7,4*np.pi*10**-7,4*np.pi*10**-7,4*np.pi*10**-7,4*np.pi*10**-7,4*np.pi*10**-7,8.575*4*np.pi*10**-7,8.575*4*np.pi*10**-7,4*np.pi*10**-7,4*np.pi*10**-7])

conductors_rho=10**-10*np.array([172.41,146.15,298.96,1666.7,1666.7,352.74,251.77,172.41,146.15,298.96,1666.7,1666.7,352.74,251.77])

rho=conductors_rho
Zc=calc_Zc(f,r,mu,rho)
print(Zc)
0.7s
Python
Python 3

关于np.frompyfunc的用法参考如下：