fortran_lorenz.f90 1.8 KB

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