Monsieur,
Je travaille sur un projet de modélisation de la montée en vitesse d'une pompe centrifuge.
J'ai trouvé un algorithme sur fortran mais je n'arrive pas à le comprendre surtout les variables.
Si quelqu'un peut m'aider et c'est urgent.
merci
Code:****************************** l'algorithme************ Methode des caractéristiques ! Exemple integer N1,N2,N3 real R,rp,kL,la,kv,L1,L3,kk,delta,L 2,kvv dimension V1(2001),Vp1(2001),H1(2001),Hp 1(2001) dimension V2(2001),Vp2(2001),H2(2001),Hp 2(2001),kL(2001) dimension X(2001),Dh(2001) dimension V3(2001),Vp3(2001),H3(2001),Hp 3(2001) open(1,file='C:\pression.dat') open(2,file='C:\vitesse.dat') ! *************** données de la pompe ***************** pi=3.14 ! N2=5 R=0.064 ! R=0.084 rp=0.019 Nm=1500 omo=Nm*pi/30 b=0.0045 c=900 kvv=8 teta=2*pi/kvv g=9.81 tp=0.15 la=0.02 !6 L2=R-rp ! ************* données des conduites amont et aval********************** D1=0.04 L1=0.5 !0.5 D3=0.04 L3=1.5 !0.5 ! write(*,*)omo,la,teta,1/2 ! stop ! *************** calcul du diametre hydraulique pompe **************************** dt=0.000005 N2=L2/dt/c dx2=(R-rp)/N2 ! write(*,*)dx2,N2 ! stop do 10 i=1,N2+1 kL(i)=(rp+(i-1)*dx2)*teta Dh(i)=4*b*kL(i)/(2*kL(i)+2*b) X(i)=rp+(i-1)*dx2 ! write(*,*) Dh(i),x(i),dh(i)/dh(1) 10 continue ! **************** calcul régime permanant pompe *************************** ! stop do 20 i=1,N2+1 H2(i)=0 V2(i)=0 20 continue ! dt=dx2/c N1=L1/(c*dt) ! L1=c*dt*N1 N3=L3/(c*dt) ! L3=c*dt*N3 ! write(*,*) N1,N3,L1,L3,N2 ! stop do 21 i=1,N1+1 H1(i)=0 V1(i)=0 21 continue do 22 i=1,N3+1 H3(i)=0 V3(i)=0 22 continue ! ****************************** *Calcul Du Régime Transitoire t=0 ! write(1,*) t,H3(N3+1) ! write(2,*) t,V2(N2+1) do 31 j=1,150000 t=t+dt if (t.ge.tp) then om=omo else om=(omo/tp)*t end if ! ******* Calcul de la conduite amont 1 ****************************** do 25 i=2, N1 Vp1(i)=0.5*((g/c)*(H1(i-1)-H1(i+1))+V1(i+1)+V1(i-1)+(g/c)*dt* 1(V1(i+1)-V1(i-1))-(la*dt/(2*D1))*(V1(i-1)*abs(V1(i-1))+V1(i+1)* 1abs(V1(i+1)))) Hp1(i)=(c/(2*g))*(V1(i-1)-V1(i+1)-(g/c)*dt*(V1(i-1)+ 1V1(i+1))+(la*dt/(2*D1))*(-V1(i-1) 1*abs(V1(i-1))+V1(i+1)*abs(V1(i+1))))+(H1 (i-1)+H1(i+1))/2 25 continue ! ******** calcul au niveau de la pompe 2 ***************************** do 26 i=2, N2 ! Vp2(i)=0.5*((g/c)*(H2(i-1)-H2(i+1))+V2(i+1)+V2(i-1)+ ! 1(g/c)*dt*(V2(i+1)-V2(i-1))-(la*dt/2)* ! 1(V2(i-1)*abs(V2(i-1))/Dh(i-1)+V2(i+1)*abs(V2(i+1))/Dh(i+1)) ! 1+(om**2)*dt*(X(i-1)+X(i+1))) ! vp2(i)=vp2(i-1)*kl(i-1)/kl(i) ! Hp2(i)=(c/(2*g))*(V2(i-1)-V2(i+1)-(g/c)*dt* ! 1(V2(i+1)+V2(i-1))+(la*dt/2)*(-V2(i-1)*abs(V2(i-1)) ! 1/Dh(i-1)+V2(i+1)*abs(V2(i+1))/Dh(i+1))+(om**2)*dt*(X(i-1)-X(i+1))) ! 1+(H2(i-1)+H2(i+1))/2 C1=g/c xl1=C1*H2(i-1)+V2(i-1)-(la*dt/2)*V2(i-1)*abs(V2(i-1))/Dh(i-1)+ 1(om**2)*dt*X(i)-C1*dt*V2(i-1) xl2=-C1*H2(i+1)+V2(i+1)-(la*dt/2)*V2(i+1)*abs(V2(i+1))/Dh(i+1)+ 1(om**2)*dt*X(i)+C1*dt*V2(i+1) Hp2(i)=(xl2-xl1)/(-2*C1) vp2(i)=xl1-C1*Hp2(i) ! write(*,*) t,xl1 ! stop 26 continue ! stop ! ************ calcul de la conduite aval 3 ****************************** * do 27 i=2,N3 Vp3(i)=0.5*((g/c)*(H3(i-1)-H3(i+1))+V3(i+1)+V3(i-1)+(g/c)*dt* 1(V3(i+1)-V3(i-1))-(la*dt/(2*D3)) 1*(V3(i-1)*abs(V3(i-1))+V3(i+1)*abs(V3(i+1)))) Hp3(i)=(c/(2*g))*(V3(i-1)-V3(i+1)-(g/c)*dt*(V3(i+1)+ 1V3(i-1))+(la*dt/(2*D3))*(-V3(i-1) 1*abs(V3(i-1))+V3(i+1)*abs(V3(i+1))))+(H3 (i-1)+H3(i+1))/2 27 continue ! ********** condition au limite à l'extrémité amont de la conduite 1 ********* Hp1(1)=0 C1=(g/c) C2=C1*H1(2)-V1(2)+(la*dt/(2*D1))*V1(2)*abs(V1(2))-C1*dt*V1(2) ! A=C1*(2.14E-8)*D1**4 ! B=1-C1*(5.26E-6)*D1**2 ! C11=C2-C1*9.2384 ! delta1=B*B-4*A*C11 ! write(*,*) t,H3(N3+1) ! stop - ! Vp1(1)=(-B+sqrt(delta1))/2/A Vp1(1)=C1*Hp1(1)-C2 ! Hp1(1)=(-Vp1(1)-C2)/C1 ! if (Hp1(1).lt.-10.33) goto 123 ! goto 124 ! 123 Hp1(1)=-10.33 !124 continue ! *********** jonction entre conduite 1 aval et pompe extrémité amont *************** s1=(pi*D1**2)/4 s3=(pi*D3**2)/4 ! sp=2*pi*rp*b ! surface latérale aspiration ! spp=2*pi*R*b ! ,, ,, refoulement Dx4=4*b*(kvv*kL(1))/(2*kvv*kL(1)+2*b) Dx5=4*b*(kvv*kL(N2+1))/(2*kvv*kL(N2+1)+2*b) sp=(pi*dx4**2)/4 spp=(pi*dx5**2)/4 ! write(*,*) dx2,dx4,dx5 ! stop kk=(sp+s1)/s1 Vp2(1)=(1/kk)*((g/c)*(H1(N1)-H2(2))+V2(2)+V1(N1)+(om**2)*dt *X(1)- 1(la*dt/2)*(V2(2)*abs(V2(2))/Dh(2)+V1(N1)*abs(V1(N1))/D1)+ 1(g/c)*dt*(-V1(N1)+V2(2))) Vp1(N1+1)=Vp2(1)*sp/s1 Hp1(N1+1)=(c/g)*(V1(N1)-Vp1(N1+1)-(la*dt/(2*D1))* 1V1(N1)*abs(V1(N1)))+H1(N1)-V1(N1)*dt Hp2(1)=Hp1(N1+1) ! write(*,*) 1/kk,Vp2(1),(spp+s3)/s3,(spp+s3)/spp ! stop ! *********** jonction entre conduite 3 amont et pompe extrémité aval *************** kk=(spp+s3)/s3 Vp2(N2+1)=(1/kk)*((g/c)*(H2(N2)-H3(2))+V2(N2)+V3(2)+(om**2)*dt * 1X(N2+1)+(g/c)*dt*(V3(2)-V2(N2))-(la*dt/2)*(V2(N2)*abs(V2(N2))/ 1Dh(N2)+V3(2)*abs(V3(2))/D3)) Vp3(1)=Vp2(N2+1)*spp/s3 Hp3(1)=(c/g)*(Vp3(1)-V3(2)+(la*dt/(2*D3))* 1V3(2)*abs(V3(2)))+H3(2)+V3(2) *dt Hp2(N2+1)=Hp3(1) ! *********** jonction entre conduite 3 aval et vanne *************** kv=0.15 !1000 ! perte de charge singulière vanne cv=1/sqrt(1+kv) cc=0.8 !0.5 ! coefficient de contaction ag=1*s3 ! ouverture de la vanne bv=cv*cc*ag ! coefficient cd=cv*cc ck=((s3/bv)**2)*(1./(2*g)) bb=g/c c3=(la*dt/(2*D3))*V3(N3)*abs(V3(N3))-bb*H3(N3)-V3(N3)+bb*V3(N3)*dt delta=1-4*bb*ck*c3 Vp3(N3+1)=(-1+sqrt(delta))/(2*ck*bb) ! Vp3(N3+1)=0 Hp3(N3+1)=(1/bb)*(-c3-Vp3(N3+1)) ! write(*,*) t,Hp3(N3+1),1/bb,delta ! stop do 28 i=1, N1+1 H1(i)=Hp1(i) V1(i)=Vp1(i) 28 continue do 29 i=1, N2+1 H2(i)=Hp2(i) V2(i)=Vp2(i) 29 continue do 30 i=1, N3+1 H3(i)=Hp3(i) V3(i)=Vp3(i) 30 continue ! if (IU/Kkk*Kkk-IU) 31,300,31 if (j/30*30-j) 39,300,39 300 write(1,*) t,H3(1)-H1(N1+1) write(2,*) t,V2(N2+1)*spp*1000 39 continue 31 continue end
-----