From db22368db2d51718e1a905595860b134f224c970 Mon Sep 17 00:00:00 2001 From: diener Date: Wed, 6 Feb 2008 15:51:04 +0000 Subject: I forgot to commit this yesterday. Changes related to the new reduction scheme and added OpenMP directives. git-svn-id: http://svn.einsteintoolkit.org/cactus/EinsteinInitialData/NoExcision/trunk@16 4ec1db94-0e4f-0410-ada3-8bed251432c9 --- src/NoExcision_mod.F90 | 193 +++++++++++++++++++++++++++---------------------- 1 file changed, 108 insertions(+), 85 deletions(-) (limited to 'src') diff --git a/src/NoExcision_mod.F90 b/src/NoExcision_mod.F90 index 26f9844..8ecf82e 100644 --- a/src/NoExcision_mod.F90 +++ b/src/NoExcision_mod.F90 @@ -15,6 +15,7 @@ module NoExcision_mod CCTK_INT :: loop_counter CCTK_REAL, dimension(16) :: infnormresid CCTK_REAL, dimension(16) :: delta_new, delta_old, alpha, beta, absres + CCTK_REAL, dimension(16) :: sumred, lsumred, infred, linfred character(len=18), dimension(16) :: red_names = (/ 'NoExcision::redgxx', & 'NoExcision::redgxy', & 'NoExcision::redgxz', & @@ -39,7 +40,7 @@ module NoExcision_mod ' kyz', ' kzz', & ' alp', 'betax', & 'betay', 'betaz' /) - integer :: nx, ny, nz + CCTK_INT :: nx, ny, nz integer :: ii logical, dimension(16) :: cont = (/ (.true., ii=1, 16) /) @@ -165,6 +166,7 @@ module NoExcision_mod CCTK_INT, dimension(:,:,:), intent(in) :: mask CCTK_REAL, intent(in) :: si CCTK_INT, intent(in) :: order + CCTK_INT :: offset if ( cont(1) ) call residual ( v1, mask, r1, si, order ) if ( cont(2) ) call residual ( v2, mask, r2, si, order ) @@ -195,7 +197,7 @@ module NoExcision_mod v9, v10, v11, v12, v13, v14, v15, v16, & r1, r2, r3, r4, r5, r6, r7, r8, & r9, r10, r11, r12, r13, r14, r15, r16, & - mask ) + mask, weight, do_inf_reduction, order ) CCTK_REAL, dimension(:,:,:), intent(in) :: u1, u2, u3, u4, u5, u6, & u7, u8, u9, u10, u11, u12, & @@ -207,106 +209,128 @@ module NoExcision_mod r7, r8, r9, r10, r11, r12, & r13, r14, r15, r16 CCTK_INT, dimension(:,:,:), intent(in) :: mask + CCTK_REAL, dimension(:,:,:), intent(in) :: weight + LOGICAL, intent(in) :: do_inf_reduction + CCTK_INT, intent(in) :: order -!$OMP PARALLEL + CCTK_INT i, j, k, offset - if ( cont(1) ) then -!$OMP WORKSHARE - where ( mask > 0 ) r1 = u1*v1 -!$OMP END WORKSHARE NOWAIT - end if + offset = order / 2 - if ( cont(2) ) then -!$OMP WORKSHARE - where ( mask > 0 ) r2 = u2*v2 -!$OMP END WORKSHARE NOWAIT - end if + lsumred = 0.0 - if ( cont(3) ) then -!$OMP WORKSHARE - where ( mask > 0 ) r3 = u3*v3 -!$OMP END WORKSHARE NOWAIT - end if + if (do_inf_reduction) linfred = 0.0 - if ( cont(4) ) then -!$OMP WORKSHARE - where ( mask > 0 ) r4 = u4*v4 -!$OMP END WORKSHARE NOWAIT - end if +!$OMP PARALLEL DO SCHEDULE(guided) PRIVATE(j,i) REDUCTION(+:lsumred) & +!$OMP REDUCTION(max:linfred) + do k = 1 + offset, nz - offset + do j = 1 + offset, ny - offset + do i = 1 + offset, nx - offset - if ( cont(5) ) then -!$OMP WORKSHARE - where ( mask > 0 ) r5 = u5*v5 -!$OMP END WORKSHARE NOWAIT - end if + if ( mask(i,j,k) > 0 ) then - if ( cont(6) ) then -!$OMP WORKSHARE - where ( mask > 0 ) r6 = u6*v6 -!$OMP END WORKSHARE NOWAIT - end if + if ( cont(1) ) then + r1(i,j,k) = u1(i,j,k)*v1(i,j,k) + lsumred(1) = lsumred(1) + weight(i,j,k)*r1(i,j,k) + if (do_inf_reduction) linfred(1) = max(linfred(1),r1(i,j,k)) + end if + + if ( cont(2) ) then + r2(i,j,k) = u2(i,j,k)*v2(i,j,k) + lsumred(2) = lsumred(2) + weight(i,j,k)*r2(i,j,k) + if (do_inf_reduction) linfred(2) = max(linfred(2),r2(i,j,k)) + end if - if ( cont(7) ) then -!$OMP WORKSHARE - where ( mask > 0 ) r7 = u7*v7 -!$OMP END WORKSHARE NOWAIT - end if + if ( cont(3) ) then + r3(i,j,k) = u3(i,j,k)*v3(i,j,k) + lsumred(3) = lsumred(3) + weight(i,j,k)*r3(i,j,k) + if (do_inf_reduction) linfred(3) = max(linfred(3),r3(i,j,k)) + end if - if ( cont(8) ) then -!$OMP WORKSHARE - where ( mask > 0 ) r8 = u8*v8 -!$OMP END WORKSHARE NOWAIT - end if + if ( cont(4) ) then + r4(i,j,k) = u4(i,j,k)*v4(i,j,k) + lsumred(4) = lsumred(4) + weight(i,j,k)*r4(i,j,k) + if (do_inf_reduction) linfred(4) = max(linfred(4),r4(i,j,k)) + end if - if ( cont(9) ) then -!$OMP WORKSHARE - where ( mask > 0 ) r9 = u9*v9 -!$OMP END WORKSHARE NOWAIT - end if + if ( cont(5) ) then + r5(i,j,k) = u5(i,j,k)*v5(i,j,k) + lsumred(5) = lsumred(5) + weight(i,j,k)*r5(i,j,k) + if (do_inf_reduction) linfred(5) = max(linfred(5),r5(i,j,k)) + end if - if ( cont(10) ) then -!$OMP WORKSHARE - where ( mask > 0 ) r10 = u10*v10 -!$OMP END WORKSHARE NOWAIT - end if + if ( cont(6) ) then + r6(i,j,k) = u6(i,j,k)*v6(i,j,k) + lsumred(6) = lsumred(6) + weight(i,j,k)*r6(i,j,k) + if (do_inf_reduction) linfred(6) = max(linfred(6),r6(i,j,k)) + end if - if ( cont(11) ) then -!$OMP WORKSHARE - where ( mask > 0 ) r11 = u11*v11 -!$OMP END WORKSHARE NOWAIT - end if + if ( cont(7) ) then + r7(i,j,k) = u7(i,j,k)*v7(i,j,k) + lsumred(7) = lsumred(7) + weight(i,j,k)*r7(i,j,k) + if (do_inf_reduction) linfred(7) = max(linfred(7),r7(i,j,k)) + end if - if ( cont(12) ) then -!$OMP WORKSHARE - where ( mask > 0 ) r12 = u12*v12 -!$OMP END WORKSHARE NOWAIT - end if + if ( cont(8) ) then + r8(i,j,k) = u8(i,j,k)*v8(i,j,k) + lsumred(8) = lsumred(8) + weight(i,j,k)*r8(i,j,k) + if (do_inf_reduction) linfred(8) = max(linfred(8),r8(i,j,k)) + end if - if ( cont(13) ) then -!$OMP WORKSHARE - where ( mask > 0 ) r13 = u13*v13 -!$OMP END WORKSHARE NOWAIT - end if + if ( cont(9) ) then + r9(i,j,k) = u9(i,j,k)*v9(i,j,k) + lsumred(9) = lsumred(9) + weight(i,j,k)*r9(i,j,k) + if (do_inf_reduction) linfred(9) = max(linfred(9),r9(i,j,k)) + end if - if ( cont(14) ) then -!$OMP WORKSHARE - where ( mask > 0 ) r14 = u14*v14 -!$OMP END WORKSHARE NOWAIT - end if + if ( cont(10) ) then + r10(i,j,k) = u10(i,j,k)*v10(i,j,k) + lsumred(10) = lsumred(10) + weight(i,j,k)*r10(i,j,k) + if (do_inf_reduction) linfred(10) = max(linfred(01),r10(i,j,k)) + end if - if ( cont(15) ) then -!$OMP WORKSHARE - where ( mask > 0 ) r15 = u15*v15 -!$OMP END WORKSHARE NOWAIT - end if + if ( cont(11) ) then + r11(i,j,k) = u11(i,j,k)*v11(i,j,k) + lsumred(11) = lsumred(11) + weight(i,j,k)*r11(i,j,k) + if (do_inf_reduction) linfred(11) = max(linfred(11),r11(i,j,k)) + end if - if ( cont(16) ) then -!$OMP WORKSHARE - where ( mask > 0 ) r16 = u16*v16 -!$OMP END WORKSHARE NOWAIT - end if + if ( cont(12) ) then + r12(i,j,k) = u12(i,j,k)*v12(i,j,k) + lsumred(12) = lsumred(12) + weight(i,j,k)*r12(i,j,k) + if (do_inf_reduction) linfred(12) = max(linfred(12),r12(i,j,k)) + end if -!$OMP END PARALLEL + if ( cont(13) ) then + r13(i,j,k) = u13(i,j,k)*v13(i,j,k) + lsumred(13) = lsumred(13) + weight(i,j,k)*r13(i,j,k) + if (do_inf_reduction) linfred(13) = max(linfred(13),r13(i,j,k)) + end if + + if ( cont(14) ) then + r14(i,j,k) = u14(i,j,k)*v14(i,j,k) + lsumred(14) = lsumred(14) + weight(i,j,k)*r14(i,j,k) + if (do_inf_reduction) linfred(14) = max(linfred(14),r14(i,j,k)) + end if + + if ( cont(15) ) then + r15(i,j,k) = u15(i,j,k)*v15(i,j,k) + lsumred(15) = lsumred(15) + weight(i,j,k)*r15(i,j,k) + if (do_inf_reduction) linfred(15) = max(linfred(15),r15(i,j,k)) + end if + + if ( cont(16) ) then + r16(i,j,k) = u16(i,j,k)*v16(i,j,k) + lsumred(16) = lsumred(16) + weight(i,j,k)*r16(i,j,k) + if (do_inf_reduction) linfred(16) = max(linfred(16),r16(i,j,k)) + end if + + end if + + end do + end do + end do +!$OMP END PARALLEL DO end subroutine multiply @@ -330,7 +354,6 @@ module NoExcision_mod v13, v14, v15, v16 CCTK_REAL, dimension(16), intent(in) :: c1, c2 CCTK_INT, dimension(:,:,:), intent(in) :: mask - !$OMP PARALLEL -- cgit v1.2.3