Serial Jacobi Iterative Scheme

# Serial Jacobi Iterative Scheme

### A single-process implementation of the Jacobi Scheme as applied to the Laplace equation is included below. Note that: F90 is used. System size, m, is determined at run time. Boundary conditions are handled by subroutine bc. Pointer arrays c, n, e, w, and s point to parts of the solution space, u. They are used to avoid unnecessary memory usage as well as to improve readability. This scheme is very slow to converge and is not used in practice. This code is shown primarily to provide a starting point for subsequent introduction of parallel concepts. It also serves as a starting point for convergence rate improvement ideas.

``````

! Solve Laplace equation using Jacobi iteration method
! Kadin Tseng, Boston University, November 1999

MODULE jacobi_module
IMPLICIT NONE
INTEGER, PARAMETER :: real4 = selected_real_kind(6,37)
INTEGER, PARAMETER :: real8 = selected_real_kind(15,307)
REAL(real8), DIMENSION(:,:), ALLOCATABLE         :: unew
REAL(real8), DIMENSION(:,:), ALLOCATABLE, TARGET :: u  ! solution array
REAL(real8) :: tol=1.d-4, gdel=1.0d0
REAL(real4) :: start_time, end_time
INTEGER :: m, iter = 0
PUBLIC

CONTAINS
SUBROUTINE bc(u, m)
! PDE: Laplacian u = 0;      0<=x<=1;  0<=y u(1:m  ,1:m  )          ! i  ,j    Current/Central for 1<=i<=m; 1<=j u(1:m  ,2:m+1)          ! i  ,j+1  North (of Current)
e => u(2:m+1,1:m  )          ! i+1,j    East  (of Current)
w => u(0:m-1,1:m  )          ! i-1,j    West  (of Current)
s => u(1:m  ,0:m-1)          ! i  ,j-1  South (of Current)

CALL bc(u, m)                ! set up boundary values

DO WHILE (gdel > tol)        ! iterate until error below threshold
iter = iter + 1            ! increment iteration counter
IF(iter > 5000) THEN
WRITE(*,*)'Iteration terminated (exceeds 5000)'
STOP                     ! nonconvergent solution
ENDIF
unew = ( n + e + w + s )*0.25 ! new solution, Eq. 3
gdel = MAXVAL(DABS(unew-c))    ! find local max error
IF(MOD(iter,10)==0) WRITE(*,"('iter,gdel:',i6,e12.4)")iter,gdel
c = unew                   ! update interior u
ENDDO

CALL CPU_TIME(end_time)      ! stop timer
PRINT *,'Total cpu time =',end_time - start_time,' x 1'
PRINT *,'Stopped at iteration =',iter
PRINT *,'The maximum error =',gdel

DEALLOCATE (unew, u)

END PROGRAM Jacobi
```
```