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))
    for i in range(N_points_1):
        for j in range(N_points_2): 
    sys = []
    ##Boundary Condition x = 0
    for i in range(N_points_1):
        BC1 = T_N[i,0] - T_x_0

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

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

    ## 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] 

    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()

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

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.

  • $\begingroup$ Do you want the transient solution or only the final temperature distribution? $\endgroup$
    – Bob
    Commented Mar 9, 2023 at 5:48
  • $\begingroup$ The entire temperature distribution which is what i get in fsolve but not on the others even with a good initial guess $\endgroup$
    – Klaus3
    Commented Mar 9, 2023 at 14:26
  • $\begingroup$ What is the physical interpretation of what you get with fsolve? are you making $\partial T / \partial t = 0$ $\endgroup$
    – Bob
    Commented Mar 9, 2023 at 14:29
  • $\begingroup$ 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. $\endgroup$
    – Klaus3
    Commented Mar 9, 2023 at 14:33


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.