Discussion :

[goiko]

bonsoir a tous
SVP j'ai besoin d'aide pour la resolution de l'equation d'ecoulement transitoire des milieu poreu.
pour la discretisation j'ai utiliser la methode des volume finie ce qui m'a donné:
H(i,j) = H(i,j ) + (T/S)*(dt/(2*dx*dx))*(H(i+1,j+1)+H(i+1,j-1)-4*H(i,j)+H(i-1,j+1)+H(i-1,j-1))
le domaine d'etude est un domaine borné un carré 0<x<xmax et 0<y<ymax.
les conditions initiales sont:
a t=0 H(i,j)=y0
pour les condition aux limites on prend sur les 4 cotés les conditions de neumann
les constantes sont:
T= 300
S= 0.002
dt = 0.01
dx=0.1
y0=100
Merci bcp .

[FR119492]

Salut!
Ton problème est posé un peu n'importe comment. En effet, si
H(i,j) = H(i,j ) + (T/S)*(dt/(2*dx*dx))*(H(i+1,j+1)+H(i+1,j-1)-4*H(i,j)+H(i-1,j+1)+H(i-1,j-1))
alors on a nécessairement
(T/S)*(dt/(2*dx*dx))*(H(i+1,j+1)+H(i+1,j-1)-4*H(i,j)+H(i-1,j+1)+H(i-1,j-1))=0
ce qui est un peu étonnant.
Alors, pose ton problème proprement et on essayera de t'aider.
Jean-Marc Blanc

[goiko]

au fait l'equation donné est l'equation discrétisé en volume finie a 2 dimensions pour l'espace et 1 dimension pour le temps.donc le H(i,j) a gauche de = est calculé pour le temps n+1 et celui a droite pour un temps n (j'ai pris le schema d'euler pour le temps).j'espere ke vous comprenez

[goiko]

l'equation en question en volume finie est
dH/dt= (T/S)*(1/(2*dx*dx))*(H(i+1,j+1)+H(i+1,j-1)-4*H(i,j)+H(i-1,j+1)+H(i-1,j-1))
donc pour le schema d'euler pour le temps l'equation devient
H(i,j,n+1) = H(i,j ,n) + (T/S)*(dt/(2*dx*dx))*(H(i+1,j+1,n)+H(i+1,j-1,n)-4*H(i,j,n)+H(i-1,j+1,n)+H(i-1,j-1,n))

[Ehouarn]

l'equation en question en volume finie est
dH/dt= (T/S)*(1/(2*dx*dx))*(H(i+1,j+1)+H(i+1,j-1)-4*H(i,j)+H(i-1,j+1)+H(i-1,j-1))
donc pour le schema d'euler pour le temps l'equation devient
H(i,j,n+1) = H(i,j ,n) + (T/S)*(dt/(2*dx*dx))*(H(i+1,j+1,n)+H(i+1,j-1,n)-4*H(i,j,n)+H(i-1,j+1,n)+H(i-1,j-1,n))

Bonjour,

Cela ressemble plutôt à une formulation en différences finies d'un laplacien.
Tu devrais préciser l'équation différentielle de départ, parce que je doute de ta version discrétisée (de plus, tu est visiblement en 2D, mais n'apparait que du dx dans ton expression. Pas de dy ? ).

De même, il faut tenir compte des conditions aux limites de ton problème (dirichlet, Neumann, mixte,... fonction du temps?).
Quant à la discrétisation temporelle par schéma d'Euler explicite, c'est hélas plutôt à proscrire (peu précis et peu stable), il faut en général utiliser des schémas plus élaborés.

[goiko]

Merci a tous j'ai pu ecrire le code et il me donne des bons resultats
je vous donne le code:
program serge

C resolution de l'equation d'ecoulement transitoire en milieu poreux
implicit double precision (a-h,o-z)
parameter(Nmax=100)
parameter( Nx=11)
parameter(Ny=11)
dimension H(Nmax,Nmax),A(Nmax,Nmax)


