forums.silverfrost.com

[curoi]

I have code trying to run multiple equations with 3-D array representing a well-mixed tank that is has some inflow = outflow where the 2-D mass transport of some contaminant mass ( c - concentration ) is defined by the advection-diffusion equation (ADE).

This tank has four boundary conditions that define the sides of the tank where

dc/dx = 0 at x = 0
dc/dx = 0 at x = L
dc/dy = 0 at y = 0
dc/dy = 0 at y = B

These boundary conditions yield 9 different finite difference approximation equations for the ADE ( 4 for the sides, 4 for the corners, and 1 for the internal cells)

The problem I am having is that I need to run all 9 equations simultaneously in an iterative way for all time ( using integer n ) so that the do loop for n = 1 to nts encompasses all 9 equations, I get this error from Silverfrost/Plato IDE at line 53:


"System. Overflow Exception

Arithmetic operation resulted in an overflow"

I had originally applied do = 1, nts loop to the internal cell equation and than another do = 1, nts loop for the other 8 equations and it ran the code even though it wasn't calculating the FDA's in the right order.

I'm wondering how to fix this error so that I can calculate all 9 equations within 1 do n = 1, nts do loop.

Fortran 95 code:

program main

! Declaration of variables

implicit none

real :: tottime, length, width, depth, dx, dy, dt, Q, u, D, mass
real, allocatable, dimension ( :: time, x, y
real, allocatable, dimension (:,:, :: conc
integer :: i, j, n, ndx, ndy, nts

! Variable values

length = 21.0 ! Length of the well-mixed tank (in meters)
width = 15.5 ! Width of the well-mixed tank (in meters)
depth = 2.0 ! Depth of the well-mixed tank (in meters)
ndx = 21 ! Number of spatial intervals in x-direction
ndy = 31 ! Number of spatial intervals in y-direction
dx = length / real( ndx ) ! Grid spacing in x-direction (in meters)
dy = width / real ( ndy ) ! Grid spacing in the y-direction (in meters)
Q = -3.0 ! Flow rate of water into and out of the tank (m^3/s)
u = Q / ( 15.0 * depth ) ! Velocity of water into and out of the tank (m/s)
D = 100.0 ! Diffusion coefficient (m^2/s) [DL = DT] - Longitudinal (x-direction) is the same as transverse (y-direction)
tottime = 25.0 ! Total time of the system
dt = 0.0005 ! Time step
nts = (tottime) / dt ! Number of time steps
mass = 100. ! Amount of mass injected into the system at x = 15, y = 7.5, and t = 0

! Allocation of time and concentration

allocate ( time( nts ), x( ndx ), y( ndy ), conc( ndy, ndx, nts ) )

conc = 0.0
time = 0.0

! Initial concentration is injected at x = 15 m (i = 16), y = 7.5 m (j = 16), and at t = 0 (n = 1)

conc( 16, 16, 1 ) = mass / ( dx * dy * depth ) * ( 1000000. / ( 100.**3 ) * 1000. )

! INTERNAL CELLS

do n = 2, nts

! INTERNAL CELLS
do i = 2, ndx - 1

do j = 2, ndy - 1

conc( j, i, n ) = conc( j, i, n - 1 ) - u * dt / dx * ( conc( j, i + 1, n - 1 ) - conc( j, i, n - 1 ) ) &
+ D * dt * ( ( 1 / dx**2 ) * ( conc( j, i + 1, n - 1 ) - 2 * conc( j, i, n - 1 ) + conc( j, i - 1, n - 1 ) ) &
+ ( 1 / dy**2 ) * ( conc( j + 1, i, n - 1 ) - 2 * conc( j, i, n - 1 ) + conc( j - 1, i, n - 1 ) ) ) <------------------- line 53

end do

end do


! BOUNDARY CONDITIONS

! CORNERS

! Bottom-left corner ! is conc supposed to be i + 1 near boundary?
conc( 1, 1, n ) = conc( 1, 1, n - 1 ) - u * dt / dx * ( conc( 1, 2, n - 1 ) - conc( 1, 1, n - 1 ) ) &
+ D * dt * ( ( 1 / dx**2 ) * ( 2 * conc( 1, 2, n - 1 ) - 2 * conc( 1, 1, n - 1 ) ) &
+ ( 1 / dy**2 ) * ( 2 * conc( 2, 1, n - 1 ) - 2 * conc( 1, 1, n - 1 ) ) )

[curoi]

! Top-left corner ! is conc supposed to be i + 1 near boundary?
conc( ndy, 1, n ) = conc( ndy, 1, n - 1 ) - u * dt / dx * ( conc( ndy, 2, n - 1 ) - conc( ndy, 1, n - 1 ) ) &
+ D * dt * ( ( 1 / dx**2 ) * ( 2 * conc( ndy, 2, n - 1 ) - 2 * conc( ndy, 1, n - 1 ) ) &
+ ( 1 / dy**2 ) * ( 2 * conc( ndy - 1, 1, n - 1 ) - 2 * conc( ndy, 1, n - 1 ) ) )

! Bottom-right corner
conc( 1, ndx, n ) = conc( 1, ndx, n - 1 ) - u * dt / dx * ( conc( 1, ndx, n - 1 ) - conc( 1, ndx - 1, n - 1 ) ) &
+ D * dt * ( ( 1 / dx**2 ) * ( 2 * conc( 1, ndx - 1, n - 1 ) - 2 * conc( 1, 1, n - 1 ) ) &
+ ( 1 / dy**2 ) * ( 2 * conc( 2, ndx, n - 1 ) - 2 * conc( 1, ndx, n - 1 ) ) )

! Top-right corner
conc( ndy, ndx, n ) = conc( ndy, ndx, n - 1 ) - u * dt / dx * ( conc( ndy, ndx, n - 1 ) - conc( ndy, ndx - 1, n - 1 ) ) &
+ D * dt * ( ( 1 / dx**2 ) * ( 2 * conc( ndy, ndx - 1, n - 1 ) - 2 * conc( ndy, ndx, n - 1 ) ) &
+ ( 1 / dy**2 ) * ( 2 * conc( ndy - 1, ndx, n - 1 ) - 2 * conc( ndy, ndx, n - 1 ) ) )

! SIDES

do i = 2, ndx - 1

! Bottom side
conc( 1, i, n ) = conc( 1, i, n - 1 ) - u * dt / dx * ( conc( 1, i, n - 1 ) - conc( 1, i - 1, n - 1 ) ) + D * dt &
* ( ( 1 / dx**2 ) * ( conc( 1, i + 1, n - 1 ) - 2 * conc( 1, i, n - 1 ) + conc( 1, i - 1, n - 1 ) ) &
+ ( 1 / dy**2 ) * ( 2 * conc( 2, i, n - 1 ) - 2 * conc( 1, 1, n - 1 ) ) )

! Top side
conc( ndy, i, n ) = conc( ndy, i, n - 1 ) - u * dt / dx * ( conc( ndy, i, n - 1 ) - conc( ndy, i - 1, n - 1 ) ) + D * dt &
* ( ( 1 / dx**2 ) * ( conc( ndy, i + 1, n - 1 ) - 2 * conc( ndy, i, n - 1 ) + conc( ndy, i - 1, n - 1 ) ) &
+ ( 1 / dy**2 ) * ( 2 * conc( ndy - 1, i, n - 1 ) - 2 * conc( ndy, i, n - 1 ) ) )

end do

do j = 2, ndy - 1

! Left side ! is conc supposed to be i + 1 near boundary?
conc( j, 1, n ) = conc( j, 1, n - 1 ) - u * dt / dx * ( conc( j, 1, n - 1 ) - conc( j, 2, n - 1 ) ) + D * dt &
* ( ( 1 / dx**2 ) * ( 2 * conc( j, 2, n - 1 ) - 2 * conc( j, 1, n - 1 ) ) &
+ ( 1 / dy**2 ) * ( conc( j + 1, 1, n - 1 ) - 2 * conc( j, 1, n - 1 ) + conc( j - 1, i, n - 1 ) ) )

! Right side
conc( j, ndx, n ) = conc( j, ndx, n - 1 ) - u * dt / dx * ( conc( j, ndx, n - 1 ) - conc( j, ndx - 1, n - 1 ) ) + D * dt &
* ( ( 1 / dx**2 ) * ( 2 * conc( j, ndx - 1, n - 1 ) - 2 * conc( j, ndx, n - 1 ) ) &
+ ( 1 / dy**2 ) * ( conc( j + 1, ndx, n - 1 ) - 2 * conc( j, ndx, n - 1 ) + conc( j - 1, ndx, n - 1 ) ) )

end do

end do


open(10, file = 'C:\g95\fortran\CEE 573\Project\diffusion25.csv', status = 'unknown')


do j = 1, ndy
y( j ) = dy * ( j - 1 )

write(10,110) ( conc( j, i, 8 ), i = 1, ndx )
110 format( 10000( f12.4, ',' ) )

end do


end program main

[PaulLaidler]

If you compile you code using CHECKMATE and then run the executable using the debugger (SDBG) then you see the line at which the program is failing and also the small (or zero) value in a denominator. It could be a programming error or simply an inappropriate numerical algorithm.

[curoi]

Okay, I tried that. I downloaded FTN95 Express in order to get the ability to use Checkmate and the debugger. But I checkmate seems spotty. I can't run Checkmate on code that I would like too fix including some of the examples that come downloaded with the software.

And now when I open up Plato IDE, there is a a dropdown menu at the top that allows me to select six different options including:

Checkmate .NET
Checkmate Win32
Debug .NET
Debug Win32

Plato IDE only lets me use the second and fourth options (Win32) and when I hit the continue debugging option (F5) it pops up another window that points to the declared real variable "length = 21."

[PaulLaidler]

You will find instructions on how to use Plato etc. in the help file FTN95.chm.