我目前正试图用Python解决我的流体力学问题,该问题导致了三个不同的方程。经过一些研究,我试图使它尽可能简单,这产生了以下代码:
from scipy.optimize import fsolve
import math
def equations(p):
d, Re, f = p
return (d**5 - 0.0165*f, 97200 - Re*d, f**(-0.5) + 2*math.log10(0.00015/(3.7*d) + 2.51/(Re*f**(0.5))))
d, Re, f = fsolve(equations, (0.22, 0.0213, 490000))
print(equations((d, Re, f)))
这给出了一个长错误:
> ValueError Traceback (most recent call
> last) <ipython-input-19-dcc7fc844bb6> in <module>()
> 6 return (d**5 - 0.0165*f, 97200 - Re*d, f**(-0.5) + 2*math.log10(0.00015/(3.7*d) + 2.51/(Re*f**(0.5))))
> 7
> ----> 8 d, Re, f = fsolve(equations, (0.22, 0.0213, 490000))
> 9
> 10 print(equations((d, Re, f)))
>
> ~Anaconda3libsite-packagesscipyoptimizeminpack.py in
> fsolve(func, x0, args, fprime, full_output, col_deriv, xtol, maxfev,
> band, epsfcn, factor, diag)
> 146 'diag': diag}
> 147
> --> 148 res = _root_hybr(func, x0, args, jac=fprime, **options)
> 149 if full_output:
> 150 x = res['x']
>
> ~Anaconda3libsite-packagesscipyoptimizeminpack.py in
> _root_hybr(func, x0, args, jac, col_deriv, xtol, maxfev, band, eps, factor, diag, **unknown_options)
> 225 with _MINPACK_LOCK:
> 226 retval = _minpack._hybrd(func, x0, args, 1, xtol, maxfev,
> --> 227 ml, mu, epsfcn, factor, diag)
> 228 else:
> 229 _check_func('fsolve', 'fprime', Dfun, x0, args, n, (n, n))
>
> <ipython-input-19-dcc7fc844bb6> in equations(p)
> 4 def equations(p):
> 5 d, Re, f = p
> ----> 6 return (d**5 - 0.0165*f, 97200 - Re*d, f**(-0.5) + 2*math.log10(0.00015/(3.7*d) + 2.51/(Re*f**(0.5))))
> 7
> 8 d, Re, f = fsolve(equations, (0.22, 0.0213, 490000))
>
> ValueError: math domain error
现在,我知道它一定与对数值有关,但我认为方程的复杂性也可能在这里发挥作用。
您可以在equations(p)
函数中添加print(d, Re, f)
,并查看值在每次迭代中的变化情况。
我对它进行了测试,这些是错误之前的结果
0.22 0.0213 490000.0
0.22 0.0213 490000.0
0.22 0.0213 490000.0
0.22000000327825547 0.0213 490000.0
0.22 0.02130000031739473 490000.0
0.22 0.0213 490000.007301569
4384243.368454826 -464.7214227491007 3096116.988633934
正如您所看到的,在计算时,最后一行值实际上对f**(-0.5) + 2*math.log10(0.00015/(3.7*d) + 2.51/(Re*f**(0.5)))
操作无效,因为0.00015/(3.7*d) + 2.51/(Re*f**(0.5))
的计算结果为-3.069524039032724e-06
,接近0,因此被识别为0。因此,由于log(0)
本身无效,它返回Domain Error
要解决此问题,您可能需要重新访问您的公式或扩展有效位数字的数量,以容纳更多数字并使其不再为0。