Help with Debugging Fortran Code | Plato Compiler

  • Context: Fortran 
  • Thread starter Thread starter 1994Bhaskar
  • Start date Start date
  • Tags Tags
    Code Fortran
Click For Summary

Discussion Overview

The discussion revolves around debugging a Fortran code related to a mechanical engineering application. Participants are addressing a run-time error encountered when executing the code, specifically focusing on the Courant subroutine and its interactions with other parts of the program.

Discussion Character

  • Technical explanation
  • Debugging
  • Homework-related

Main Points Raised

  • One participant reports encountering a run-time error indicating a reference to an undefined variable or array element when executing the Fortran code.
  • The same participant seeks assistance in understanding whether the declaration of the Courant subroutine is correct and if the variables inside it were properly declared.
  • Another participant notes that they cannot access the provided link to the code, suggesting that the code should be posted directly in the thread.
  • A third participant echoes the inability to access the link and requests the code to be shared directly, indicating a potential barrier to effective troubleshooting.
  • The original poster provides a substantial portion of the Fortran code, including the main program and several subroutines, to facilitate debugging efforts.

Areas of Agreement / Disagreement

There is no consensus on the cause of the run-time error, and multiple views exist regarding the accessibility of the code. The discussion remains unresolved as participants have not yet identified a solution to the debugging issue.

Contextual Notes

The discussion highlights potential limitations in understanding the Fortran programming structure, as the original poster is new to the language. There may also be assumptions about variable declarations and subroutine calls that have not been fully clarified.

1994Bhaskar
Messages
134
Reaction score
0
Hey everybody, I am a mechanical engineer and new to Fortran(I have decent experience in programming in C, Matlab and understand languages similar to them). However the Fortran programming structure is new to me. Hence I am facing problems debugging a Fortran code.
Here is the link of the code which I want to run.

http://www2.meteo.uni-bonn.de/forsch...rtran/coad1d.f

I am using a Plato Fortran Compiler.

It runs without showing error, but when I enter the values initially asked, the program stops mid-way showing an error:
Run-time Error
*** Error 112, Reference to undefined variable, array element or function result (/UNDEF)

COURANT - in file help.for at line 149 [+01de]

main - in file help.for at line 45 [+05dc]

So there is an error in the Courant subroutine or error in calling it. Again i am a little unfamiliar to fortran subroutines but is it possible for any of you to also check what the problem here is ?
As far as i could see the variables inside the subroutine were declared beforehand. And the logic looks correct.
Also is the declaration line 45 calling courant subroutine done correctly ?

Any help is appreciated. I have already received awesome help from you all here.
 
Last edited by a moderator:
Technology news on Phys.org
I can't get the link to open (HTTP 404 Not Found error). Please post the code directly in this page.
 
Mark44 said:
I can't get the link to open (HTTP 404 Not Found error). Please post the code directly in this page.

Mark44, here it is:(it is huge)

program coad1d
parameter (n=400)
implicit double precision (a-h,o-z)
common /const/ rq0,dlnr,scal,ax
common /cour/ c(n,n),ima(n,n)
common /grid/ g(n),r(n),e(n)
common /kern/ ck(n,n),ec(n,n)
dimension rri(n),eei(n)
data emin,pi,tmax/1.e-9,3.141592654,3600./
c dt : time step (input in sec)
c g : spectral mass distribution (mg/cm**3)
c e : droplet mass grid (mg)
c r : droplet radius grid (um)
c dlnr: constant grid distance of logarithmic grid
c rq0 : mode radius of initial distribution (input in um)
c xmw : total water content (input in g/m**3)
c xn0 : mean initial droplet mass (mg)
c xn1 : total initial droplet number concentration (1/cm^3)
c ax : growth factor for consecutive masses
c scal: scaling factor for calculation of ax, see below
c isw : collision kernel: 0 (long), 1 (hall), 2 (golovin)
c read input variables
print *,'enter dt,rq0,xmw,scal,isw'
read *,dt,rq0,xmw,scal,isw
rq0=rq0*1.e-04
xmw=xmw*1.d-3
xn0=4./3.*pi*1000.*exp(log(rq0)*3.)
xn1=xmw/xn0
dlnr=dlog(2.d0)/(3.*scal)
ax=2.d0**(1.0/scal)
c mass and radius grid
e(1)=emin*0.5*(ax+1.)
r(1)=1000.*dexp(dlog(3.*e(1)/(4.*pi))/3.)
do i=2,n
e(i)=ax*e(i-1)
r(i)=1000.*dexp(dlog(3.*e(i)/(4.*pi))/3.)
enddo
c initial mass distribution
x0=xn1/xn0
do i=1,n
x1=e(i)
g(i)=3.*x1*x1*x0*dexp(-x1/xn0)
enddo
c courant numbers
call courant
c kernel
call trkern (isw)
c output for plotting
do i=1,n
rri(i)=min(1.d38,r(i))
eei(i)=min(1.d38,e(i))
enddo
open (16,file='boplot00.out',status='old')
write (16,6100) rri,eei
6100 format (5e16.8)
c nt: number of iterations
nt=int(tmax/dt)
c multiply kernel with constant timestep and logarithmic grid distance
do i=1,n
do j=1,n
ck(i,j)=ck(i,j)*dt*dlnr
enddo
enddo
c time integration
tlmin=1.d-6
t=1.d-6
lmin=0
do ij=1,nt
t=t+dt
tlmin=tlmin+dt
c collision
call coad
c output for plotting
if (tlmin.ge.60.) then
tlmin=tlmin-60.
lmin=lmin+1
print *,'time in minutes:',lmin
write (16,6100) t,g
c mass balance
x0=0.
x1=0.
do i=1,n
x0=x0+g(i)*dlnr
x1=max(x1,g(i))
if (dabs(x1-g(i)).lt.1.d-9) imax=i
enddo
print *,'mass ',x0,'max ',x1,'imax ',imax
endif
enddo
close (16)
stop 'stop coad1d'
end

