c--------------------------------------------------------------------- c--------------------------------------------------------------------- subroutine ssor(niter) c--------------------------------------------------------------------- c--------------------------------------------------------------------- c--------------------------------------------------------------------- c to perform pseudo-time stepping SSOR iterations c for five nonlinear pde's. c--------------------------------------------------------------------- implicit none integer niter include 'mpinpb.h' include 'applu.incl' c--------------------------------------------------------------------- c local variables c--------------------------------------------------------------------- integer i, j, k, m integer istep double precision tmp double precision delunm(5), tv(5,isiz1,isiz2) external timer_read double precision wtime, timer_read integer IERROR ROOT = 0 c--------------------------------------------------------------------- c begin pseudo-time stepping iterations c--------------------------------------------------------------------- tmp = 1.0d+00 / ( omega * ( 2.0d+00 - omega ) ) c--------------------------------------------------------------------- c initialize a,b,c,d to zero (guarantees that page tables have been c formed, if applicable on given architecture, before timestepping). c--------------------------------------------------------------------- do m=1,isiz2 do k=1,isiz1 do j=1,5 do i=1,5 a(i,j,k,m) = 0.d0 b(i,j,k,m) = 0.d0 c(i,j,k,m) = 0.d0 d(i,j,k,m) = 0.d0 enddo enddo enddo enddo c--------------------------------------------------------------------- c compute the steady-state residuals c--------------------------------------------------------------------- call rhs c--------------------------------------------------------------------- c compute the L2 norms of newton iteration residuals c--------------------------------------------------------------------- call l2norm( isiz1, isiz2, isiz3, nx0, ny0, nz0, > ist, iend, jst, jend, > rsd, rsdnm ) call MPI_BARRIER( MPI_COMM_WORLD, IERROR ) call timer_clear(1) call timer_start(1) c--------------------------------------------------------------------- c the timestep loop c--------------------------------------------------------------------- do istep = 1, niter if (id .eq. 0) then if (mod ( istep, 20) .eq. 0 .or. > istep .eq. itmax .or. > istep .eq. 1) then if (niter .gt. 1) write( *, 200) istep 200 format(' Time step ', i4) endif endif c--------------------------------------------------------------------- c perform SSOR iteration c--------------------------------------------------------------------- do k = 2, nz - 1 do j = jst, jend do i = ist, iend do m = 1, 5 rsd(m,i,j,k) = dt * rsd(m,i,j,k) end do end do end do end do DO k = 2, nz -1 c--------------------------------------------------------------------- c form the lower triangular part of the jacobian matrix c--------------------------------------------------------------------- call jacld(k) c--------------------------------------------------------------------- c perform the lower triangular solution c--------------------------------------------------------------------- call blts( isiz1, isiz2, isiz3, > nx, ny, nz, k, > omega, > rsd, > a, b, c, d, > ist, iend, jst, jend, > nx0, ny0, ipt, jpt) END DO DO k = nz - 1, 2, -1 c--------------------------------------------------------------------- c form the strictly upper triangular part of the jacobian matrix c--------------------------------------------------------------------- call jacu(k) c--------------------------------------------------------------------- c perform the upper triangular solution c--------------------------------------------------------------------- call buts( isiz1, isiz2, isiz3, > nx, ny, nz, k, > omega, > rsd, tv, > d, a, b, c, > ist, iend, jst, jend, > nx0, ny0, ipt, jpt) END DO c--------------------------------------------------------------------- c update the variables c--------------------------------------------------------------------- do k = 2, nz-1 do j = jst, jend do i = ist, iend do m = 1, 5 u( m, i, j, k ) = u( m, i, j, k ) > + tmp * rsd( m, i, j, k ) end do end do end do end do c--------------------------------------------------------------------- c compute the max-norms of newton iteration corrections c--------------------------------------------------------------------- if ( mod ( istep, inorm ) .eq. 0 ) then call l2norm( isiz1, isiz2, isiz3, nx0, ny0, nz0, > ist, iend, jst, jend, > rsd, delunm ) c if ( ipr .eq. 1 .and. id .eq. 0 ) then c write (*,1006) ( delunm(m), m = 1, 5 ) c else if ( ipr .eq. 2 .and. id .eq. 0 ) then c write (*,'(i5,f15.6)') istep,delunm(5) c end if end if c--------------------------------------------------------------------- c compute the steady-state residuals c--------------------------------------------------------------------- call rhs c--------------------------------------------------------------------- c compute the max-norms of newton iteration residuals c--------------------------------------------------------------------- if ( ( mod ( istep, inorm ) .eq. 0 ) .or. > ( istep .eq. itmax ) ) then call l2norm( isiz1, isiz2, isiz3, nx0, ny0, nz0, > ist, iend, jst, jend, > rsd, rsdnm ) c if ( ipr .eq. 1.and.id.eq.0 ) then c write (*,1007) ( rsdnm(m), m = 1, 5 ) c end if end if c--------------------------------------------------------------------- c check the newton-iteration residuals against the tolerance levels c--------------------------------------------------------------------- if ( ( rsdnm(1) .lt. tolrsd(1) ) .and. > ( rsdnm(2) .lt. tolrsd(2) ) .and. > ( rsdnm(3) .lt. tolrsd(3) ) .and. > ( rsdnm(4) .lt. tolrsd(4) ) .and. > ( rsdnm(5) .lt. tolrsd(5) ) ) then c if (ipr .eq. 1 .and. id.eq.0) then c write (*,1004) istep c end if return end if end do call timer_stop(1) wtime = timer_read(1) call MPI_ALLREDUCE( wtime, > maxtime, > 1, > MPI_DOUBLE_PRECISION, > MPI_MAX, > MPI_COMM_WORLD, > IERROR ) return 1001 format (1x/5x,'pseudo-time SSOR iteration no.=',i4/) 1004 format (1x/1x,'convergence was achieved after ',i4, > ' pseudo-time steps' ) 1006 format (1x/1x,'RMS-norm of SSOR-iteration correction ', > 'for first pde = ',1pe12.5/, > 1x,'RMS-norm of SSOR-iteration correction ', > 'for second pde = ',1pe12.5/, > 1x,'RMS-norm of SSOR-iteration correction ', > 'for third pde = ',1pe12.5/, > 1x,'RMS-norm of SSOR-iteration correction ', > 'for fourth pde = ',1pe12.5/, > 1x,'RMS-norm of SSOR-iteration correction ', > 'for fifth pde = ',1pe12.5) 1007 format (1x/1x,'RMS-norm of steady-state residual for ', > 'first pde = ',1pe12.5/, > 1x,'RMS-norm of steady-state residual for ', > 'second pde = ',1pe12.5/, > 1x,'RMS-norm of steady-state residual for ', > 'third pde = ',1pe12.5/, > 1x,'RMS-norm of steady-state residual for ', > 'fourth pde = ',1pe12.5/, > 1x,'RMS-norm of steady-state residual for ', > 'fifth pde = ',1pe12.5) end