Dpcpp ports IDR

[yhmtsai]

This PR adds IDR support
Note. there is a weird barrier in IDR kernel. I do not think it need the barrier but it is always failed when I add barrier in other place without it. The issue is only happened in CPU

TODO:

• use oneDPL random number generator to generate complex
• refine atomic_max such that be similar impl as cuda? we need the memoryspace template, so we can not do similar with cuda

[codecov]

Codecov Report

Merging #849 (aad7412) into develop (ebba714) will decrease coverage by 0.00%.
The diff coverage is 100.00%.

@@             Coverage Diff             @@
##           develop     #849      +/-   ##
===========================================
- Coverage    94.54%   94.54%   -0.01%     
===========================================
  Files          411      411              
  Lines        33134    33104      -30     
===========================================
- Hits         31326    31297      -29     
+ Misses        1808     1807       -1     

Impacted Files	Coverage Δ
omp/test/solver/idr_kernels.cpp	100.00% <100.00%> (ø)
omp/reorder/rcm_kernels.cpp	98.13% <0.00%> (+0.60%)	⬆️

---

Continue to review full report at Codecov.

Legend - Click here to learn more
Δ = absolute <relative> (impact), ø = not affected, ? = missing data
Powered by Codecov. Last update ebba714...aad7412. Read the comment docs.

[tcojean]

LGTM!

[tcojean]

Should it be like MKL?

Suggested change

	find_package(oneDPL REQUIRED HINTS "$ENV{DPL_ROOT}")
	find_package(oneDPL CONFIG REQUIRED HINTS "$ENV{DPL_ROOT}")

[yhmtsai]

I use the example from https://github.com/oneapi-src/oneDPL/tree/main/cmake#how-to-use-onedpl-package-from-cmake

[upsj]

They are not necessarily consistent inside Intel's documentation, and since there is no FindoneDPL.cmake module inside CMake (and I don't think there will be in the forseeable future), it doesn't really matter, except for slightly changing the error message.

[tcojean]

Is the oneMKL dot better on one RHS as well like in CUDA? Maybe our implementation is better?

[yhmtsai]

I do not test it actually

[tcojean]

Same question.

[upsj]

LGTM! I actually hadn't looked at the IDR kernels yet. Most of them look pretty well-suited for common kernels.

[upsj]

I am certain that we don't need this. We do not read or write any memory before this that would be accessed afterwards.

[yhmtsai]

I also think it should not need this one, but it is the place I can make it correct on cpu.
I do not have time to investigate it yet

[yhmtsai]

I am thinking it requires the global space sync because we update the value in global space later and use it in next loop.
Originally, I thought it always accessed by the same thread (with same tid), so it should be okay even it is not visible from global. However, it may not be true in CPU?

[upsj]

We have two other barriers before we access values again, and the column is being used exclusively by that thread. Does the code break if you remove this synchronization? If so, that sounds like a compiler issue.

[yhmtsai]

Yes, it will break if I remove the barrier.
I hope it's a compile issue. In my understanding, the value should be update because it is on the same thread.
It is only happened in cpu and we also need this weird kind barrier in coo or csr, too.

[upsj]

We only seem to be using one of them at a time. Can we avoid the duplicate allocation by storing everything in complex and reinterpret-casting to real for the norm? (that is safe by the standard, the other way round not necessarily).

[yhmtsai]

I guess I need to add sync before reinterpret-cast?
It may not be needed due to data dependency on values though

[upsj]

We only need a barrier between the dot and norm computation, or make sure that we don't accidentally claim they are non-aliasing with __restrict__

[upsj]

Suggested change

	UninitializedArray<ValueType, default_dot_dim *(default_dot_dim + 1)>
	*reduction_helper_array)
	UninitializedArray<ValueType, default_dot_dim *(default_dot_dim + 1)>
	&reduction_helper_array)

[yhmtsai]

I think it is the pointer not reference?

[upsj]

If you passed it by reference, you wouldn't need to dereference it later.

[yhmtsai]

@upsj I will not move some of them into common in this pr because some kernel in the function can not be in common.
For example, initialize contain initialize_m and initialize_subspace_vectors (able to common), but orthonormalize_subspace_vectors not. I think I can not put individual kernel into common?

[upsj]

@yhmtsai I see, that's fine by me then :)

[yhmtsai]

@upsj is it what you means?
also the calling,

cgh.parallel_for(
            sycl_nd_range(grid, block), [=](sycl::nd_item<3> item_ct1) {
                orthonormalize_subspace_vectors_kernel<block_size>(
                    num_rows, num_cols, values, stride, item_ct1,
                    *reduction_helper_array_acc_ct1.get_pointer(),
                    reduction_helper_real_acc_ct1.get_pointer().get());
            });

but I am a little unsure.
Is there any reason such that dpct do not use this way?

[upsj]

If you pass only the complex reduction helper, you can use

ValueType *__restrict__ reduction_helper = reduction_helper_array;
auto real_reduction_helper = reinterpret_cast<remove_complex<ValueType>*>(reduction_helper);

That should be fine.

[yhmtsai]

I mean the reference part

[upsj]

yes

[yhmtsai]

does real_reduction_helper not lead issue with reduction_helper?
reduction_helper thinks no one overlap on it, but real_reduction_helper does.

[upsj]

That's why it's important to set the __restrict__ on the first variable, but not the second! The compiler can't assume that pointers derived from the same pointer don't overlap.

[yhmtsai]

the second is overlapped on the first, so I thought it destroys the first __restrict__ property?
Might the first be arranged after the second (if we do not consider other dependency) because the first is __restrict__?

[upsj]

Now we are going deep into standard-ese, but from my reading of the C standard, pointers that are derived from a __restrict__ pointer may alias with them. Details: http://www.open-std.org/jtc1/sc22/wg14/www/docs/n1124.pdf 6.7.3.1 Formal definition of restrict. Anyways, I don't think using restrict here gives us any benefits, since we have memory barriers between writes and reads to the values. I am mostly concerned about unnecessary shared memory storage.

[tcojean]

rebase!

[yhmtsai]

rebase!

[upsj]

I think we can remove this and move the first __syncthreads() before reduction_helper[tidx] = dot, because reduce already synchronizes before accessing memory. Not sure whether the CUDA compiler removes duplicate barrier.sync instructions.

[yhmtsai]

I guess you mean removing the second one?
I will need to add another one before reduction_helper_real[tidx] = norm; to avoid reduction_helper_real[tidx] issue because they are on the same memory.

[upsj]

I mean

for (...) {
__syncthreads();
reduction_helper[tidx] = dot;
reduce(...);
__syncthreads();
dot = reduction_helper[0];
}
...
__syncthreads();
reduction_helper_real[tidx] = norm;
reduce(...);
__syncthreads();
norm = reduction_helper_real[0];

That makes sure that every time we write a reduction helper, we first make sure that nobody still needs to read it, and every time we read it, we first make sure that the reduction on it has already finished.

[yhmtsai]

I see. I am fine with both version

[sonarcloud]

Kudos, SonarCloud Quality Gate passed!   

0 Bugs
0 Vulnerabilities
0 Security Hotspots
1 Code Smell

No Coverage information
3.2% Duplication