Bod Blog

Equation de la chaleur en coordonnées polaires


        # -*- coding: utf-8 -*-
        """
        Created on Mon Nov 20 20:37:28 2023

        @author: THOMAS-SCBoduch
        """

        # plot a nice polar function
        # Try with animation


        import math
        import matplotlib.pyplot as plt
        import numpy as np
        from matplotlib.animation import FuncAnimation

        num_r=120
        num_theta=25
        max_iter=2000
        radius=1
        delta_r = radius / num_r
        delta_theta = 2*np.pi / num_theta

        alpha=1.1
        delta_t = (delta_r ** 2)/(4 * alpha)
        gamma = (alpha * delta_t) / (delta_r ** 2)


        ##########################Def of the grid#############################

        r=np.linspace(0, radius, num_r)
        theta=np.linspace(0,2*np.pi, num_theta)
        R,Theta = np.meshgrid(r,theta)

        #T=np.ones((num_theta,num_r))
        T=np.ones((num_theta,num_r))

        ######################################################################
        ####Def of the plot
        fig, ax1  = plt.subplots(subplot_kw=dict(projection='polar'))
        #plt.rcParams['pcolor.shading'] ='nearest'
        quad = ax1.pcolormesh(Theta,R,T,shading='gouraud')

        ######################################################################

        #other setups,
        #ax1.plot(Theta, R , T, color='r', ls='none', marker='.')

        #def des fonctions initiales


        #conditions aux limites
        T[:,-1]=20.0
        T[:,0]=25.0
        T[:,1]=23.0
        #T[0,0]=25.0
        #T[0,1]=23.0
        #T[-1,0]=25.0
        #T[-1,1]=23.0
        #Conditions initiales

        def animate(frame):
            global T, quad
            T_temp = T.copy()  # Initialisez T_temp en dehors de la boucle

            for i in range(1, num_r - 1):
                r_i = i * delta_r
                T_temp[:, i] = T[:, i] + delta_t * ((alpha / delta_r**2) * (T[:, i + 1] - 2 * T[:, i] + T[:, i - 1]) + (alpha / r_i / delta_r) * (T[:, i + 1] - T[:, i]))
            #print(T_temp[:,1])
            # Copiez les conditions initiales de la frame précédente
            T_temp[:, 0] = T[:, 0]
            T_temp[:, 1] = T[:, 1]
            T_temp[:, -1] = T[:, -1]

            #T_temp[0, 0] = T[0,0]
            #T_temp[0, 1] = T[0,1]
            #T_temp[-1, 0] = T[-1,0]
            #T_temp[-1, 1] = T[-1,1]
            T = T_temp.copy()

            #quad.set_array(T.flatten())

            # Mise à jour du plot
            ax1.clear()  # Effacer l'axe pour éviter la superposition des images
            quad = ax1.pcolormesh(Theta, R, T, shading='gouraud', cmap=plt.cm.jet, vmin=-10 , vmax=50)
            ax1.set_title('{:.1f} ms'.format(frame*delta_t*1000))

            return quad,

        """ Laplacien steady state R
        #champ scalaire de temperature
        T=np.ones((num_theta,num_r))

        #conditions aux limites
        T[:,-1]=20.0
        #Conditions initiales
        T[:,0]=25.0
        T[:,1]=23.0


        for i in range(1,num_r-1):
            r_i = i * delta_r
            T[:,i+1]= (T[:,i] * (delta_r + 2*r_i) - r_i * T[:,i-1] ) / ( r_i + delta_r )
            print(T[0,i+1])
        """

        """ Laplacien steady state R and Theta
        #champ scalaire de temperature
        T=np.ones((num_theta,num_r))

        #conditions aux limites
        T[:,-1]=20.0
        #Conditions initiales
        T[:,0]=25.0
        T[:,1]=23.0


        for i in range(1,num_r-1):
            r_i = i * delta_r
            for j in range(1,num_theta-1):
                T[:,i+1]= (T[:,i] * (delta_r + 2*r_i) - r_i * T[:,i-1] ) / ( r_i + delta_r )
                T[j+1,:]= 2*T[j,:] - T[j-1,:]



        #"""

        anim=FuncAnimation(fig,animate,interval=50, frames=max_iter,repeat=False)
        anim.save("heat_equation_solution.gif")
        plt.colorbar(quad)
        plt.show()


      

Plot de l'animation

Initial condition are as follow :

The center of the cylinder is at 25°C (this is where the gas filling is going on). The surface of the cylinder is at 20°C (ambiant temperature)

Bod Blog