1program main 2 implicit none 3 4 integer, parameter :: dp = 8 5 real(dp), dimension(1:3) :: x 6 integer, parameter :: nstep = 20000000 7 real(dp) :: t = 0.0_dp 8 real(dp) :: h = 1.0e-10_dp 9 integer, parameter :: nb_loops = 21 10 integer, parameter :: n = 3 11 integer :: k 12 integer :: time_begin 13 integer :: time_end 14 integer :: count_rate 15 real(dp) :: time 16 real(dp) :: min_time = 100.0 17 18 do k = 1, nb_loops 19 x = [ 8.5_dp, 3.1_dp, 1.2_dp ] 20 call system_clock(time_begin, count_rate) 21 call rk4sys(n, t, x, h, nstep) 22 call system_clock(time_end, count_rate) 23 time = real(time_end - time_begin, dp) / real(count_rate, dp) 24 min_time = min(time, min_time) 25 write (*,*) time, x(1) 26 end do 27 write (*,*) "Minimal Runtime:", min_time 28contains 29 subroutine xpsys(x,f) 30 real(dp), dimension(1:3), intent(in) :: x 31 real(dp), dimension(1:3), intent(out) :: f 32 f(1) = 10.0_dp * ( x(2) - x(1) ) 33 f(2) = 28.0_dp * x(1) - x(2) - x(1) * x(3) 34 f(3) = x(1) * x(2) - (8.0_dp / 3.0_dp) * x(3) 35 end subroutine xpsys 36 37 subroutine rk4sys(n, t, x, h, nstep) 38 integer, intent(in) :: n 39 real(dp), intent(in) :: t 40 real(dp), dimension(1:n), intent(inout) :: x 41 real(dp), intent(in) :: h 42 integer, intent(in) :: nstep 43 ! Local variables 44 real(dp) :: h2 45 real(dp), dimension(1:n) :: y, f1, f2, f3, f4 46 integer :: i, k 47 48 h2 = 0.5_dp * h 49 do k = 1, nstep 50 call xpsys(x, f1) 51 y = x + h2 * f1 52 call xpsys(y, f2) 53 y = x + h2 * f2 54 call xpsys(y, f3) 55 y = x + h * f3 56 call xpsys(y, f4) 57 x = x + h * (f1 + 2.0_dp * (f2 + f3) + f4) / 6.0_dp 58 end do 59 end subroutine rk4sys 60end program main 61