subroutine coad
c collision subroutine, exponential approach
parameter (n=400)
implicit double precision (a-h,o-z)
common /cour/ c(n,n),ima(n,n)
common /grid/ g(n),r(n),e(n)
common /kern/ ck(n,n),ec(n,n)
data gmin /1.d-60/
c lower and upper integration limit i0,i1
do i=1,n-1
i0=i
if (g(i).gt.gmin) go to 2000
enddo
2000 continue
do i=n-1,1,-1
i1=i
if (g(i).gt.gmin) go to 2010
enddo
2010 continue
do i=i0,i1
do j=i,i1
k=ima(i,j)
kp=k+1
x0=ck(i,j)*g(i)*g(j)
x0=min(x0,g(i)*e(j))
if (j.ne.k) x0=min(x0,g(j)*e(i))
gsi=x0/e(j)
gsj=x0/e(i)
gsk=gsi+gsj
g(i)=g(i)-gsi
g(j)=g(j)-gsj
gk=g(k)+gsk
if (gk.gt.gmin) then
x1=dlog(g(kp)/gk+1.d-60)
flux=gsk/x1*(dexp(0.5*x1)-dexp(x1*(0.5-c(i,j))))
flux=min(flux,gk)
g(k)=gk-flux
g(kp)=g(kp)+flux
endif
enddo
enddo
return
end

subroutine courant
parameter (n=400)
implicit double precision (a-h,o-z)
common /const/ rq0,dlnr,scal,ax
common /cour/ c(n,n),ima(n,n)
common /grid/ g(n),r(n),e(n)
do i=1,n
do j=i,n
x0=e(i)+e(j)
do k=j,n
if (e(k).ge.x0.and.e(k-1).lt.x0) then
if (c(i,j).lt.1.-1.d-08) then
kk=k-1
c(i,j)=dlog(x0/e(k-1))/(3.d0*dlnr)
else
c(i,j)=0.
kk=k
endif
ima(i,j)=min(n-1,kk)
go to 2000
endif
enddo
2000 continue
c(j,i)=c(i,j)
ima(j,i)=ima(i,j)
enddo
enddo
return
end

subroutine trkern (isw)
parameter (n=400)
implicit double precision (a-h,o-z)
common /grid/ g(n),r(n),e(n)
common /kern/ ck(n,n),ec(n,n)
common /veloc/ winf(n),rr(n)
dimension cck(n,n)
data pi/3.141592654/
c terminal velocity
call fallg
if (isw.eq.0) then
c long kernel
do j=1,n
do i=1,j
if(r(j).le.50.) then
effi=4.5d-4*r(j)*r(j)*
& (1.d0-3.d0/(max(3.d0,dble(r(i)))+1.d-2))
else
effi=1.d0
endif
cck(j,i)=pi*(rr(j)+rr(i))*(rr(j)+rr(i))*effi*
& abs(winf(j)-winf(i))
cck(i,j)=cck(j,i)
enddo
enddo
elseif (isw.eq.1) then
c hall kernel
call effic
do j=1,n
do i=1,j
cck(j,i)=pi*(rr(j)+rr(i))*(rr(j)+rr(i))*ec(j,i)*
& abs(winf(j)-winf(i))
cck(i,j)=cck(j,i)
enddo
enddo
else
c golovin kernel
do j=1,n
do i=1,j
cck(j,i)=1.5*(e(j)+e(i))
cck(i,j)=cck(j,i)
enddo
enddo
endif
c two-dimensional linear interpolation of kernel
do i=1,n
do j=1,n
jm=max0(j-1,1)
im=max0(i-1,1)
jp=min0(j+1,n)
ip=min0(i+1,n)
ck(i,j)=0.125*(cck(i,jm)+cck(im,j)+cck(ip,j)+cck(i,jp))
& +.5*cck(i,j)
if (i.eq.j) ck(i,j)=0.5*ck(i,j)
enddo
enddo
return
end

