Fortran90 wrong output

Question

나는 대학에서 수업을 받기 위해 작은 프로그램을 공부하고 있는데, 거의 끝났지 만 어떻게 든 작동하기를 원하는 것처럼 작동하지 않는다.

이제 출력 파일 gravity1.dat는 나에게 0이 아닌 값을 주어야합니다.하지만 그렇지 않습니다 ... 마지막 공식의 어딘가에서 g (서핑)을 계산하면 변수 중 하나가 0입니다. 거의 모든 것을 시도한 경우 그것을 수정하는 나의 힘에서 그러나 나는 그것을 고치는 것처럼 보이지 않는다.

!=========================================== 

subroutine write_two_columns (nnn,xxx,yyy,filename) 
implicit none 
integer i,nnn 
real(8) xxx(nnn),yyy(nnn) 
character(LEN=*) filename 

open(unit=123,file=filename,action='write') 
do i=1,nnn 
write(123,*) xxx(i),yyy(i) 
end do 
close(123) 
end subroutine 

다른 서브 루틴 :

program gravity 

implicit none 
real(8) Lx,Ly,sx,sy,xsphere,ysphere,r,A,rho1,rho2,dx,G 
integer np,nel,nelx,nely,i,nnx,nny,j,counter,nsurf 
real(8),dimension(:),allocatable :: xcgrid 
real(8),dimension(:),allocatable :: ycgrid 
real(8),dimension(:),allocatable :: xgrid 
real(8),dimension(:),allocatable :: ygrid 
real(8),dimension(:),allocatable :: rho 
real(8),dimension(:),allocatable :: xsurf, gsurf 
real(8),dimension(:),allocatable :: ysurf 

nnx=11. 
nny=11. 
Lx=10. 
Ly=10. 
nelx=nnx-1. 
nely=nny-1. 
nel=nelx*nely 
np=nnx*nny 
sx=Lx/nelx 
sy=Ly/nely 
xsphere=5. 
ysphere=5. 
r=3. 
nsurf=7 !number of gravimeters 

dx=Lx/(nsurf-1.) 

!========================================================== 

allocate(xgrid(np)) 
allocate(ygrid(np)) 

counter=0 
do i=1,nnx 
    do j=1,nny 
    counter=counter+1 
    xgrid(counter)=dble(i-1)*sx 
    ygrid(counter)=dble(j-1)*sy 
    end do 
end do 

call write_two_columns(np,xgrid,ygrid,'grid_init1.dat') 
!========================================================== 

allocate(xcgrid(np)) 
allocate(ycgrid(np)) 


counter=0 
do i=1,nnx-1 
    do j=1,nny-1 
    counter=counter+1 
    xcgrid(counter)=dble(i-1)*sx+0.5*sx 
    ycgrid(counter)=dble(j-1)*sy+0.5*sy 
    end do 
end do 

call write_two_columns(np,xcgrid,ycgrid,'gridc_init1.dat') 
!========================================================== 

allocate(rho(nel)) 

rho1=3000. !kg/m^3 
rho2=3200. !kg/m^3 

do i=1,nel 
    if (sqrt((xsphere-xcgrid(i))**2)+((ysphere-ycgrid(i))**2)<r) then 
    rho(i)=3200. 
    else 
    rho(i)=3000. 
    end if 
end do 

call write_three_columns(nel,xcgrid,ycgrid,rho,'inclusion1.dat') 
!========================================================== 

allocate(xsurf(nsurf)) 
allocate(ysurf(nsurf)) 

do i=1,nsurf 
xsurf(i)=(i-1)*dx 
ysurf(i)=ly 
end do 

call write_two_columns(nsurf,xsurf,ysurf,'surf_init1.dat') 
!========================================================== 

allocate(gsurf(nsurf)) 

G=0.000000000066738480 !m^3 kg^-1 s^-2 

do i=1,nsurf 
    do j=1,nel 
    gsurf(i)=gsurf(i)+(-2.*G*(((rho(i)-rho1)*(ycgrid(counter)-ysurf(i)))/((xcgrid(counter)-xsurf(i))**2.+(ycgrid(counter)-ysurf(i))**2.))*sx*sy) 

    end do 
end do 

call write_two_columns(nsurf,xsurf,gsurf,'gravity1.dat') 

deallocate(xgrid) 
deallocate(ygrid) 
deallocate(xcgrid) 
deallocate(ycgrid) 
deallocate(xsurf) 
deallocate(ysurf) 
deallocate(gsurf) 
end program 

서브 루틴을 사용

!=========================================== 

subroutine write_three_columns (nnn,xxx,yyy,zzz,filename) 
implicit none 
integer i,nnn 
real(8) xxx(nnn),yyy(nnn),zzz(nnn) 
character(LEN=*) filename 
open(unit=123,file=filename,action='write') 
do i=1,nnn 
write(123,*) xxx(i),yyy(i),zzz(i) 
end do 
close(123) 
end subroutine 

!=========================================== 

Answer 1

이 (rho(j)-rho1)을해야하지? 만 사용하면 rho(1:nel)을 채울 수 있습니다.

덧붙여서 변수 초기화에주의하십시오. real을 integer에 할당하고 혼합형 연산을 수행하십시오. 예기치 않은 결과가 발생할 수 있으므로주의하십시오. 컴파일러를 사용하여 이러한 문제점을 감지하십시오.