대기오염농도를 구하는 포트란 프로그램입니다...ㅜㅠ도와주세요

급합니다. 제발...고수님들 도와주세요.. 열심히 포트란을 하려하지만,,,역시나 프로그램쪽으로는 타고난 머리가 필요한듯..ㅋ 가우시안 식을 이용하여 굴뚝에서 나온 연기의 농도를 계산하는 프로그램입니다. 바람의 세기의 차이에 따라 농도(c)의 값의 변화를 알려고 하는데요.. 배열이 들어가니깐 도무지...어려워서.. 다행인지 불행인지 문법에러는 안뜹니다. 그런데 문제는 결과가 아예나오지 않습니다. 실행파일을 시키면 커서만 깜박일뿐 runtime erro도 나지 않는데.. 제 생각에는 무한 do문이 돌아가는 것 같기도 하고...에러가 없다고 나오니깐 어느부분이 틀렸는지 도저히 감이 안옵니다. 아래의 것이 제가 짠 프로그램입니다. program change_windspeed implicit none integer::i integer,dimension(0:15)::u real::x real,parameter :: d=10., v=10., stack=10. , pi=3.14, a=61.141, b=0.91465, Q=500. real,dimension(0:15) :: Hs, h real,dimension(100)::yth,zth, T real,dimension(16,100)::c !Hs : 굴뚝의 유효 높이 h : 연기상승높이, v : 배출속도(m/s), D: 연돌내경(m), U :연돌주위 풍속 stack: 굴뚝의 높이 Q:배출량(kg) open(11,file="c:/windspeed.txt")! 바람의 세기를 0-16m/s가 들어있는 파일 굳이 필요하지 않을 것 같기도 합니다.<-있어야 하는지 없어도 되는지...도 갈쳐주세요.. open(12,file="c:/result.txt") !********************유효굴뚝의 높이구하기******************** do 3 read(11,*,end=3) u 3 continue do i=0, 15 !아래의 식들은 오염농도를 구하기 위해 필요한 식입니다... !배열에 문제가 없는듯 한데... h(i)=1.5*v*d/u(i) Hs(i)=stack+h(i) do x=1,100 T(x)=12.5-1.0857*2.7*log(x) ! pasquill에 의한 안정도 c를 선택하였다.저번 숙제에 부산의 안정도는 보통 c로 나타났기 때문에... yth(x)=1000*x*tan(T(x))/2.15 zth(x)=a*x**b c(i,x)=(Q/(2*pi*u(i)*yth(x)*zth(x)))*(1+exp(-2*Hs(i)**2/2*zth(x)**2)) write(12,*)u(i), Hs(i), c(i,x)!<-write 문을 쓰는데 결과 파일에는 아무것도 적히지 않습니다....ㅜㅠ enddo enddo end program change_windspeed