subroutine fallg
c terminal velocity of falling drops
parameter (n=400)
implicit double precision (a-h,o-z)
common /grid/ g(n),r(n),e(n)
common /kern/ ck(n,n),ec(n,n)
common /veloc/ winf(n),rr(n)
dimension b(7),c(6),rat(20),r0(15),ecoll(15,20)
data b /-0.318657e1,0.992696,-0.153193e-2,-0.987059e-3,
& -0.578878e-3,0.855176e-4,-0.327815e-5/
data c /-0.500015e1,0.523778e1,-0.204914e1,0.475294,-0.542819e-1,
& 0.238449e-2/
data pi /3.141592654/
eta=1.818e-4
xlamb=6.62e-6
rhow=1.
rhoa=1.225e-3
grav=980.665
cunh=1.257*xlamb
t0=273.15
sigma=76.1-0.155*(293.15-t0)
stok=2.*grav*(rhow-rhoa)/(9.*eta)
stb=32.*rhoa*(rhow-rhoa)*grav/(3.*eta*eta)
phy=sigma*sigma*sigma*rhoa*rhoa/((eta**4)*grav*(rhow-rhoa))
py=phy**(1./6.)
c rr: radius in cm-units
do j=1,n
rr(j)=r(j)*1.e-4
enddo
do j=1,n
if (rr(j).le.1.e-3) then
winf(j)=stok*(rr(j)*rr(j)+cunh*rr(j))
elseif (rr(j).gt.1.e-3.and.rr(j).le.5.35e-2) then
x=log(stb*rr(j)*rr(j)*rr(j))
y=0.
do i=1,7
y=y+b(i)*(x**(i-1))
enddo
xrey=(1.+cunh/rr(j))*exp(y)
winf(j)=xrey*eta/(2.*rhoa*rr(j))
elseif (rr(j).gt.5.35e-2) then
bond=grav*(rhow-rhoa)*rr(j)*rr(j)/sigma
if (rr(j).gt.0.35) bond=grav*(rhow-rhoa)*0.35*0.35/sigma
x=log(16.*bond*py/3.)
y=0.
do i=1,6
y=y+c(i)*(x**(i-1))
enddo
xrey=py*exp(y)
winf(j)=xrey*eta/(2.*rhoa*rr(j))
if (rr(j).gt.0.35) winf(j)=xrey*eta/(2.*rhoa*0.35)
endif
enddo
return
end

