From: Felipe Abrão on 5 Jul 2010 15:37 Please, could somebody explain me what is happening. I simply allocated a matrix called "a", but figured out that the program is inserting some strange values to some of its elements (like 1.057229231618489E-307 or 6.941570408926931E+228), even before I set any value to it! Shouldn't all the values be zero before any desired value setting? Below is my code, but you don't have to read it all. You can check the problem just after line 65, where I wrote a command in order to read the values of matrix "a". The same thing is happening with matrixs "b" and "c". In advance I thank you all! PS: I use Fortran 90 in Windows XP. ! EXERCICIO.f90 ! ! FUNCTIONS: ! EXERCICIO - Entry point of console application. ! ! **************************************************************************** ! ! PROGRAM: EXERCICIO ! ! PURPOSE: Entry point for the console application. ! ! **************************************************************************** program EXERCICIO use IMSL implicit none ! Variables integer(4)::i,j,k,l,m,Nx,Ny integer(4),parameter::ipath=1 integer(4),dimension(6)::iparam integer(4),allocatable,dimension(:)::irow_u,jcol_u,irow_v,jcol_v,irow_p,jcol_p real(8)::dx,dy,D_w,D_e,D_s,D_n,D_EE,D_F,D_G,D_H,D_I,D_J,D_K,D_L real(8),parameter::Re=1.0d2,gama=1.0d0/ Re,u_E=1.0d0,u_G=0.0d0,u_H=0.0d0,v_I=0.0d0,v_J=0.0d0,& v_K=0.0d0,v_L=0.0d0,error=1.0d-6,alpha_u=5.0d-1,alpha_v=5.0d-1,alpha_p=5.0d-1 !real(8),dimension(N*N)::theta_old real(8),dimension(5)::rparam real(8),allocatable,dimension(:)::u,v,p,theta,S,S_v,S_p,X,Y,u_F,AA,BB,CC !real(8),dimension(N,N)::theta_fake real(8),allocatable,dimension(:,:)::a,a_inv,b,c,p_fake,p_old,u_fake,u_old,v_fake,v_old,a_P,& a_W,a_E,a_S,a_N,b_P,b_W,b_E,b_S,b_N,c_P,c_W,c_E,c_S,c_N,F_w_u,F_e_u,F_s_u,F_n_u,S_fake,& S_v_fake,S_p_fake,F_w_v,F_e_v,F_s_v,F_n_v ! Body of EXERCICIO_1 iparam(1)=0 Nx=5 Ny=5 allocate(a(Ny*Nx,Ny*Nx),a_inv(Ny*Nx,Ny*Nx),b(Ny*Nx,Ny*Nx),c(Ny*Nx,Ny*Nx)) allocate(a_P(Ny,Nx),a_W(Ny,Nx),a_E(Ny,Nx),a_S(Ny,Nx),a_N(Ny,Nx),b_P(Ny,Nx),b_W(Ny,Nx),& b_E(Ny,Nx),b_S(Ny,Nx),b_N(Ny,Nx),c_P(Ny,Nx),c_W(Ny,Nx),c_E(Ny,Nx),c_S(Ny,Nx),c_N(Ny,Nx),& p_fake(Ny,Nx),p_old(Ny,Nx),u_fake(Ny,Nx),u_old(Ny,Nx),v_fake(Ny,Nx),v_old(Ny,Nx),F_w_u(Ny,Nx),& F_e_u(Ny,Nx),F_s_u(Ny,Nx),F_n_u(Ny,Nx),S_fake(Ny,Nx),S_v_fake(Ny,Nx),S_p_fake(Ny,Nx),& F_w_v(Ny,Nx),F_e_v(Ny,Nx),F_s_v(Ny,Nx),F_n_v(Ny,Nx)) allocate(u(Ny*Nx),v(Ny*Nx),p(Ny*Nx),theta(Ny*Nx),S(Ny*Nx),S_v(Ny*Nx),S_p(Ny*Nx)) allocate(Y(Ny),u_F(Ny)) allocate(X(Nx)) open(7,file='pre.dat') do i=1,Ny*Nx do j=1,Ny*Nx write(7,*) i,j,a(i,j) end do end do dx=1.0d0/dble(Nx) dy=1.0d0/dble(Ny) D_w=gama/dx D_e=gama/dx D_s=gama/dy D_n=gama/dy D_EE=2.0d0*D_w D_F=2.0d0*D_e D_G=2.0d0*D_s D_H=2.0d0*D_n D_I=D_EE D_J=D_F D_K=D_G D_L=D_H p_fake=5.0d-1 u_fake=5.0d-1 v_fake=5.0d-1 ! do i=1,Ny ! do j=1,Nx ! u_fake(i,j)=dble(i+j) ! v_fake(i,j)=dble(i+j) ! p_fake(i,j)=dble(i+j) ! end do ! end do theta=5.0d-1 !theta_old=6.0d-1 do i=1,Nx X(i)=1.0d0/dble((Nx-1))*dble((i-1)) end do do i=1,Ny Y(i)=1.0d0/dble((Ny-1))*dble((i-1)) u_F(i)=4.0d0*(Y(i)-Y(i)**2.0d0) end do p_old=1.0d0 u_old=1.0d0 v_old=1.0d0 ! do while(maxval(dabs(p_fake-p_old)).gt.error) !do while((maxval(dabs(p_fake-p_old)).gt.error).or.(maxval(dabs(u_fake- u_old)).gt.error).or.& !(maxval(dabs(v_fake-v_old)).gt.error)) !do while((maxval(dabs(p_fake(2:Ny-1,2:Nx-1)- p_old(2:Ny-1,2:Nx-1))).gt.error).or.(maxval(dabs(& !u_fake(2:Ny-1,2:Nx-1)-u_old(2:Ny-1,2:Nx-1))).gt.error).or. (maxval(dabs(v_fake(2:Ny-1,2:Nx-1)-& !v_old(2:Ny-1,2:Nx-1))).gt.error)) ! write(*,*) 'diferença de pressão =',maxval(dabs(p_fake-p_old)) ! write(*,*) 'diferença de velocidade u =',maxval(dabs(u_fake-u_old)) ! write(*,*) 'diferença de velocidade v =',maxval(dabs(v_fake-v_old)) ! p_old=p_fake ! u_old=u_fake ! v_old=v_fake do i=1,Ny do j=1,Nx if (j.eq.1) then F_w_u(2:Ny-1,j)=1.0d0 else F_w_u(i,j)=1.0d0/2.0d0*(u_fake(i,j)+u_fake(i,j-1)) end if if (i.eq.1) then F_w_v(i,j)=0.0d0 !Tenho dúvida. else F_w_v(i,j)=1.0d0/2.0d0*(u_fake(i,j)+u_fake(i-1,j)) end if if (j.eq.Nx) then F_e_u(i,j)=4.0d0*(Y(i)-Y(i)**2.0d0) else F_e_u(i,j)=1.0d0/2.0d0*(u_fake(i,j+1)+u_fake(i,j)) end if if (i.eq.1) then F_e_v(i,j)=0.0d0 !Tenho dúvida. else if (j.eq.Nx) then F_e_v(i,j)=4.0d0*(Y(i)-Y(i)**2.0d0) !Tenho dúvida. else F_e_v(i,j)=1.0d0/2.0d0*(u_fake(i,j+1)+u_fake(i-1,j+1)) end if if (i.eq.1) then F_s_u(i,j)=0.0d0 else if (j.eq.1) then F_s_u(i,j)=0.0d0 else F_s_u(i,j)=1.0d0/2.0d0*(v_fake(i,j)+v_fake(i,j-1)) end if if (i.eq.1) then F_s_v(i,j)=0.0d0 else F_s_v(i,j)=1.0d0/2.0d0*(v_fake(i-1,j)+v_fake(i,j)) end if if (i.eq.Ny) then F_n_u(i,j)=0.0d0 else if (j.eq.1) then F_n_u(i,j)=0.0d0 else F_n_u(i,j)=1.0d0/2.0d0*(v_fake(i+1,j)+v_fake(i+1,j-1)) end if if (i.eq.Ny) then F_n_v(i,j)=0.0d0 else F_n_v(i,j)=1.0d0/2.0d0*(v_fake(i,j)+v_fake(i+1,j)) end if end do end do !interior do i=3,Ny-1 do j=3,Nx-1 !coeficientes de u a_P(i,j)=(F_e_u(i,j)*dy+F_n_u(i,j)*dx+D_w*dy+D_e*dy+D_s*dx+D_n*dx)/ alpha_u a_W(i,j)=F_w_u(i,j)*dy+F_e_u(i,j)*1.0d0/8.0d0*dy+D_w*dy a_E(i,j)=D_e*dy a_S(i,j)=F_s_u(i,j)*dx+F_n_u(i,j)*1.0d0/8.0d0*dx+D_s*dx a_N(i,j)=D_n*dx S_fake(i,j)=(3.0d0/8.0d0*u_fake(i,j)*dy-2.0d0/8.0d0*u_fake(i,j-1)*dy-1.0d0/8.0d0& *u_fake(i,j-2)*dy)*F_w_u(i,j)+(-3.0d0/8.0d0*u_fake(i,j+1)*dy +2.0d0/8.0d0*& u_fake(i,j)*dy)*F_e_u(i,j) +(3.0d0/8.0d0*u_fake(i,j)*dx-2.0d0/8.0d0*u_fake(i-1,j)& *dx-1.0d0/8.0d0*u_fake(i-2,j)*dx)*F_s_u(i,j)+(-3.0d0/8.0d0*u_fake(i +1,j)*dx+& 2.0d0/8.0d0*u_fake(i,j)*dx)*F_n_u(i,j)+(p_fake(i,j-1)- p_fake(i,j))*dy+((1.0d0-& alpha_u)*a_P(i,j))*u_old(i,j) !coeficientes de v b_P(i,j)=(F_e_v(i,j)*dy+F_n_v(i,j)*dx+D_w*dy+D_e*dy+D_s*dx+D_n*dx)/ alpha_v b_W(i,j)=F_w_v(i,j)*dy+F_e_v(i,j)*1.0d0/8.0d0*dy+D_w*dy b_E(i,j)=D_e*dy b_S(i,j)=F_s_v(i,j)*dx+F_n_v(i,j)*1.0d0/8.0d0*dx+D_s*dx b_N(i,j)=D_n*dx S_v_fake(i,j)=(3.0d0/8.0d0*v_fake(i,j)*dy-2.0d0/8.0d0*v_fake(i,j-1)*dy-1.0d0/ & 8.0d0*v_fake(i,j-2)*dy)*F_w_v(i,j)+(-3.0d0/8.0d0*v_fake(i,j+1)*dy +2.0d0/8.0d0*& v_fake(i,j)*dy)*F_e_v(i,j) +(3.0d0/8.0d0*v_fake(i,j)*dx-2.0d0/8.0d0*v_fake(i-1,j)& *dx-1.0d0/8.0d0*v_fake(i-2,j)*dx)*F_s_v(i,j)+(-3.0d0/8.0d0*v_fake(i +1,j)*dx+& 2.0d0/8.0d0*v_fake(i,j)*dx)*F_n_v(i,j)+(p_fake(i-1,j)- p_fake(i,j))*dx+((1.0d0-& alpha_v)*b_P(i,j))*v_old(i,j) !coeficientes de p c_W(i,j)=dy/a_P(i,j)*dy c_E(i,j)=dy/a_E(i,j)*alpha_u*dy c_S(i,j)=dx/b_P(i,j)*dx c_N(i,j)=dx/b_N(i,j)*alpha_v*dx c_P(i,j)=c_W(i,j)+c_E(i,j)+c_S(i,j)+c_N(i,j) S_p_fake(i,j)=u_fake(i,j)*dy-u_fake(i,j+1)*dy+v_fake(i,j)*dx- v_fake(i+1,j)*dx end do end do !coluna 2 - com exceção da linha 2 do i=3,Ny-1 !coeficientes de u a_P(i,2)=(F_e_u(i,2)*dy+F_n_u(i,2)*dx+D_w*dy+D_e*dy+D_s*dx+D_n*dx)/ alpha_u a_W(i,2)=F_w_u(i,2)*7.0d0/8.0d0*dy+F_e_u(i,2)*1.0d0/8.0d0*dy+D_w*dy a_E(i,2)=D_e*dy a_S(i,2)=F_s_u(i,2)*dx+F_n_u(i,2)*1.0d0/8.0d0*dx+D_s*dx a_N(i,2)=D_n*dx S_fake(i,2)=(3.0d0/8.0d0*u_fake(i,2)*dy-2.0d0/8.0d0*u_E*dy)*F_w_u(i, 2)+(-3.0d0/8.0d0& *u_fake(i,3)*dy+2.0d0/8.0d0*u_fake(i,2)*dy)*F_e_u(i, 2)+(3.0d0/8.0d0*u_fake(i,2)*dx-& 2.0d0/8.0d0*u_fake(i-1,2)*dx-1.0d0/8.0d0*u_fake(i-2,2)*dx)*F_s_u(i, 2)+(-3.0d0/8.0d0*& u_fake(i+1,2)*dx+2.0d0/8.0d0*u_fake(i,2)*dx)*F_n_u(i,2)+(p_fake(i, 1)-p_fake(i,2))*dy& +((1.0d0-alpha_u)*a_P(i,2))*u_old(i,2) !coeficientes de v b_P(i,2)=(F_e_v(i,2)*dy+F_n_v(i,2)*dx+D_w*dy+D_e*dy+D_s*dx+D_n*dx)/ alpha_v b_W(i,2)=F_w_v(i,2)*7.0d0/8.0d0*dy+F_e_v(i,2)*1.0d0/8.0d0*dy+D_w*dy b_E(i,2)=D_e*dy b_S(i,2)=F_s_v(i,2)*dx+F_n_v(i,2)*1.0d0/8.0d0*dx+D_s*dx b_N(i,2)=D_n*dx S_v_fake(i,2)=(3.0d0/8.0d0*v_fake(i, 2)*dy-2.0d0/8.0d0*v_I*dy)*F_w_v(i,2)+(-3.0d0/& 8.0d0*v_fake(i,3)*dy+2.0d0/8.0d0*v_fake(i,2)*dy)*F_e_v(i, 2)+(3.0d0/8.0d0*v_fake(i,2)& *dx-2.0d0/8.0d0*v_fake(i-1,2)*dx-1.0d0/8.0d0*v_fake(i-2,2)*dx)*F_s_v(i, 2)+(-3.0d0/& 8.0d0*v_fake(i+1,2)*dx+2.0d0/8.0d0*v_fake(i,2)*dx)*F_n_v(i,2)+ (p_fake(i-1,2)-& p_fake(i,2))*dx+((1.0d0-alpha_v)*b_P(i,2))*v_old(i,2) !coeficientes de p c_W(i,2)=dy/a_P(i,2)*dy c_E(i,2)=dy/a_E(i,2)*alpha_u*dy c_S(i,2)=dx/b_P(i,2)*dx c_N(i,2)=dx/b_N(i,2)*alpha_v*dx c_P(i,2)=c_W(i,2)+c_E(i,2)+c_S(i,2)+c_N(i,2) S_p_fake(i,2)=u_fake(i,2)*dy-u_fake(i,3)*dy+v_fake(i,2)*dx-v_fake(i +1,2)*dx end do !linha 2 - com exceção da coluna 2 do j=3,Nx-1 !coeficientes de u a_P(2,j)=(F_e_u(2,j)*dy+F_n_u(2,j)*dx+D_w*dy+D_e*dy+D_s*dx+D_n*dx)/ alpha_u a_W(2,j)=F_w_u(2,j)*dy+F_e_u(2,j)*1.0d0/8.0d0*dy+D_w*dy a_E(2,j)=D_e*dy a_S(2,j)=F_s_u(2,j)*7.0d0/8.0d0*dx+F_n_u(2,j)*1.0d0/8.0d0*dx+D_s*dx a_N(2,j)=D_n*dx S_fake(2,j)=(3.0d0/8.0d0*u_fake(2,j)*dy-2.0d0/8.0d0*u_fake(2,j-1)*dy-1.0d0/8.0d0*& u_fake(2,j-2)*dy)*F_w_u(2,j)+(-3.0d0/8.0d0*u_fake(2,j+1)*dy +2.0d0/8.0d0*u_fake(2,j)*& dy)*F_e_u(2,j) +(3.0d0/8.0d0*u_fake(2,j)*dx-2.0d0/8.0d0*u_G*dx)*F_s_u(2,j)+(-3.0d0/& 8.0d0*u_fake(3,j)*dx+2.0d0/8.0d0*u_fake(2,j)*dx)*F_n_u(2,j)+ (p_fake(2,j-1)-& p_fake(2,j))*dy+((1.0d0-alpha_u)*a_P(2,j))*u_old(2,j) !coeficientes de v b_P(2,j)=(F_e_v(2,j)*dy+F_n_v(2,j)*dx+D_w*dy+D_e*dy+D_s*dx+D_n*dx)/ alpha_v b_W(2,j)=F_w_v(2,j)*dy+F_e_v(2,j)*1.0d0/8.0d0*dy+D_w*dy b_E(2,j)=D_e*dy b_S(2,j)=F_s_v(2,j)*7.0d0/8.0d0*dx+F_n_v(2,j)*1.0d0/8.0d0*dx+D_s*dx b_N(2,j)=D_n*dx S_v_fake(2,j)=(3.0d0/8.0d0*v_fake(2,j)*dy-2.0d0/8.0d0*v_fake(2,j-1)*dy-1.0d0/8.0d0*& v_fake(2,j-2)*dy)*F_w_v(2,j)+(-3.0d0/8.0d0*v_fake(2,j+1)*dy +2.0d0/8.0d0*v_fake(2,j)*& dy)*F_e_v(2,j) +(3.0d0/8.0d0*v_fake(2,j)*dx-2.0d0/8.0d0*v_K*dx)*F_s_v(2,j)+(-3.0d0/& 8.0d0*v_fake(3,j)*dx+2.0d0/8.0d0*v_fake(2,j)*dx)*F_n_v(2,j)+ (p_fake(1,j)-& p_fake(2,j))*dx+((1.0d0-alpha_v)*b_P(2,j))*v_old(2,j) !coeficientes de p c_W(2,j)=dy/a_P(2,j)*dy c_E(2,j)=dy/a_E(2,j)*alpha_u*dy c_S(2,j)=dx/b_P(2,j)*dx c_N(2,j)=dx/b_N(2,j)*alpha_v*dx c_P(2,j)=c_W(2,j)+c_E(2,j)+c_S(2,j)+c_N(2,j) S_p_fake(2,j)=u_fake(2,j)*dy-u_fake(2,j+1)*dy+v_fake(2,j)*dx- v_fake(3,j)*dx end do !coluna 2 - linha 2 !coeficientes de u a_P(2,2)=(F_e_u(2,2)*dy+F_n_u(2,2)*dx+D_w*dy+D_e*dy+D_s*dx+D_n*dx)/ alpha_u a_W(2,2)=F_w_u(2,2)*7.0d0/8.0d0*dy+F_e_u(2,2)*1.0d0/8.0d0*dy +D_w*dy a_E(2,2)=D_e*dy a_S(2,2)=F_s_u(2,2)*7.0d0/8.0d0*dx+F_n_u(2,2)*1.0d0/8.0d0*dx +D_s*dx a_N(2,2)=D_n*dx S_fake(2,2)=(3.0d0/8.0d0*u_fake(2,2)*dy-2.0d0/8.0d0*u_E*dy)*F_w_u(2,2)+ (-3.0d0/& 8.0d0*u_fake(2,3)*dy +2.0d0/8.0d0*u_fake(2,2)*dy)*F_e_u(2,2)+(3.0d0/8.0d0*& u_fake(2,2)*dx-2.0d0/8.0d0*u_G*dx)*F_s_u(2,2)+ (-3.0d0/8.0d0*u_fake(3,2)*dx+2.0d0& /8.0d0*u_fake(2,2)*dx)*F_n_u(2,2)+(p_fake(2,1)-p_fake(2,2))*dy+ ((1.0d0-alpha_u)*& a_P(2,2))*u_old(2,2) !coeficientes de v b_P(2,2)=(F_e_v(2,2)*dy+F_n_v(2,2)*dx+D_w*dy+D_e*dy+D_s*dx+D_n*dx)/ alpha_v b_W(2,2)=F_w_v(2,2)*7.0d0/8.0d0*dy+F_e_v(2,2)*1.0d0/8.0d0*dy +D_w*dy b_E(2,2)=D_e*dy b_S(2,2)=F_s_v(2,2)*7.0d0/8.0d0*dx+F_n_v(2,2)*1.0d0/8.0d0*dx +D_s*dx b_N(2,2)=D_n*dx S_v_fake(2,2)=(3.0d0/8.0d0*v_fake(2,2)*dy-2.0d0/8.0d0*v_I*dy)*F_w_v(2,2)+ (-3.0d0& /8.0d0*v_fake(2,3)*dy +2.0d0/8.0d0*v_fake(2,2)*dy)*F_e_v(2,2)+(3.0d0/8.0d0*& v_fake(2,2)*dx-2.0d0/8.0d0*v_K*dx)*F_s_v(2,2)+ (-3.0d0/8.0d0*v_fake(3,2)*dx+2.0d0& /8.0d0*v_fake(2,2)*dx)*F_n_v(2,2)+(p_fake(1,2)-p_fake(2,2))*dx+ ((1.0d0-alpha_v)*& b_P(2,2))*v_old(2,2) !coeficientes de p c_W(2,2)=dy/a_P(2,2)*dy c_E(2,2)=dy/a_E(2,2)*alpha_u*dy c_S(2,2)=dx/b_P(2,2)*dx c_N(2,2)=dx/b_N(2,2)*alpha_v*dx c_P(2,2)=c_W(2,2)+c_E(2,2)+c_S(2,2)+c_N(2,2) S_p_fake(2,2)=u_fake(2,2)*dy-u_fake(2,3)*dy+v_fake(2,2)*dx- v_fake(3,2)*dx !fronteira W - com exceção da linha 2 do i=3,Ny-1 !coeficientes de u a_P(i,1)=(F_e_u(i,1)*7.0d0/8.0d0*dy+F_n_u(i,1)*dx+D_e*dy+D_s*dx +D_n*dx+D_EE*9.0d0/& 3.0d0*dy)/alpha_u a_W(i,1)=1.0d-20 !0.0d0 a_E(i,1)=D_e*dy+D_EE*1.0d0/3.0d0*dy a_S(i,1)=F_s_u(i,1)*dx+F_n_u(i,1)*1.0d0/8.0d0*dx+D_s*dx a_N(i,1)=D_n*dx S_fake(i,1)=(u_E*dy)*F_w_u(i,1)+(-3.0d0/8.0d0*u_fake(i,2)*dy +2.0d0/8.0d0*u_E*dy)*& F_e_u(i,1)+(3.0d0/8.0d0*u_fake(i, 1)*dx-2.0d0/8.0d0*u_fake(i-1,1)*dx-1.0d0/8.0d0*& u_fake(i-2,1)*dx)*F_s_u(i,1)+(-3.0d0/8.0d0*u_fake(i+1,1)*dx +2.0d0/8.0d0*u_fake(i,1)*& dx)*F_n_u(i,1)+D_EE*8.0d0/3.0d0*u_E*dy+((1.0d0-alpha_u)*a_P(i, 1))*u_old(i,1) !Falta o termo de pressão. ! a_P(i,1)=1.0d0 ! a_W(i,1)=1.0d-20 ! a_E(i,1)=1.0d-20 ! a_S(i,1)=1.0d-20 ! a_N(i,1)=1.0d-20 ! S_fake(i,1)=u_E !coeficientes de v b_P(i,1)=(F_e_v(i,1)*7.0d0/8.0d0*dy+F_n_v(i,1)*dx+D_e*dy+D_s*dx +D_n*dx+D_I*9.0d0/& 3.0d0*dy)/alpha_v b_W(i,1)=1.0d-20 !0.0d0 b_E(i,1)=D_e*dy+D_I*1.0d0/3.0d0*dy b_S(i,1)=F_s_v(i,1)*dx+F_n_v(i,1)*1.0d0/8.0d0*dx+D_s*dx b_N(i,1)=D_n*dx S_v_fake(i,1)=(v_I*dy)*F_w_v(i,1)+(-3.0d0/8.0d0*v_fake(i,2)*dy +2.0d0/8.0d0*v_I*dy)*& F_e_v(i,1)+(3.0d0/8.0d0*v_fake(i, 1)*dx-2.0d0/8.0d0*v_fake(i-1,1)*dx-1.0d0/8.0d0*& v_fake(i-2,1)*dx)*F_s_v(i,1)+(-3.0d0/8.0d0*v_fake(i+1,1)*dx +2.0d0/8.0d0*v_fake(i,1)*& dx)*F_n_v(i,1)+(p_fake(i-1,1)-p_fake(i,1))*dx+D_I*8.0d0/3.0d0*v_I*dy +((1.0d0-alpha_v& )*b_P(i,1))*v_old(i,1) ! b_P(i,1)=1.0d0 ! b_W(i,1)=1.0d-20 ! b_E(i,1)=1.0d-20 ! b_S(i,1)=1.0d-20 ! b_N(i,1)=1.0d-20 ! S_v_fake(i,1)=v_I !coeficientes de p c_W(i,1)=dy/a_P(i,1)*dy c_E(i,1)=dy/a_E(i,1)*alpha_u*dy c_S(i,1)=dx/b_P(i,1)*dx c_N(i,1)=dx/b_N(i,1)*alpha_v*dx c_P(i,1)=c_W(i,1)+c_E(i,1)+c_S(i,1)+c_N(i,1) S_p_fake(i,1)=u_fake(i,1)*dy-u_fake(i,2)*dy+v_fake(i,1)*dx-v_fake(i +1,1)*dx end do !fronteira W - linha 2 !coeficientes de u a_P(2,1)=(F_e_u(2,1)*7.0d0/8.0d0*dy+F_n_u(2,1)*dx+D_e*dy+D_s*dx +D_n*dx+D_EE*9.0d0/& 3.0d0*dy)/alpha_u a_W(2,1)=1.0d-20 !0.0d0 a_E(2,1)=D_e*dy+D_EE*1.0d0/3.0d0*dy a_S(2,1)=F_s_u(2,1)*7.0d0/8.0d0*dx+F_n_u(2,1)*1.0d0/8.0d0*dx+D_s*dx a_N(2,1)=D_n*dx S_fake(2,1)=(u_E*dy)*F_w_u(2,1)+(-3.0d0/8.0d0*u_fake(2,2)*dy +2.0d0/8.0d0*u_E*dy)*& F_e_u(2,1)+(3.0d0/8.0d0*u_fake(2,1)*dx-2.0d0/8.0d0*u_G*dx)*F_s_u(2,1)+ (-3.0d0/8.0d0*& u_fake(3,1)*dx +2.0d0/8.0d0*u_fake(2,1)*dx)*F_n_u(2,1)+D_EE*8.0d0/3.0d0*u_E*dy+((& 1.0d0-alpha_u)*a_P(2,1))*u_old(2,1) !Falta o termo de pressão. ! a_P(2,1)=1.0d0 ! a_W(2,1)=1.0d-20 ! a_E(2,1)=1.0d-20 ! a_S(2,1)=1.0d-20 ! a_N(2,1)=1.0d-20 ! S_fake(2,1)=u_E !coeficientes de v b_P(2,1)=(F_e_v(2,1)*7.0d0/8.0d0*dy+F_n_v(2,1)*dx+D_e*dy+D_s*dx +D_n*dx+D_I*9.0d0/& 3.0d0*dy)/alpha_v b_W(2,1)=1.0d-20 !0.0d0 b_E(2,1)=D_e*dy+D_I*1.0d0/3.0d0*dy b_S(2,1)=F_s_v(2,1)*7.0d0/8.0d0*dx+F_n_v(2,1)*1.0d0/8.0d0*dx+D_s*dx b_N(2,1)=D_n*dx S_v_fake(2,1)=(v_I*dy)*F_w_v(2,1)+(-3.0d0/8.0d0*v_fake(2,2)*dy +2.0d0/8.0d0*v_I*dy)*& F_e_v(2,1)+(3.0d0/8.0d0*v_fake(2,1)*dx-2.0d0/8.0d0*v_K*dx)*F_s_v(2,1)+ (-3.0d0/8.0d0*& v_fake(3,1)*dx+2.0d0/8.0d0*v_fake(2,1)*dx)*F_n_v(2,1)+(p_fake(1,1)- p_fake(2,1))*dx+& D_I*8.0d0/3.0d0*v_I*dy+((1.0d0-alpha_v)*b_P(2,1))*v_old(2,1) ! b_P(2,1)=1.0d0 ! b_W(2,1)=1.0d-20 ! b_E(2,1)=1.0d-20 ! b_S(2,1)=1.0d-20 ! b_N(2,1)=1.0d-20 ! S_v_fake(2,1)=v_I !coeficientes de p c_W(2,1)=dy/a_P(2,1)*dy c_E(2,1)=dy/a_E(2,1)*alpha_u*dy c_S(2,1)=dx/b_P(2,1)*dx c_N(2,1)=dx/b_N(2,1)*alpha_v*dx c_P(2,1)=c_W(2,1)+c_E(2,1)+c_S(2,1)+c_N(2,1) S_p_fake(2,1)=u_fake(2,1)*dy-u_fake(2,2)*dy+v_fake(2,1)*dx- v_fake(3,1)*dx !fronteira E - com exceção da linha 2 do i=3,Ny-1 !coeficientes de u a_P(i,Nx)=(F_n_u(i,Nx)*dx+D_w*dy+D_s*dx+D_n*dx+D_F*9.0d0/3.0d0*dy)/ alpha_u a_W(i,Nx)=F_w_u(i,Nx)*dy+D_w*dy+D_F*1.0d0/3.0d0*dy a_E(i,Nx)=1.0d-20 !0.0d0 a_S(i,Nx)=F_s_u(i,Nx)*dx+F_n_u(i,Nx)*1.0d0/8.0d0*dx+D_s*dx a_N(i,Nx)=D_n*dx S_fake(i,Nx)=(3.0d0/8.0d0*u_fake(i,Nx)*dy-2.0d0/8.0d0*u_fake(i,Nx-1)*dy-1..0d0/8.0d0*& u_fake(i,Nx-2)*dy)*F_w_u(i,Nx)+(-u_F(i)*dy)*F_e_u(i,Nx) +(3.0d0/8.0d0*u_fake(i,Nx)*dx& -2.0d0/8.0d0*u_fake(i-1,Nx)*dx-1.0d0/8.0d0*u_fake(i-2,Nx)*dx)*F_s_u(i,Nx) +(-3.0d0/& 8.0d0*u_fake(i+1,Nx)*dx+2.0d0/8.0d0*u_fake(i,Nx)*dx)*F_n_u(i,Nx)+ (p_fake(i,Nx-1)-& p_fake(i,Nx))*dy+D_F*8.0d0/3.0d0*u_F(i)*dy+((1.0d0- alpha_u)*a_P(i,Nx))*u_old(i,Nx) ! a_P(i,Nx)=1.0d0 ! a_W(i,Nx)=1.0d-20 ! a_E(i,Nx)=1.0d-20 ! a_S(i,Nx)=1.0d-20 ! a_N(i,Nx)=1.0d-20 ! S_fake(i,Nx)=u_F(i) !coeficientes de v b_P(i,Nx)=(F_n_v(i,Nx)*dx+D_w*dy+D_s*dx+D_n*dx+D_J*9.0d0/3.0d0*dy)/ alpha_v b_W(i,Nx)=F_w_v(i,Nx)*dy+D_w*dy+D_J*1.0d0/3.0d0*dy b_E(i,Nx)=1.0d-20 !0.0d0 b_S(i,Nx)=F_s_v(i,Nx)*dx+F_n_v(i,Nx)*1.0d0/8.0d0*dx+D_s*dx b_N(i,Nx)=D_n*dx S_v_fake(i,Nx)=(3.0d0/8.0d0*v_fake(i,Nx)*dy-2.0d0/8.0d0*v_fake(i,Nx-1)*dy-1.0d0/ & 8.0d0*v_fake(i,Nx-2)*dy)*F_w_v(i,Nx)+(-v_J*dy)*F_e_v(i,Nx) +(3.0d0/8.0d0*v_fake(i,Nx)& *dx-2.0d0/8.0d0*v_fake(i-1,Nx)*dx-1.0d0/8.0d0*v_fake(i-2,Nx)*dx)*F_s_v(i,Nx) +(-3.0d0& /8.0d0*v_fake(i+1,Nx)*dx+2.0d0/8.0d0*v_fake(i,Nx)*dx)*F_n_v(i,Nx)+ (p_fake(i-1,Nx)-& p_fake(i,Nx))*dx+D_J*8.0d0/3.0d0*v_J*dy+((1.0d0- alpha_v)*b_P(i,Nx))*v_old(i,Nx) ! b_P(i,Nx)=1.0d0 ! b_W(i,Nx)=1.0d-20 ! b_E(i,Nx)=1.0d-20 ! b_S(i,Nx)=1.0d-20 ! b_N(i,Nx)=1.0d-20 ! S_v_fake(i,Nx)=v_J !coeficientes de p c_W(i,Nx)=dy/a_P(i,Nx)*dy c_E(i,Nx)=dy/a_E(i,Nx)*alpha_u*dy c_S(i,Nx)=dx/b_P(i,Nx)*dx c_N(i,Nx)=dx/b_N(i,Nx)*alpha_v*dx c_P(i,Nx)=c_W(i,Nx)+c_E(i,Nx)+c_S(i,Nx)+c_N(i,Nx) S_p_fake(i,Nx)=v_fake(i,Nx)*dx-v_fake(i+1,Nx)*dx !Faltam termos de u. end do !fronteira E - linha 2 !coeficientes de u a_P(2,Nx)=(F_n_u(2,Nx)*dx+D_w*dy+D_s*dx+D_n*dx+D_F*9.0d0/3.0d0*dy)/ alpha_u a_W(2,Nx)=F_w_u(2,Nx)*dy+D_w*dy+D_F*1.0d0/3.0d0*dy a_E(2,Nx)=1.0d-20 !0.0d0 a_S(2,Nx)=F_s_u(2,Nx)*7.0d0/8.0d0*dx+F_n_u(2,Nx)*1.0d0/8.0d0*dx +D_s*dx a_N(2,Nx)=D_n*dx S_fake(2,Nx)=(3.0d0/8.0d0*u_fake(2,Nx)*dy-2.0d0/8.0d0*u_fake(2,Nx-1)*dy-1..0d0/8.0d0*& u_fake(2,Nx-2)*dy)*F_w_u(2,Nx)+(-u_F(2)*dy)*F_e_u(2,Nx) +(3.0d0/8.0d0*u_fake(2,Nx)*dx& -2.0d0/8.0d0*u_G*dx)*F_s_u(2,Nx)+(-3.0d0/8.0d0*u_fake(3,Nx)*dx +2.0d0/8.0d0*& u_fake(2,Nx)*dx)*F_n_u(2,Nx)+(p_fake(2,Nx-1)-p_fake(2,Nx))*dy +D_F*8.0d0/3.0d0*u_F(2)& *dy+((1.0d0-alpha_u)*a_P(2,Nx))*u_old(2,Nx) ! a_P(2,Nx)=1.0d0 ! a_W(2,Nx)=1.0d-20 ! a_E(2,Nx)=1.0d-20 ! a_S(2,Nx)=1.0d-20 ! a_N(2,Nx)=1.0d-20 ! S_fake(2,Nx)=u_F(2) !coeficientes de v b_P(2,Nx)=(F_n_v(2,Nx)*dx+D_w*dy+D_s*dx+D_n*dx+D_J*9.0d0/3.0d0*dy)/ alpha_v b_W(2,Nx)=F_w_v(2,Nx)*dy+D_w*dy+D_J*1.0d0/3.0d0*dy b_E(2,Nx)=1.0d-20 !0.0d0 b_S(2,Nx)=F_s_v(2,Nx)*7.0d0/8.0d0*dx+F_n_v(2,Nx)*1.0d0/8.0d0*dx +D_s*dx b_N(2,Nx)=D_n*dx S_v_fake(2,Nx)=(3.0d0/8.0d0*v_fake(2,Nx)*dy-2.0d0/8.0d0*v_fake(2,Nx-1)*dy-1.0d0/ & 8.0d0*v_fake(2,Nx-2)*dy)*F_w_v(2,Nx)+(-v_J*dy)*F_e_v(2,Nx) +(3.0d0/8.0d0*v_fake(2,Nx)& *dx-2.0d0/8.0d0*v_K*dx)*F_s_v(2,Nx)+(-3.0d0/8.0d0*v_fake(3,Nx)*dx +2.0d0/8.0d0*& v_fake(2,Nx)*dx)*F_n_v(2,Nx)+(p_fake(1,Nx)-p_fake(2,Nx))*dx +D_J*8.0d0/3.0d0*v_J*dy+(& (1.0d0-alpha_v)*b_P(2,Nx))*v_old(2,Nx) ! b_P(2,Nx)=1.0d0 ! b_W(2,Nx)=1.0d-20 ! b_E(2,Nx)=1.0d-20 ! b_S(2,Nx)=1.0d-20 ! b_N(2,Nx)=1.0d-20 ! S_v_fake(2,Nx)=v_J !coeficientes de p c_W(2,Nx)=dy/a_P(2,Nx)*dy c_E(2,Nx)=dy/a_E(2,Nx)*alpha_u*dy c_S(2,Nx)=dx/b_P(2,Nx)*dx c_N(2,Nx)=dx/b_N(2,Nx)*alpha_v*dx c_P(2,Nx)=c_W(2,Nx)+c_E(2,Nx)+c_S(2,Nx)+c_N(2,Nx) S_p_fake(2,Nx)=v_fake(2,Nx)*dx-v_fake(3,Nx)*dx !Faltam termos de u. !fronteira S - com exceção da coluna 2 do j=3,Nx-1 !coeficientes de u a_P(1,j)=(F_e_u(1,j)*dy+F_n_u(1,j)*7.0d0/8.0d0*dx+D_w*dy+D_e*dy +D_n*dx+D_G*9.0d0/& 3.0d0*dx)/alpha_u a_W(1,j)=F_w_u(1,j)*dy+F_e_u(1,j)*1.0d0/8.0d0*dy+D_w*dy a_E(1,j)=D_e*dy a_S(1,j)=1.0d-20 !0.0d0 a_N(1,j)=D_n*dx+D_G*1.0d0/3.0d0*dx S_fake(1,j)=(3.0d0/8.0d0*u_fake(1,j)*dy-2.0d0/8.0d0*u_fake(1,j-1)*dy-1.0d0/8.0d0*& u_fake(1,j-2)*dy)*F_w_u(1,j)+(-3.0d0/8.0d0*u_fake(1,j+1)*dy +2.0d0/8.0d0*u_fake(1,j)*& dy)*F_e_u(1,j)+(u_G*dx)*F_s_u(1,j)+(-3.0d0/8.0d0*u_fake(2,j)*dx +2.0d0/8.0d0*u_G*dx)*& F_n_u(1,j)+(p_fake(1,j-1)-p_fake(1,j))*dy+D_G*8.0d0/3.0d0*u_G*dx+ ((1.0d0-alpha_u)*& a_P(1,j))*u_old(1,j) ! a_P(1,j)=1.0d0 ! a_W(1,j)=1.0d-20 ! a_E(1,j)=1.0d-20 ! a_S(1,j)=1.0d-20 ! a_N(1,j)=1.0d-20 ! S_fake(1,j)=u_G !coeficientes de v b_P(1,j)=(F_e_v(1,j)*dy+F_n_v(1,j)*7.0d0/8.0d0*dx+D_w*dy+D_e*dy +D_n*dx+D_K*9.0d0/& 3.0d0*dx)/alpha_v b_W(1,j)=F_w_v(1,j)*dy+F_e_v(1,j)*1.0d0/8.0d0*dy+D_w*dy b_E(1,j)=D_e*dy b_S(1,j)=1.0d-20 !0.0d0 b_N(1,j)=D_n*dx+D_K*1.0d0/3.0d0*dx S_v_fake(1,j)=(3.0d0/8.0d0*v_fake(1,j)*dy-2.0d0/8.0d0*v_fake(1,j-1)*dy-1.0d0/8.0d0*& v_fake(1,j-2)*dy)*F_w_v(1,j)+(-3.0d0/8.0d0*v_fake(1,j+1)*dy +2.0d0/8.0d0*v_fake(1,j)*& dy)*F_e_v(1,j)+(v_K*dx)*F_s_v(1,j)+(-3.0d0/8.0d0*v_fake(2,j)*dx +2.0d0/8.0d0*v_K*dx)*& F_n_v(1,j)+D_K*8.0d0/3.0d0*v_K*dx+((1.0d0- alpha_v)*b_P(1,j))*v_old(1,j) !Falta o termo de pressão. ! b_P(1,j)=1.0d0 ! b_W(1,j)=1.0d-20 ! b_E(1,j)=1.0d-20 ! b_S(1,j)=1.0d-20 ! b_N(1,j)=1.0d-20 ! S_v_fake(1,j)=v_K !coeficientes de p c_W(1,j)=dy/a_P(1,j)*dy c_E(1,j)=dy/a_E(1,j)*alpha_u*dy c_S(1,j)=dx/b_P(1,j)*dx c_N(1,j)=dx/b_N(1,j)*alpha_v*dx c_P(1,j)=c_W(1,j)+c_E(1,j)+c_S(1,j)+c_N(1,j) S_p_fake(1,j)=u_fake(1,j)*dy-u_fake(1,j+1)*dy+v_fake(1,j)*dx- v_fake(2,j)*dx end do !fronteira S - coluna 2 !coeficientes de u a_P(1,2)=(F_e_u(1,2)*dy+F_n_u(1,2)*7.0d0/8.0d0*dx+D_w*dy+D_e*dy+D_n*dx +D_G*9.0d0/3.0d0*dx)/& alpha_u a_W(1,2)=F_w_u(1,2)*7.0d0/8.0d0*dy+F_e_u(1,2)*1.0d0/8.0d0*dy+D_w*dy a_E(1,2)=D_e*dy a_S(1,2)=1.0d-20 !0.0d0 a_N(1,2)=D_n*dx+D_G*1.0d0/3.0d0*dx S_fake(1,2)=(3.0d0/8.0d0*u_fake(1,2)*dy-2.0d0/8.0d0*u_E*dy)*F_w_u(1,2)+ (-3.0d0/8.0d0*& u_fake(1,3)*dy+2.0d0/8.0d0*u_fake(1,2)*dy)*F_e_u(1,2)+ (u_G*dx)*F_s_u(1,2)+(-3.0d0/8.0d0*& u_fake(2,2)*dx+2.0d0/8.0d0*u_G*dx)*F_n_u(1,2)+(p_fake(1,1)- p_fake(1,2))*dy+D_G*8.0d0/3.0d0*& u_G*dx+((1.0d0-alpha_u)*a_P(1,2))*u_old(1,2) ! a_P(1,2)=1.0d0 ! a_W(1,2)=1.0d-20 ! a_E(1,2)=1.0d-20 ! a_S(1,2)=1.0d-20 ! a_N(1,2)=1.0d-20 ! S_fake(1,2)=u_G !coeficientes de v b_P(1,2)=(F_e_v(1,2)*dy+F_n_v(1,2)*7.0d0/8.0d0*dx+D_w*dy+D_e*dy+D_n*dx +D_K*9.0d0/3.0d0*dx)/& alpha_v b_W(1,2)=F_w_v(1,2)*7.0d0/8.0d0*dy+F_e_v(1,2)*1.0d0/8.0d0*dy+D_w*dy b_E(1,2)=D_e*dy b_S(1,2)=1.0d-20 !0.0d0 b_N(1,2)=D_n*dx+D_K*1.0d0/3.0d0*dx S_v_fake(1,2)=(3.0d0/8.0d0*v_fake(1,2)*dy-2.0d0/8.0d0*v_I*dy)*F_w_v(1,2)+ (-3.0d0/8.0d0*& v_fake(1,3)*dy+2.0d0/8.0d0*v_fake(1,2)*dy)*F_e_v(1,2)+ (v_K*dx)*F_s_v(1,2)+(-3.0d0/8.0d0*& v_fake(2,2)*dx+2.0d0/8.0d0*v_K*dx)*F_n_v(1,2)+D_K*8.0d0/3.0d0*v_K*dx+ ((1.0d0-alpha_v)*& b_P(1,2))*v_old(1,2) !Falta o termo de pressão. ! b_P(1,2)=1.0d0 ! b_W(1,2)=1.0d-20 ! b_E(1,2)=1.0d-20 ! b_S(1,2)=1.0d-20 ! b_N(1,2)=1.0d-20 ! S_v_fake(1,2)=v_K !coeficientes de p c_W(1,2)=dy/a_P(1,2)*dy c_E(1,2)=dy/a_E(1,2)*alpha_u*dy c_S(1,2)=dx/b_P(1,2)*dx c_N(1,2)=dx/b_N(1,2)*alpha_v*dx c_P(1,2)=c_W(1,2)+c_E(1,2)+c_S(1,2)+c_N(1,2) S_p_fake(1,2)=u_fake(1,2)*dy-u_fake(1,3)*dy+v_fake(1,2)*dx- v_fake(2,2)*dx !fronteira N - com exceção da coluna 2 do j=3,Nx-1 !coeficientes de u a_P(Ny,j)=(F_e_u(Ny,j)*dy+D_w*dy+D_e*dy+D_s*dx+D_H*9.0d0/3.0d0*dx)/ alpha_u a_W(Ny,j)=F_w_u(Ny,j)*dy+F_e_u(Ny,j)*1.0d0/8.0d0*dy+D_w*dy a_E(Ny,j)=D_e*dy a_S(Ny,j)=F_s_u(Ny,j)*dx+D_s*dx+D_H*1.0d0/3.0d0*dx a_N(Ny,j)=1.0d-20 !0.0d0 S_fake(Ny,j)=(3.0d0/8.0d0*u_fake(Ny,j)*dy-2.0d0/8.0d0*u_fake(Ny,j-1)*dy-1..0d0/8.0d0*& u_fake(Ny,j-2)*dy)*F_w_u(Ny,j)+(-3.0d0/8.0d0*u_fake(Ny,j+1)*dy +2.0d0/8.0d0*& u_fake(Ny,j)*dy)*F_e_u(Ny,j) +(3.0d0/8.0d0*u_fake(Ny,j)*dx-2.0d0/8.0d0*u_fake(Ny-1,j)& *dx-1.0d0/8.0d0*u_fake(Ny-2,j)*dx)*F_s_u(Ny,j)+(-u_H*dx)*F_n_u(Ny,j) +(p_fake(Ny,j-1)& -p_fake(Ny,j))*dy+D_H*8.0d0/3.0d0*u_H*dx+((1.0d0- alpha_u)*a_P(Ny,j))*u_old(Ny,j) ! a_P(Ny,j)=1.0d0 ! a_W(Ny,j)=1.0d-20 ! a_E(Ny,j)=1.0d-20 ! a_S(Ny,j)=1.0d-20 ! a_N(Ny,j)=1.0d-20 ! S_fake(Ny,j)=u_H !coeficientes de v b_P(Ny,j)=(F_e_v(Ny,j)*dy+D_w*dy+D_e*dy+D_s*dx+D_L*9.0d0/3.0d0*dx)/ alpha_v b_W(Ny,j)=F_w_v(Ny,j)*dy+F_e_v(Ny,j)*1.0d0/8.0d0*dy+D_w*dy b_E(Ny,j)=D_e*dy b_S(Ny,j)=F_s_v(Ny,j)*dx+D_s*dx+D_L*1.0d0/3.0d0*dx b_N(Ny,j)=1.0d-20 !0.0d0 S_v_fake(Ny,j)=(3.0d0/8.0d0*v_fake(Ny,j)*dy-2.0d0/8.0d0*v_fake(Ny,j-1)*dy-1.0d0/ & 8.0d0*v_fake(Ny,j-2)*dy)*F_w_v(Ny,j)+(-3.0d0/8.0d0*v_fake(Ny,j+1)*dy +2.0d0/8.0d0*& v_fake(Ny,j)*dy)*F_e_v(Ny,j) +(3.0d0/8.0d0*v_fake(Ny,j)*dx-2.0d0/8.0d0*v_fake(Ny-1,j)& *dx-1.0d0/8.0d0*v_fake(Ny-2,j)*dx)*F_s_v(Ny,j)+(-v_L*dx)*F_n_v(Ny,j) +(p_fake(Ny-1,j)& -p_fake(Ny,j))*dx+D_L*8.0d0/3.0d0*v_L*dx+((1.0d0- alpha_v)*b_P(Ny,j))*v_old(Ny,j) ! b_P(Ny,j)=1.0d0 ! b_W(Ny,j)=1.0d-20 ! b_E(Ny,j)=1.0d-20 ! b_S(Ny,j)=1.0d-20 ! b_N(Ny,j)=1.0d-20 ! S_v_fake(Ny,j)=v_L !coeficientes de p c_W(Ny,j)=dy/a_P(Ny,j)*dy c_E(Ny,j)=dy/a_E(Ny,j)*alpha_u*dy c_S(Ny,j)=dx/b_P(Ny,j)*dx c_N(Ny,j)=dx/b_N(Ny,j)*alpha_v*dx c_P(Ny,j)=c_W(Ny,j)+c_E(Ny,j)+c_S(Ny,j)+c_N(Ny,j) S_p_fake(Ny,j)=u_fake(Ny,j)*dy-u_fake(Ny,j+1)*dy !Faltam termos de v. end do !fronteira N - coluna 2 !coeficientes de u a_P(Ny,2)=(F_e_u(Ny,2)*dy+D_w*dy+D_e*dy+D_s*dx+D_H*9.0d0/3.0d0*dx)/ alpha_u a_W(Ny,2)=F_w_u(Ny,2)*7.0d0/8.0d0*dy+F_e_u(Ny,2)*1.0d0/8.0d0*dy +D_w*dy a_E(Ny,2)=D_e*dy a_S(Ny,2)=F_s_u(Ny,2)*dx+D_s*dx+D_H*1.0d0/3.0d0*dx a_N(Ny,2)=1.0d-20 !0.0d0 S_fake(Ny,2)=(3.0d0/8.0d0*u_fake(Ny, 2)*dy-2.0d0/8.0d0*u_E*dy)*F_w_u(Ny,2)+(-3.0d0/& 8.0d0*u_fake(Ny,3)*dy+2.0d0/8.0d0*u_fake(Ny,2)*dy)*F_e_u(Ny, 2)+(3.0d0/8.0d0*& u_fake(Ny, 2)*dx-2.0d0/8.0d0*u_fake(Ny-1,2)*dx-1.0d0/8.0d0*u_fake(Ny-2,2)*dx)*& F_s_u(Ny,2)+(-u_H*dx)*F_n_u(Ny,2)+(p_fake(Ny,1)-p_fake(Ny,2))*dy +D_H*8.0d0/3.0d0*u_H& *dx+((1.0d0-alpha_u)*a_P(Ny,2))*u_old(Ny,2) ! a_P(Ny,2)=1.0d0 ! a_W(Ny,2)=1.0d-20 ! a_E(Ny,2)=1.0d-20 ! a_S(Ny,2)=1.0d-20 ! a_N(Ny,2)=1.0d-20 ! S_fake(Ny,2)=u_H !coeficientes de v b_P(Ny,2)=(F_e_v(Ny,2)*dy+D_w*dy+D_e*dy+D_s*dx+D_L*9.0d0/3.0d0*dx)/ alpha_v b_W(Ny,2)=F_w_v(Ny,2)*7.0d0/8.0d0*dy+F_e_v(Ny,2)*1.0d0/8.0d0*dy +D_w*dy b_E(Ny,2)=D_e*dy b_S(Ny,2)=F_s_v(Ny,2)*dx+D_s*dx+D_L*1.0d0/3.0d0*dx b_N(Ny,2)=1.0d-20 !0.0d0 S_v_fake(Ny,2)=(3.0d0/8.0d0*v_fake(Ny, 2)*dy-2.0d0/8.0d0*v_I*dy)*F_w_v(Ny,2)+(-3.0d0/& 8.0d0*v_fake(Ny,3)*dy+2.0d0/8.0d0*v_fake(Ny,2)*dy)*F_e_v(Ny, 2)+(3.0d0/8.0d0*& v_fake(Ny, 2)*dx-2.0d0/8.0d0*v_fake(Ny-1,2)*dx-1.0d0/8.0d0*v_fake(Ny-2,2)*dx)*& F_s_v(Ny,2)+(-v_L*dx)*F_n_v(Ny,2)+(p_fake(Ny-1,2)-p_fake(Ny,2))*dx +D_L*8.0d0/3.0d0*& v_L*dx+((1.0d0-alpha_v)*b_P(Ny,2))*v_old(Ny,2) ! b_P(Ny,2)=1.0d0 ! b_W(Ny,2)=1.0d-20 ! b_E(Ny,2)=1.0d-20 ! b_S(Ny,2)=1.0d-20 ! b_N(Ny,2)=1.0d-20 ! S_v_fake(Ny,2)=v_L !coeficientes de p c_W(Ny,2)=dy/a_P(Ny,2)*dy c_E(Ny,2)=dy/a_E(Ny,2)*alpha_u*dy c_S(Ny,2)=dx/b_P(Ny,2)*dx c_N(Ny,2)=dx/b_N(Ny,2)*alpha_v*dx c_P(Ny,2)=c_W(Ny,2)+c_E(Ny,2)+c_S(Ny,2)+c_N(Ny,2) S_p_fake(Ny,2)=u_fake(Ny,2)*dy-u_fake(Ny,3)*dy !Faltam termos de v. !fronteira SW !coeficientes de u a_P(1,1)=(F_e_u(1,1)*7.0d0/8.0d0*dy+F_n_u(1,1)*7.0d0/8.0d0*dx+D_e*dy +D_n*dx+D_EE*& 9.0d0/3.0d0*dy+D_G*9.0d0/3.0d0*dx)/alpha_u a_W(1,1)=1.0d-20 !0.0d0 a_E(1,1)=D_e*dy+D_EE*1.0d0/3.0d0*dy a_S(1,1)=1.0d-20 !0.0d0 a_N(1,1)=D_n*dx+D_G*1.0d0/3.0d0*dx S_fake(1,1)=(u_E*dy)*F_w_u(1,1)+(-3.0d0/8.0d0*u_fake(1,2)*dy +2.0d0/8.0d0*u_E*dy)*& F_e_u(1,1)+(u_G*dx)*F_s_u(1,1)+(-3.0d0/8.0d0*u_fake(2,1)*dx +2.0d0/8.0d0*u_G*dx)*& F_n_u(1,1)+D_EE*8.0d0/3.0d0*u_E*dy+D_G*8.0d0/3.0d0*u_G*dx+((1.0d0- alpha_u)*a_P(1,1))& *u_old(1,1) !Falta o termo de pressão. ! a_P(1,1)=1.0d0 ! a_W(1,1)=1.0d-20 ! a_E(1,1)=1.0d-20 ! a_S(1,1)=1.0d-20 ! a_N(1,1)=1.0d-20 ! S_fake(1,1)=1.0d-30 !coeficientes de v b_P(1,1)=(F_e_v(1,1)*7.0d0/8.0d0*dy+F_n_v(1,1)*7.0d0/8.0d0*dx+D_e*dy +D_n*dx+D_I*& 9.0d0/3.0d0*dy+D_K*9.0d0/3.0d0*dx)/alpha_v b_W(1,1)=1.0d-20 !0.0d0 b_E(1,1)=D_e*dy+D_I*1.0d0/3.0d0*dy b_S(1,1)=1.0d-20 !0.0d0 b_N(1,1)=D_n*dx+D_K*1.0d0/3.0d0*dx S_v_fake(1,1)=(v_I*dy)*F_w_v(1,1)+(-3.0d0/8.0d0*v_fake(1,2)*dy +2.0d0/8.0d0*v_I*dy)*& F_e_v(1,1)+(v_K*dx)*F_s_v(1,1)+(-3.0d0/8.0d0*v_fake(2,1)*dx +2.0d0/8.0d0*v_K*dx)*& F_n_v(1,1)+D_I*8.0d0/3.0d0*v_I*dy+D_K*8.0d0/3.0d0*v_K*dx+((1.0d0- alpha_v)*b_P(1,1))& *v_old(1,1) !Falta o termo de pressão. ! b_P(1,1)=1.0d0 ! b_W(1,1)=1.0d-20 ! b_E(1,1)=1.0d-20 ! b_S(1,1)=1.0d-20 ! b_N(1,1)=1.0d-20 ! S_v_fake(1,1)=1.0d-30 !coeficientes de p c_W(1,1)=dy/a_P(1,1)*dy c_E(1,1)=dy/a_E(1,1)*alpha_u*dy c_S(1,1)=dx/b_P(1,1)*dx c_N(1,1)=dx/b_N(1,1)*alpha_v*dx c_P(1,1)=c_W(1,1)+c_E(1,1)+c_S(1,1)+c_N(1,1) S_p_fake(1,1)=u_fake(1,1)*dy-u_fake(1,2)*dy+v_fake(1,1)*dx- v_fake(2,1)*dx !fronteira NW !coeficientes de u a_P(Ny,1)=(F_e_u(Ny,1)*7.0d0/8.0d0*dy+D_e*dy+D_s*dx +D_EE*9.0d0/3.0d0*dy+D_H*9.0d0/& 3.0d0*dx)/alpha_u a_W(Ny,1)=1.0d-20 !0.0d0 a_E(Ny,1)=D_e*dy+D_EE*1.0d0/3.0d0*dy a_S(Ny,1)=F_s_u(Ny,1)*dx+D_s*dx+D_H*1.0d0/3.0d0*dx a_N(Ny,1)=1.0d-20 !0.0d0 S_fake(Ny,1)=(u_E*dy)*F_w_u(Ny,1)+(-3.0d0/8.0d0*u_fake(Ny,2)*dy +2.0d0/8.0d0*u_E*dy)*& F_e_u(Ny,1)+(3.0d0/8.0d0*u_fake(Ny, 1)*dx-2.0d0/8.0d0*u_fake(Ny-1,1)*dx-1.0d0/8.0d0*& u_fake(Ny-2,1)*dx)*F_s_u(Ny,1)+(-u_H*dx)*F_n_u(Ny, 1)+D_EE*8.0d0/3.0d0*u_E*dy+D_H*& 8.0d0/3.0d0*u_H*dx+((1.0d0-alpha_u)*a_P(Ny,1))*u_old(Ny,1) ! Falta o termo de pressão. ! a_P(Ny,1)=1.0d0 ! a_W(Ny,1)=1.0d-20 ! a_E(Ny,1)=1.0d-20 ! a_S(Ny,1)=1.0d-20 ! a_N(Ny,1)=1.0d-20 ! S_fake(Ny,1)=1.0d-30 !coeficientes de v b_P(Ny,1)=(F_e_v(Ny,1)*7.0d0/8.0d0*dy+D_e*dy+D_s*dx +D_I*9.0d0/3.0d0*dy+D_L*9.0d0/& 3.0d0*dx)/alpha_v b_W(Ny,1)=1.0d-20 !0.0d0 b_E(Ny,1)=D_e*dy+D_I*1.0d0/3.0d0*dy b_S(Ny,1)=F_s_v(Ny,1)*dx+D_s*dx+D_L*1.0d0/3.0d0*dx b_N(Ny,1)=1.0d-20 !0.0d0 S_v_fake(Ny,1)=(v_I*dy)*F_w_v(Ny,1)+(-3.0d0/8.0d0*v_fake(Ny,2)*dy +2.0d0/8.0d0*v_I*dy& )*F_e_v(Ny,1)+(3.0d0/8.0d0*v_fake(Ny, 1)*dx-2.0d0/8.0d0*v_fake(Ny-1,1)*dx-1.0d0/8.0d0& *v_fake(Ny-2,1)*dx)*F_s_v(Ny,1)+(-v_L*dx)*F_n_v(Ny,1)+ (p_fake(Ny-1,1)-p_fake(Ny,1))*& dx+D_I*8.0d0/3.0d0*v_I*dy+D_L*8.0d0/3.0d0*v_L*dx+((1.0d0- alpha_v)*b_P(Ny,1))*& v_old(Ny,1) ! b_P(Ny,1)=1.0d0 ! b_W(Ny,1)=1.0d-20 ! b_E(Ny,1)=1.0d-20 ! b_S(Ny,1)=1.0d-20 ! b_N(Ny,1)=1.0d-20 ! S_v_fake(Ny,1)=1.0d-30 !coeficientes de p c_W(Ny,1)=dy/a_P(Ny,1)*dy c_E(Ny,1)=dy/a_E(Ny,1)*alpha_u*dy c_S(Ny,1)=dx/b_P(Ny,1)*dx c_N(Ny,1)=dx/b_N(Ny,1)*alpha_v*dx c_P(Ny,1)=c_W(Ny,1)+c_E(Ny,1)+c_S(Ny,1)+c_N(Ny,1) S_p_fake(Ny,1)=u_fake(Ny,1)*dy-u_fake(Ny,2)*dy !Faltam termos de v. !fronteira SE !coeficientes de u a_P(1,Nx)=(F_n_u(1,Nx)*7.0d0/8.0d0*dx+D_w*dy+D_n*dx +D_F*9.0d0/3.0d0*dy+D_G*9.0d0/& 3.0d0*dx)/alpha_u a_W(1,Nx)=F_w_u(1,Nx)*dy+D_w*dy+D_F*1.0d0/3.0d0*dy a_E(1,Nx)=1.0d-20 !0.0d0 a_S(1,Nx)=1.0d-20 !0.0d0 a_N(1,Nx)=D_n*dx+D_G*1.0d0/3.0d0*dx S_fake(1,Nx)=(3.0d0/8.0d0*u_fake(1,Nx)*dy-2.0d0/8.0d0*u_fake(1,Nx-1)*dy-1..0d0/8.0d0*& u_fake(1,Nx-2)*dy)*F_w_u(1,Nx)+(-u_F(1)*dy)*F_e_u(1,Nx)+ (u_G*dx)*F_s_u(1,Nx)+(-3.0d0& /8.0d0*u_fake(2,Nx)*dx+2.0d0/8.0d0*u_G*dx)*F_n_u(1,Nx)+ (p_fake(1,Nx-1)-p_fake(1,Nx))& *dy+D_F*8.0d0/3.0d0*u_F(1)*dy+D_G*8.0d0/3.0d0*u_G*dx+((1.0d0- alpha_u)*a_P(1,Nx))*& u_old(1,Nx) ! a_P(1,Nx)=1.0d0 ! a_W(1,Nx)=1.0d-20 ! a_E(1,Nx)=1.0d-20 ! a_S(1,Nx)=1.0d-20 ! a_N(1,Nx)=1.0d-20 ! S_fake(1,Nx)=1.0d-30 !coeficientes de v b_P(1,Nx)=(F_n_v(1,Nx)*7.0d0/8.0d0*dx+D_w*dy+D_n*dx +D_J*9.0d0/3.0d0*dy+D_K*9.0d0/& 3.0d0*dx)/alpha_v b_W(1,Nx)=F_w_v(1,Nx)*dy+D_w*dy+D_J*1.0d0/3.0d0*dy b_E(1,Nx)=1.0d-20 !0.0d0 b_S(1,Nx)=1.0d-20 !0.0d0 b_N(1,Nx)=D_n*dx+D_K*1.0d0/3.0d0*dx S_v_fake(1,Nx)=(3.0d0/8.0d0*v_fake(1,Nx)*dy-2.0d0/8.0d0*v_fake(1,Nx-1)*dy-1.0d0/ & 8.0d0*v_fake(1,Nx-2)*dy)*F_w_v(1,Nx)+(-v_J*dy)*F_e_v(1,Nx)+ (v_K*dx)*F_s_v(1,Nx)+(& -3.0d0/8.0d0*v_fake(2,Nx)*dx+2.0d0/8.0d0*v_K*dx)*F_n_v(1,Nx) +D_J*8.0d0/3.0d0*v_J*dy+& D_K*8.0d0/3.0d0*v_K*dx+((1.0d0-alpha_v)*b_P(1,Nx))*v_old(1,Nx) ! Falta o termo de pressão. ! b_P(1,Nx)=1.0d0 ! b_W(1,Nx)=1.0d-20 ! b_E(1,Nx)=1.0d-20 ! b_S(1,Nx)=1.0d-20 ! b_N(1,Nx)=1.0d-20 ! S_v_fake(1,Nx)=1.0d-30 !coeficientes de p c_W(1,Nx)=dy/a_P(1,Nx)*dy c_E(1,Nx)=dy/a_E(1,Nx)*alpha_u*dy c_S(1,Nx)=dx/b_P(1,Nx)*dx c_N(1,Nx)=dx/b_N(1,Nx)*alpha_v*dx c_P(1,Nx)=c_W(1,Nx)+c_E(1,Nx)+c_S(1,Nx)+c_N(1,Nx) S_p_fake(1,Nx)=v_fake(1,Nx)*dx-v_fake(2,Nx)*dx !Faltam termos de u. !fronteira NE !coeficientes de u a_P(Ny,Nx)=(D_w*dy+D_s*dx+D_F*9.0d0/3.0d0*dy+D_H*9.0d0/3.0d0*dx)/ alpha_u a_W(Ny,Nx)=F_w_u(Ny,Nx)*dy+D_w*dy+D_F*1.0d0/3.0d0*dy a_E(Ny,Nx)=1.0d-20 !0.0d0 a_S(Ny,Nx)=F_s_u(Ny,Nx)*dx+D_s*dx+D_H*1.0d0/3.0d0*dx a_N(Ny,Nx)=1.0d-20 !0.0d0 S_fake(Ny,Nx)=(3.0d0/8.0d0*u_fake(Ny,Nx)*dy-2.0d0/8.0d0*u_fake(Ny,Nx-1)*dy-1.0d0/ & 8.0d0*u_fake(Ny,Nx-2)*dy)*F_w_u(Ny,Nx)+(-u_F(Ny)*dy)*F_e_u(Ny,Nx) +(3.0d0/8.0d0*& u_fake(Ny,Nx)*dx-2.0d0/8.0d0*u_fake(Ny-1,Nx)*dx-1.0d0/8.0d0*u_fake(Ny-2,Nx)*dx)*& F_s_u(Ny,Nx)+(-u_H*dx)*F_n_u(Ny,Nx)+(p_fake(Ny,Nx-1)- p_fake(Ny,Nx))*dy+D_F*8.0d0/& 3.0d0*u_F(Ny)*dy+D_H*8.0d0/3.0d0*u_H*dx+((1.0d0- alpha_u)*a_P(Ny,Nx))*u_old(Ny,Nx) ! a_P(Ny,Nx)=1.0d0 ! a_W(Ny,Nx)=1.0d-20 ! a_E(Ny,Nx)=1.0d-20 ! a_S(Ny,Nx)=1.0d-20 ! a_N(Ny,Nx)=1.0d-20 ! S_fake(Ny,Nx)=1.0d-30 !coeficientes de v b_P(Ny,Nx)=(D_w*dy+D_s*dx+D_J*9.0d0/3.0d0*dy+D_L*9.0d0/3.0d0*dx)/ alpha_v b_W(Ny,Nx)=F_w_v(Ny,Nx)*dy+D_w*dy+D_J*1.0d0/3.0d0*dy b_E(Ny,Nx)=1.0d-20 !0.0d0 b_S(Ny,Nx)=F_s_v(Ny,Nx)*dx+D_s*dx+D_L*1.0d0/3.0d0*dx b_N(Ny,Nx)=1.0d-20 !0.0d0 S_v_fake(Ny,Nx)=(3.0d0/8.0d0*v_fake(Ny,Nx)*dy-2.0d0/8.0d0*v_fake(Ny,Nx-1)*dy-1.0d0/ & 8.0d0*v_fake(Ny,Nx-2)*dy)*F_w_v(Ny,Nx)+(-v_J*dy)*F_e_v(Ny,Nx) +(3.0d0/8.0d0*& v_fake(Ny,Nx)*dx-2.0d0/8.0d0*v_fake(Ny-1,Nx)*dx-1.0d0/8.0d0*v_fake(Ny-2,Nx)*dx)*& F_s_v(Ny,Nx)+(-v_L*dx)*F_n_v(Ny,Nx)+(p_fake(Ny-1,Nx)- p_fake(Ny,Nx))*dx+D_J*8.0d0/& 3.0d0*v_J*dy+D_L*8.0d0/3.0d0*v_L*dx+((1.0d0- alpha_v)*b_P(Ny,Nx))*v_old(Ny,Nx) ! b_P(Ny,Nx)=1.0d0 ! b_W(Ny,Nx)=1.0d-20 ! b_E(Ny,Nx)=1.0d-20 ! b_S(Ny,Nx)=1.0d-20 ! b_N(Ny,Nx)=1.0d-20 ! S_v_fake(Ny,Nx)=1.0d-30 !coeficientes de p c_W(Ny,Nx)=dy/a_P(Ny,Nx)*dy c_E(Ny,Nx)=dy/a_E(Ny,Nx)*alpha_u*dy c_S(Ny,Nx)=dx/b_P(Ny,Nx)*dx c_N(Ny,Nx)=dx/b_N(Ny,Nx)*alpha_v*dx c_P(Ny,Nx)=c_W(Ny,Nx)+c_E(Ny,Nx)+c_S(Ny,Nx)+c_N(Ny,Nx) S_p_fake(Ny,Nx)=0.0d0 !Faltam termos de u e de v. !Fonte k=1 do i=1,Ny do j=1,Nx S(k)=S_fake(i,j) S_v(k)=S_v_fake(i,j) S_p(k)=S_p_fake(i,j) k=k+1 end do end do !Matriz dos Coeficientes k=1 do i=1,Ny do j=1,Nx a(k,k)=a_P(i,j) b(k,k)=b_P(i,j) c(k,k)=c_P(i,j) k=k+1 end do end do k=2 do i=1,Ny if (i.eq.1) then do j=2,Nx a(k,k-1)=-a_W(i,j) b(k,k-1)=-b_W(i,j) c(k,k-1)=-c_W(i,j) k=k+1 end do else do j=1,Nx a(k,k-1)=-a_W(i,j) b(k,k-1)=-b_W(i,j) c(k,k-1)=-c_W(i,j) k=k+1 end do end if end do k=2 do i=1,Ny if (i.ne.Ny) then do j=1,Nx a(k-1,k)=-a_E(i,j) b(k-1,k)=-b_E(i,j) c(k-1,k)=-c_E(i,j) k=k+1 end do else do j=1,Nx-1 a(k-1,k)=-a_E(i,j) b(k-1,k)=-b_E(i,j) c(k-1,k)=-c_E(i,j) k=k+1 end do end if end do k=Nx+1 do i=2,Ny do j=1,Nx a(k,k-Nx)=-a_S(i,j) b(k,k-Nx)=-b_S(i,j) c(k,k-Nx)=-c_S(i,j) k=k+1 end do end do k=1 do i=1,Ny-1 do j=1,Nx a(k,k+Nx)=-a_N(i,j) b(k,k+Nx)=-b_N(i,j) c(k,k+Nx)=-c_N(i,j) k=k+1 end do end do k=0 l=0 m=0 do i=1,Ny*Nx do j=1,Ny*Nx if (dabs(a(i,j)).ne.0.0d0) then k=k+1 end if if (dabs(b(i,j)).ne.0.0d0) then l=l+1 end if if (dabs(c(i,j)).ne.0.0d0) then m=m+1 end if end do end do allocate(AA(k),irow_u(k),jcol_u(k)) allocate(BB(l),irow_v(l),jcol_v(l)) allocate(CC(m),irow_p(m),jcol_p(m)) k=1 l=1 m=1 do i=1,Ny*Nx do j=1,Ny*Nx if (dabs(a(i,j)).ne.0.0d0) then AA(k)=a(i,j) irow_u(k)=i jcol_u(k)=j k=k+1 end if if (dabs(b(i,j)).ne.0.0d0) then BB(l)=b(i,j) irow_v(l)=i jcol_v(l)=j l=l+1 end if if (dabs(c(i,j)).ne.0.0d0) then CC(m)=c(i,j) irow_p(m)=i jcol_p(m)=j m=m+1 end if end do end do call DL4LXG(iparam,rparam) iparam(5)=1000000 !call DLINRG(Ny*Nx,a,Ny*Nx,a_inv,Ny*Nx) !CALL DLSARG (Ny*Nx,a,Ny*Nx,S,ipath,u) call DLSLXG(Ny*Nx,k-1,AA,irow_u,jcol_u,S,ipath,iparam,rparam,u) call DLSLXG(Ny*Nx,l-1,BB,irow_v,jcol_v,S_v,ipath,iparam,rparam,v) call DLSLXG(Ny*Nx,m-1,CC,irow_p,jcol_p,S_p,ipath,iparam,rparam,p) ! k=1 ! do i=1,Ny ! do j=1,Nx ! p_fake(i,j)=p(k) ! if (j.ne.1) then ! u_fake(i,j)=u(k)+dy/(a_P(i,j)*alpha_u)*(p_fake(i,j-1)- p_fake(i,j)) ! end if ! if (i.ne.1) then ! v_fake(i,j)=v(k)+dx/(b_P(i,j)*alpha_v)*(p_fake(i-1,j)- p_fake(i,j)) ! end if ! p_fake(i,j)=p_old(i,j)+alpha_p*p(k) ! k=k+1 ! end do ! end do p_old=p_fake u_old=u_fake v_old=v_fake k=1 do i=1,Ny do j=1,Nx p_fake(i,j)=p(k) !p_fake(i,j)=p_old(i,j)+alpha_p*p(k) !p_fake(i,j)=alpha_p*(k)+(1.0d0-alpha_p)*p_old(i,j) !if (i.ne.1) then ! if (i.ne.Ny) then ! if (j.ne.1) then ! if (j.ne.Nx) then !u_fake(i,j)=u(k)+dy/(a_P(i,j)*alpha_u)*(p_fake(i,j-1)- p_fake(i,j)) !v_fake(i,j)=v(k)+dx/(b_P(i,j)*alpha_v)*(p_fake(i-1,j)- p_fake(i,j)) ! u_fake(i,j)=u(k)+dy/a_P(i,j)*(p_fake(i,j-1)-p_fake(i,j)) ! v_fake(i,j)=v(k)+dx/b_P(i,j)*(p_fake(i-1,j)-p_fake(i,j)) ! u_fake(i,j)=alpha_u*u(k)+(1.0d0-alpha_u)*u_old(i,j) ! v_fake(i,j)=alpha_v*v(k)+(1.0d0-alpha_v)*v_old(i,j) ! end if ! end if ! end if !end if if (j.ne.1) then u_fake(i,j)=u(k)+dy/(a_P(i,j)*alpha_u)*(p_fake(i,j-1)-p_fake(i,j)) !u_fake(i,j)=u(k)+dy/a_P(i,j)*(p_fake(i,j-1)-p_fake(i,j)) !u_fake(i,j)=alpha_u*u(k) !u_fake(i,j)=alpha_u*u(k)+(1.0d0-alpha_u)*u_old(i,j) ! write(*,*) u_fake(i,j) end if if (i.ne.1) then v_fake(i,j)=v(k)+dx/(b_P(i,j)*alpha_v)*(p_fake(i-1,j)-p_fake(i,j)) !v_fake(i,j)=v(k)+dx/b_P(i,j)*(p_fake(i-1,j)-p_fake(i,j)) !v_fake(i,j)=alpha_v*v(k) !v_fake(i,j)=alpha_v*v(k)+(1.0d0-alpha_v)*v_old(i,j) end if p_fake(i,j)=p_old(i,j)+alpha_p*p(k) k=k+1 end do end do deallocate(AA,irow_u,jcol_u,BB,irow_v,jcol_v,CC,irow_p,jcol_p) !u=matmul(a_inv,S) !theta=matmul(a_inv,S) ! end do open(1,file='u.dat') write(1,11) Ny,Nx 11 format('variables=X,Y,u',/,'zone t="time','" i=',i4,' j=',i4,' f=point') k=1 do i=1,Ny do j=1,Nx write(1,*) Y(i),X(j),u(k) k=k+1 end do end do open(2,file='v.dat') write(2,12) Ny,Nx 12 format('variables=X,Y,v',/,'zone t="time','" i=',i4,' j=',i4,' f=point') k=1 do i=1,Ny do j=1,Nx write(2,*) Y(i),X(j),v(k) k=k+1 end do end do open(3,file='p.dat') write(3,13) Ny,Nx 13 format('variables=X,Y,p',/,'zone t="time','" i=',i4,' j=',i4,' f=point') k=1 do i=1,Ny do j=1,Nx write(3,*) Y(i),X(j),p(k) k=k+1 end do end do open(4,file='S.dat') write(4,14) Ny,Nx 14 format('variables=X,Y,p',/,'zone t="time','" i=',i4,' j=',i4,' f=point') k=1 do i=1,Ny do j=1,Nx write(4,*) S_fake(i,j),k,S(k) k=k+1 end do end do open(5,file='a.dat') !write(5,15) Ny,Nx !15 format('variables=X,Y,p',/,'zone t="time','" i=',i4,' j=',i4,' f=point') do i=1,Ny*Nx do j=1,Ny*Nx write(5,*) i,j,a(i,j) end do end do open(6,file='aa.dat') !write(6,16) Ny,Nx !16 format('variables=X,Y,p',/,'zone t="time','" i=',i4,' j=',i4,' f=point') do i=1,Ny do j=1,Nx write(6,*) i,j,a_P(i,j) end do end do end program EXERCICIO
From: Gordon Sande on 5 Jul 2010 15:59 On 2010-07-05 16:37:27 -0300, Felipe Abr�o <felipeabrao(a)gmail.com> said: > Please, could somebody explain me what is happening. I simply > allocated a matrix called "a", but figured out that the program is > inserting some strange values to some of its elements (like > 1.057229231618489E-307 or 6.941570408926931E+228), even before I set > any value to it! Shouldn't all the values be zero before any desired > value setting? No! Before you set the value it has NO value! The "history" is that back when computers were $1000 per hour and programmers were paid $3.50 per hour it was more efficient to just use whatever value was left over from the previous program that was executing. It is called an undefined value. Fortran starts things with an undefined value. This is not as silly as it sounds as you might have a subroutine that gets used twice. What values should be used the second time? If you had been using a good debugging compiler, and had set the suitable options to have it check, then you would have gotten a diagnostic about an undefined value. Since you are using Windows you might try the Silverfrost Fortran 90 (it used to be Salford who developed it for student use) which is free for personal use. Be sure to set the /undef option! > Below is my code, but you don't have to read it all. > You can check the problem just after line 65, where I wrote a command > in order to read the values of matrix "a". The same thing is happening > with matrixs "b" and "c". In advance I thank you all! > > PS: I use Fortran 90 in Windows XP. > > > ! EXERCICIO.f90 > ! > ! FUNCTIONS: > ! EXERCICIO - Entry point of console application. > ! > > ! > *************************************************************************** > * > ! > ! PROGRAM: EXERCICIO > ! > ! PURPOSE: Entry point for the console application. > ! > ! > *************************************************************************** > * > > program EXERCICIO > > use IMSL > > implicit none > > ! Variables > > integer(4)::i,j,k,l,m,Nx,Ny > integer(4),parameter::ipath=1 > integer(4),dimension(6)::iparam > > integer(4),allocatable,dimension(:)::irow_u,jcol_u,irow_v,jcol_v,irow_p,jco > l_p > real(8)::dx,dy,D_w,D_e,D_s,D_n,D_EE,D_F,D_G,D_H,D_I,D_J,D_K,D_L > real(8),parameter::Re=1.0d2,gama=1.0d0/ > Re,u_E=1.0d0,u_G=0.0d0,u_H=0.0d0,v_I=0.0d0,v_J=0.0d0,& > > v_K=0.0d0,v_L=0.0d0,error=1.0d-6,alpha_u=5.0d-1,alpha_v=5.0d-1,al > pha_p=5.0d-1 > !real(8),dimension(N*N)::theta_old > real(8),dimension(5)::rparam > > real(8),allocatable,dimension(:)::u,v,p,theta,S,S_v,S_p,X,Y,u_F,AA,BB,CC > !real(8),dimension(N,N)::theta_fake > > real(8),allocatable,dimension(:,:)::a,a_inv,b,c,p_fake,p_old,u_fake,u_old,v > _fake,v_old,a_P,& > > a_W,a_E,a_S,a_N,b_P,b_W,b_E,b_S,b_N,c_P,c_W,c_E,c_S,c_N,F_w_u,F_e_u,F_s_u,F > _n_u,S_fake,& > S_v_fake,S_p_fake,F_w_v,F_e_v,F_s_v,F_n_v > > ! Body of EXERCICIO_1 > > iparam(1)=0 > > Nx=5 > Ny=5 > > > > > allocate(a(Ny*Nx,Ny*Nx),a_inv(Ny*Nx,Ny*Nx),b(Ny*Nx,Ny*Nx),c(Ny*Nx,Ny*Nx)) > > allocate(a_P(Ny,Nx),a_W(Ny,Nx),a_E(Ny,Nx),a_S(Ny,Nx),a_N(Ny,Nx),b_P(Ny,Nx), > b_W(Ny,Nx),& > > b_E(Ny,Nx),b_S(Ny,Nx),b_N(Ny,Nx),c_P(Ny,Nx),c_W(Ny,Nx),c_E(Ny,Nx),c_S(Ny,Nx > ),c_N(Ny,Nx),& > > p_fake(Ny,Nx),p_old(Ny,Nx),u_fake(Ny,Nx),u_old(Ny,Nx),v_fake(Ny,Nx),v_old(N > y,Nx),F_w_u(Ny,Nx),& > > F_e_u(Ny,Nx),F_s_u(Ny,Nx),F_n_u(Ny,Nx),S_fake(Ny,Nx),S_v_fake(Ny,Nx),S_p_fa > ke(Ny,Nx),& > F_w_v(Ny,Nx),F_e_v(Ny,Nx),F_s_v(Ny,Nx),F_n_v(Ny,Nx)) > > allocate(u(Ny*Nx),v(Ny*Nx),p(Ny*Nx),theta(Ny*Nx),S(Ny*Nx),S_v(Ny*Nx),S_p(Ny > *Nx)) > allocate(Y(Ny),u_F(Ny)) > allocate(X(Nx)) > > > open(7,file='pre.dat') > > do i=1,Ny*Nx > do j=1,Ny*Nx > write(7,*) i,j,a(i,j) > end do > end do > > > > dx=1.0d0/dble(Nx) > dy=1.0d0/dble(Ny) > D_w=gama/dx > D_e=gama/dx > D_s=gama/dy > D_n=gama/dy > D_EE=2.0d0*D_w > D_F=2.0d0*D_e > D_G=2.0d0*D_s > D_H=2.0d0*D_n > D_I=D_EE > D_J=D_F > D_K=D_G > D_L=D_H > > p_fake=5.0d-1 > u_fake=5.0d-1 > v_fake=5.0d-1 > > ! do i=1,Ny > ! do j=1,Nx > ! u_fake(i,j)=dble(i+j) > ! v_fake(i,j)=dble(i+j) > ! p_fake(i,j)=dble(i+j) > ! end do > ! end do > > > > theta=5.0d-1 > !theta_old=6.0d-1 > > do i=1,Nx > X(i)=1.0d0/dble((Nx-1))*dble((i-1)) > end do > > do i=1,Ny > Y(i)=1.0d0/dble((Ny-1))*dble((i-1)) > u_F(i)=4.0d0*(Y(i)-Y(i)**2.0d0) > end do > > p_old=1.0d0 > u_old=1.0d0 > v_old=1.0d0 > > ! do while(maxval(dabs(p_fake-p_old)).gt.error) > > !do while((maxval(dabs(p_fake-p_old)).gt.error).or.(maxval(dabs(u_fake- > u_old)).gt.error).or.& > !(maxval(dabs(v_fake-v_old)).gt.error)) > > !do while((maxval(dabs(p_fake(2:Ny-1,2:Nx-1)- > p_old(2:Ny-1,2:Nx-1))).gt.error).or.(maxval(dabs(& > !u_fake(2:Ny-1,2:Nx-1)-u_old(2:Ny-1,2:Nx-1))).gt.error).or. > (maxval(dabs(v_fake(2:Ny-1,2:Nx-1)-& > !v_old(2:Ny-1,2:Nx-1))).gt.error)) > > > > ! write(*,*) 'diferen�a de press�o =',maxval(dabs(p_fake-p_old)) > ! write(*,*) 'diferen�a de velocidade u =',maxval(dabs(u_fake-u_old)) > ! write(*,*) 'diferen�a de velocidade v =',maxval(dabs(v_fake-v_old)) > > > ! p_old=p_fake > ! u_old=u_fake > ! v_old=v_fake > > > do i=1,Ny > do j=1,Nx > if (j.eq.1) then > F_w_u(2:Ny-1,j)=1.0d0 > else > F_w_u(i,j)=1.0d0/2.0d0*(u_fake(i,j)+u_fake(i,j-1)) > end if > if (i.eq.1) then > F_w_v(i,j)=0.0d0 !Tenho d�vida. > else > F_w_v(i,j)=1.0d0/2.0d0*(u_fake(i,j)+u_fake(i-1,j)) > end if > if (j.eq.Nx) then > F_e_u(i,j)=4.0d0*(Y(i)-Y(i)**2.0d0) > else > F_e_u(i,j)=1.0d0/2.0d0*(u_fake(i,j+1)+u_fake(i,j)) > end if > if (i.eq.1) then > F_e_v(i,j)=0.0d0 !Tenho d�vida. > else if (j.eq.Nx) then > F_e_v(i,j)=4.0d0*(Y(i)-Y(i)**2.0d0) !Tenho d�vida. > else > F_e_v(i,j)=1.0d0/2.0d0*(u_fake(i,j+1)+u_fake(i-1,j+1)) > end if > if (i.eq.1) then > F_s_u(i,j)=0.0d0 > else if (j.eq.1) then > F_s_u(i,j)=0.0d0 > else > F_s_u(i,j)=1.0d0/2.0d0*(v_fake(i,j)+v_fake(i,j-1)) > end if > if (i.eq.1) then > F_s_v(i,j)=0.0d0 > else > F_s_v(i,j)=1.0d0/2.0d0*(v_fake(i-1,j)+v_fake(i,j)) > end if > if (i.eq.Ny) then > F_n_u(i,j)=0.0d0 > else if (j.eq.1) then > F_n_u(i,j)=0.0d0 > else > F_n_u(i,j)=1.0d0/2.0d0*(v_fake(i+1,j)+v_fake(i+1,j-1)) > end if > if (i.eq.Ny) then > F_n_v(i,j)=0.0d0 > else > F_n_v(i,j)=1.0d0/2.0d0*(v_fake(i,j)+v_fake(i+1,j)) > end if > end do > end do > > > > > > > > !interior > > do i=3,Ny-1 > do j=3,Nx-1 > !coeficientes de u > a_P(i,j)=(F_e_u(i,j)*dy+F_n_u(i,j)*dx+D_w*dy+D_e*dy+D_s*dx+D_n*dx)/ > alpha_u > a_W(i,j)=F_w_u(i,j)*dy+F_e_u(i,j)*1.0d0/8.0d0*dy+D_w*dy > a_E(i,j)=D_e*dy > a_S(i,j)=F_s_u(i,j)*dx+F_n_u(i,j)*1.0d0/8.0d0*dx+D_s*dx > a_N(i,j)=D_n*dx > > S_fake(i,j)=(3.0d0/8.0d0*u_fake(i,j)*dy-2.0d0/8.0d0*u_fake(i,j-1)*dy-1.0d > 0/8.0d0& > *u_fake(i,j-2)*dy)*F_w_u(i,j)+(-3.0d0/8.0d0*u_fake(i,j+1)*dy > +2.0d0/8.0d0*& > u_fake(i,j)*dy)*F_e_u(i,j) > +(3.0d0/8.0d0*u_fake(i,j)*dx-2.0d0/8.0d0*u_fake(i-1,j)& > *dx-1.0d0/8.0d0*u_fake(i-2,j)*dx)*F_s_u(i,j)+(-3.0d0/8.0d0*u_fake(i > +1,j)*dx+& > 2.0d0/8.0d0*u_fake(i,j)*dx)*F_n_u(i,j)+(p_fake(i,j-1)- > p_fake(i,j))*dy+((1.0d0-& > alpha_u)*a_P(i,j))*u_old(i,j) > > !coeficientes de v > b_P(i,j)=(F_e_v(i,j)*dy+F_n_v(i,j)*dx+D_w*dy+D_e*dy+D_s*dx+D_n*dx)/ > alpha_v > b_W(i,j)=F_w_v(i,j)*dy+F_e_v(i,j)*1.0d0/8.0d0*dy+D_w*dy > b_E(i,j)=D_e*dy > b_S(i,j)=F_s_v(i,j)*dx+F_n_v(i,j)*1.0d0/8.0d0*dx+D_s*dx > b_N(i,j)=D_n*dx > > S_v_fake(i,j)=(3.0d0/8.0d0*v_fake(i,j)*dy-2.0d0/8.0d0*v_fake(i,j-1)*dy-1. > 0d0/ > & > 8.0d0*v_fake(i,j-2)*dy)*F_w_v(i,j)+(-3.0d0/8.0d0*v_fake(i,j+1)*dy > +2.0d0/8.0d0*& > v_fake(i,j)*dy)*F_e_v(i,j) > +(3.0d0/8.0d0*v_fake(i,j)*dx-2.0d0/8.0d0*v_fake(i-1,j)& > *dx-1.0d0/8.0d0*v_fake(i-2,j)*dx)*F_s_v(i,j)+(-3.0d0/8.0d0*v_fake(i > +1,j)*dx+& > 2.0d0/8.0d0*v_fake(i,j)*dx)*F_n_v(i,j)+(p_fake(i-1,j)- > p_fake(i,j))*dx+((1.0d0-& > alpha_v)*b_P(i,j))*v_old(i,j) > > !coeficientes de p > c_W(i,j)=dy/a_P(i,j)*dy > c_E(i,j)=dy/a_E(i,j)*alpha_u*dy > c_S(i,j)=dx/b_P(i,j)*dx > c_N(i,j)=dx/b_N(i,j)*alpha_v*dx > c_P(i,j)=c_W(i,j)+c_E(i,j)+c_S(i,j)+c_N(i,j) > S_p_fake(i,j)=u_fake(i,j)*dy-u_fake(i,j+1)*dy+v_fake(i,j)*dx- > v_fake(i+1,j)*dx > > end do > end do > > !coluna 2 - com exce��o da linha 2 > > do i=3,Ny-1 > !coeficientes de u > a_P(i,2)=(F_e_u(i,2)*dy+F_n_u(i,2)*dx+D_w*dy+D_e*dy+D_s*dx+D_n*dx)/ > alpha_u > a_W(i,2)=F_w_u(i,2)*7.0d0/8.0d0*dy+F_e_u(i,2)*1.0d0/8.0d0*dy+D_w*dy > a_E(i,2)=D_e*dy > a_S(i,2)=F_s_u(i,2)*dx+F_n_u(i,2)*1.0d0/8.0d0*dx+D_s*dx > a_N(i,2)=D_n*dx > S_fake(i,2)=(3.0d0/8.0d0*u_fake(i,2)*dy-2.0d0/8.0d0*u_E*dy)*F_w_u(i, > 2)+(-3.0d0/8.0d0& > *u_fake(i,3)*dy+2.0d0/8.0d0*u_fake(i,2)*dy)*F_e_u(i, > 2)+(3.0d0/8.0d0*u_fake(i,2)*dx-& > 2.0d0/8.0d0*u_fake(i-1,2)*dx-1.0d0/8.0d0*u_fake(i-2,2)*dx)*F_s_u(i, > 2)+(-3.0d0/8.0d0*& > u_fake(i+1,2)*dx+2.0d0/8.0d0*u_fake(i,2)*dx)*F_n_u(i,2)+(p_fake(i, > 1)-p_fake(i,2))*dy& > +((1.0d0-alpha_u)*a_P(i,2))*u_old(i,2) > > !coeficientes de v > b_P(i,2)=(F_e_v(i,2)*dy+F_n_v(i,2)*dx+D_w*dy+D_e*dy+D_s*dx+D_n*dx)/ > alpha_v > b_W(i,2)=F_w_v(i,2)*7.0d0/8.0d0*dy+F_e_v(i,2)*1.0d0/8.0d0*dy+D_w*dy > b_E(i,2)=D_e*dy > b_S(i,2)=F_s_v(i,2)*dx+F_n_v(i,2)*1.0d0/8.0d0*dx+D_s*dx > b_N(i,2)=D_n*dx > S_v_fake(i,2)=(3.0d0/8.0d0*v_fake(i, > 2)*dy-2.0d0/8.0d0*v_I*dy)*F_w_v(i,2)+(-3.0d0/& > 8.0d0*v_fake(i,3)*dy+2.0d0/8.0d0*v_fake(i,2)*dy)*F_e_v(i, > 2)+(3.0d0/8.0d0*v_fake(i,2)& > > *dx-2.0d0/8.0d0*v_fake(i-1,2)*dx-1.0d0/8.0d0*v_fake(i-2,2)*dx)*F_s_v(i, > 2)+(-3.0d0/& > 8.0d0*v_fake(i+1,2)*dx+2.0d0/8.0d0*v_fake(i,2)*dx)*F_n_v(i,2)+ > (p_fake(i-1,2)-& > p_fake(i,2))*dx+((1.0d0-alpha_v)*b_P(i,2))*v_old(i,2) > > !coeficientes de p > c_W(i,2)=dy/a_P(i,2)*dy > c_E(i,2)=dy/a_E(i,2)*alpha_u*dy > c_S(i,2)=dx/b_P(i,2)*dx > c_N(i,2)=dx/b_N(i,2)*alpha_v*dx > c_P(i,2)=c_W(i,2)+c_E(i,2)+c_S(i,2)+c_N(i,2) > S_p_fake(i,2)=u_fake(i,2)*dy-u_fake(i,3)*dy+v_fake(i,2)*dx-v_fake(i > +1,2)*dx > > end do > > !linha 2 - com exce��o da coluna 2 > > do j=3,Nx-1 > !coeficientes de u > a_P(2,j)=(F_e_u(2,j)*dy+F_n_u(2,j)*dx+D_w*dy+D_e*dy+D_s*dx+D_n*dx)/ > alpha_u > a_W(2,j)=F_w_u(2,j)*dy+F_e_u(2,j)*1.0d0/8.0d0*dy+D_w*dy > a_E(2,j)=D_e*dy > a_S(2,j)=F_s_u(2,j)*7.0d0/8.0d0*dx+F_n_u(2,j)*1.0d0/8.0d0*dx+D_s*dx > a_N(2,j)=D_n*dx > > S_fake(2,j)=(3.0d0/8.0d0*u_fake(2,j)*dy-2.0d0/8.0d0*u_fake(2,j-1)*dy-1.0d > 0/8.0d0*& > u_fake(2,j-2)*dy)*F_w_u(2,j)+(-3.0d0/8.0d0*u_fake(2,j+1)*dy > +2.0d0/8.0d0*u_fake(2,j)*& > dy)*F_e_u(2,j) > +(3.0d0/8.0d0*u_fake(2,j)*dx-2.0d0/8.0d0*u_G*dx)*F_s_u(2,j)+(-3.0d0/& > 8.0d0*u_fake(3,j)*dx+2.0d0/8.0d0*u_fake(2,j)*dx)*F_n_u(2,j)+ > (p_fake(2,j-1)-& > p_fake(2,j))*dy+((1.0d0-alpha_u)*a_P(2,j))*u_old(2,j) > > !coeficientes de v > b_P(2,j)=(F_e_v(2,j)*dy+F_n_v(2,j)*dx+D_w*dy+D_e*dy+D_s*dx+D_n*dx)/ > alpha_v > b_W(2,j)=F_w_v(2,j)*dy+F_e_v(2,j)*1.0d0/8.0d0*dy+D_w*dy > b_E(2,j)=D_e*dy > b_S(2,j)=F_s_v(2,j)*7.0d0/8.0d0*dx+F_n_v(2,j)*1.0d0/8.0d0*dx+D_s*dx > b_N(2,j)=D_n*dx > > S_v_fake(2,j)=(3.0d0/8.0d0*v_fake(2,j)*dy-2.0d0/8.0d0*v_fake(2,j-1)*dy-1. > 0d0/8.0d0*& > v_fake(2,j-2)*dy)*F_w_v(2,j)+(-3.0d0/8.0d0*v_fake(2,j+1)*dy > +2.0d0/8.0d0*v_fake(2,j)*& > dy)*F_e_v(2,j) > +(3.0d0/8.0d0*v_fake(2,j)*dx-2.0d0/8.0d0*v_K*dx)*F_s_v(2,j)+(-3.0d0/& > 8.0d0*v_fake(3,j)*dx+2.0d0/8.0d0*v_fake(2,j)*dx)*F_n_v(2,j)+ > (p_fake(1,j)-& > p_fake(2,j))*dx+((1.0d0-alpha_v)*b_P(2,j))*v_old(2,j) > > !coeficientes de p > c_W(2,j)=dy/a_P(2,j)*dy > c_E(2,j)=dy/a_E(2,j)*alpha_u*dy > c_S(2,j)=dx/b_P(2,j)*dx > c_N(2,j)=dx/b_N(2,j)*alpha_v*dx > c_P(2,j)=c_W(2,j)+c_E(2,j)+c_S(2,j)+c_N(2,j) > S_p_fake(2,j)=u_fake(2,j)*dy-u_fake(2,j+1)*dy+v_fake(2,j)*dx- > v_fake(3,j)*dx > > end do > > !coluna 2 - linha 2 > > !coeficientes de u > a_P(2,2)=(F_e_u(2,2)*dy+F_n_u(2,2)*dx+D_w*dy+D_e*dy+D_s*dx+D_n*dx)/ > alpha_u > a_W(2,2)=F_w_u(2,2)*7.0d0/8.0d0*dy+F_e_u(2,2)*1.0d0/8.0d0*dy > +D_w*dy > a_E(2,2)=D_e*dy > a_S(2,2)=F_s_u(2,2)*7.0d0/8.0d0*dx+F_n_u(2,2)*1.0d0/8.0d0*dx > +D_s*dx > a_N(2,2)=D_n*dx > > S_fake(2,2)=(3.0d0/8.0d0*u_fake(2,2)*dy-2.0d0/8.0d0*u_E*dy)*F_w_u(2,2)+ > (-3.0d0/& > 8.0d0*u_fake(2,3)*dy > +2.0d0/8.0d0*u_fake(2,2)*dy)*F_e_u(2,2)+(3.0d0/8.0d0*& > u_fake(2,2)*dx-2.0d0/8.0d0*u_G*dx)*F_s_u(2,2)+ > (-3.0d0/8.0d0*u_fake(3,2)*dx+2.0d0& > /8.0d0*u_fake(2,2)*dx)*F_n_u(2,2)+(p_fake(2,1)-p_fake(2,2))*dy+ > ((1.0d0-alpha_u)*& > a_P(2,2))*u_old(2,2) > > !coeficientes de v > b_P(2,2)=(F_e_v(2,2)*dy+F_n_v(2,2)*dx+D_w*dy+D_e*dy+D_s*dx+D_n*dx)/ > alpha_v > b_W(2,2)=F_w_v(2,2)*7.0d0/8.0d0*dy+F_e_v(2,2)*1.0d0/8.0d0*dy > +D_w*dy > b_E(2,2)=D_e*dy > b_S(2,2)=F_s_v(2,2)*7.0d0/8.0d0*dx+F_n_v(2,2)*1.0d0/8.0d0*dx > +D_s*dx > b_N(2,2)=D_n*dx > > S_v_fake(2,2)=(3.0d0/8.0d0*v_fake(2,2)*dy-2.0d0/8.0d0*v_I*dy)*F_w_v(2,2)+ > (-3.0d0& > /8.0d0*v_fake(2,3)*dy > +2.0d0/8.0d0*v_fake(2,2)*dy)*F_e_v(2,2)+(3.0d0/8.0d0*& > v_fake(2,2)*dx-2.0d0/8.0d0*v_K*dx)*F_s_v(2,2)+ > (-3.0d0/8.0d0*v_fake(3,2)*dx+2.0d0& > /8.0d0*v_fake(2,2)*dx)*F_n_v(2,2)+(p_fake(1,2)-p_fake(2,2))*dx+ > ((1.0d0-alpha_v)*& > b_P(2,2))*v_old(2,2) > > !coeficientes de p > c_W(2,2)=dy/a_P(2,2)*dy > c_E(2,2)=dy/a_E(2,2)*alpha_u*dy > c_S(2,2)=dx/b_P(2,2)*dx > c_N(2,2)=dx/b_N(2,2)*alpha_v*dx > c_P(2,2)=c_W(2,2)+c_E(2,2)+c_S(2,2)+c_N(2,2) > S_p_fake(2,2)=u_fake(2,2)*dy-u_fake(2,3)*dy+v_fake(2,2)*dx- > v_fake(3,2)*dx > > > !fronteira W - com exce��o da linha 2 > > do i=3,Ny-1 > !coeficientes de u > a_P(i,1)=(F_e_u(i,1)*7.0d0/8.0d0*dy+F_n_u(i,1)*dx+D_e*dy+D_s*dx > +D_n*dx+D_EE*9.0d0/& > 3.0d0*dy)/alpha_u > a_W(i,1)=1.0d-20 !0.0d0 > a_E(i,1)=D_e*dy+D_EE*1.0d0/3.0d0*dy > a_S(i,1)=F_s_u(i,1)*dx+F_n_u(i,1)*1.0d0/8.0d0*dx+D_s*dx > a_N(i,1)=D_n*dx > S_fake(i,1)=(u_E*dy)*F_w_u(i,1)+(-3.0d0/8.0d0*u_fake(i,2)*dy > +2.0d0/8.0d0*u_E*dy)*& > F_e_u(i,1)+(3.0d0/8.0d0*u_fake(i, > 1)*dx-2.0d0/8.0d0*u_fake(i-1,1)*dx-1.0d0/8.0d0*& > u_fake(i-2,1)*dx)*F_s_u(i,1)+(-3.0d0/8.0d0*u_fake(i+1,1)*dx > +2.0d0/8.0d0*u_fake(i,1)*& > dx)*F_n_u(i,1)+D_EE*8.0d0/3.0d0*u_E*dy+((1.0d0-alpha_u)*a_P(i, > 1))*u_old(i,1) !Falta o termo de press�o. > > > ! a_P(i,1)=1.0d0 > ! a_W(i,1)=1.0d-20 > ! a_E(i,1)=1.0d-20 > ! a_S(i,1)=1.0d-20 > ! a_N(i,1)=1.0d-20 > ! S_fake(i,1)=u_E > > > > > > > !coeficientes de v > b_P(i,1)=(F_e_v(i,1)*7.0d0/8.0d0*dy+F_n_v(i,1)*dx+D_e*dy+D_s*dx > +D_n*dx+D_I*9.0d0/& > 3.0d0*dy)/alpha_v > b_W(i,1)=1.0d-20 !0.0d0 > b_E(i,1)=D_e*dy+D_I*1.0d0/3.0d0*dy > b_S(i,1)=F_s_v(i,1)*dx+F_n_v(i,1)*1.0d0/8.0d0*dx+D_s*dx > b_N(i,1)=D_n*dx > S_v_fake(i,1)=(v_I*dy)*F_w_v(i,1)+(-3.0d0/8.0d0*v_fake(i,2)*dy > +2.0d0/8.0d0*v_I*dy)*& > F_e_v(i,1)+(3.0d0/8.0d0*v_fake(i, > 1)*dx-2.0d0/8.0d0*v_fake(i-1,1)*dx-1.0d0/8.0d0*& > v_fake(i-2,1)*dx)*F_s_v(i,1)+(-3.0d0/8.0d0*v_fake(i+1,1)*dx > +2.0d0/8.0d0*v_fake(i,1)*& > dx)*F_n_v(i,1)+(p_fake(i-1,1)-p_fake(i,1))*dx+D_I*8.0d0/3.0d0*v_I*dy > +((1.0d0-alpha_v& > )*b_P(i,1))*v_old(i,1) > > ! b_P(i,1)=1.0d0 > ! b_W(i,1)=1.0d-20 > ! b_E(i,1)=1.0d-20 > ! b_S(i,1)=1.0d-20 > ! b_N(i,1)=1.0d-20 > ! S_v_fake(i,1)=v_I > > > > > !coeficientes de p > > c_W(i,1)=dy/a_P(i,1)*dy > c_E(i,1)=dy/a_E(i,1)*alpha_u*dy > c_S(i,1)=dx/b_P(i,1)*dx > c_N(i,1)=dx/b_N(i,1)*alpha_v*dx > c_P(i,1)=c_W(i,1)+c_E(i,1)+c_S(i,1)+c_N(i,1) > S_p_fake(i,1)=u_fake(i,1)*dy-u_fake(i,2)*dy+v_fake(i,1)*dx-v_fake(i > +1,1)*dx > > end do > > !fronteira W - linha 2 > > !coeficientes de u > a_P(2,1)=(F_e_u(2,1)*7.0d0/8.0d0*dy+F_n_u(2,1)*dx+D_e*dy+D_s*dx > +D_n*dx+D_EE*9.0d0/& > 3.0d0*dy)/alpha_u > a_W(2,1)=1.0d-20 !0.0d0 > a_E(2,1)=D_e*dy+D_EE*1.0d0/3.0d0*dy > a_S(2,1)=F_s_u(2,1)*7.0d0/8.0d0*dx+F_n_u(2,1)*1.0d0/8.0d0*dx+D_s*dx > a_N(2,1)=D_n*dx > S_fake(2,1)=(u_E*dy)*F_w_u(2,1)+(-3.0d0/8.0d0*u_fake(2,2)*dy > +2.0d0/8.0d0*u_E*dy)*& > > F_e_u(2,1)+(3.0d0/8.0d0*u_fake(2,1)*dx-2.0d0/8.0d0*u_G*dx)*F_s_u(2,1)+ > (-3.0d0/8.0d0*& > u_fake(3,1)*dx > +2.0d0/8.0d0*u_fake(2,1)*dx)*F_n_u(2,1)+D_EE*8.0d0/3.0d0*u_E*dy+((& > 1.0d0-alpha_u)*a_P(2,1))*u_old(2,1) !Falta o termo de press�o. > > > ! a_P(2,1)=1.0d0 > ! a_W(2,1)=1.0d-20 > ! a_E(2,1)=1.0d-20 > ! a_S(2,1)=1.0d-20 > ! a_N(2,1)=1.0d-20 > ! S_fake(2,1)=u_E > > > > > !coeficientes de v > b_P(2,1)=(F_e_v(2,1)*7.0d0/8.0d0*dy+F_n_v(2,1)*dx+D_e*dy+D_s*dx > +D_n*dx+D_I*9.0d0/& > 3.0d0*dy)/alpha_v > b_W(2,1)=1.0d-20 !0.0d0 > b_E(2,1)=D_e*dy+D_I*1.0d0/3.0d0*dy > b_S(2,1)=F_s_v(2,1)*7.0d0/8.0d0*dx+F_n_v(2,1)*1.0d0/8.0d0*dx+D_s*dx > b_N(2,1)=D_n*dx > S_v_fake(2,1)=(v_I*dy)*F_w_v(2,1)+(-3.0d0/8.0d0*v_fake(2,2)*dy > +2.0d0/8.0d0*v_I*dy)*& > > F_e_v(2,1)+(3.0d0/8.0d0*v_fake(2,1)*dx-2.0d0/8.0d0*v_K*dx)*F_s_v(2,1)+ > (-3.0d0/8.0d0*& > v_fake(3,1)*dx+2.0d0/8.0d0*v_fake(2,1)*dx)*F_n_v(2,1)+(p_fake(1,1)- > p_fake(2,1))*dx+& > D_I*8.0d0/3.0d0*v_I*dy+((1.0d0-alpha_v)*b_P(2,1))*v_old(2,1) > > ! b_P(2,1)=1.0d0 > ! b_W(2,1)=1.0d-20 > ! b_E(2,1)=1.0d-20 > ! b_S(2,1)=1.0d-20 > ! b_N(2,1)=1.0d-20 > ! S_v_fake(2,1)=v_I > > > > > !coeficientes de p > c_W(2,1)=dy/a_P(2,1)*dy > c_E(2,1)=dy/a_E(2,1)*alpha_u*dy > c_S(2,1)=dx/b_P(2,1)*dx > c_N(2,1)=dx/b_N(2,1)*alpha_v*dx > c_P(2,1)=c_W(2,1)+c_E(2,1)+c_S(2,1)+c_N(2,1) > S_p_fake(2,1)=u_fake(2,1)*dy-u_fake(2,2)*dy+v_fake(2,1)*dx- > v_fake(3,1)*dx > > > > !fronteira E - com exce��o da linha 2 > > do i=3,Ny-1 > !coeficientes de u > a_P(i,Nx)=(F_n_u(i,Nx)*dx+D_w*dy+D_s*dx+D_n*dx+D_F*9.0d0/3.0d0*dy)/ > alpha_u > a_W(i,Nx)=F_w_u(i,Nx)*dy+D_w*dy+D_F*1.0d0/3.0d0*dy > a_E(i,Nx)=1.0d-20 !0.0d0 > a_S(i,Nx)=F_s_u(i,Nx)*dx+F_n_u(i,Nx)*1.0d0/8.0d0*dx+D_s*dx > a_N(i,Nx)=D_n*dx > > S_fake(i,Nx)=(3.0d0/8.0d0*u_fake(i,Nx)*dy-2.0d0/8.0d0*u_fake(i,Nx-1)*dy-1 > .0d0/8.0d0*& > u_fake(i,Nx-2)*dy)*F_w_u(i,Nx)+(-u_F(i)*dy)*F_e_u(i,Nx) > +(3.0d0/8.0d0*u_fake(i,Nx)*dx& > -2.0d0/8.0d0*u_fake(i-1,Nx)*dx-1.0d0/8.0d0*u_fake(i-2,Nx)*dx)*F_s_u(i,Nx > ) > +(-3.0d0/& > 8.0d0*u_fake(i+1,Nx)*dx+2.0d0/8.0d0*u_fake(i,Nx)*dx)*F_n_u(i,Nx)+ > (p_fake(i,Nx-1)-& > p_fake(i,Nx))*dy+D_F*8.0d0/3.0d0*u_F(i)*dy+((1.0d0- > alpha_u)*a_P(i,Nx))*u_old(i,Nx) > > > ! a_P(i,Nx)=1.0d0 > ! a_W(i,Nx)=1.0d-20 > ! a_E(i,Nx)=1.0d-20 > ! a_S(i,Nx)=1.0d-20 > ! a_N(i,Nx)=1.0d-20 > ! S_fake(i,Nx)=u_F(i) > > > > > > !coeficientes de v > b_P(i,Nx)=(F_n_v(i,Nx)*dx+D_w*dy+D_s*dx+D_n*dx+D_J*9.0d0/3.0d0*dy)/ > alpha_v > b_W(i,Nx)=F_w_v(i,Nx)*dy+D_w*dy+D_J*1.0d0/3.0d0*dy > b_E(i,Nx)=1.0d-20 !0.0d0 > b_S(i,Nx)=F_s_v(i,Nx)*dx+F_n_v(i,Nx)*1.0d0/8.0d0*dx+D_s*dx > b_N(i,Nx)=D_n*dx > > S_v_fake(i,Nx)=(3.0d0/8.0d0*v_fake(i,Nx)*dy-2.0d0/8.0d0*v_fake(i,Nx-1)*dy > -1.0d0/ > & > 8.0d0*v_fake(i,Nx-2)*dy)*F_w_v(i,Nx)+(-v_J*dy)*F_e_v(i,Nx) > +(3.0d0/8.0d0*v_fake(i,Nx)& > > *dx-2.0d0/8.0d0*v_fake(i-1,Nx)*dx-1.0d0/8.0d0*v_fake(i-2,Nx)*dx)*F_s_v(i,Nx > ) > +(-3.0d0& > /8.0d0*v_fake(i+1,Nx)*dx+2.0d0/8.0d0*v_fake(i,Nx)*dx)*F_n_v(i,Nx)+ > (p_fake(i-1,Nx)-& > p_fake(i,Nx))*dx+D_J*8.0d0/3.0d0*v_J*dy+((1.0d0- > alpha_v)*b_P(i,Nx))*v_old(i,Nx) > > > ! b_P(i,Nx)=1.0d0 > ! b_W(i,Nx)=1.0d-20 > ! b_E(i,Nx)=1.0d-20 > ! b_S(i,Nx)=1.0d-20 > ! b_N(i,Nx)=1.0d-20 > ! S_v_fake(i,Nx)=v_J > > > > > !coeficientes de p > c_W(i,Nx)=dy/a_P(i,Nx)*dy > c_E(i,Nx)=dy/a_E(i,Nx)*alpha_u*dy > c_S(i,Nx)=dx/b_P(i,Nx)*dx > c_N(i,Nx)=dx/b_N(i,Nx)*alpha_v*dx > c_P(i,Nx)=c_W(i,Nx)+c_E(i,Nx)+c_S(i,Nx)+c_N(i,Nx) > S_p_fake(i,Nx)=v_fake(i,Nx)*dx-v_fake(i+1,Nx)*dx !Faltam termos de > u. > > end do > > !fronteira E - linha 2 > > !coeficientes de u > a_P(2,Nx)=(F_n_u(2,Nx)*dx+D_w*dy+D_s*dx+D_n*dx+D_F*9.0d0/3.0d0*dy)/ > alpha_u > a_W(2,Nx)=F_w_u(2,Nx)*dy+D_w*dy+D_F*1.0d0/3.0d0*dy > a_E(2,Nx)=1.0d-20 !0.0d0 > a_S(2,Nx)=F_s_u(2,Nx)*7.0d0/8.0d0*dx+F_n_u(2,Nx)*1.0d0/8.0d0*dx > +D_s*dx > a_N(2,Nx)=D_n*dx > > S_fake(2,Nx)=(3.0d0/8.0d0*u_fake(2,Nx)*dy-2.0d0/8.0d0*u_fake(2,Nx-1)*dy-1 > .0d0/8.0d0*& > u_fake(2,Nx-2)*dy)*F_w_u(2,Nx)+(-u_F(2)*dy)*F_e_u(2,Nx) > +(3.0d0/8.0d0*u_fake(2,Nx)*dx& > -2.0d0/8.0d0*u_G*dx)*F_s_u(2,Nx)+(-3.0d0/8.0d0*u_fake(3,Nx)*dx > +2.0d0/8.0d0*& > u_fake(2,Nx)*dx)*F_n_u(2,Nx)+(p_fake(2,Nx-1)-p_fake(2,Nx))*dy > +D_F*8.0d0/3.0d0*u_F(2)& > *dy+((1.0d0-alpha_u)*a_P(2,Nx))*u_old(2,Nx) > > > ! a_P(2,Nx)=1.0d0 > ! a_W(2,Nx)=1.0d-20 > ! a_E(2,Nx)=1.0d-20 > ! a_S(2,Nx)=1.0d-20 > ! a_N(2,Nx)=1.0d-20 > ! S_fake(2,Nx)=u_F(2) > > > > !coeficientes de v > b_P(2,Nx)=(F_n_v(2,Nx)*dx+D_w*dy+D_s*dx+D_n*dx+D_J*9.0d0/3.0d0*dy)/ > alpha_v > b_W(2,Nx)=F_w_v(2,Nx)*dy+D_w*dy+D_J*1.0d0/3.0d0*dy > b_E(2,Nx)=1.0d-20 !0.0d0 > b_S(2,Nx)=F_s_v(2,Nx)*7.0d0/8.0d0*dx+F_n_v(2,Nx)*1.0d0/8.0d0*dx > +D_s*dx > b_N(2,Nx)=D_n*dx > > S_v_fake(2,Nx)=(3.0d0/8.0d0*v_fake(2,Nx)*dy-2.0d0/8.0d0*v_fake(2,Nx-1)*dy > -1.0d0/ > & > 8.0d0*v_fake(2,Nx-2)*dy)*F_w_v(2,Nx)+(-v_J*dy)*F_e_v(2,Nx) > +(3.0d0/8.0d0*v_fake(2,Nx)& > *dx-2.0d0/8.0d0*v_K*dx)*F_s_v(2,Nx)+(-3.0d0/8.0d0*v_fake(3,Nx)*dx > +2.0d0/8.0d0*& > v_fake(2,Nx)*dx)*F_n_v(2,Nx)+(p_fake(1,Nx)-p_fake(2,Nx))*dx > +D_J*8.0d0/3.0d0*v_J*dy+(& > (1.0d0-alpha_v)*b_P(2,Nx))*v_old(2,Nx) > > > ! b_P(2,Nx)=1.0d0 > ! b_W(2,Nx)=1.0d-20 > ! b_E(2,Nx)=1.0d-20 > ! b_S(2,Nx)=1.0d-20 > ! b_N(2,Nx)=1.0d-20 > ! S_v_fake(2,Nx)=v_J > > !coeficientes de p > c_W(2,Nx)=dy/a_P(2,Nx)*dy > c_E(2,Nx)=dy/a_E(2,Nx)*alpha_u*dy > c_S(2,Nx)=dx/b_P(2,Nx)*dx > c_N(2,Nx)=dx/b_N(2,Nx)*alpha_v*dx > c_P(2,Nx)=c_W(2,Nx)+c_E(2,Nx)+c_S(2,Nx)+c_N(2,Nx) > S_p_fake(2,Nx)=v_fake(2,Nx)*dx-v_fake(3,Nx)*dx !Faltam termos de u. > > > !fronteira S - com exce��o da coluna 2 > > do j=3,Nx-1 > !coeficientes de u > a_P(1,j)=(F_e_u(1,j)*dy+F_n_u(1,j)*7.0d0/8.0d0*dx+D_w*dy+D_e*dy > +D_n*dx+D_G*9.0d0/& > 3.0d0*dx)/alpha_u > a_W(1,j)=F_w_u(1,j)*dy+F_e_u(1,j)*1.0d0/8.0d0*dy+D_w*dy > a_E(1,j)=D_e*dy > a_S(1,j)=1.0d-20 !0.0d0 > a_N(1,j)=D_n*dx+D_G*1.0d0/3.0d0*dx > > S_fake(1,j)=(3.0d0/8.0d0*u_fake(1,j)*dy-2.0d0/8.0d0*u_fake(1,j-1)*dy-1.0d > 0/8.0d0*& > u_fake(1,j-2)*dy)*F_w_u(1,j)+(-3.0d0/8.0d0*u_fake(1,j+1)*dy > +2.0d0/8.0d0*u_fake(1,j)*& > dy)*F_e_u(1,j)+(u_G*dx)*F_s_u(1,j)+(-3.0d0/8.0d0*u_fake(2,j)*dx > +2.0d0/8.0d0*u_G*dx)*& > F_n_u(1,j)+(p_fake(1,j-1)-p_fake(1,j))*dy+D_G*8.0d0/3.0d0*u_G*dx+ > ((1.0d0-alpha_u)*& > a_P(1,j))*u_old(1,j) > > > ! a_P(1,j)=1.0d0 > ! a_W(1,j)=1.0d-20 > ! a_E(1,j)=1.0d-20 > ! a_S(1,j)=1.0d-20 > ! a_N(1,j)=1.0d-20 > ! S_fake(1,j)=u_G > > > > > > > !coeficientes de v > b_P(1,j)=(F_e_v(1,j)*dy+F_n_v(1,j)*7.0d0/8.0d0*dx+D_w*dy+D_e*dy > +D_n*dx+D_K*9.0d0/& > 3.0d0*dx)/alpha_v > b_W(1,j)=F_w_v(1,j)*dy+F_e_v(1,j)*1.0d0/8.0d0*dy+D_w*dy > b_E(1,j)=D_e*dy > b_S(1,j)=1.0d-20 !0.0d0 > b_N(1,j)=D_n*dx+D_K*1.0d0/3.0d0*dx > > S_v_fake(1,j)=(3.0d0/8.0d0*v_fake(1,j)*dy-2.0d0/8.0d0*v_fake(1,j-1)*dy-1. > 0d0/8.0d0*& > v_fake(1,j-2)*dy)*F_w_v(1,j)+(-3.0d0/8.0d0*v_fake(1,j+1)*dy > +2.0d0/8.0d0*v_fake(1,j)*& > dy)*F_e_v(1,j)+(v_K*dx)*F_s_v(1,j)+(-3.0d0/8.0d0*v_fake(2,j)*dx > +2.0d0/8.0d0*v_K*dx)*& > F_n_v(1,j)+D_K*8.0d0/3.0d0*v_K*dx+((1.0d0- > alpha_v)*b_P(1,j))*v_old(1,j) !Falta o termo de press�o. > > ! b_P(1,j)=1.0d0 > ! b_W(1,j)=1.0d-20 > ! b_E(1,j)=1.0d-20 > ! b_S(1,j)=1.0d-20 > ! b_N(1,j)=1.0d-20 > ! S_v_fake(1,j)=v_K > > > > > > !coeficientes de p > c_W(1,j)=dy/a_P(1,j)*dy > c_E(1,j)=dy/a_E(1,j)*alpha_u*dy > c_S(1,j)=dx/b_P(1,j)*dx > c_N(1,j)=dx/b_N(1,j)*alpha_v*dx > c_P(1,j)=c_W(1,j)+c_E(1,j)+c_S(1,j)+c_N(1,j) > S_p_fake(1,j)=u_fake(1,j)*dy-u_fake(1,j+1)*dy+v_fake(1,j)*dx- > v_fake(2,j)*dx > > end do > > !fronteira S - coluna 2 > > !coeficientes de u > a_P(1,2)=(F_e_u(1,2)*dy+F_n_u(1,2)*7.0d0/8.0d0*dx+D_w*dy+D_e*dy+D_n*dx > +D_G*9.0d0/3.0d0*dx)/& > alpha_u > a_W(1,2)=F_w_u(1,2)*7.0d0/8.0d0*dy+F_e_u(1,2)*1.0d0/8.0d0*dy+D_w*dy > a_E(1,2)=D_e*dy > a_S(1,2)=1.0d-20 !0.0d0 > a_N(1,2)=D_n*dx+D_G*1.0d0/3.0d0*dx > > S_fake(1,2)=(3.0d0/8.0d0*u_fake(1,2)*dy-2.0d0/8.0d0*u_E*dy)*F_w_u(1,2)+ > (-3.0d0/8.0d0*& > u_fake(1,3)*dy+2.0d0/8.0d0*u_fake(1,2)*dy)*F_e_u(1,2)+ > (u_G*dx)*F_s_u(1,2)+(-3.0d0/8.0d0*& > u_fake(2,2)*dx+2.0d0/8.0d0*u_G*dx)*F_n_u(1,2)+(p_fake(1,1)- > p_fake(1,2))*dy+D_G*8.0d0/3.0d0*& > u_G*dx+((1.0d0-alpha_u)*a_P(1,2))*u_old(1,2) > > ! a_P(1,2)=1.0d0 > ! a_W(1,2)=1.0d-20 > ! a_E(1,2)=1.0d-20 > ! a_S(1,2)=1.0d-20 > ! a_N(1,2)=1.0d-20 > ! S_fake(1,2)=u_G > > > > > > > > > !coeficientes de v > b_P(1,2)=(F_e_v(1,2)*dy+F_n_v(1,2)*7.0d0/8.0d0*dx+D_w*dy+D_e*dy+D_n*dx > +D_K*9.0d0/3.0d0*dx)/& > alpha_v > b_W(1,2)=F_w_v(1,2)*7.0d0/8.0d0*dy+F_e_v(1,2)*1.0d0/8.0d0*dy+D_w*dy > b_E(1,2)=D_e*dy > b_S(1,2)=1.0d-20 !0.0d0 > b_N(1,2)=D_n*dx+D_K*1.0d0/3.0d0*dx > > S_v_fake(1,2)=(3.0d0/8.0d0*v_fake(1,2)*dy-2.0d0/8.0d0*v_I*dy)*F_w_v(1,2)+ > (-3.0d0/8.0d0*& > v_fake(1,3)*dy+2.0d0/8.0d0*v_fake(1,2)*dy)*F_e_v(1,2)+ > (v_K*dx)*F_s_v(1,2)+(-3.0d0/8.0d0*& > v_fake(2,2)*dx+2.0d0/8.0d0*v_K*dx)*F_n_v(1,2)+D_K*8.0d0/3.0d0*v_K*dx+ > ((1.0d0-alpha_v)*& > b_P(1,2))*v_old(1,2) !Falta o termo de press�o. > > ! b_P(1,2)=1.0d0 > ! b_W(1,2)=1.0d-20 > ! b_E(1,2)=1.0d-20 > ! b_S(1,2)=1.0d-20 > ! b_N(1,2)=1.0d-20 > ! S_v_fake(1,2)=v_K > > > > > > !coeficientes de p > c_W(1,2)=dy/a_P(1,2)*dy > c_E(1,2)=dy/a_E(1,2)*alpha_u*dy > c_S(1,2)=dx/b_P(1,2)*dx > c_N(1,2)=dx/b_N(1,2)*alpha_v*dx > c_P(1,2)=c_W(1,2)+c_E(1,2)+c_S(1,2)+c_N(1,2) > S_p_fake(1,2)=u_fake(1,2)*dy-u_fake(1,3)*dy+v_fake(1,2)*dx- > v_fake(2,2)*dx > > > !fronteira N - com exce��o da coluna 2 > > do j=3,Nx-1 > !coeficientes de u > a_P(Ny,j)=(F_e_u(Ny,j)*dy+D_w*dy+D_e*dy+D_s*dx+D_H*9.0d0/3.0d0*dx)/ > alpha_u > a_W(Ny,j)=F_w_u(Ny,j)*dy+F_e_u(Ny,j)*1.0d0/8.0d0*dy+D_w*dy > a_E(Ny,j)=D_e*dy > a_S(Ny,j)=F_s_u(Ny,j)*dx+D_s*dx+D_H*1.0d0/3.0d0*dx > a_N(Ny,j)=1.0d-20 !0.0d0 > > S_fake(Ny,j)=(3.0d0/8.0d0*u_fake(Ny,j)*dy-2.0d0/8.0d0*u_fake(Ny,j-1)*dy-1 > .0d0/8.0d0*& > u_fake(Ny,j-2)*dy)*F_w_u(Ny,j)+(-3.0d0/8.0d0*u_fake(Ny,j+1)*dy > +2.0d0/8.0d0*& > u_fake(Ny,j)*dy)*F_e_u(Ny,j) > +(3.0d0/8.0d0*u_fake(Ny,j)*dx-2.0d0/8.0d0*u_fake(Ny-1,j)& > *dx-1.0d0/8.0d0*u_fake(Ny-2,j)*dx)*F_s_u(Ny,j)+(-u_H*dx)*F_n_u(Ny,j) > +(p_fake(Ny,j-1)& > -p_fake(Ny,j))*dy+D_H*8.0d0/3.0d0*u_H*dx+((1.0d0- > alpha_u)*a_P(Ny,j))*u_old(Ny,j) > > > ! a_P(Ny,j)=1.0d0 > ! a_W(Ny,j)=1.0d-20 > ! a_E(Ny,j)=1.0d-20 > ! a_S(Ny,j)=1.0d-20 > ! a_N(Ny,j)=1.0d-20 > ! S_fake(Ny,j)=u_H > > > > > !coeficientes de v > b_P(Ny,j)=(F_e_v(Ny,j)*dy+D_w*dy+D_e*dy+D_s*dx+D_L*9.0d0/3.0d0*dx)/ > alpha_v > b_W(Ny,j)=F_w_v(Ny,j)*dy+F_e_v(Ny,j)*1.0d0/8.0d0*dy+D_w*dy > b_E(Ny,j)=D_e*dy > b_S(Ny,j)=F_s_v(Ny,j)*dx+D_s*dx+D_L*1.0d0/3.0d0*dx > b_N(Ny,j)=1.0d-20 !0.0d0 > > S_v_fake(Ny,j)=(3.0d0/8.0d0*v_fake(Ny,j)*dy-2.0d0/8.0d0*v_fake(Ny,j-1)*dy > -1.0d0/ > & > 8.0d0*v_fake(Ny,j-2)*dy)*F_w_v(Ny,j)+(-3.0d0/8.0d0*v_fake(Ny,j+1)*dy > +2.0d0/8.0d0*& > v_fake(Ny,j)*dy)*F_e_v(Ny,j) > +(3.0d0/8.0d0*v_fake(Ny,j)*dx-2.0d0/8.0d0*v_fake(Ny-1,j)& > *dx-1.0d0/8.0d0*v_fake(Ny-2,j)*dx)*F_s_v(Ny,j)+(-v_L*dx)*F_n_v(Ny,j) > +(p_fake(Ny-1,j)& > -p_fake(Ny,j))*dx+D_L*8.0d0/3.0d0*v_L*dx+((1.0d0- > alpha_v)*b_P(Ny,j))*v_old(Ny,j) > > > ! b_P(Ny,j)=1.0d0 > ! b_W(Ny,j)=1.0d-20 > ! b_E(Ny,j)=1.0d-20 > ! b_S(Ny,j)=1.0d-20 > ! b_N(Ny,j)=1.0d-20 > ! S_v_fake(Ny,j)=v_L > > > > > > !coeficientes de p > c_W(Ny,j)=dy/a_P(Ny,j)*dy > c_E(Ny,j)=dy/a_E(Ny,j)*alpha_u*dy > c_S(Ny,j)=dx/b_P(Ny,j)*dx > c_N(Ny,j)=dx/b_N(Ny,j)*alpha_v*dx > c_P(Ny,j)=c_W(Ny,j)+c_E(Ny,j)+c_S(Ny,j)+c_N(Ny,j) > S_p_fake(Ny,j)=u_fake(Ny,j)*dy-u_fake(Ny,j+1)*dy !Faltam termos de > v. > > > end do > > !fronteira N - coluna 2 > > !coeficientes de u > a_P(Ny,2)=(F_e_u(Ny,2)*dy+D_w*dy+D_e*dy+D_s*dx+D_H*9.0d0/3.0d0*dx)/ > alpha_u > a_W(Ny,2)=F_w_u(Ny,2)*7.0d0/8.0d0*dy+F_e_u(Ny,2)*1.0d0/8.0d0*dy > +D_w*dy > a_E(Ny,2)=D_e*dy > a_S(Ny,2)=F_s_u(Ny,2)*dx+D_s*dx+D_H*1.0d0/3.0d0*dx > a_N(Ny,2)=1.0d-20 !0.0d0 > S_fake(Ny,2)=(3.0d0/8.0d0*u_fake(Ny, > 2)*dy-2.0d0/8.0d0*u_E*dy)*F_w_u(Ny,2)+(-3.0d0/& > 8.0d0*u_fake(Ny,3)*dy+2.0d0/8.0d0*u_fake(Ny,2)*dy)*F_e_u(Ny, > 2)+(3.0d0/8.0d0*& > u_fake(Ny, > 2)*dx-2.0d0/8.0d0*u_fake(Ny-1,2)*dx-1.0d0/8.0d0*u_fake(Ny-2,2)*dx)*& > F_s_u(Ny,2)+(-u_H*dx)*F_n_u(Ny,2)+(p_fake(Ny,1)-p_fake(Ny,2))*dy > +D_H*8.0d0/3.0d0*u_H& > *dx+((1.0d0-alpha_u)*a_P(Ny,2))*u_old(Ny,2) > > ! a_P(Ny,2)=1.0d0 > ! a_W(Ny,2)=1.0d-20 > ! a_E(Ny,2)=1.0d-20 > ! a_S(Ny,2)=1.0d-20 > ! a_N(Ny,2)=1.0d-20 > ! S_fake(Ny,2)=u_H > > > > > !coeficientes de v > b_P(Ny,2)=(F_e_v(Ny,2)*dy+D_w*dy+D_e*dy+D_s*dx+D_L*9.0d0/3.0d0*dx)/ > alpha_v > b_W(Ny,2)=F_w_v(Ny,2)*7.0d0/8.0d0*dy+F_e_v(Ny,2)*1.0d0/8.0d0*dy > +D_w*dy > b_E(Ny,2)=D_e*dy > b_S(Ny,2)=F_s_v(Ny,2)*dx+D_s*dx+D_L*1.0d0/3.0d0*dx > b_N(Ny,2)=1.0d-20 !0.0d0 > S_v_fake(Ny,2)=(3.0d0/8.0d0*v_fake(Ny, > 2)*dy-2.0d0/8.0d0*v_I*dy)*F_w_v(Ny,2)+(-3.0d0/& > 8.0d0*v_fake(Ny,3)*dy+2.0d0/8.0d0*v_fake(Ny,2)*dy)*F_e_v(Ny, > 2)+(3.0d0/8.0d0*& > v_fake(Ny, > 2)*dx-2.0d0/8.0d0*v_fake(Ny-1,2)*dx-1.0d0/8.0d0*v_fake(Ny-2,2)*dx)*& > F_s_v(Ny,2)+(-v_L*dx)*F_n_v(Ny,2)+(p_fake(Ny-1,2)-p_fake(Ny,2))*dx > +D_L*8.0d0/3.0d0*& > v_L*dx+((1.0d0-alpha_v)*b_P(Ny,2))*v_old(Ny,2) > > ! b_P(Ny,2)=1.0d0 > ! b_W(Ny,2)=1.0d-20 > ! b_E(Ny,2)=1.0d-20 > ! b_S(Ny,2)=1.0d-20 > ! b_N(Ny,2)=1.0d-20 > ! S_v_fake(Ny,2)=v_L > > > > > > > > !coeficientes de p > c_W(Ny,2)=dy/a_P(Ny,2)*dy > c_E(Ny,2)=dy/a_E(Ny,2)*alpha_u*dy > c_S(Ny,2)=dx/b_P(Ny,2)*dx > c_N(Ny,2)=dx/b_N(Ny,2)*alpha_v*dx > c_P(Ny,2)=c_W(Ny,2)+c_E(Ny,2)+c_S(Ny,2)+c_N(Ny,2) > S_p_fake(Ny,2)=u_fake(Ny,2)*dy-u_fake(Ny,3)*dy !Faltam termos de v. > > > !fronteira SW > > !coeficientes de u > a_P(1,1)=(F_e_u(1,1)*7.0d0/8.0d0*dy+F_n_u(1,1)*7.0d0/8.0d0*dx+D_e*dy > +D_n*dx+D_EE*& > 9.0d0/3.0d0*dy+D_G*9.0d0/3.0d0*dx)/alpha_u > a_W(1,1)=1.0d-20 !0.0d0 > a_E(1,1)=D_e*dy+D_EE*1.0d0/3.0d0*dy > a_S(1,1)=1.0d-20 !0.0d0 > a_N(1,1)=D_n*dx+D_G*1.0d0/3.0d0*dx > S_fake(1,1)=(u_E*dy)*F_w_u(1,1)+(-3.0d0/8.0d0*u_fake(1,2)*dy > +2.0d0/8.0d0*u_E*dy)*& > F_e_u(1,1)+(u_G*dx)*F_s_u(1,1)+(-3.0d0/8.0d0*u_fake(2,1)*dx > +2.0d0/8.0d0*u_G*dx)*& > F_n_u(1,1)+D_EE*8.0d0/3.0d0*u_E*dy+D_G*8.0d0/3.0d0*u_G*dx+((1.0d0- > alpha_u)*a_P(1,1))& > *u_old(1,1) !Falta o termo de press�o. > > > ! a_P(1,1)=1.0d0 > ! a_W(1,1)=1.0d-20 > ! a_E(1,1)=1.0d-20 > ! a_S(1,1)=1.0d-20 > ! a_N(1,1)=1.0d-20 > ! S_fake(1,1)=1.0d-30 > > > > > > > > > !coeficientes de v > b_P(1,1)=(F_e_v(1,1)*7.0d0/8.0d0*dy+F_n_v(1,1)*7.0d0/8.0d0*dx+D_e*dy > +D_n*dx+D_I*& > 9.0d0/3.0d0*dy+D_K*9.0d0/3.0d0*dx)/alpha_v > b_W(1,1)=1.0d-20 !0.0d0 > b_E(1,1)=D_e*dy+D_I*1.0d0/3.0d0*dy > b_S(1,1)=1.0d-20 !0.0d0 > b_N(1,1)=D_n*dx+D_K*1.0d0/3.0d0*dx > S_v_fake(1,1)=(v_I*dy)*F_w_v(1,1)+(-3.0d0/8.0d0*v_fake(1,2)*dy > +2.0d0/8.0d0*v_I*dy)*& > F_e_v(1,1)+(v_K*dx)*F_s_v(1,1)+(-3.0d0/8.0d0*v_fake(2,1)*dx > +2.0d0/8.0d0*v_K*dx)*& > F_n_v(1,1)+D_I*8.0d0/3.0d0*v_I*dy+D_K*8.0d0/3.0d0*v_K*dx+((1.0d0- > alpha_v)*b_P(1,1))& > *v_old(1,1) !Falta o termo de press�o. > > ! b_P(1,1)=1.0d0 > ! b_W(1,1)=1.0d-20 > ! b_E(1,1)=1.0d-20 > ! b_S(1,1)=1.0d-20 > ! b_N(1,1)=1.0d-20 > ! S_v_fake(1,1)=1.0d-30 > > > > > !coeficientes de p > c_W(1,1)=dy/a_P(1,1)*dy > c_E(1,1)=dy/a_E(1,1)*alpha_u*dy > c_S(1,1)=dx/b_P(1,1)*dx > c_N(1,1)=dx/b_N(1,1)*alpha_v*dx > c_P(1,1)=c_W(1,1)+c_E(1,1)+c_S(1,1)+c_N(1,1) > S_p_fake(1,1)=u_fake(1,1)*dy-u_fake(1,2)*dy+v_fake(1,1)*dx- > v_fake(2,1)*dx > > > !fronteira NW > > !coeficientes de u > a_P(Ny,1)=(F_e_u(Ny,1)*7.0d0/8.0d0*dy+D_e*dy+D_s*dx > +D_EE*9.0d0/3.0d0*dy+D_H*9.0d0/& > 3.0d0*dx)/alpha_u > a_W(Ny,1)=1.0d-20 !0.0d0 > a_E(Ny,1)=D_e*dy+D_EE*1.0d0/3.0d0*dy > a_S(Ny,1)=F_s_u(Ny,1)*dx+D_s*dx+D_H*1.0d0/3.0d0*dx > a_N(Ny,1)=1.0d-20 !0.0d0 > S_fake(Ny,1)=(u_E*dy)*F_w_u(Ny,1)+(-3.0d0/8.0d0*u_fake(Ny,2)*dy > +2.0d0/8.0d0*u_E*dy)*& > F_e_u(Ny,1)+(3.0d0/8.0d0*u_fake(Ny, > 1)*dx-2.0d0/8.0d0*u_fake(Ny-1,1)*dx-1.0d0/8.0d0*& > u_fake(Ny-2,1)*dx)*F_s_u(Ny,1)+(-u_H*dx)*F_n_u(Ny, > 1)+D_EE*8.0d0/3.0d0*u_E*dy+D_H*& > 8.0d0/3.0d0*u_H*dx+((1.0d0-alpha_u)*a_P(Ny,1))*u_old(Ny,1) ! Falta > o termo de press�o. > > > ! a_P(Ny,1)=1.0d0 > ! a_W(Ny,1)=1.0d-20 > ! a_E(Ny,1)=1.0d-20 > ! a_S(Ny,1)=1.0d-20 > ! a_N(Ny,1)=1.0d-20 > ! S_fake(Ny,1)=1.0d-30 > > > > > > !coeficientes de v > b_P(Ny,1)=(F_e_v(Ny,1)*7.0d0/8.0d0*dy+D_e*dy+D_s*dx > +D_I*9.0d0/3.0d0*dy+D_L*9.0d0/& > 3.0d0*dx)/alpha_v > b_W(Ny,1)=1.0d-20 !0.0d0 > b_E(Ny,1)=D_e*dy+D_I*1.0d0/3.0d0*dy > b_S(Ny,1)=F_s_v(Ny,1)*dx+D_s*dx+D_L*1.0d0/3.0d0*dx > b_N(Ny,1)=1.0d-20 !0.0d0 > S_v_fake(Ny,1)=(v_I*dy)*F_w_v(Ny,1)+(-3.0d0/8.0d0*v_fake(Ny,2)*dy > +2.0d0/8.0d0*v_I*dy& > )*F_e_v(Ny,1)+(3.0d0/8.0d0*v_fake(Ny, > 1)*dx-2.0d0/8.0d0*v_fake(Ny-1,1)*dx-1.0d0/8.0d0& > *v_fake(Ny-2,1)*dx)*F_s_v(Ny,1)+(-v_L*dx)*F_n_v(Ny,1)+ > (p_fake(Ny-1,1)-p_fake(Ny,1))*& > dx+D_I*8.0d0/3.0d0*v_I*dy+D_L*8.0d0/3.0d0*v_L*dx+((1.0d0- > alpha_v)*b_P(Ny,1))*& > v_old(Ny,1) > > > ! b_P(Ny,1)=1.0d0 > ! b_W(Ny,1)=1.0d-20 > ! b_E(Ny,1)=1.0d-20 > ! b_S(Ny,1)=1.0d-20 > ! b_N(Ny,1)=1.0d-20 > ! S_v_fake(Ny,1)=1.0d-30 > > > !coeficientes de p > c_W(Ny,1)=dy/a_P(Ny,1)*dy > c_E(Ny,1)=dy/a_E(Ny,1)*alpha_u*dy > c_S(Ny,1)=dx/b_P(Ny,1)*dx > c_N(Ny,1)=dx/b_N(Ny,1)*alpha_v*dx > c_P(Ny,1)=c_W(Ny,1)+c_E(Ny,1)+c_S(Ny,1)+c_N(Ny,1) > S_p_fake(Ny,1)=u_fake(Ny,1)*dy-u_fake(Ny,2)*dy !Faltam termos de v. > > > !fronteira SE > > !coeficientes de u > a_P(1,Nx)=(F_n_u(1,Nx)*7.0d0/8.0d0*dx+D_w*dy+D_n*dx > +D_F*9.0d0/3.0d0*dy+D_G*9.0d0/& > 3.0d0*dx)/alpha_u > a_W(1,Nx)=F_w_u(1,Nx)*dy+D_w*dy+D_F*1.0d0/3.0d0*dy > a_E(1,Nx)=1.0d-20 !0.0d0 > a_S(1,Nx)=1.0d-20 !0.0d0 > a_N(1,Nx)=D_n*dx+D_G*1.0d0/3.0d0*dx > > S_fake(1,Nx)=(3.0d0/8.0d0*u_fake(1,Nx)*dy-2.0d0/8.0d0*u_fake(1,Nx-1)*dy-1 > .0d0/8.0d0*& > u_fake(1,Nx-2)*dy)*F_w_u(1,Nx)+(-u_F(1)*dy)*F_e_u(1,Nx)+ > (u_G*dx)*F_s_u(1,Nx)+(-3.0d0& > /8.0d0*u_fake(2,Nx)*dx+2.0d0/8.0d0*u_G*dx)*F_n_u(1,Nx)+ > (p_fake(1,Nx-1)-p_fake(1,Nx))& > *dy+D_F*8.0d0/3.0d0*u_F(1)*dy+D_G*8.0d0/3.0d0*u_G*dx+((1.0d0- > alpha_u)*a_P(1,Nx))*& > u_old(1,Nx) > > ! a_P(1,Nx)=1.0d0 > ! a_W(1,Nx)=1.0d-20 > ! a_E(1,Nx)=1.0d-20 > ! a_S(1,Nx)=1.0d-20 > ! a_N(1,Nx)=1.0d-20 > ! S_fake(1,Nx)=1.0d-30 > > > > > > > > > !coeficientes de v > b_P(1,Nx)=(F_n_v(1,Nx)*7.0d0/8.0d0*dx+D_w*dy+D_n*dx > +D_J*9.0d0/3.0d0*dy+D_K*9.0d0/& > 3.0d0*dx)/alpha_v > b_W(1,Nx)=F_w_v(1,Nx)*dy+D_w*dy+D_J*1.0d0/3.0d0*dy > b_E(1,Nx)=1.0d-20 !0.0d0 > b_S(1,Nx)=1.0d-20 !0.0d0 > b_N(1,Nx)=D_n*dx+D_K*1.0d0/3.0d0*dx > > S_v_fake(1,Nx)=(3.0d0/8.0d0*v_fake(1,Nx)*dy-2.0d0/8.0d0*v_fake(1,Nx-1)*dy > -1.0d0/ > & > 8.0d0*v_fake(1,Nx-2)*dy)*F_w_v(1,Nx)+(-v_J*dy)*F_e_v(1,Nx)+ > (v_K*dx)*F_s_v(1,Nx)+(& > -3.0d0/8.0d0*v_fake(2,Nx)*dx+2.0d0/8.0d0*v_K*dx)*F_n_v(1,Nx) > +D_J*8.0d0/3.0d0*v_J*dy+& > D_K*8.0d0/3.0d0*v_K*dx+((1.0d0-alpha_v)*b_P(1,Nx))*v_old(1,Nx) ! > Falta o termo de press�o. > > > ! b_P(1,Nx)=1.0d0 > ! b_W(1,Nx)=1.0d-20 > ! b_E(1,Nx)=1.0d-20 > ! b_S(1,Nx)=1.0d-20 > ! b_N(1,Nx)=1.0d-20 > ! S_v_fake(1,Nx)=1.0d-30 > > > > > !coeficientes de p > c_W(1,Nx)=dy/a_P(1,Nx)*dy > c_E(1,Nx)=dy/a_E(1,Nx)*alpha_u*dy > c_S(1,Nx)=dx/b_P(1,Nx)*dx > c_N(1,Nx)=dx/b_N(1,Nx)*alpha_v*dx > c_P(1,Nx)=c_W(1,Nx)+c_E(1,Nx)+c_S(1,Nx)+c_N(1,Nx) > S_p_fake(1,Nx)=v_fake(1,Nx)*dx-v_fake(2,Nx)*dx !Faltam termos de u. > > > !fronteira NE > > !coeficientes de u > a_P(Ny,Nx)=(D_w*dy+D_s*dx+D_F*9.0d0/3.0d0*dy+D_H*9.0d0/3.0d0*dx)/ > alpha_u > a_W(Ny,Nx)=F_w_u(Ny,Nx)*dy+D_w*dy+D_F*1.0d0/3.0d0*dy > a_E(Ny,Nx)=1.0d-20 !0.0d0 > a_S(Ny,Nx)=F_s_u(Ny,Nx)*dx+D_s*dx+D_H*1.0d0/3.0d0*dx > a_N(Ny,Nx)=1.0d-20 !0.0d0 > > S_fake(Ny,Nx)=(3.0d0/8.0d0*u_fake(Ny,Nx)*dy-2.0d0/8.0d0*u_fake(Ny,Nx-1)*d > y-1.0d0/ > & > 8.0d0*u_fake(Ny,Nx-2)*dy)*F_w_u(Ny,Nx)+(-u_F(Ny)*dy)*F_e_u(Ny,Nx) > +(3.0d0/8.0d0*& > > u_fake(Ny,Nx)*dx-2.0d0/8.0d0*u_fake(Ny-1,Nx)*dx-1.0d0/8.0d0*u_fake(Ny-2,Nx) > *dx)*& > F_s_u(Ny,Nx)+(-u_H*dx)*F_n_u(Ny,Nx)+(p_fake(Ny,Nx-1)- > p_fake(Ny,Nx))*dy+D_F*8.0d0/& > 3.0d0*u_F(Ny)*dy+D_H*8.0d0/3.0d0*u_H*dx+((1.0d0- > alpha_u)*a_P(Ny,Nx))*u_old(Ny,Nx) > > ! a_P(Ny,Nx)=1.0d0 > ! a_W(Ny,Nx)=1.0d-20 > ! a_E(Ny,Nx)=1.0d-20 > ! a_S(Ny,Nx)=1.0d-20 > ! a_N(Ny,Nx)=1.0d-20 > ! S_fake(Ny,Nx)=1.0d-30 > > > > > !coeficientes de v > b_P(Ny,Nx)=(D_w*dy+D_s*dx+D_J*9.0d0/3.0d0*dy+D_L*9.0d0/3.0d0*dx)/ > alpha_v > b_W(Ny,Nx)=F_w_v(Ny,Nx)*dy+D_w*dy+D_J*1.0d0/3.0d0*dy > b_E(Ny,Nx)=1.0d-20 !0.0d0 > b_S(Ny,Nx)=F_s_v(Ny,Nx)*dx+D_s*dx+D_L*1.0d0/3.0d0*dx > b_N(Ny,Nx)=1.0d-20 !0.0d0 > > S_v_fake(Ny,Nx)=(3.0d0/8.0d0*v_fake(Ny,Nx)*dy-2.0d0/8.0d0*v_fake(Ny,Nx-1) > *dy-1.0d0/ > & > 8.0d0*v_fake(Ny,Nx-2)*dy)*F_w_v(Ny,Nx)+(-v_J*dy)*F_e_v(Ny,Nx) > +(3.0d0/8.0d0*& > > v_fake(Ny,Nx)*dx-2.0d0/8.0d0*v_fake(Ny-1,Nx)*dx-1.0d0/8.0d0*v_fake(Ny-2,Nx) > *dx)*& > F_s_v(Ny,Nx)+(-v_L*dx)*F_n_v(Ny,Nx)+(p_fake(Ny-1,Nx)- > p_fake(Ny,Nx))*dx+D_J*8.0d0/& > 3.0d0*v_J*dy+D_L*8.0d0/3.0d0*v_L*dx+((1.0d0- > alpha_v)*b_P(Ny,Nx))*v_old(Ny,Nx) > > ! b_P(Ny,Nx)=1.0d0 > ! b_W(Ny,Nx)=1.0d-20 > ! b_E(Ny,Nx)=1.0d-20 > ! b_S(Ny,Nx)=1.0d-20 > ! b_N(Ny,Nx)=1.0d-20 > ! S_v_fake(Ny,Nx)=1.0d-30 > > > > > > !coeficientes de p > c_W(Ny,Nx)=dy/a_P(Ny,Nx)*dy > c_E(Ny,Nx)=dy/a_E(Ny,Nx)*alpha_u*dy > c_S(Ny,Nx)=dx/b_P(Ny,Nx)*dx > c_N(Ny,Nx)=dx/b_N(Ny,Nx)*alpha_v*dx > c_P(Ny,Nx)=c_W(Ny,Nx)+c_E(Ny,Nx)+c_S(Ny,Nx)+c_N(Ny,Nx) > S_p_fake(Ny,Nx)=0.0d0 !Faltam termos de u e de v. > > > !Fonte > > k=1 > > do i=1,Ny > do j=1,Nx > S(k)=S_fake(i,j) > S_v(k)=S_v_fake(i,j) > S_p(k)=S_p_fake(i,j) > k=k+1 > end do > end do > > > !Matriz dos Coeficientes > > k=1 > > do i=1,Ny > do j=1,Nx > a(k,k)=a_P(i,j) > b(k,k)=b_P(i,j) > c(k,k)=c_P(i,j) > k=k+1 > end do > end do > > > k=2 > > do i=1,Ny > if (i.eq.1) then > do j=2,Nx > a(k,k-1)=-a_W(i,j) > b(k,k-1)=-b_W(i,j) > c(k,k-1)=-c_W(i,j) > k=k+1 > end do > else > do j=1,Nx > a(k,k-1)=-a_W(i,j) > b(k,k-1)=-b_W(i,j) > c(k,k-1)=-c_W(i,j) > k=k+1 > end do > end if > end do > > k=2 > > do i=1,Ny > if (i.ne.Ny) then > do j=1,Nx > a(k-1,k)=-a_E(i,j) > b(k-1,k)=-b_E(i,j) > c(k-1,k)=-c_E(i,j) > k=k+1 > end do > else > do j=1,Nx-1 > a(k-1,k)=-a_E(i,j) > b(k-1,k)=-b_E(i,j) > c(k-1,k)=-c_E(i,j) > k=k+1 > end do > end if > end do > > k=Nx+1 > > do i=2,Ny > do j=1,Nx > a(k,k-Nx)=-a_S(i,j) > b(k,k-Nx)=-b_S(i,j) > c(k,k-Nx)=-c_S(i,j) > k=k+1 > end do > end do > > k=1 > > do i=1,Ny-1 > do j=1,Nx > a(k,k+Nx)=-a_N(i,j) > b(k,k+Nx)=-b_N(i,j) > c(k,k+Nx)=-c_N(i,j) > k=k+1 > end do > end do > > > > k=0 > l=0 > m=0 > > do i=1,Ny*Nx > do j=1,Ny*Nx > if (dabs(a(i,j)).ne.0.0d0) then > k=k+1 > end if > if (dabs(b(i,j)).ne.0.0d0) then > l=l+1 > end if > if (dabs(c(i,j)).ne.0.0d0) then > m=m+1 > end if > end do > end do > > allocate(AA(k),irow_u(k),jcol_u(k)) > allocate(BB(l),irow_v(l),jcol_v(l)) > allocate(CC(m),irow_p(m),jcol_p(m)) > > k=1 > l=1 > m=1 > > do i=1,Ny*Nx > do j=1,Ny*Nx > if (dabs(a(i,j)).ne.0.0d0) then > AA(k)=a(i,j) > irow_u(k)=i > jcol_u(k)=j > k=k+1 > end if > if (dabs(b(i,j)).ne.0.0d0) then > BB(l)=b(i,j) > irow_v(l)=i > jcol_v(l)=j > l=l+1 > end if > if (dabs(c(i,j)).ne.0.0d0) then > CC(m)=c(i,j) > irow_p(m)=i > jcol_p(m)=j > m=m+1 > end if > end do > end do > > call DL4LXG(iparam,rparam) > > iparam(5)=1000000 > > > !call DLINRG(Ny*Nx,a,Ny*Nx,a_inv,Ny*Nx) > !CALL DLSARG (Ny*Nx,a,Ny*Nx,S,ipath,u) > call DLSLXG(Ny*Nx,k-1,AA,irow_u,jcol_u,S,ipath,iparam,rparam,u) > call DLSLXG(Ny*Nx,l-1,BB,irow_v,jcol_v,S_v,ipath,iparam,rparam,v) > call DLSLXG(Ny*Nx,m-1,CC,irow_p,jcol_p,S_p,ipath,iparam,rparam,p) > > ! k=1 > > > ! do i=1,Ny > ! do j=1,Nx > ! p_fake(i,j)=p(k) > ! if (j.ne.1) then > ! u_fake(i,j)=u(k)+dy/(a_P(i,j)*alpha_u)*(p_fake(i,j-1)- > p_fake(i,j)) > ! end if > ! if (i.ne.1) then > ! v_fake(i,j)=v(k)+dx/(b_P(i,j)*alpha_v)*(p_fake(i-1,j)- > p_fake(i,j)) > ! end if > ! p_fake(i,j)=p_old(i,j)+alpha_p*p(k) > ! k=k+1 > ! end do > ! end do > > > p_old=p_fake > u_old=u_fake > v_old=v_fake > > > k=1 > > do i=1,Ny > do j=1,Nx > p_fake(i,j)=p(k) > !p_fake(i,j)=p_old(i,j)+alpha_p*p(k) > !p_fake(i,j)=alpha_p*(k)+(1.0d0-alpha_p)*p_old(i,j) > !if (i.ne.1) then > ! if (i.ne.Ny) then > ! if (j.ne.1) then > ! if (j.ne.Nx) then > !u_fake(i,j)=u(k)+dy/(a_P(i,j)*alpha_u)*(p_fake(i,j-1)- > p_fake(i,j)) > !v_fake(i,j)=v(k)+dx/(b_P(i,j)*alpha_v)*(p_fake(i-1,j)- > p_fake(i,j)) > ! u_fake(i,j)=u(k)+dy/a_P(i,j)*(p_fake(i,j-1)-p_fake(i,j)) > ! v_fake(i,j)=v(k)+dx/b_P(i,j)*(p_fake(i-1,j)-p_fake(i,j)) > ! u_fake(i,j)=alpha_u*u(k)+(1.0d0-alpha_u)*u_old(i,j) > ! v_fake(i,j)=alpha_v*v(k)+(1.0d0-alpha_v)*v_old(i,j) > ! end if > ! end if > ! end if > !end if > if (j.ne.1) then > u_fake(i,j)=u(k)+dy/(a_P(i,j)*alpha_u)*(p_fake(i,j-1)-p_fake(i,j)) > !u_fake(i,j)=u(k)+dy/a_P(i,j)*(p_fake(i,j-1)-p_fake(i,j)) > !u_fake(i,j)=alpha_u*u(k) > !u_fake(i,j)=alpha_u*u(k)+(1.0d0-alpha_u)*u_old(i,j) > ! write(*,*) u_fake(i,j) > end if > if (i.ne.1) then > v_fake(i,j)=v(k)+dx/(b_P(i,j)*alpha_v)*(p_fake(i-1,j)-p_fake(i,j)) > !v_fake(i,j)=v(k)+dx/b_P(i,j)*(p_fake(i-1,j)-p_fake(i,j)) > !v_fake(i,j)=alpha_v*v(k) > !v_fake(i,j)=alpha_v*v(k)+(1.0d0-alpha_v)*v_old(i,j) > end if > p_fake(i,j)=p_old(i,j)+alpha_p*p(k) > k=k+1 > end do > end do > > > > > deallocate(AA,irow_u,jcol_u,BB,irow_v,jcol_v,CC,irow_p,jcol_p) > > !u=matmul(a_inv,S) > !theta=matmul(a_inv,S) > > > > ! end do > > open(1,file='u.dat') > write(1,11) Ny,Nx > 11 format('variables=X,Y,u',/,'zone t="time','" i=',i4,' j=',i4,' > f=point') > > k=1 > > do i=1,Ny > do j=1,Nx > write(1,*) Y(i),X(j),u(k) > k=k+1 > end do > end do > > open(2,file='v.dat') > write(2,12) Ny,Nx > 12 format('variables=X,Y,v',/,'zone t="time','" i=',i4,' j=',i4,' > f=point') > > k=1 > > do i=1,Ny > do j=1,Nx > write(2,*) Y(i),X(j),v(k) > k=k+1 > end do > end do > > open(3,file='p.dat') > write(3,13) Ny,Nx > 13 format('variables=X,Y,p',/,'zone t="time','" i=',i4,' j=',i4,' > f=point') > > k=1 > > do i=1,Ny > do j=1,Nx > write(3,*) Y(i),X(j),p(k) > k=k+1 > end do > end do > > open(4,file='S.dat') > write(4,14) Ny,Nx > 14 format('variables=X,Y,p',/,'zone t="time','" i=',i4,' j=',i4,' > f=point') > > k=1 > > do i=1,Ny > do j=1,Nx > write(4,*) S_fake(i,j),k,S(k) > k=k+1 > end do > end do > > open(5,file='a.dat') > !write(5,15) Ny,Nx > !15 format('variables=X,Y,p',/,'zone t="time','" i=',i4,' j=',i4,' > f=point') > > > do i=1,Ny*Nx > do j=1,Ny*Nx > write(5,*) i,j,a(i,j) > end do > end do > > open(6,file='aa.dat') > !write(6,16) Ny,Nx > !16 format('variables=X,Y,p',/,'zone t="time','" i=',i4,' j=',i4,' > f=point') > > > do i=1,Ny > do j=1,Nx > write(6,*) i,j,a_P(i,j) > end do > end do > > > > end program EXERCICIO
From: mecej4 on 5 Jul 2010 15:59 Felipe Abrão wrote: > Please, could somebody explain me what is happening. I simply > allocated a matrix called "a", but figured out that the program is > inserting some strange values to some of its elements (like > 1.057229231618489E-307 or 6.941570408926931E+228), even before I set > any value to it! Shouldn't all the values be zero before any desired > value setting? Below is my code, but you don't have to read it all. > You can check the problem just after line 65, where I wrote a command > in order to read the values of matrix "a". The same thing is happening > with matrixs "b" and "c". In advance I thank you all! > > PS: I use Fortran 90 in Windows XP. <-- long program listing --> It is your job to initialize variables before using them. Some compilers have a switch to initialize variables by default to zero. This, however, does not work with allocatable variables -- you cannot initialize variables that have not yet been allocated. -- mecej4
From: e p chandler on 5 Jul 2010 16:01 "Felipe Abr�o" <felipeabrao(a)gmail.com> wrote in message news:bf5a29d3-df12-47e5-82bf-951daffefbea(a)u7g2000yqm.googlegroups.com... Please, could somebody explain me what is happening. I simply allocated a matrix called "a", but figured out that the program is inserting some strange values to some of its elements (like 1.057229231618489E-307 or 6.941570408926931E+228), even before I set any value to it! Shouldn't all the values be zero before any desired value setting? Below is my code, but you don't have to read it all. You can check the problem just after line 65, where I wrote a command in order to read the values of matrix "a". The same thing is happening with matrixs "b" and "c". In advance I thank you all! PS: I use Fortran 90 in Windows XP. [really long code snipped] ------> That is what is supposed to happen. If you don't assign a value, then you should get whatever happens to be in memory. These small and large numbers are typical "garbage" values. If you want zero, assign zero to it. The good news is that you can assign a value to an entire array with one statement. [rant mode on] While there are some compilers that will zero memory first and there are some compiler options that may zero memory first, you should not count on that taking place. Un-initialized variables are a programming error. Consider the case of a recent poster in this newsgroup. The poor soul has inherited a code base of approximately 600,000 lines of code. His code requires that un-initialized variables be set to zero (along with something else not quite relevant to this discussion.). His compiler has an option to do that. But he wants to port his code to a different compiler and possibly different platforms. He has a monumental task which I am extremely happy is NOT mine! [rant mode off] --- Elliot
From: Louis Krupp on 5 Jul 2010 16:04
On 7/5/2010 1:37 PM, Felipe Abr�o wrote: > Please, could somebody explain me what is happening. I simply > allocated a matrix called "a", but figured out that the program is > inserting some strange values to some of its elements (like > 1.057229231618489E-307 or 6.941570408926931E+228), even before I set > any value to it! Shouldn't all the values be zero before any desired > value setting? Below is my code, but you don't have to read it all. You don't have to post it all. Write a small program to demonstrate the problem you're having, and post it. People will thank you for doing that. Louis |