@@ -147,34 +147,49 @@ subroutine der_univ_subs_omp(du, recv_u_s, recv_u_e, n, dist_sa, dist_sc)
147147
148148 ! Local variables
149149 integer :: i, j! , b
150- real (dp) :: ur, bl, recp, du_1, du_n
151-
152- bl = dist_sa(1 )
153- ur = dist_sc(n)
154- recp = 1._dp / (1._dp - ur* bl)
150+ real (dp) :: ur, bl, recp, du_s, du_e
155151
156152 ! $omp simd
157153 do i = 1 , SZ
158- du_1 = recp* (du(i, 1 ) - bl* recv_u_s(i, 1 ))
159- du_n = recp* (du(i, n) - ur* recv_u_e(i, 1 ))
154+ ! A small trick we do here is valid for symmetric Toeplitz matrices.
155+ ! In our case our matrices satisfy this criteria in the (5:n-4) region
156+ ! and as long as a rank has around at least 20 entries the assumptions
157+ ! we make here are perfectly valid.
158+
159+ ! bl is the bottom left entry in the 2x2 matrix
160+ ! ur is the upper right entry in the 2x2 matrix
161+
162+ ! Start
163+ ! At the start we have the 'bl', and assume 'ur'
164+ bl = dist_sa(1 )
165+ ur = dist_sa(1 )
166+ recp = 1._dp / (1._dp - ur* bl)
167+ du_s = recp* (du(i, 1 ) - bl* recv_u_s(i, 1 ))
168+
169+ ! End
170+ ! At the end we have the 'ur', and assume 'bl'
171+ bl = dist_sc(n)
172+ ur = dist_sc(n)
173+ recp = 1._dp / (1._dp - ur* bl)
174+ du_e = recp* (du(i, n) - ur* recv_u_e(i, 1 ))
160175 end do
161176 ! $omp end simd
162177
163178 ! $omp simd
164179 do i = 1 , SZ
165- du(i, 1 ) = du_1
180+ du(i, 1 ) = du_s
166181 end do
167182 ! $omp end simd
168183 do j = 2 , n - 1
169184 ! $omp simd
170185 do i = 1 , SZ
171- du(i, j) = (du(i, j) - dist_sa(j)* du_1 - dist_sc(j)* du_n )
186+ du(i, j) = (du(i, j) - dist_sa(j)* du_s - dist_sc(j)* du_e )
172187 end do
173188 ! $omp end simd
174189 end do
175190 ! $omp simd
176191 do i = 1 , SZ
177- du(i, n) = du_n
192+ du(i, n) = du_e
178193 end do
179194 ! $omp end simd
180195
0 commit comments