subroutine effic
c collision efficiencies of hall kernel
parameter (n=400)
implicit double precision (a-h,o-z)
common /grid/ g(n),r(n),e(n)
common /kern/ ck(n,n),ec(n,n)
dimension rat(21),r0(15),ecoll(15,21)
data r0 /6.,8.,10.,15.,20.,25.,30.,40.,50.,
& 60.,70.,100.,150.,200.,300./
data rat /0.,0.05,0.1,0.15,0.2,0.25,0.3,0.35,0.4,0.45,0.5,
& 0.55,0.6,0.65,0.7,0.75,0.8,0.85,0.9,0.95,1.0/
data ecoll /
& 0.001,0.001,0.001,0.001,0.001,0.001,0.001,0.001,0.001,0.001
& ,0.001,0.001,0.001,0.001,0.001,0.003,0.003,0.003,0.004,0.005
& ,0.005,0.005,0.010,0.100,0.050,0.200,0.500,0.770,0.870,0.970
& ,0.007,0.007,0.007,0.008,0.009,0.010,0.010,0.070,0.400,0.430
& ,0.580,0.790,0.930,0.960,1.000,0.009,0.009,0.009,0.012,0.015
& ,0.010,0.020,0.280,0.600,0.640,0.750,0.910,0.970,0.980,1.000
& ,0.014,0.014,0.014,0.015,0.016,0.030,0.060,0.500,0.700,0.770
& ,0.840,0.950,0.970,1.000,1.000,0.017,0.017,0.017,0.020,0.022
& ,0.060,0.100,0.620,0.780,0.840,0.880,0.950,1.000,1.000,1.000
& ,0.030,0.030,0.024,0.022,0.032,0.062,0.200,0.680,0.830,0.870
& ,0.900,0.950,1.000,1.000,1.000,0.025,0.025,0.025,0.036,0.043
& ,0.130,0.270,0.740,0.860,0.890,0.920,1.000,1.000,1.000,1.000
& ,0.027,0.027,0.027,0.040,0.052,0.200,0.400,0.780,0.880,0.900
& ,0.940,1.000,1.000,1.000,1.000,0.030,0.030,0.030,0.047,0.064
& ,0.250,0.500,0.800,0.900,0.910,0.950,1.000,1.000,1.000,1.000
& ,0.040,0.040,0.033,0.037,0.068,0.240,0.550,0.800,0.900,0.910
& ,0.950,1.000,1.000,1.000,1.000,0.035,0.035,0.035,0.055,0.079
& ,0.290,0.580,0.800,0.900,0.910,0.950,1.000,1.000,1.000,1.000
& ,0.037,0.037,0.037,0.062,0.082,0.290,0.590,0.780,0.900,0.910
& ,0.950,1.000,1.000,1.000,1.000,0.037,0.037,0.037,0.060,0.080
& ,0.290,0.580,0.770,0.890,0.910,0.950,1.000,1.000,1.000,1.000
& ,0.037,0.037,0.037,0.041,0.075,0.250,0.540,0.760,0.880,0.920
& ,0.950,1.000,1.000,1.000,1.000,0.037,0.037,0.037,0.052,0.067
& ,0.250,0.510,0.770,0.880,0.930,0.970,1.000,1.000,1.000,1.000
& ,0.037,0.037,0.037,0.047,0.057,0.250,0.490,0.770,0.890,0.950
& ,1.000,1.000,1.000,1.000,1.000,0.036,0.036,0.036,0.042,0.048
& ,0.230,0.470,0.780,0.920,1.000,1.020,1.020,1.020,1.020,1.020
& ,0.040,0.040,0.035,0.033,0.040,0.112,0.450,0.790,1.010,1.030
& ,1.040,1.040,1.040,1.040,1.040,0.033,0.033,0.033,0.033,0.033
& ,0.119,0.470,0.950,1.300,1.700,2.300,2.300,2.300,2.300,2.300
& ,0.027,0.027,0.027,0.027,0.027,0.125,0.520,1.400,2.300,3.000
& ,4.000,4.000,4.000,4.000,4.000/
c two-dimensional linear interpolation of the collision efficiency
do j=1,n
do i=1,j
do k=2,15
if (r(j).le.r0(k).and.r(j).ge.r0(k-1)) then
ir=k
elseif (r(j).gt.r0(15)) then
ir=16
elseif (r(j).lt.r0(1)) then
ir=1
endif
enddo
rq=r(i)/r(j)
do kk=2,21
if (rq.le.rat(kk).and.rq.gt.rat(kk-1)) iq=kk
enddo
if (ir.lt.16) then
if (ir.ge.2) then
p=(r(j)-r0(ir-1))/(r0(ir)-r0(ir-1))
q=(rq-rat(iq-1))/(rat(iq)-rat(iq-1))
ec(j,i)=(1.-p)*(1.-q)*ecoll(ir-1,iq-1)+
& p*(1.-q)*ecoll(ir,iq-1)+
& q*(1.-p)*ecoll(ir-1,iq)+
& p*q*ecoll(ir,iq)
else
q=(rq-rat(iq-1))/(rat(iq)-rat(iq-1))
ec(j,i)=(1.-q)*ecoll(1,iq-1)+q*ecoll(1,iq)
endif
else
q=(rq-rat(iq-1))/(rat(iq)-rat(iq-1))
ek=(1.-q)*ecoll(15,iq-1)+q*ecoll(15,iq)
ec(j,i)=min(ek,1.d0)
endif
ec(i,j)=ec(j,i)
if (ec(i,j).lt.1.e-20) stop 99
enddo
enddo
return
end
 
Isn't this the same question you asked at the end of your other thread? If so, I posted what I think might be wrong in that thread.
 

Similar threads

Fortran Passing variables in FORTRAN
  • · Replies 25 ·
Replies
25
Views
4K
Fortran Fortran Code Help for Mechanical Engineers | Solving Debugging Issues
  • · Replies 17 ·
Replies
17
Views
4K
Fortran Errors in Fortran code for solving Laplace's equation in 3D
  • · Replies 4 ·
Replies
4
Views
3K
Fortran Fortran subroutine solving for the velocity and displacement of a spring
  • · Replies 8 ·
Replies
8
Views
2K
Fortran Fixing Fortran Runtime Error for Beginner
  • · Replies 17 ·
Replies
17
Views
7K
Fortran What are common errors when compiling a simple Fortran code?
  • · Replies 2 ·
Replies
2
Views
3K
Fortran Fortran external functions vs subroutines
  • · Replies 8 ·
Replies
8
Views
4K
Fortran How does Fortran know the Inputs and Outputs of Subroutines
  • · Replies 59 ·
2
Replies
59
Views
12K
Fortran Assign 2D Array with Reals in Fortran - Help Needed
  • · Replies 4 ·
Replies
4
Views
2K
Fortran Runtime error while reading file in fortran90
  • · Replies 5 ·
Replies
5
Views
5K