The Mandelbrot Set (Widget Program, Fortran)

      module myglob
        implicit none
        integer, parameter :: N=600
        integer :: id_draw,id_xmin,id_xmax,id_ymin,id_ymax,id_zoom,id_reset
        integer :: id_iter,id_scale,id_progress,id_cancel,id_undo
        integer :: iplot,izoom,iundo,izlog,iold,niter
        real, dimension (2) :: xscl,yscl,zscl,xold,yold
        real, dimension (N,N) :: zmat,zold 
      end module myglob

      program mandel
      use dislin
      use myglob
      implicit none

      interface
        subroutine myplot (id)
          implicit none
          integer, intent (in) :: id
        end subroutine myplot

        subroutine resetplot (id)
          implicit none
          integer, intent (in) :: id
        end subroutine resetplot

        subroutine zoomplot (id)
          implicit none
          integer, intent (in) :: id
        end subroutine zoomplot

        subroutine undoplot (id)
          implicit none
          integer, intent (in) :: id
        end subroutine undoplot
      end interface

      integer :: ip,ip1,ip2,id,id_but

      iplot = 0
      izoom = 0
      iundo = 0
      niter = 100
      izlog = 0
      iold = 0

      call swgtit ('DISLIN Mandelbrot Plot')

      call wgini ('hori', ip)
      call swgwth (-70)
      call wgbas (ip, 'vert', ip1)
      call swgwth (-15)
      call wgbas (ip, 'vert', ip2)

      call swgdrw (2100.0/2970.0)
      call wgdraw (ip1, id_draw)

      call wgltxt (ip2, 'xmin:', '-2.000000', 50, id_xmin) 
      call wgltxt (ip2, 'xmax:', '1.000000', 50, id_xmax) 
      call wgltxt (ip2, 'ymin:', '-1.000000', 50, id_ymin) 
      call wgltxt (ip2, 'ymax:', '1.000000', 50, id_ymax) 

      call wglab (ip2, ' ', id)
      call wgbut (ip2, 'Log. Colous', 0, id_scale) 
      call swgopt ('smooth', 'pbar')

      call wglab (ip2, ' ', id)
      call wglab (ip2, 'Progress:', id)
      call wgpbar (ip2, 0.0, 100.0, 5.0, id_progress) 

      call wglab (ip2, ' ', id)
      call wgltxt (ip2, 'Iterations:', '100', 40, id_iter) 

      call wglab (ip2, ' ', id)
      call wgpbut (ip2, 'Zoom', id_zoom)
      call swgcbk (id_zoom, zoomplot) 

      call wgpbut (ip2, 'Undo Zoom', id_undo)
      call swgcbk (id_undo, undoplot) 

      call wgpbut (ip2, 'Cancel', id_cancel)
      call wgpbut (ip2, 'Reset', id_reset)
      call swgcbk (id_reset, resetplot) 

      call wglab (ip2, ' ', id)
      call wgquit (ip2, id)
      call wgpbut (ip2, 'Plot', id_but)
      call swgcbk (id_but, myplot) 
      call wgfin ()
      end program mandel

      subroutine myplot (id)
      use dislin
      use myglob
      implicit none

      interface
        function iterate (cx, cy) result (m)
          implicit none
          double precision, intent (in) :: cx, cy
          integer :: m
        end function iterate
      end interface

      integer, intent (in) :: id
      integer :: i, j, nclr, isel, nx1, ny1, nx2, ny2,ibut
      double precision :: cx, cy, xd, yd
      real :: xa, xe, xor, xstp, ya, ye, yor, ystp, za, ze, zor, zstp,xv

      xa = -1.0
      xe =  1.0
      xor = -1.0
      xstp = 0.2

      ya = -1.0
      ye =  1.0
      yor = -1.0
      ystp = 0.2
 
      za = 1.0
      ze =  100.0
      zor = 10.0
      zstp = 10.0

      call setxid (id_draw, 'widget') 

      call gwgflt (id_xmin, xscl(1))
      call gwgflt (id_xmax, xscl(2))
      call gwgflt (id_ymin, yscl(1))
      call gwgflt (id_ymax, yscl(2))
      call gwgint (id_iter, niter)
      call gwgbut (id_scale, izlog)  

      xd = (xscl(2) - xscl(1)) / (N - 1)
      yd = (yscl(2) - yscl(1)) / (N - 1)
      call metafl ('cons')
      call scrmod ('revers')

      call disini ()
      if (izoom .eq. 0) call erase ()

      zscl(1) = 1.0
      zscl(2) = real (niter)

      call setscl (xscl, 2, 'x')
      call setscl (yscl, 2, 'y')
      call setscl (zscl, 2, 'z')

      call nochek ()
      call axspos (300, 1900)
      call ax3len (2200,1700,1700)

      if (izlog .eq. 1) then
        call axsscl ('log', 'z')
        call labels ('log', 'z')
      else
        call axsscl ('lin', 'z')
        call labels ('float', 'z')
      end if

      call graf3 (xa, xe, xor, xstp, ya, ye, yor, ystp, za, ze, zor, zstp) 
      call sendbf ()

      if ((izoom .eq. 0) .and. (iundo .eq. 0)) then  
        do i=1,N
          cx = xscl(1) + (i - 1) * xd
          call doevnt ()
          call gwgbut (id_cancel, ibut) 
          if (ibut .eq. 1) then 
            call swgbut (id_cancel, 0) 
            call erase ()
            iplot = 0
            call disfin ()
            return
          end if

          xv = (i - 1) * 100
          xv = xv / N 
          call swgval (id_progress, xv)
          do j=1,N
            cy = yscl(1) + (j - 1) * yd
            zmat(i,j) = iterate (cx, cy)
          end do
        end do
      end if

      call crvmat (zmat, N, N, 1, 1)

      if (izoom .eq. 1) then
        call csrrec (nx1, ny1, nx2, ny2)

        xold(1) = xscl(1)
        xold(2) = xscl(2)
        yold(1) = yscl(1)
        yold(2) = yscl(2)

        do i=1,N
          do j=1,N
  	    zold(i,j) = zmat(i,j)
          end do
        end do

        iold = 1
        xscl(1) = xinvrs (nx1)
        xscl(2) = xinvrs (nx2)
        yscl(1) = yinvrs (ny1)
        yscl(2) = yinvrs (ny2)

        call swgflt (id_xmin, xscl(1), 6)
        call swgflt (id_xmax, xscl(2), 6)
        call swgflt (id_ymin, yscl(1), 6)
        call swgflt (id_ymax, yscl(2), 6)
      end if

      call disfin ()
      iplot = 1
      end subroutine myplot

      subroutine zoomplot (id)
      use dislin
      use myglob
      implicit none
      integer, intent (in) :: id 

      if (iplot .eq. 0) return  

      izoom = 1
      call myplot (id)                
      if (iplot .eq. 0) return        

      izoom = 0
      call myplot (id)                ! replot with new scaling  
      end subroutine zoomplot

      subroutine undoplot (id)
      use dislin
      use myglob
      implicit none

      integer, intent (in) :: id   
      integer :: i,j   

      if (iold .eq. 0) return

      xscl(1) = xold(1)
      xscl(2) = xold(2)
      yscl(1) = yold(1)
      yscl(2) = yold(2)

      call swgflt (id_xmin, xscl(1), 6)
      call swgflt (id_xmax, xscl(2), 6)
      call swgflt (id_ymin, yscl(1), 6)
      call swgflt (id_ymax, yscl(2), 6)

      do i=1,N
        do j=1,N
          zmat(i,j) = zold(i,j)
        end do
      end do

      iundo = 1
      call myplot (id)
      iundo = 0
      end subroutine undoplot

      subroutine resetplot (id)
      use dislin
      use myglob
      implicit none

      integer, intent (in) :: id
       
      xscl(1) = -2.0
      xscl(2) = 1.0
      yscl(1) = -1.0
      yscl(2) = 1.0

      call swgflt (id_xmin, xscl(1), 6)
      call swgflt (id_xmax, xscl(2), 6)
      call swgflt (id_ymin, yscl(1), 6)
      call swgflt (id_ymax, yscl(2), 6)

      call myplot (id)
      end subroutine resetplot

      function iterate (cx, cy) result (m)
      use dislin
      use myglob
      implicit none

      double precision, intent (in) :: cx, cy
      integer :: m
      double precision :: x, y, x2, y2, zb, zbmax

      m = 0
      x = 0.0 
      y = 0.0
      zb = 0.0
      zbmax = 4.0

      do while ((zb .lt. zbmax) .and. (m .lt. niter)) 
        x2 = x * x - y * y + cx
        y2 = 2 * x * y + cy
        x = x2
        y = y2
        m = m + 1
        zb = x * x + y * y
      end do
      return
      end function iterate

News