Errors in Fortran code for solving Laplace's equation in 3D

  • #1
physicsuser2023
3
0
TL;DR Summary
I encountered errors while running a Fortran code to calculate potential values inside a cube, with the first error being a reference to an undefined variable and the second error being insufficient storage allocation.
hey everyone. I wrote a code in Fortran to calculate potential values or solve the Laplace equation inside a cube according to the boundary conditions mentioned in the code.
this is my code:
Fortran:
program laplace_cubic

  implicit none

  REAL*8 :: LX, LY, LZ, DELTA, MAX_ERR, ERR
  INTEGER :: NX, NY, NZ, I, J, K, M
  INTEGER*8 :: MAX_COUNT
  REAL*8, DIMENSION(:,:,:), ALLOCATABLE :: PHI, PHIO
 
  write(*, *) 'LX='
  read(*, *) LX
  write(*, *) 'LY='
  read(*, *) LY
  write(*, *) 'LZ='
  read(*, *) LZ
  write(*, *) 'ENTER MAX COUNTER='
  read(*, *) MAX_COUNT
  write(*, *) 'ENTER MAX ERROR='
  read(*, *) MAX_ERR

  DELTA = 0.01
  NX = INT(LX/DELTA)
  NY = INT(LY/DELTA)
  NZ = INT(LZ/DELTA)
 
  ALLOCATE(PHI(0:NX, 0:NY, 0:NZ))
  ALLOCATE(PHIO(0:NX, 0:NY, 0:NZ))

  CALL INITIALIZE(PHI, NX, NY, NZ, LX, LY, LZ, DELTA)

DO M = 1, MAX_COUNT
    PHIO=PHI
    ERR=0.0
    DO I = 1, NX - 1
      DO J = 1, NY - 1
        DO K = 1, NZ - 1
          PHI(I, J, K) = (PHIO(I + 1, J, K) + PHIO(I - 1, J, K) + PHIO(I, J + 1, K) + &
               PHIO(I, J - 1, K) + PHIO(I, J, K + 1) + PHIO(I, J, K - 1)) / 6.0
          ERR = ERR + (PHI(I, J, K) - PHIO(I, J, K))**2
        END DO
      END DO
    END DO
    IF (ERR < MAX_ERR) EXIT
END DO

OPEN(100, FILE = 'laplaceoutput.TXT')
  DO I = 0, NX
    DO J = 0, NY
      DO K = 0, NZ
        WRITE(100, *) I * DELTA, J * DELTA, K * DELTA, PHI(I, J, K)
      END DO
    END DO
  END DO

end program laplace_cubic

SUBROUTINE INITIALIZE(PHI,NX,NY,NZ,LX,LY,LZ,DELTA)
   
    IMPLICIT NONE

    REAL,PARAMETER::PI=3.1415
    INTEGER::I,J,K
    INTEGER,INTENT(IN)::NX,NY,NZ
    REAL*8::DELTA,X,Y,Z,LX,LY,LZ
    REAL*8,DIMENSION(0:NX,0:NY,0:NZ)::PHI

DO J = 1, NY - 1
  DO K = 1, NZ - 1
    PHI(0, J, K) = 1.0  ! potential inside the back face (zy plane)
    PHI(NX,J,K) = 0.0
  END DO
END DO

DO J = 0, NY
  Y = J * DELTA
  PHI(0, J, 0) = SIN(PI*Y/LY)  !sinusoidal potential along y-axis
END DO

DO K = 0, NZ
  Z = K * DELTA
  PHI(0, 0, K) = SIN(PI*Z/LZ)  ! sinusoidal potential along z-axis
END DO

DO J = 1, NY-1
  DO I = 1, NX-1
    PHI(I, J, 0) = 1.0  ! potential inside the bottom face (xy plane)
    PHI(I,J,NZ) = 0.0
  END DO
END DO

DO I = 0, NX
   X= I * DELTA
  PHI(I, 0, 0) = SIN(PI*X/LX)  ! sinusoidal potential along x-axis
END DO

DO I = 1, NX-1
  DO K = 1, NZ-1
    PHI(I, 0, K) = 1.0  ! potential inside the left face (XZ plane)
    PHI(I,NY,K) = 0.0
  END DO
END DO

DO K = 1, NZ-1
  DO J = 1, NY-1
    DO I = 1, NX-1
      PHI(I, J, K) = 0.5  ! potential inside the cube
    END DO
  END DO
END DO

END SUBROUTINE INITIALIZE
I first entered the dimensions of the cube as 3x3x3, Iteration or MAX COUNTER as 1000, and MAX ERROR as 0.000001. After these inputs, I enter to make the program perform the operation, but it showed me an error: error 112, reference to undefined variable, array element or function result.
I then set the dimensions of the cube to 10x10x10 and the rest of the inputs the same as before. But this time I received another error: ALLOCATE was unable to obtain sufficient storage.
I tested a code similar to the above code for the two-dimensional case and it worked.
what's the solution?
 
Technology news on Phys.org
  • #2
With DELTA at .01, 10x10x10 results in both PHI and PHI0 arrays being dimensioned 1000x1000x1000.
So you were trying to allocate memory for 2 billion real values (16,000,000,000 bytes). So there is no mystery about the allocation error.

For the other error, it will be a lot easier to track down if you note where this error occurred in the code. I don't know what debug tools you have, but at a minimum, you can run the code after putting some "WRITE(*,*) 'at line #999' statements in.
 
  • Like
Likes Vanadium 50 and physicsuser2023
  • #3
Thank you for your guidance. Now I understand that the error related to 10*10*10 is a normal and reasonable thing. So can you suggest me some suitable input values to fix this error? Values for delta and other inputs. thanks
But the error related to 3*3*3. its said there is in line 32
 
  • #4
OK. So you have not fully initialized PHI. Some of the 301x301x301 entries in that array were not set in function INITIALIZE. On the face of it, it appears that you are never setting array elements with indexes of 300.
In Fortran, an array of 0:N has N+1 elements.
(So my calculation above should have been (1001x1001x1001).)
 
  • Like
Likes physicsuser2023
  • #5
Thank you for guidance, i check the code again
 

Similar threads

  • Programming and Computer Science
Replies
1
Views
907
Replies
1
Views
1K
  • Programming and Computer Science
Replies
3
Views
1K
  • Programming and Computer Science
Replies
7
Views
3K
  • Programming and Computer Science
Replies
12
Views
3K
  • Programming and Computer Science
Replies
5
Views
1K
  • Atomic and Condensed Matter
Replies
3
Views
737
  • Programming and Computer Science
Replies
4
Views
1K
  • Programming and Computer Science
Replies
8
Views
1K
  • Programming and Computer Science
Replies
5
Views
1K
Back
Top