dx=20
dy=10
dt =0.001
c=0.02
y0=100
T=300
S=0.002
mt=1/(dt)
c création des fichiers
open(UNIT=10,FILE='temps0,01j.DAT')
write(10,*)'variable:"X" , "Y" , "H"'
open(UNIT=9,FILE='temps0,025j.DAT')
write(9,*)'variable:"X" , "Y" , "H"'
open(UNIT=8,FILE='temps0,05j.DAT')
write(8,*)'variable:"X" , "Y" , "H"'
open(UNIT=7,FILE='temps0,1j.DAT')
write(7,*)'variable:"X" , "Y" , "H"'
open(UNIT=6,FILE='temps0,5j.DAT')
write(6,*)'variable:"X" , "Y" , "H"'
open(UNIT=5,FILE='temps1j.DAT')
write(5,*)'variable:"X" , "Y" , "H"'
c condition initiale
do 10 i=1,Nx
do 10 j=1,Ny
H(i,j)=y0
10 continue

C boucle de calcul de la solution approché

do 30 k=1,mt

C conditions aux limites
do i=1,Nx
H(i,1)= 0.5*H(i,2)
enddo
do j=1,Ny
H(Nx,j)=0.5*H(Nx-1,j)
enddo
do j=1,Ny
H(1,j)= 0.5*H(2,j)
enddo
do i=2,Nx
x=(i-1)*dx
H(i,Ny)= c*x+y0
enddo
do i=2,Nx
x= (i-1)*dx
do j=2,Ny
y= (j-1)*dy
A(i,j) = H(i,j )+ (T/S)*(dt/(dx*dx))*(H(i+1,j)+H(i-1,j)-4*H(i,j)
+H(i,j+1)+H(i,j-1))

if(k==900)then
print *,x,y,H(i,j)
write(10,*)x,y,H(i,j)
end if
if(k==960)then
print *,x,y,H(i,j)
write(9,*)x,y,H(i,j)
end if
if(k==980)then
print *,x,y,H(i,j)
write(8,*)x,y,H(i,j)
end if
if(k==990)then
print *,x,y,H(i,j)
write(7,*)x,y,H(i,j)
end if
if(k==998)then
print *,x,y,H(i,j)
write(6,*)x,y,H(i,j)
end if
if(k==999)then
print *,x,y,H(i,j)
write(5,*)x,y,H(i,j)
end if

enddo
enddo

C affectation
do i= 1,Nx
do j=1,Ny
H(i,j)=A(i,j)
enddo
enddo
30 continue


end program serge



Maintenant j'ai un nouveau problème je suis débutant et je dois tracer les isocharges pour un certain temps donné et je ne maitrise pas bien tecplot alors si vous vous y connaissez aidez moi a trouvez comment faire pour tracer les isocharges pour un temps par exemple t=0.01
Merci d'avance

[Leflix]

Salut!
Ton problème est posé un peu n'importe comment. En effet, si
H(i,j) = H(i,j ) + (T/S)*(dt/(2*dx*dx))*(H(i+1,j+1)+H(i+1,j-1)-4*H(i,j)+H(i-1,j+1)+H(i-1,j-1))
alors on a nécessairement
(T/S)*(dt/(2*dx*dx))*(H(i+1,j+1)+H(i+1,j-1)-4*H(i,j)+H(i-1,j+1)+H(i-1,j-1))=0
ce qui est un peu étonnant.
Alors, pose ton problème proprement et on essayera de t'aider.
Jean-Marc Blanc

Jean Marc soit un peu moins acerbe avec les débutants ou ceux qui posent des questions...
Son problème n'est absolument pas posé n'importe comment comme tu le dis si brutalement...
au vu de son équation discrète il est évident qu'il cherchait à résoudre dH/dt=(T/S)*Lap(H) ou Lap est l'operateur laplacien , avec un schema d'Euler d'ordre 1. D'où le H(i,j)=H(i,j)+......
je suis surpris qu'un expert comme toi ne s'en soit pas rendu compte au premier coup d'oeil...

