Since you helped me so much, you deserve my full code! haha
(I'm not sure that it's fully correct (that is, whether it actually solve the exercise I was trying to do), but it does work and produces something!)
Enjoy:
#Fichier TP3
#Méthode des différences finies
#Equation de la chaleur
from numpy import *
from numpy.linalg import *
from matplotlib.pyplot import *
from mpl_toolkits.mplot3d import Axes3D
#Définition des fonctions
def uex(x,y):
z=sin(2*pi*x)*sinh(2*pi*y)
return z
def ud(x):
y=sin(2*pi*x)*sinh(2*pi)
return y
#Maillage
print('Choix du nombre Ns de points intérieurs du maillage')
Ns = int(input('Ns = '))
h = 1/(Ns+2)
X=linspace(0,1,Ns+2) #maillage complet
#Définition de la matrice A et du vecteur B
u=-1*ones(Ns-1)
D=4*eye(Ns)+diag(u,1)+ diag(u,-1)
E=-eye(Ns)
N=Ns*Ns
A=array([0]*N*N).reshape(N,N)
for i in range (0,N,Ns):
for j in range (i,N-Ns,Ns):
A[i:i+Ns,i:i+Ns]=D
A[i:i+Ns,i+Ns:i+2*Ns]=E
A[i+Ns:i+2*Ns,i:i+Ns]=E
A[N-Ns:N,N-Ns:N]=D
print('A = ', A)
B=zeros((N,1))
C=zeros((Ns,1))
for i in range (0,Ns):
C[i]=ud(X[i])
B[N-Ns:N]=C
#Définition de Uh
Uh=dot(inv(A),B)
#Définition de Zh
Zh=Uh.reshape(Ns,Ns)
print("Zh = ",Zh)
#Rectification des bords de Zh
for i in range (0,Ns):
for j in range (0,Ns):
Zh[i,Ns-1]=ud(X[i])
Zh[0,j]=0
Zh[Ns-1,j]=0
Zh[i,0]=0
print ('Zh = ',Zh)
#Définition de la solution exacte
Uex=zeros((Ns,Ns))
for i in range (0,Ns):
for j in range (0,Ns):
Uex[i,j]=uex(X[i],X[j])
print('Uex = ', Uex)
#Calcul de l'erreur
err=amax(Zh-Uex)
print("Erreur = ", err)
#Affichage des résultats
coordX=zeros(Ns)
coordY=zeros(Ns)
for k in range (0,Ns):
for i in range (0,Ns):
coordX[k]=X[k]
coordY[i]=X[i]
fig = figure()
ax = Axes3D(fig)
coordX, coordY = np.meshgrid(coordX,coordY)
#ax.plot_surface(coordX, coordY, abs(Zh-Uex), rstride=1, cstride=1, cmap='hot')
ax.plot_surface(coordX, coordY, Zh, rstride=1, cstride=1, cmap='hot')
#ax.plot_surface(coordX, coordY, Uex,rstride=1, cstride=1, cmap='hot')
show()
I let you see that output yourself!