1+ ! ==============================================================
2+ !
3+ ! SAMPLE SOURCE CODE - SUBJECT TO THE TERMS OF SAMPLE CODE LICENSE AGREEMENT,
4+ ! http://software.intel.com/en-us/articles/intel-sample-source-code-license-agreement/
5+ !
6+ ! Copyright 2016 Intel Corporation
7+ !
8+ ! THIS FILE IS PROVIDED "AS IS" WITH NO WARRANTIES, EXPRESS OR IMPLIED, INCLUDING BUT
9+ ! NOT LIMITED TO ANY IMPLIED WARRANTY OF MERCHANTABILITY, FITNESS FOR A PARTICULAR
10+ ! PURPOSE, NON-INFRINGEMENT OF INTELLECTUAL PROPERTY RIGHTS.
11+ !
12+ ! =============================================================
13+ !
14+ ! Part of the Coarray Tutorial. For information, please read
15+ ! Tutorial: Using Fortran Coarrays
16+ ! Getting Started Tutorials document
17+ program mcpi
18+ ! This program demonstrates using Fortran coarrays to implement the classic
19+ ! method of computing the mathematical value pi using a Monte Carlo technique.
20+ ! A good explanation of this method can be found at:
21+ ! http://www.mathcs.emory.edu/~cheung/Courses/170/Syllabus/07/compute-pi.html
22+ !
23+ ! Compiler options: /Qcoarray
24+ ! -coarray
25+ !
26+ ! Note for Visual Studio users - this source is excluded from building in the
27+ ! tutorial project. If you wish to change that, right click on the file in
28+ ! Solution Explorer and select Properties. Then go to Configuration Properties >
29+ ! General. Make sure that All Configurations and All Platforms are selected, then
30+ ! change Exclude File From Build to No. Be sure to remove or exclude the sequential
31+ ! source file from the build.
32+ implicit none
33+ ! Declare kind values for large integers, single and double precision
34+ integer , parameter :: K_BIGINT = selected_int_kind (15 )
35+ integer , parameter :: K_DOUBLE = selected_real_kind (15 ,300 )
36+ ! Number of trials per image. The bigger this is, the better the result
37+ ! This value must be evenly divisible by the number of images.
38+ integer (K_BIGINT), parameter :: num_trials = 600000000_K_BIGINT
39+ ! Actual value of PI to 18 digits for comparison
40+ real (K_DOUBLE), parameter :: actual_pi = 3.141592653589793238_K_DOUBLE
41+ ! Declare scalar coarray that will exist on each image
42+ integer (K_BIGINT) :: total[* ] ! Per-image subtotal
43+ ! Local variables
44+ real (K_DOUBLE) :: x,y
45+ real (K_DOUBLE) :: computed_pi
46+ integer :: i
47+ integer (K_BIGINT) :: bigi
48+ integer (K_BIGINT) :: clock_start,clock_end,clock_rate
49+ integer , allocatable :: seed_array(:)
50+ integer :: seed_size
51+ ! Image 1 initialization
52+ if (THIS_IMAGE() == 1 ) then
53+ ! Make sure that num_trials is divisible by the number of images
54+ if (MOD (num_trials,INT (NUM_IMAGES(),K_BIGINT)) /= 0_K_BIGINT ) &
55+ error stop " Number of trials not evenly divisible by number of images!"
56+ print ' (A,I0,A,I0,A)' , " Computing pi using " ,num_trials," trials across " ,NUM_IMAGES()," images"
57+ call SYSTEM_CLOCK (clock_start)
58+ end if
59+ ! Set the initial random number seed to an unpredictable value, with a different
60+ ! sequence on each image. The Fortran 2015 standard specifies a new RANDOM_INIT
61+ ! intrinsic subroutine that does this, but Intel Fortran doesn't yet support it.
62+ !
63+ ! What we do here is first call RANDOM_SEED with no arguments. The standard doesn't
64+ ! specify that behavior, but Intel Fortran sets the seed to a value based on the
65+ ! system clock. Then, because it's likely that multiple threads get the same seed,
66+ ! we modify the seed based on the image number.
67+ call RANDOM_SEED () ! Initialize based on time
68+ ! Alter the seed values per-image
69+ call RANDOM_SEED (seed_size) ! Get size of seed array
70+ allocate (seed_array(seed_size))
71+ call RANDOM_SEED (GET= seed_array) ! Get the current seed
72+ seed_array(1 ) = seed_array(1 ) + (37 * THIS_IMAGE()) ! Ignore potential overflow
73+ call RANDOM_SEED (PUT= seed_array) ! Set the new seed
74+ ! Initialize our subtotal
75+ total = 0_K_BIGINT
76+ ! Run the trials, with each image doing its share of the trials.
77+ !
78+ ! Get a random X and Y and see if the position
79+ ! is within a circle of radius 1. If it is, add one to the subtotal
80+ do bigi= 1_K_BIGINT ,num_trials/ int (NUM_IMAGES(),K_BIGINT)
81+ call RANDOM_NUMBER (x)
82+ call RANDOM_NUMBER (y)
83+ if ((x* x)+ (y* y) <= 1.0_K_DOUBLE ) total = total + 1_K_BIGINT
84+ end do
85+ ! Wait for everyone
86+ sync all
87+ ! Image 1 end processing
88+ if (this_image() == 1 ) then
89+ ! Sum all of the images' subtotals
90+ do i= 2 ,num_images()
91+ total = total + total[i]
92+ end do
93+ ! total/num_trials is an approximation of pi/4
94+ computed_pi = 4.0_K_DOUBLE * (REAL (total,K_DOUBLE)/ REAL (num_trials,K_DOUBLE))
95+ print ' (A,G0.8,A,G0.3)' , " Computed value of pi is " , computed_pi, &
96+ " , Relative Error: " ,ABS ((computed_pi- actual_pi)/ actual_pi)
97+ ! Show elapsed time
98+ call SYSTEM_CLOCK (clock_end,clock_rate)
99+ print ' (A,G0.3,A)' , " Elapsed time is " , &
100+ REAL (clock_end- clock_start)/ REAL (clock_rate)," seconds"
101+ end if
102+ end program mcpi
0 commit comments