program laplace
	implicit double precision (a-h,o-z)

	parameter (mn=201)
      common /t/t(mn,mn),told(mn,mn)

      real e,etime,et(2)

      open(unit=1,file='t.dat',form='formatted')

c         Initialize

      error=1.0d-06
      nc=0

	dx=1.0d0/dble(mn-1)
	do 10 i=1,mn
	do 10 j=1,mn
	t(i,j)=0.0d0
10    continue

c         Boundary Conditions

      ny=mn/10
	do 20 j=2,ny
      t(1,j)=dble(j-1)/dble(ny)
20    continue

      do 30 j=ny+1,mn
      t(1,j)=1.0d0
30    continue

      do 40 i=2,mn
      t(i,mn)=1.0d0
40    continue

c         to(i,j)

      do 45 i=1,mn
      do 45 j=1,mn
      told(i,j)=t(i,j)
45    continue

c         Laplace calculation

      e=etime(et)

1     nc=nc+1
      do 50 i=2,mn-1
      do 50 j=2,mn-1
          te=t(i+1,j)
          tw=t(i-1,j)
          tn=t(i,j+1)
          ts=t(i,j-1)
          t(i,j)=0.25d0*(te+tw+tn+ts)
50    continue

      do 60 j=2,mn-1
          t(mn,j)=t(mn-1,j)
60    continue

      res=0.0d0
      do 70 i=2,mn-1
      do 70 j=2,mn-1
      tnewv=t(i,j)
      toldv=told(i,j)
      res=res+dabs(tnewv-toldv)
      told(i,j)=t(i,j)
70    continue

      write(*,71) nc,dlog10(res)
71    format(2x,i7,2x,f6.3)

      if (res.le.error) go to 100
      go to 1

100   e=etime(et)

c         Output

      write(*,*) e,et(1),et(2)

      WRITE(1,*) mn,mn
      Y=-DX
      DO 110 J=1,MN
          Y=Y+DX
      X=-DX
      DO 110 I=1,MN
          X=X+DX
      WRITE(1,102) X,Y,T(I,J)
110   CONTINUE

102   FORMAT(3(2X,D15.8))

	end program laplace