program Mandelbrot
      parameter( nmax=100 )
      parameter( mmax=100 )
      parameter( kmax=56, vmin=0, vmax=5 )
      real x(mmax+1), y(mmax+1), u(mmax+1,mmax+1)
      real dy, dx
      real xmin, xmax, ymin, ymax
      complex z, z0, zn, ima
      real pi(2,kmax+1)
      ima=(0.0,1.0)
*-- 出力ファイルをダンプ (保存) するための命令 ---
*      CALL SWpSET( 'LDUMP', .TRUE. )
*      CALL SWpSET( 'LWAIT', .FALSE. )
*      CALL SWpSET( 'LWAIT1', .FALSE. )
*-- 描画領域の設定 ---
      xmin=-2.0
      xmax=1.0
      ymin=-1.5
      ymax=1.5
*-- 座標変数の設定 ---
      dx=(xmax-xmin)/mmax
      dy=(ymax-ymin)/mmax
      do 10 i=1,mmax+1
         x(i)=xmin+dx*(i-1)
         y(i)=ymin+dy*(i-1)
  10  continue
*-- 漸化式の z0 の設定 (マンデルブロ集合なら z0=0) ---
      z0=(0.0,0.0)
*-- 漸化式の計算 ---
      do 20 j=1,mmax+1
         do 30 i=1,mmax+1
            z=x(i)+ima*y(j)
            zn=z0**2+z
            do 40 k=1,nmax
               zn=zn**2+z
               if(abs(zn).gt.2.0)then
                  u(i,j)=5.0
                  go to 30
               end if
  40  continue
            u(i,j)=abs(zn)            
  30  continue
  20  continue


*-- グラフの出力装置の指定 (IWS=4 は GTK 指定) ----
      WRITE(*,*) ' WORKSTATION ID (I)  ? ;'
      CALL SGPWSN
*      IWS=4
      READ (*,*) IWS

*-- 出力装置の呼び出し
      CALL GROPN( IWS )

*-- トーン塗り分けの設定 1 ---
*-- UESTLV は 100 以下しか呼び出せないので, 苦肉の策 ---
         dp=(vmax-vmin)/kmax
         do 110 k=1,kmax
            TLEV1=vmin+(k-1)*dp
            TLEV2=TLEV1+dp
            IPAT=(29+k)*1000+999
            CALL UESTLV( TLEV1, TLEV2, IPAT )
  110  continue
         do 120 k=1,kmax+1
            pi(1,k)=vmin+(k-1)*dp
            pi(2,k)=vmin+(k-1)*dp
  120  continue

      CALL GRFRM
      CALL GRSWND(  xmin,  xmax,  ymin,  ymax )
      CALL GRSVPT(  0.2,  0.8,  0.2,  0.8 )
*      CALL USPFIT
      CALL GRSTRN( 1 )
      CALL GRSTRF

         CALL UETONE( u, mmax+1, mmax+1, mmax+1 )  ! 格子点の個数を与える
         CALL USDAXS
         CALL UDCNTR( u, mmax+1, mmax+1, mmax+1 )
*         CALL UDCNTZ( u, mmax+1, mmax+1, mmax+1, u, 10000 )  ! 配列が100 を越えたときに使用する
*-- トーンバー（凡例）の設定 2 ---
         CALL GRSWND( 0.0, 1.0, vmin, vmax )
         CALL GRSVPT( 0.85, 0.9, 0.2, 0.8 )
         CALL GRSTRN( 1 )
         CALL GRSTRF

      CALL UETONE( pi, 2, 2, kmax+1 )

         CALL SLPVPR( 3 )            ! これは GRSTRN より後ろで宣言される
         CALL UZLSET( 'LABELYR', .TRUE. )
*         CALL UZFACT( 0.8 )
         CALL UYSFMT( '(F4.1)' )
         CALL UYAXDV( 'R', 0.1, 0.2 )

      CALL USSTTL( 'TIME', '', 'DEPARTMENT', '' )

      CALL GRCLS

      stop
      end