# Scipy.root not converging even when provided with initial guesses very close to solution

I've made a previous question here and also in SO wondering why only the fsolve solver converges for the simple one dimensional unsteady conduction problem

$$\frac{\partial T}{\partial t} = \alpha \frac{\partial^2T}{\partial x^2}$$

$$\alpha = 0.001$$

$$T(x=0,t) = 60°C$$

$$T(x=L,t) = 25°C$$

$$T(x,t=0) = 25°C$$

One of the answers alluded to initial values, from which i started to test and, to make things easier, i changed the value of alpha to

$$\alpha = 0.000001$$

Which means that the result profile shall be very close to simply

$$T(x=0,t) = 60°C$$

$$T(x=L,t) = 25°C$$

$$T(x,t) = 25°C$$

Which is what i used as initial condition for this new code

import matplotlib.pyplot as plt
import numpy as np
import sympy as sp
from scipy.optimize import fsolve
from scipy import optimize
import pandas as pd

t_max = 25
x_max = 1
N_points_1 = 100
N_points_2 = 11
dt = t_max/(N_points_1-1)
dx = x_max/(N_points_2-1)
N_variables = N_points_1 * N_points_2
alpha = 0.000001
T_x_0 = 60
T_x_L = 25
T_t_0 = 25

def EDP (v):
T_N =  np.ones((N_points_1, N_points_2))

#Variables
k=0
for i in range(N_points_1):
for j in range(N_points_2):
T_N[i,j]=v[k]
k=k+1

#Model
sys = []
##Boundary Condition x = 0
for i in range(N_points_1):
BC1 = T_N[i,0] - T_x_0
sys.append(BC1)

##Boundary Condition x = L
for i in range(N_points_1):
BC2 = T_N[i,-1] - T_x_L
sys.append(BC2)

##Initial Condition
for j in range(1,N_points_2-1):
IC1 = T_N[0,j] - T_t_0
sys.append(IC1)

## Energy Balance
for i in range(N_points_1-1):
for j in range (1,N_points_2 -1):
EB = alpha*(T_N[i,j+1]-2*T_N[i,j]+T_N[i,j-1])/dx**2 - (T_N[i+1,j]-T_N[i,j])/dt  #partial(T_N,t_1,1)[i,j] - alpha*partial_2(T_N,x,2)[i,j]
sys.append(EB)

sys=np.array(sys)

return sys

initial_guess = 25*np.ones((N_points_1,N_points_2))
initial_guess[:,0] = T_x_0
initial_guess[:,-1] = T_x_L
initial_guess = initial_guess.flatten()
print(np.sum(EDP(initial_guess)))
0.3464999999999999


As you can see the value of the sum of the values of each equation is already close to 0

fsolve does converge, but even so, it takes a long number of iterations. Others either do not converge or straigth up diverge

Solution_T = fsolve(EDP,initial_guess,maxfev = 1200)

Solution_T2 = optimize.anderson(EDP,initial_guess, verbose = True,maxiter = 1000)
104:  |F(x)| = 2.36365e+93; step 1
OverflowError: (34, 'Result too large')

Solution_T3 = optimize.newton_krylov(EDP,initial_guess,verbose=True,maxiter = 1000)
999:  |F(x)| = 0.00372199; step 1
NoConvergence: [60.         24.99999991 24.99999972 ... 24.99990242 24.99977307
24.99968681]

Solution_T4 = optimize.broyden1(EDP,initial_guess,verbose=True,maxiter = 1000)
659:  |F(x)| = 4.37351e+151; step 1
OverflowError: (34, 'Result too large')

Solution_T5 = optimize.broyden2(EDP,initial_guess,verbose=True,maxiter = 1000)
999:  |F(x)| = 2.81888e+10; step 1.74723e-09
NoConvergence: [ 6.00000000e+01 -2.63924280e+06 -2.98161396e+04 ...  3.39690111e+08
2.14121070e+07 -2.54528649e+08]


So, what actually is happening here? I can only converge these algorithms when i feed the solution of fsolve into them, from which the objective function is already below the tolerance for convergence. I know that fsolve did converge, but i am just running tests for much larger system of equations, from which the large scale solvers, those above besides fsolve, are required.

• Do you want the transient solution or only the final temperature distribution?
– Bob
Commented Mar 9, 2023 at 5:48
• The entire temperature distribution which is what i get in fsolve but not on the others even with a good initial guess Commented Mar 9, 2023 at 14:26
• What is the physical interpretation of what you get with fsolve? are you making $\partial T / \partial t = 0$
– Bob
Commented Mar 9, 2023 at 14:29
• So basically by putting alpha to be a very small value i am getting pretty close to a situation with $\frac{\partial T}{\partial t} = 0$ which means that the initial profile should be very close to the final profile(not the steady state one, but on each dynamic step). It means that the effective thermal conductivity is so low that there is virtually no heat flux between the high 60 degree temperature end to the 25 temperature end. Commented Mar 9, 2023 at 14:33