fortran - Division by zero in geany -
i trying run sub routine in geany, bur keeps me giving following warning
npj(i,j) = dble(((2)/((x(i+1))-(x(i-1))))*dble(-((1)/(x(i+1)-x(i)))-((1)/(x(i)- x(i-1)))))+ & 1 warning: possible change of value in conversion real(8) integer(4) @ (1) the code follows by
c(i,j) = -(re(i,j))/npj(i,j) on next line.
everytime run program gives getting divisions zero.
the code here:
! projeto 1 - mÉtodos experimentais em hidrodinÂmica !****************************************** ! * ! programa principal * ! * !****************************************** implicit double precision (a-h,o-z) parameter (nx = 100, ny = 100) ! declaraÇÃo de variÁveis common/grid/x(nx),y(ny) common/resi/re(nx,ny) common/ptn/ pot(nx,ny) character*30 bobo ! definiÇÃo dos arquivos de entrada e saÍda open(3,file = 'placa.dat') open(4,file = 'output.dat') ! entrada de dados read(3,*) bobo,imax write(*,'(a30,i10)')bobo,imax read(3,*) bobo,jmax write(*,'(a30,i10)')bobo,jmax read(3,*) bobo,t write(*,'(a30,f10.3)')bobo,t read(3,*) bobo,nmax write(*,'(a30,i10)')bobo,nmax read(3,*) bobo,ua write(*,'(a30,d10.3)')bobo,ua read(3,*) bobo,ub write(*,'(a30,d10.3)')bobo,ub read(3,*) bobo,uc write(*,'(a30,d10.3)')bobo,uc read(3,*) bobo,ud write(*,'(a30,d10.3)')bobo,ud read(3,*) bobo,up write(*,'(a30,d10.3)')bobo,up read(3,*) bobo,prec write(*,'(a30,d10.3)')bobo,prec read(3,*) bobo,npr write(*,'(a30,i10)')bobo,npr read(3,*) bobo,ite write(*,'(a30,i10)')bobo,ite read(3,*) bobo,ile write(*,'(a30,i10)')bobo,ile read(3,*) bobo,xsf write(*,'(a30,f10.3)')bobo,xsf read(3,*) bobo,ysf write(*,'(a30,f10.3)')bobo,ysf !$$$$$$ !$$$$$$ write(*,*)"os dados de entrada estao corretos?" !$$$$$$ write(*,*)"1--------sim" !$$$$$$ write(*,*)"2--------nao" !$$$$$$ read(*,*)inf !$$$$$$ if(inf.eq.2) stop ! geraÇÃo da malha computacional call malha(imax,jmax,dx,ite,ile,xsf,ysf,dy) ! ! condiÇÃo inicial call inicial(imax,jmax,up) !!! inÍcio das iteraÇÕes call solver(imax,jmax,nmax,prec,n,npr,dy,ua,ub,ud,ile,ite,t) !! fim da execuÇÃo stop end !****************************************** ! * ! subrotina malha * ! * !****************************************** subroutine malha(imax,jmax,dx,ite,ile,xsf,ysf,dy) implicit double precision (a-h,o-z) parameter (nx=100,ny=100) common/grid/ x(nx),y(ny) dx=1.0d0/dble(ite-ile) i=ile,ite x(i)=dx*dble(i-ile) end i=ite,imax x(i)=x(i-1)+((x(i-1)-x(i-2))*xsf) end i=ile-1,1,-1 x(i)=x(i+1)+((x(i+1)-x(i+2))*xsf) end y(1) = (-dx)/2.0d0 y(2) = dx/2.0d0 dy = y(2)-y(1) j=3, jmax y(j)=y(j-1)+((y(j-1)-y(j-2))*ysf) end return end !****************************************** ! * ! subrotina inicial * ! * !****************************************** subroutine inicial(imax,jmax,up) implicit double precision (a-h,o-z) parameter (nx=100,ny=100) common/ptn/ pot(nx,ny) common/grid/ x(nx),y(ny) = 1.0d0 j=1,jmax = 1,imax pot(i,j)=up*x(i) end end return end !****************************************** ! * ! subrotina contorno * ! * !****************************************** subroutine contorno(imax,jmax,ua,ub,ud,ile,ite,dy,t) implicit double precision (a-h,o-z) parameter (nx=100, ny=100) common/grid/ x(nx),y(ny) common/ptn/ pot(nx,ny) ! entrada e saída j=1,jmax pot(1,j)=ua*x(1) pot(imax,j)=ub*x(imax) end ! fronteira superior i=1, imax pot(i,jmax)=ud*x(i) end !simetria i=1, imax pot(i,1) = pot(i,2) end ! sobre o corpo i=ile,ite pot(i,1) = pot(i,2) - 2.0d0*dy*ua*t*(1-(2.0d0*x(i))) end return end !****************************************** ! * ! subrotina residuo * ! * !****************************************** subroutine residuo(imax,jmax,teste) implicit double precision (a-h,o-z) parameter (nx=100, ny=100) common/resi/ re(nx,ny) common/grid/ x(nx),y(ny) common/ptn/ pot(nx,ny) ! calculo residuo j=2,jmax-1 i=2,imax-1 re(i,j) = ((2.0d0/((x(i+1))-(x(i-1))))*(((pot(i+1,j)-pot(i,j))/(x(i+1)-x(i)))-((pot(i,j)-pot(i-1,j))/(x(i)-x(i-1)))))+ & &((2.0d0/((y(j+1))-(y(j-1))))*(((pot(i,j+1)-pot(i,j))/(y(j+1)-y(j)))-((pot(i,j)-pot(i,j-1))/(y(j)-y(j-1))))) end end ! calculo residuo maximo teste=0.0d0 j=2,jmax-1 i=2,imax-1 if(dabs(re(i,j)).gt.teste) teste=dabs(re(i,j)) end if end end return end !****************************************** ! * ! subrotina solver * ! * !****************************************** subroutine solver(imax,jmax,nmax,prec,n,npr,dy,ua,ub,ud,ile,ite,t) implicit double precision (a-h,o-z) parameter (nx=100,ny=100) common/grid/ x(nx),y(ny) common/resi/ re(nx,ny) common/npjaco/ npj(nx,ny) common/correc/ c(nx,ny) common/ptn/ pot(nx,ny) teste=100.0d0 n=1 imp=npr-1 open(10,file='remax.dat') ! inicio das iteraÇÕes call contorno(imax,jmax,ua,ub,ud,ile,ite,dy,t) while((n.ne.(nmax+1)).and.(teste.gt.prec)) call residuo(imax,jmax,teste) write(*,*) n,teste write(10,*) n,teste j=2,jmax-1 i=2,imax-1 npj(i,j) = (((2.0d0)/((x(i+1))-(x(i-1))))*(-((1.0d0)/(x(i+1)-x(i)))-((1.0d0)/(x(i)-x(i-1)))))+ & & (((2.0d0)/((y(j+1))-(y(j-1))))*(-((1.0d0)/(y(j+1)-y(j)))-((1.0d0)/(y(j)-y(j-1))))) c(i,j) = -(re(i,j))/npj(i,j) pot(i,j) = pot(i,j)+c(i,j) end end n=n+1 imp=imp+1 !========================================================================================================== ! saÍda de resultados if(imp.eq.npr) write(4,10) imax,jmax 10 format('title = " malha cartesiana "',/,& & 'variables = x, y, pot, npj, re',/,& & 'zone t ="zone-one", i=',i5,'j=',i5,',f=point') j=1,jmax i=1,imax write(4,*) 'x',i,':', x(i),y(j), pot(i,j), npj(i,j), re(i,j) end end imp=0 end if ! end ! fim das iteraÇÕes return end
you not declaring variables in code, counting on implicit type. hence npj integer array. bad habit. declare types, , put implicit none in each program unit. may require more coding, it's worth it.
in assignment npj(i,j), if right-hand side small, truncated 0 [1]. since line followed c(i,j) = -(re(i,j))/npj(i,j), division zero.
in case, npj should declared double precision. didn't investigate it's supposed anyway.
[1] more precisely, double precision value (call a) converted if there int(a,k) in code, k integer kind of npj (should default integer here). see section 7.2.1.3 #8 of fortran 2008 standard. find in section 13.7.81 #5 int(a) 0 when |a|<1.
Comments
Post a Comment