Ceci dit goiko, je suis d'accord avec Ehouarn, ton programme donne certainement maintenant des resultats mais je doute que ce soit les bons. la discretisation de ton equation aux dérivees partielles ne sera valable que si dx=dy.
Comme le fait remarquer aussi justement Ehouarn, tu dois ecrire pour tes conditions aux limites de Neumann...
H(1,j)=H(2,j)
H(Nx,j)=H(Nx-1,j)
H(i,1)=H(i,2)

de plus....

si ton domaine est [0,xmax]x[0,ymax] et que tu commences à numéroter tes noeuds de 1 jusqu'à NX et 1 jusqu'à NY
alors dx=xmax/(NX-1) et dy=ymax/(NY-1)
ensuite tu fixes si tu veux NX et NY et le programme calculera les bons dx et dy. Dans ton cas la valeur que tu donnes à dx devrait être égal à celle de dy et ces valeurs laissent supposer que xmax=200 et ymax=100..
Il faudra veiller à ce qu'en fonction de xmax et ymax qui définissent ton pb et en fonction de NX et NY que tu donnes, tu aies dx=dy sinon ton equation discrete:
A(i,j) = H(i,j )+ (T/S)*(dt/(dx*dx))*(H(i+1,j)+H(i-1,j)-4*H(i,j)
+H(i,j+1)+H(i,j-1))

est fausse...

enfin à cause de ton implicit double precision (a-h,o-m) mt sera considéré comme un entier car ça commence par "m"
fait attention à ce que tu as écrit: mt=1/dt avec dt=0,001
écris plutôt mt= nint (1/dt) nint= "l'entier le plus proche" car sinon tu peux avoir de mauvaises surprises avec certains compilateurs... J'imagine que tu t'attends à faire 1000 itérations en temps sur un temps de simulation qui vaut 1s...mais là tu ne feras que 999 itérations....mais bon ça ce n'est pas très grave...

[FR119492]

Salut!

Jean Marc soit un peu moins acerbe avec les débutants ou ceux qui posent des questions...

Mon but n'était pas de massacrer un débutant, mais de contribuer à sa formation. Je constate en effet de plus en plus des lacunes en ce qui concerne la méthodologie chez les débutants (et parfois chez leurs professeurs). Cette méthodologie a été définie par G. Polya dans son excellent livre: "Comment poser et résoudre un problème" (angl.: "How to solve it").

1. Comprendre le problème.
2. Concevoir un plan.
3. Mettre le plan à exécution.
4. Examiner la solution obtenue.

Dans le premier message que goiko a envoyé, il ne commence pas comme il aurait dû, c'est-à-dire par les équations aux dérivées partielles qu'il doit intégrer qui font intervenir dH/dt, d2H/dx2 et d2H/dy2, mais ce qu'il a obtenu en appliquant la méthode des différences finies et la méthode d'Euler. Or la cause de ses difficultés est dû au choix de ces méthodes et non à des erreurs dans la suite des calculs.
Jean-Marc Blanc

[Leflix]

Salut!

Dans le premier message que goiko a envoyé, il ne commence pas comme il aurait dû, c'est-à-dire par les équations aux dérivées partielles qu'il doit intégrer qui font intervenir dH/dt, d2H/dx2 et d2H/dy2, mais ce qu'il a obtenu en appliquant la méthode des différences finies et la méthode d'Euler. Or la cause de ses difficultés est dû au choix de ces méthodes et non à des erreurs dans la suite des calculs.
Jean-Marc Blanc

Non je ne suis pas d'accord avec toi, les méthodes choisies ne sont pas en cause ici. La methode des differences finies et le schema d'euler en particulier permettent tres bien de resoudre son EDP. Les erreurs provenaient comme l'a pointé Ehouarn de l'implementation des conditions aux limite et quelques inconsistances au niveau de dx et dy et des indices de boucles...

Son EDP qui est une equation parabolique se resoud tres bien avec ce type de methode et de schema. Absolument aucun probleme ! Stabillité, précision aucun soucis !!!

Et j'emploie deliberement le mot "resoudre" au lieu "d'integrer" car peu importe ça c'est de la semantique Jean Marc
Et je pense que ce n'est pas ce genre de remarque que les gens qui posent des questions ici ont envie d'entendre...ça ne fait pas avancer leur problème.