見出し画像

【Fortranコードあり】美しいマンデルブロー集合をparaviewで可視化

美しいマンデルブロー集合をFortranで解くプログラムです。

Fortran90:ソルバ
paraview:可視化

Fortranでのソースコード(paraviewでの可視化用ファイル出力あり)を添付しておきます。
Fortranを学びたい方よろしくお願いいたします。
※Fortran90/95の書き方です。

paraviewで可視化すると美しいマンデルブロー集合が現れます。

美しいジュリア集合もお試しあれ。

Fortranのソースコード

program main
implicit none

integer, parameter :: Nmax=50, Nx=1024, Ny=1024
double precision, parameter :: Lx=4.d0, Ly=3.d0, dx=Lx/dble(Nx), dy=Ly/dble(Ny)
complex(kind(0d0)), parameter :: uso=(0.d0,1.d0)
integer :: i, j, n, m
double precision, allocatable, dimension(:) :: x, y
complex(kind(0d0)) , allocatable, dimension(:,:) ::z, c
character(len=40) filename

allocate( x(0:Nx-1) , y(0:Nx-1), z(0:Nx-1,0:Ny-1), c(0:Nx-1,0:Ny-1))

call system("mkdir data")
m = 0

!-----------------------------------------------------------------------------------
z = 0
do j=0,Ny-1
   do i=0,Nx-1
      x(i)     = -Lx/2.d0+dx*dble(i)
      y(j)     = -Ly/2.d0+dy*dble(j)
      c(i,j) = x(i) + uso*y(j)

      do n=0,Nmax
         if(cdabs(z(i,j))>2.0d0 .or. n == Nmax) exit
          z(i,j) = z(i,j)**2 + c(i,j)
      end do

      ! 発散しなかったCの集合(Cを上書き)
      if (n == Nmax) then
         c(i,j) = 0.0
      else
      ! 発散したCの集合(Cを上書き)
         c(i,j) = dble(n)
      end if
   end do
end do

deallocate(z)

!-----------------------------------------------------------------------------------
!output of movie file
write(filename,"(a,i5.5,a)") "data/julia_",int(m),".vtk"
open(10,file=filename)
write(10,"('# vtk DataFile Version 3.0')")
write(10,"('test')")
write(10,"('ASCII ')")

write(10,"('DATASET STRUCTURED_GRID')")
write(10,"('DIMENSIONS ',3(1x,i4))") Nx, Ny, 1

write(10,"('POINTS ',i9,' float')") Nx*Ny*1
do j=0,Ny-1
do i=0,Nx-1
   write(10,"(3(f9.4,1x))") x(i), y(j), 0.d0
end do
end do

write(10,"('POINT_DATA ',i9)") Nx*Ny*1

!date input
write(10,"('SCALARS julia float')")
write(10,"('LOOKUP_TABLE default')")
do j=0,Ny-1
   do i=0,Nx-1
      write(10,*)  cdabs(c(i,j))
   end do
end do

close(10)


end program main

Twitter➡@t_kun_kamakiri
Instagram➡kamakiri1225
youtube➡https://www.youtube.com/channel/UCbG6_Q9ZRqqVT6YZOpcjDlQ
ブログ➡宇宙に入ったカマキリ(物理ブログ)
ココナラ➡物理の質問サポートサービス
コミュニティ➡製造業ブロガー


この記事が気に入ったらサポートをしてみませんか?