The first Mandelbrot set plot was published in a paper by Brooks and Matelski in 1978. It is unclear what programming language they used, but it might have been F77.
This article has a short Fortran program trying to reproduce it:
But does not quite get it exactly. I was wondering if anybody can figure out exactly the parameters to get every star in the Mandelbrot “ascii-art” correctly.
Interesting. There is also fig. 1 probably the first plot of a Julia set, with c = 0.1+0.6i.
They thank Henry Laufer, a Princeton university mathematician, for the use of the computer.
The UNIVAC mainframe at Stony Brook had FORTRAN supplied by the vendor. It was almost certainly not F77 in Fall 1978.
The values of parameters are analyzed in fair detail in math.stachexchange.com
Just search for brooks matelski on that site to get
“Grid Spacing …” The cells are rectangular with aspect 5/3
Prof Laufer was a senior prof in Fall 1978. If someone was a gofer, was it a green postdoc
or a senior prof? He had the key to the Math dept
dumb terminal and acoustic modem which he let
me use. The department had no output device so
I had to hike to the computer building to retrieve
every printout.
Regards JPM
Thanks for the post @JPM and welcome to the forum!
Strange questions can arise in a bored brain. What would have they seen on a color device? In those days, computer monitors were most of the time monochrome (black & white, or green, or ochre…) and printers were noisy dot matrix printers. What would have bring a good color device? (supposing they would have think to the color code typically used today for Mandelbrot sets)
I did not use the exact algorithme of the 1978 paper but the mandelbrot_pixbuf.f90 gtk-fortran example. I set the size of the picture to 70x70 pixels (instead of 700x700). And I first used a maximum of 200 iterations for (z_n):
It is approximately the size of an icon, so I zoomed in Gimp:
With a maximum of 10000 iterations, the frontier is more precise and therefore more white pixels appear:
Finally it is interesting to consider that Gaston Julia worked on that kind of sequences around 1918 without seeing the esthetic beauty of that kind of images, not even the complexity of the object on a matrix printer. But he saw the intrinsic mathematical beauty of those sequences.
Nice work!
If you are wondering what the Mandelbrot set sounds like as it is printed, you don’t need to look further than this video of a Honeywell H316 running a Mandelbrot program in Fortran IV (printing starts at time 8:00 in the video):
The Fortran program, and other H316-related material can be found at the page of Philipp Hachtmann: http://h316.org/
In 1969, those flashing lights were almost certainly low-voltage incandescent bulbs, not LEDs which would become popular a few years later. And the paper tape was probably read mechanically, with electrical contacts meeting through the holes, and not optically. I transferred data and source code using paper tape a few years later in the mid 1970s, and by then there were optical readers, although I usually used a nearby slower mechanical reader.
My father brought us some used paper rolls like those in the video to play with a couple of times. I was about 6. Now I see how they where punched.
Nice video.
Thanks Ivan, it is a very interesting video about the history of technology. It shows the lot of time necessary to load, compile and run a program. And it shows the old input/output devices in action.
At the end of the video, they show a Mandelbrot set printed with half-tones. I have found that online generator:
https://pauln.github.io/asciibrot/
Finally, it is also important to preserve sounds of the past. We have images drawn thousands of years ago, but no sounds recording before the end of the 19th century. It reminds me the book The Tuning of the World by R. Murray Schafer, about soundscapes.
In the early 80’s, 8-bits home computers were often using audio tapes to load/store programs (floppy drivers were still expensive). So you could hear the sound if you unplugged the jack of the tape recorder. In the 90’s, with 56kbits modems you could also hear the same kind of sounds when connecting to your Internet service provider. The computer soundscape has changed.
For reference, I adapted the code from math history - Grid spacing, iterations used in the 1978 first published rendering of the Mandelbrot set? - Mathematics Stack Exchange to work for me (I only changed the “write” statements):
C *M*-SET PLOT FROM BROOKS AND MATELSKI PAPER
C OS/8 FORTRAN IV FOR DEC PDP-8 - SCRUSS, 2023-04
C THANKS TO JPM FOR PARAMETERS
C HTTPS://MATH.STACKEXCHANGE.COM/A/3510304/578897
C LPT IS ON LOGICAL UNIT 3
INTEGER IX,IY,K
DOUBLE PRECISION CR,CI,ZR,ZI,TI,D
D=.035D0
DO 40 IY=-15,15
C FORTRAN IV CANNOT COUNT BACKWARDS, SO NEGATE Y
CI=-1.66D0*D*IY
DO 30 IX=-35,35
CR=-.75D0+D*IX
ZR=0D0
ZI=0D0
DO 10 K=1,200
TI=2D0*ZR*ZI
ZR=ZR*ZR-ZI*ZI+CR
ZI=TI+CI
IF (ZR*ZR+ZI*ZI-4D0) 10,10,20
10 CONTINUE
C PRINT '*' IF C IS IN SET
WRITE (*,50,advance="no")
GO TO 30
C PRINT SPACE IF C IS NOT IN SET
20 WRITE (*,60,advance="no")
30 CONTINUE
C PRINT NEW LINE
WRITE (*,70)
40 CONTINUE
C FORTRAN CARRIAGE CONTROL: '+' SUPPRESSES NEW LINE, ' ' DOES NOT
C 'H' DENOTES HOLLERITH CONSTANTS: OBSOLETE STRING HANDLING
50 FORMAT (1H*)
60 FORMAT (1H )
70 FORMAT (1H )
END
This prints:
$ gfortran -std=legacy a.f && ./a.out
*
****
******
*****
* *********
*** ****************
******************** **
*************************
****************************
*******************************
******************************
* ***** ********************************
*********** ********************************
************* *******************************
** ************** ******************************
****************************************************************
** ************** ******************************
************* *******************************
*********** ********************************
* ***** ********************************
******************************
*******************************
****************************
*************************
******************** **
*** ****************
* *********
*****
******
****
*
Compared to the image from the paper:
I think the match is exact! Thank you @JPM, @vmagnin and @scruss for reproducing it.
FYI, R.A. Vowels in his book “Algorithms and Data Structures in F and Fortran” has a program for displaying Mandelbrot sets using the Salford graphics package that I assume is part of the Salford compiler. He also has a program to generate the Julia (as in Gaston Julia who discovered fractals not the language) fractal set. This is a nice book but I wish someone would take it upon themselves to write a proper algorithm and data structures book thats targeted exclusively at Modern Fortran. I think this would be an important contribution to the current literature and could go a long way to demonstrate to the world just how far Fortran has moved from F77. I doubt however you would be able to get one of the major publishing houses to publish it.
Asking ChatGPT 4 to “Translate this fixed format FORTRAN code to modern Fortran.” gives
! THANKS TO JPM FOR PARAMETERS
! HTTPS://MATH.STACKEXCHANGE.COM/A/3510304/578897
! LPT IS ON LOGICAL UNIT 3
program mandelbrot
implicit none
integer :: ix, iy, k
real(kind=8) :: cr, ci, zr, zi, ti, d
d = 0.035d0
do iy = -15, 15
! FORTRAN IV CANNOT COUNT BACKWARDS, SO NEGATE Y
ci = -1.66d0 * d * real(iy, kind=8)
do ix = -35, 35
cr = -0.75d0 + d * real(ix, kind=8)
zr = 0d0
zi = 0d0
do k = 1, 200
ti = 2d0 * zr * zi
zr = zr * zr - zi * zi + cr
zi = ti + ci
if (zr * zr + zi * zi > 4d0) exit
end do
! PRINT '*' IF C IS IN SET
if (k == 201) then
write(*, '(A)', advance="no") '*'
else
! PRINT SPACE IF C IS NOT IN SET
write(*, '(A)', advance="no") ' '
endif
end do
! PRINT NEW LINE
write(*, *)
end do
end program mandelbrot
which compiles with gfortran -std=f2018 -Wall -Wextra. When I run the executables of the original FORTRAN and translated Fortran code, the output looks the same to me in Notepad, but the Windows CMD fc (file compare) command says all the lines are different. I wonder what the difference is and how the Fortran 2018 code could be modified to exactly replicate the output of the original code.
Here is an equivalent modern Fortran version that I did by hand:
program mandelbrot
implicit none
integer, parameter :: dp = kind(0.d0)
real(dp), parameter :: D = 0.035_dp
integer :: ix, iy, k
complex(dp) :: c, z
logical escaped
do iy=-15,15
do ix=-35,35
c = cmplx(-0.75_dp+D*ix, -1.66_dp*D*iy, dp)
z = 0
escaped = .false.
do k=1,200
z = z**2 + c
if (z%re**2 + z%im**2 > 4) then
escaped = .true.
exit
end if
end do
if (escaped) then
write(*,"(a)",advance="no") " "
else
write(*,"(a)",advance="no") "*"
end if
end do
write(*,*)
end do
end program
Just a matter of taste, but this can be written as:
do k=1,200
z = z**2 + c
escaped = z%re**2 + z%im**2 > 4
if (escaped) exit
end do
or even more concisely as
do k=1,200
z = z**2 + c
escaped = abs(z) > 2
if (escaped) exit
end do