1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
| module charchal
use, intrinsic :: iso_fortran_env, only : wp=>real64
implicit none
private
real, parameter :: pi = 4.*atan(1.)
public :: plot_charchal
contains
subroutine plot_charchal
integer(wp), parameter :: m = 2e5
real(wp) :: a(2,2), b(2), c(2,2), d(2), e(2,2), f(2), x(2), n, data(2, m)
integer(wp) :: i, clr(m), fn
call randinit(a)
call randinit(c)
call randinit(e)
b = [0.0, 2*rand()]
d = [0.0, rand()]
f = [0.0, rand()]
x = 0
clr = 0
do i = 1, m
call random_number(n)
if (n <= 0.85) then
x = matmul(a, x) + b
clr(i) = randint(3)
else if (n <= 0.94) then
x = matmul(c, x) + d
clr(i) = randint(3) + 1
else if (n <= 1) then
x = matmul(e, x) + f
clr(i) = 2 +randint(4)
end if
data(:, i) = x
! clr(i) = randint(3)
end do
open(newunit=fn, file='charchal.dat')
do i= 1, m
write(fn, *) data(:,i), clr(i)
end do
close(fn)
call execute_command_line('gnuplot -p charchal.plt')
end subroutine plot_charchal
subroutine randinit(a)
real(wp), intent(inout) :: a(2,2)
call random_number(a)
a = 2*a - 1
end subroutine randinit
function randint(seed) result(val)
integer, intent(in) :: seed
integer :: val
val = floor(seed*(1+rand()))
end function randint
end module charchal
program main
use charchal, only : plot_charchal
implicit none
call plot_charchal()
end program main
|
