Решите нелинейное уравнение numpy.

Вот решение, использующее интерполяцию 1D для вычисления обратной функции F(h). Поскольку используется не стандартный метод поиска корня, ошибка не контролируется, и сетка дискретизации должна выбираться с осторожностью. Однако интерполированная обратная функция может быть непосредственно вычислена по массиву.

примечание: определение F изменено, проблема теперь Solve h for F(h) = TR

import numpy as np
from scipy.interpolate import interp1d
import matplotlib.pylab as plt

# The function to inverse:
t = 10
alfa = 1.1*10**(-7)
k = 0.18
T1 = 20
Tpow = 100

def F(h):
    A = np.exp(h**2*alfa*t/k**2)
    B = h**3*2/(3*np.sqrt(3))*(alfa*t)**(3/2)/k**3
    return -(Tpow-T1)*( 1 - A + B )

# Interpolation 
h_eval = np.linspace(40, 100, 50)   # critical step: define the discretization grid
F_inverse = interp1d( F(h_eval), h_eval, kind='cubic', bounds_error=True )


# Some random data:
TR = np.array([[13, 10, 13, 13, 13],
       [ 9, 11, 10, 12,  9],
       [13, 13, 10, 10, 13],
       [12, 11, 12,  9, 11],
       [11, 13, 13, 11, 13]])

# Compute the array h for a given array TR
h = F_inverse(TR)
print(h)

# Graph to verify the interpolation 
plt.plot(h_eval, F(h_eval), '.-', label='discretized F(h)');
plt.plot(h.ravel(), TR.ravel(), 'or', label='interpolated values')
plt.xlabel('h'); plt.ylabel('F(h) or TR'); plt.legend();

С помощью другой функции изменяются следующие строки:

from scipy.special import erf
def F(h):
    return (Tpow-T1)*(1-np.exp((h**2*alfa*t)/k**2)*(1.0-erf(h*np.sqrt(alfa*t)/k)))

# Interpolation 
h_eval = np.linspace(15, 35, 50)   # the range is changed