1. Лабораторная работа 01 Вводные замечания
 Скачать 104.44 Kb.
|
|
1. Лабораторная работа №0 1.1. Вводные замечания Каждая лабораторная состоит из нескольких заданий. Выполнение каждого задания заключается в напи- сании одной или нескольких функций или подпрограмм на языке Fortran. Перед началом выполнения заданий прочтите теоретическую часть, где приводятся необходимые формулы и даются дополнительные пояснения. Всю информацию по языку Fortran и OpenMP можно найти в лекциях и методическом пособии. При написании программы следует организовать код определенным образом. Для примера смотрите тре- нировочную лабораторную № 0. В каждом задании необходимо провести замеры времени работы функций и добиться того, чтобы парал- лельная версия работала быстрее последовательной. Кроме того, при увеличении количества потоков время выполнения программы должно уменьшаться. Как производить замеры и как нарисовать графики времени выполнения, ускорения и эффективности смотрите в нижеследующем описании лабораторной №0. Продемон- стрируем образец выполнения лабораторной работы на примере двух простейших заданий. Основное внимание уделим тому, как правильно оформлять код и тестировать созданную программу. Принципы оформления и тестирования в своей основе не зависят от выбранного языка программирования и их понимание пригодится при написании программ на любом другом языке программирования. 1.2. Задание №1 Написать функцию или подпрограмму, которая решает квадратное уравнение 𝑎𝑥 2 + 𝑏𝑥 + 𝑐 = 0 с действи- тельными коэффициентами 𝑎, 𝑏 и 𝑐. Протестировать программу на корректность работы. 1.3. Задание №2 Написать функцию, находящую сумму элементов одномерного массива с использованием OpenMP. Про- тестировать функцию на корректность и на производительность, то есть удостоверится, что при увеличении количества потоков время вычислений уменьшается. 1.4. Общие пояснения Вышеприведенные задания крайне просты и обычно даются учащимся в самом начале обучения програм- мирования. Однако мы подойдем к выполнению данных заданий используя методику написания крупных программ. • В первую очередь мы организуем программу следуя принципам структурного программирования: разо- бьем ее на функции, подпрограммы и исходные файлы. • Во вторую очередь мы напишем тесты, автоматически проверяющие наш код, что позволит отловить большинство потенциальных ошибок. Начнем с того, что создадим каталог lab00 и все исходные файлы, касающиеся лабораторной №0 будем хранить исключительно в нем. Далее мысленно представим из каких частей будет состоять наша первая про- грамма, решающая квадратное уравнение. • Непосредственное решение уравнения. Можно оформить его в виде главной программы, но так как мы хотим далее протестировать этот код, указывая разные коэффициенты 𝑎, 𝑏 и 𝑐, то выделить его в виде отдельной подпрограммы или функции. • Проверка работы нашей программы (тестирование). На первый взгляд проверять в такой простой про- грамме нечего, но далее мы увидим, что это не так. Для второго задания программа кроме вышеперечисленных двух частей будет содержать еще две. • Проверка эффективности работы параллельной версии. 1 • Визуализация результатов замеров быстродействия параллельной версии программы, то есть построение графиков времени выполнения, ускорения и эффективности. Кроме исходного кода в процессе сборки и компиляции программы могут появится временные служебные файлы .mod и .o, непосредственно исполняемые файлы программ, а также файлы с изображениями графиков. Даже в такой небольшой программе как наша, будет неудобно, если все эти файлы будут находиться в одном каталоге. Поэтому создадим следующую иерархию каталогов. • bin — в данный каталог будем сохранять созданные исполняемые файлы; изначально пуст. • imgs— в данный каталог будем сохранять графики времени выполнения, ускорения и эффективности; изначально пуст. • mod — в данный каталог сохраняются созданные компилятором mod файлы; изначально пуст. • src — в данном каталоге сохраним файлы с исходным кодом, который непосредственно ответственен за выполнение заданий. • test — в данном каталоге создадим файлы с исходным кодом тестирующих программ. Так как наши программы будут сопровождаться проверками на корректность выполнения и на произво- дительность, разумно команды компиляции и запуска не вводить каждый раз с консоли, а автоматизировать. Существует масса утилит для такой автоматизации (системы сборки). Для Unix (GNU Linux, macOS) простей- шей является утилита Make, а для Windows можно обойтись bat или powershell скриптом. Поэтому у нас в каталоге будут присутствовать еще два файла. • Makefile.bat — сценарий для запуска компиляции и тестирующих программ под Windows. • Makefile — make файл для сборки программ и запуска тестирования под GNU Linux и macOS. Незабудем также, что для второго задания следует нарисовать ряд графиков. Будем использовать язык Python с библиотекой Matplotlib. Скрипт назовем plot.py и напишем так, чтобы с минимальными изменениями он работал и для других заданий. Скриптами для сборки можно и не пользоваться, однако они существенно облегчат работу. 1.5. Выполнение задания №1 Начнем выполнение задания №1. Для вычисления корней квадратного уравнения 𝑎𝑥 2 +𝑏𝑥+𝑐 = 0 необходимо вычислить дискриминант 𝐷 = 𝑏 2 − 4𝑎𝑐 и далее корни по формуле 𝑥 1,2 = −𝑏 ± √ 𝐷 2𝑎 Следует отдельно учесть случаи 𝐷 > 0, 𝐷 < 0 и 𝐷 = 0. Так как Fortran поддерживает комплексные числа, то наша программа может выдать корректный результат и при 𝐷 < 0, а не сообщение об отсутствии дей- ствительных корней. Кроме того, с вычислительной точки зрения эффективней будет поделить уравнение на коэффициент 𝑎 и решать его в виде 𝑥 2 + 𝑝𝑥 + 𝑞 = 0 , где 𝑝 = 𝑏/𝑎 и 𝑞 = 𝑐/𝑎. В результате подпрограмма решающая квадратное уравнение может выглядить следующим образом. 31 subroutine quadroots(a, b, c, real_roots, complex_roots, is_real) 32 implicit none 33 real (real64), intent (in) :: a, b, c 34 real (real64), intent (out), dimension ( 1 : 2 ) :: real_roots 35 complex (real64), intent (out), dimension ( 1 : 2 ) :: complex_roots 36 logical , intent (out) :: is_real 37 38 real (real64) :: p, q 39 real (real64) :: D 40 2 41 p = b / a 42 q = c / a 43 D = p * p - 4.0 * q 44 is_real = .true. 45 if (D > 0 ) then 46 real_roots( 1 ) = 0.5 * ( - p - sqrt (D)) 47 reaL_roots( 2 ) = 0.5 * ( - p + sqrt (D)) 48 else if (D < 0 ) then 49 is_real = .false. 50 complex_roots( 1 ) = 0.5 * ( - p - sqrt ( cmplx (D, kind = 8 ))) 51 complex_roots( 2 ) = 0.5 * ( - p + sqrt ( cmplx (D, kind = 8 ))) 52 else if ( abs (D) < epsilon ( 1.0d0 )) then 53 real_roots( 1 ) = - 0.5 * p 54 real_roots( 2 ) = - 0.5 * p 55 end if 56 57 end subroutine quadroots Рассмотрим какие аргументы принимает данная подпрограмма. • Аргументы a, b, c соответствуют коэффициентами уравнения 𝑎, 𝑏 и 𝑐 и имеют тип real(real64) или иначе double precision. Подпрограмма не должна в процессе своей работы каким либо образом моди- фицировать аргументы a, b, c, на что указывает атрибут intent(in). • Аргумент real_roots является одномерным массивом из двух элементов типа double precision. В этот аргумент подпрограмма должна записать результат вычисления корней уравнения, в случае если корни действительные. • Аргумент complex_roots является одномерным массивом из двух элементов типа complex(real64). В этот аргумент подпрограмма должна записать результат вычисления корней уравнения, в случае если корни комплексные или чисто мнимые. • Аргумент is_real имеет тип logical(kind=1). Ему присваивается значение .true. в случае если корни уравнения действительные и .false. в случае комплексных корней. Обратите внимание, что проверка на равенство нулю дискриминанта сделана как 52 else if ( abs (D) < epsilon ( 1.0d0 )) then так как для действительных чисел необходимо учитывать погрешность представления с помощью чисел с плавающей запятой. Дополним нашу подпрограмму еще двумя функциями, вычисляющими значение квадратного трехчлена. Эти функции понадобятся при тестировании результатов вычисления. 10 ! Вычисление квадратного трехчлена 11 pure function real_p2(a, b, c, x) 12 implicit none 13 real (real64), intent (in) :: a, b, c 14 real (real64), intent (in) :: x 15 real (real64) :: real_p2 16 ! Оптимизируем число операций 17 real_p2 = (x + (b / a)) * x + (c / a) 18 end function real_p2 19 20 ! Вычисление квадратного трехчлена от комплексных чисел 21 pure function complex_p2(a, b, c, x) 22 implicit none 3 23 real (real64), intent (in) :: a, b, c 24 complex (real64), intent (in) :: x 25 real (real64) :: complex_p2 26 ! оптимизируем число операций 27 complex_p2 = dble ((x + (b / a)) * x + (c / a)) 28 end function complex_p2 Сохраним все процедуры в файле quadequation.f90, заключив их в модуль quadequation для последую- щего удобного вызова. В результате получим следующий исходный код. 1 module quadequation 2 use iso_fortran_env, only : int32, int64, real32, real64 3 implicit none 4 private 5 6 ! Доступные вовне функции и подпрограммы 7 public :: real_p2, complex_p2, quadroots 8 9 contains 10 ! Вычисление квадратного трехчлена 11 pure function real_p2(a, b, c, x) 12 implicit none 13 real (real64), intent (in) :: a, b, c 14 real (real64), intent (in) :: x 15 real (real64) :: real_p2 16 ! Оптимизируем число операций 17 real_p2 = (x + (b / a)) * x + (c / a) 18 end function real_p2 19 20 ! Вычисление квадратного трехчлена от комплексных чисел 21 pure function complex_p2(a, b, c, x) 22 implicit none 23 real (real64), intent (in) :: a, b, c 24 complex (real64), intent (in) :: x 25 real (real64) :: complex_p2 26 ! оптимизируем число операций 27 complex_p2 = dble ((x + (b / a)) * x + (c / a)) 28 end function complex_p2 29 30 ! Вычисление корней квадратного уравнения 31 subroutine quadroots(a, b, c, real_roots, complex_roots, is_real) 32 implicit none 33 real (real64), intent (in) :: a, b, c 34 real (real64), intent (out), dimension ( 1 : 2 ) :: real_roots 35 complex (real64), intent (out), dimension ( 1 : 2 ) :: complex_roots 36 logical , intent (out) :: is_real 37 38 real (real64) :: p, q 39 real (real64) :: D 40 41 p = b / a 42 q = c / a 43 D = p * p - 4.0 * q 44 is_real = .true. 45 if (D > 0 ) then 4 46 real_roots( 1 ) = 0.5 * ( - p - sqrt (D)) 47 reaL_roots( 2 ) = 0.5 * ( - p + sqrt (D)) 48 else if (D < 0 ) then 49 is_real = .false. 50 complex_roots( 1 ) = 0.5 * ( - p - sqrt ( cmplx (D, kind = 8 ))) 51 complex_roots( 2 ) = 0.5 * ( - p + sqrt ( cmplx (D, kind = 8 ))) 52 else if ( abs (D) < epsilon ( 1.0d0 )) then 53 real_roots( 1 ) = - 0.5 * p 54 real_roots( 2 ) = - 0.5 * p 55 end if 56 57 end subroutine quadroots 58 end module quadequation После того, как написали процедуру можно приступать к ее проверке. Необходима подавать на вход процедуре quadroots разные коэффициенты 𝑎, 𝑏 и 𝑐 и проверять правильность вычисленных корней. Даже для такой про- стой программы как наша можно придумать большое количество проверок. Проверяющая программа может выглядеть следующим образом. 1 include "../src/quadequation.f90" 2 3 ! Тестирование решений квадратного уравнения 4 subroutine test(a, b, c) 5 use iso_fortran_env, only : int32, int64, real32, real64 6 use quadequation 7 implicit none 8 9 real (real64), intent (in) :: a, b, c 10 11 real (real64), dimension ( 1 : 2 ) :: x 12 complex (real64), dimension ( 1 : 2 ) :: cx 13 logical :: is_real 14 15 is_real = .true. 16 17 call quadroots(a, b, c, x, cx, is_real) 18 19 ! действительные корни 20 if (is_real) then 21 ! проверяем корни 22 * , "test 01: " , abs (real_p2(a, b, c, x( 1 ))) <= epsilon ( 1.0d0 ) 23 * , "test 02: " , abs (real_p2(a, b, c, x( 2 ))) <= epsilon ( 1.0d0 ) 24 ! теорема Виета 25 * , "test 03: " , abs (x( 1 ) + x( 2 ) + b / a) <= epsilon ( 1.0d0 ) 26 * , "test 04: " , abs (x( 1 ) * x( 2 ) - c / a) <= epsilon ( 1.0d0 ) 27 ! комплексные корни 28 else if (. not . is_real) then 29 ! проверяем корни 30 * , "test 01: " , abs (complex_p2(a, b, c, cx( 1 ))) <= epsilon ( 1.0d0 ) 31 * , "test 02: " , abs (complex_p2(a, b, c, cx( 2 ))) <= epsilon ( 1.0d0 ) 32 ! теорема Виета 33 * , "test 03: " , abs (cx( 1 ) + cx( 2 ) + b / a) <= epsilon ( 1.0d0 ) 34 * , "test 04: " , abs (cx( 1 ) * cx( 2 ) - c / a) <= epsilon ( 1.0d0 ) 35 end if 36 end subroutine test 5 37 38 program test_quad 39 use iso_fortran_env, only : int32, int64, real32, real64 40 implicit none 41 real (real64) :: a, b, c 42 43 * , "Целые коэффициенты" 44 a = 1 ; b = 2 ; c = 1 45 call test(a, b, c) 46 47 * , "Рациональные коэффициенты" 48 a = 1.2 ; b = 2.5 ; c = 0.5 49 call test(a, b, c) 50 51 * , "Иррациональные коэффициенты" 52 a = 1 ; b = sqrt ( 2.0d0 ); c = 0.5 53 call test(a, b, c) 54 55 * , "Комплексные корни" 56 a = 4.0 ; b = 2.0 ; c = 1.0 57 call test(a, b, c) 58 59 * , "Большие числа" 60 a = 45678907.0 ; b = 20987654.0 ; c = 109876.0 61 call test(a, b, c) 62 63 * , "Большие отрицательные числа" 64 a = - 45678907.0 ; b = 20987654.0 ; c = - 109876.0 65 call test(a, b, c) 66 end program test_quad Мы едва ли перебрали все возможные проверки. Например, наша программа даст сбой, при 𝑎 = 0, так как в теле самой программы происходит деление на 𝑎. Для более крупных программ перебрать все возможные случаи и протестировать все ветки бывает зачастую невозможно. Для запуска программы выполним следующую команду 1 gfortran -Wall ./test/test_quad.f90 -o ./bin/test_quad.exe -J mod Так как в файле test_quad.f90 в первой строке содержится директива 1 include "../src/quadequation.f90" то компилятор автоматически подключит текст модуля quadequation и все переменные и процедуры до- ступные вовне модуля станет возможно вызывать и использовать в основной программе. Опция -J указывает, что файл mod, который компилятор генерирует автоматически, надо сохранить в каталог mod, а не создавать непосредственно в корневом каталоге. 1.6. Выполнение задания №2 Во втором задании следует используя директивы OpenMP реализовать итеративное суммирование эле- ментов массива. Организуем код точно также в виде файла с кодом модуля и заключенных в нем процедур. Хотя в данном случае можно обойтись всего одной функцией, но модуль все равно нужен, чтобы компилятор автоматически создал интерфейс для вызова функции. 1 module summation 2 use iso_fortran_env, only : int32, int64, real32, real64 3 implicit none 6 4 private 5 public :: iterative_sum 6 contains 7 8 function iterative_sum( array , threads_num) result (res) 9 implicit none 10 real (real64), dimension ( 1 :), intent (in) :: array 11 integer (int32), intent (in) :: threads_num 12 real (real64) :: res 13 integer (int32) :: length, i 14 15 length = size( array , 1 ) 16 res = 0.0 17 !$omp parallel shared(array) num_threads(threads_num) 18 !$omp do reduction(+ : res) 19 do i = 1 ,length, 1 20 res = res + array (i) 21 end do 22 !$omp end do 23 !$omp end parallel 24 end function iterative_sum 25 26 end module summation • Аргумент array — массив типа double precision, сумму элементов которого необходимо найти. • Аргумент threads_num — число потоков, которые функция может создать для проведения вычислений. Имеет тип integer. В отличие от предыдущего задания, в данном задании есть две тестирующие программы. Первая про- грамма проверяет функцию на корректность, а вторая замеряет время работы, ускорение и эффективность в зависимости от числа потоков. Приведем исходный код этих программ. 1 include "../src/summation.f90" 2 3 program test_sum 4 use iso_fortran_env, only : int32, int64, real32, real64 5 use summation 6 implicit none 7 8 ! Количество испытаний 9 integer (int32), parameter :: iterations_num = 100 10 ! Длина массива 11 integer (int64) :: array_length 12 ! Количество нитей 13 ! integer(int32) :: threads_num 14 15 ! Динамический массив для проверки функции 16 real (real64), dimension (:), allocatable :: test_array 17 real (real64) :: correct_sum 18 real (real64) :: test_res 19 integer (int64) :: i 20 ! Вначале проверяем на корректность вычислений 21 22 7 23 ! -------- Тест 1 --------- 24 array_length = 10 25 allocate (test_array( 1 :array_length)) 26 27 ! Массив натуральных чисел от 1 до конца массива 28 test_array = [(i, i = 1 ,array_length)] 29 correct_sum = 0.5 * (test_array( 1 ) + test_array(array_length)) * array_length 30 test_res = iterative_sum(test_array, 1 ) 31 32 * , test_res == correct_sum 33 34 deallocate (test_array) 35 ! ----------------- 36 37 ! -------- Тест 2 --------- 38 array_length = 1000 39 allocate (test_array( 1 :array_length)) 40 41 ! Массив отрицательных целых чисел 42 test_array = [( - i, i = 1 ,array_length)] 43 correct_sum = 0.5 * (test_array( 1 ) + test_array(array_length)) * array_length 44 test_res = iterative_sum(test_array, 1 ) 45 * , test_res == correct_sum 46 47 deallocate (test_array) 48 ! ----------------- 49 50 ! -------- Тест 3 --------- 51 array_length = 1000 52 allocate (test_array( 1 :array_length)) 53 54 ! Наполняем случайными числами 55 call random_number (test_array) 56 57 correct_sum = sum (test_array) 58 test_res = iterative_sum(test_array, 1 ) 59 * , abs (test_res - correct_sum) <= epsilon ( 1.0d0 ) 60 61 deallocate (test_array) 62 ! ----------------- 63 end program test_sum 1 include "../src/summation.f90" 2 3 program test_sum 4 use iso_fortran_env, only : int32, int64, real32, real64 5 use summation 6 implicit none 7 include "omp_lib.h" 8 9 ! Количество испытаний 10 integer (int32), parameter :: iterations_num = 1000 11 ! Максимальное количество нитей 12 integer (int32), parameter :: max_threads_num = 8 8 13 ! Длина массива 14 integer (int64) :: array_length 15 ! Количество нитей 16 integer (int32) :: threads_num 17 ! Для замера времени 18 real (real64) :: t1, t2 19 real (real64), dimension (iterations_num) :: T 20 21 ! Динамический массив для проверки функции 22 real (real64), dimension (:), allocatable :: test_array 23 real (real64) :: test_res 24 integer (int64) :: i 25 26 ! Создаем большой массив для тестирования производительности 27 array_length = 1000 * 1000 28 allocate (test_array( 1 :array_length)) 29 30 call random_number (test_array) 31 32 do threads_num = 1 ,max_threads_num, 1 33 do i = 1 ,iterations_num, 1 34 t1 = omp_get_wtime() 35 test_res = iterative_sum(test_array, threads_num) 36 t2 = omp_get_wtime() 37 T(i) = t2 - t1 38 end do 39 (i2,",",G0) , threads_num, sum (T) / dble (iterations_num) 40 end do 41 42 deallocate (test_array) 43 end program test_sum Для компиляции проверяющих программ и сборки исполняемых файлов следует выполнить следующие команды: 1 gfortran -fopenmp -Wall test_sum.f90 -o test_sum -J mod 2 gfortran -fopenmp -Wall test_sum_omp.f90 -o test_sum_omp -J mod В результате компиляции будет созданы два исполняемых файла: test_sum и test_sum_omp. Для визуализации результатов замеров времени выполняем следующую команду 1 ./test_sum_omp | python plot.py Пример получившихся графиков можно видеть на рисунках 1, 2 и 3. Обратите внимание, что вам необходимо выяснить сколько параллельных потоков поддерживает ваш процессор и указать в исходном коде программы test_sum_omp соответствующее число. Код программы plot.py специально разбирать не будем. При желании за пояснениями можете обратиться к преподавателю. 9 1 2 3 4 5 6 7 8 Количество потоков 𝑝 0.0005 0.0010 0.0015 0.0020 0.0025 0.0030 0.0035 Время выполнения [сек.] Рис. 1: Время, затраченное на вычисление суммы элементов массива в зависимости от количества порожденных потоков 1 2 3 4 5 6 7 8 Количество потоков 𝑝 1 2 3 4 5 6 7 Ускорение 𝑆 𝑝 = 𝑇 1 /𝑇 𝑝 Рис. 2: Время, затраченное на вычисление суммы элементов массива в зависимости от количества порожденных потоков 10 1 2 3 4 5 6 7 8 Количество потоков 𝑝 0.88 0.90 0.92 0.94 0.96 0.98 1.00 Эффективность 𝑆 𝑝 /𝑝 Рис. 3: Время, затраченное на вычисление суммы элементов массива в зависимости от количества порожденных потоков 11 |