Skip to content

Commit ecaffbe

Browse files
committed
krr
1 parent 5f7b4f6 commit ecaffbe

File tree

2 files changed

+247
-0
lines changed

2 files changed

+247
-0
lines changed

example/KRR.jl

Lines changed: 66 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,66 @@
1+
using Dagger
2+
using Base: read, reinterpret
3+
using LinearAlgebra
4+
#using Distances
5+
6+
function read_binary_to_matrix(filename::String, nrows::Int, ncols::Int)
7+
# Open the file in read mode
8+
open(filename, "r") do io
9+
# Read the entire file content as raw bytes
10+
bytes = read(io)
11+
# Reinterpret the raw bytes as an array of Float32
12+
data = reinterpret(Float32, bytes)
13+
# Reshape the array into the desired matrix dimensions
14+
new_data = reshape(data, nrows, ncols)
15+
16+
return new_data#[1:10, 1:20]
17+
end
18+
end
19+
20+
function all_elements_less_than(matrix::Array{Float32,2}, threshold::Float32)
21+
return all(matrix .< threshold)
22+
end
23+
24+
function gaussian_kernel(matrix::Array{Float32,2}, sigma::Float32)
25+
return exp.(-matrix.^2 ./ (2 * sigma^2))
26+
end
27+
28+
# Example usage:
29+
filename = "/ibex/ai/home/omairyrm/gbgwas/hicmamain/scripts/genotype09.bin" # Replace with your actual file path
30+
nrows, ncols = 1024, 30720 # Replace with the desired matrix dimensions
31+
32+
A = read_binary_to_matrix(filename, nrows, ncols)
33+
println(size(A))
34+
A = Matrix(A)
35+
36+
#nrows, ncols = 10, 20
37+
#A = rand(0:2, nrows, ncols)
38+
#display(A)
39+
println(size(A))
40+
Acpy = A
41+
AA= A'*A #syrk
42+
#display(AA)
43+
44+
C = diag(AA) .* ones(ncols, ncols)
45+
#display(C)
46+
D = C'
47+
#display(D)
48+
49+
C = ((-2) .* AA) + (C) + (D)
50+
#C = LowerTriangular(C) + LowerTriangular(D)
51+
#display(C)
52+
C= sqrt.(C)
53+
C[diagind(C)] .+= 10
54+
C = Float32.(C)
55+
sigma = Float32(1e1)
56+
println(typeof(C))
57+
58+
#Dl = pairwise(Euclidean(), Acpy)
59+
60+
gaussian_matrix = gaussian_kernel(C, sigma)
61+
display(gaussian_matrix)
62+
63+
DA = view(gaussian_matrix, Blocks(1024, 1024));
64+
65+
threshold = Float32(10^-3)
66+
MP = Dagger.adapt_precision(DA, threshold)

src/array/adapt_precision2.jl

Lines changed: 181 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,181 @@
1+
"""
2+
tile_precision(uplo, global_norm, scalar_factor, tolerance, A)
3+
4+
it receives tile and it compute required precision per tile
5+
6+
### Input
7+
- `A` -- tile of size m x n
8+
- `global_norm` -- global norm of the whole matrix
9+
- `scalar_factor` -- scale tile by this value which is the number of tiles
10+
- `tolerance` -- user defined tolerance as required aby the application
11+
12+
### Output
13+
The required precision of the tile
14+
15+
"""
16+
function tile_precision(A, global_norm, scalar_factor, tolerance)
17+
18+
tile_sqr = mapreduce(LinearAlgebra.norm_sqr, +, A)
19+
20+
tile_norm = sqrt(tile_sqr)
21+
22+
cal = tile_norm * scalar_factor / global_norm
23+
decision_hp = tile_norm * scalar_factor / global_norm < tolerance / eps(Float16)
24+
decision_sp = tile_norm * scalar_factor / global_norm < tolerance / eps(Float32)
25+
26+
#We are planning in near future to support fp8 E4M3 and E5M2
27+
decision_fp8 = tile_norm * scalar_factor / global_norm < tolerance / 0.0625
28+
if decision_fp8
29+
return "Float8"
30+
elseif decision_hp
31+
return "Float16"
32+
elseif decision_sp
33+
#@show m, n
34+
return "Float32"
35+
else
36+
return "Float64"
37+
end
38+
end
39+
40+
"""
41+
function adapt_precision( A::UpperTriangular{T,<:DArray{T,2}},
42+
MP::UpperTriangular{String,<:DArray{String,2}}, tolerance::Float64) where {T}
43+
44+
it iterates over all tiles and calculates the required precision per tile based on formulation from Nicholas J. Higham
45+
46+
### Input
47+
- `A` -- Dagger UpperTriangular array of tiles with real values
48+
- `MP` -- Dagger UpperTriangular array to associate precision with each tile
49+
- `tolerance` -- User defined tolerance as required aby the application
50+
51+
### Output
52+
The Dagger array shows the required precision of each tile
53+
54+
"""
55+
"""
56+
function adapt_precision(A::UpperTriangular{T,<:DArray{T,2}}, tolerance::Float64) where {T}
57+
58+
Ac = parent(A).chunks
59+
mt, nt = size(Ac)
60+
61+
global_norm = LinearAlgebra.norm2(A)
62+
63+
MP = fill("T", mt, nt)
64+
DMP = view(MP, Blocks(1, 1))
65+
MPc = parent(DMP).chunks
66+
67+
for n in range(1, nt)
68+
for m in range(1, n)
69+
if m == n
70+
MPc[m, n] = Dagger.@spawn tile_precision(
71+
UpperTriangular(Ac[m, n]),
72+
global_norm,
73+
max(mt, nt),
74+
tolerance)
75+
else
76+
MPc[m, n] = Dagger.@spawn tile_precision(
77+
Ac[m, n],
78+
global_norm,
79+
max(mt, nt),
80+
tolerance)
81+
end
82+
83+
end
84+
end
85+
86+
return UpperTriangular(collect(DMP))
87+
end
88+
"""
89+
90+
"""
91+
adapt_precision( A::LowerTriangular{T,<:DArray{T,2}},
92+
MP::LowerTriangular{String,<:DArray{String,2}}, tolerance::Float64) where {T}
93+
94+
it iterates over all tiles and calculates the required precision per tile based on formulation from Nicholas J. Higham
95+
96+
### Input
97+
- `A` -- Dagger LowerTriangular array of tiles with real values
98+
- `MP` -- Dagger LowerTriangular array to associate precision with each tile
99+
- `tolerance` -- User defined tolerance as required aby the application
100+
101+
### Output
102+
The Dagger array shows the required precision of each tile
103+
104+
"""
105+
106+
"""
107+
function adapt_precision(A::LowerTriangular{T,<:DArray{T,2}}, tolerance::T) where {T}
108+
109+
Ac = parent(A).chunks
110+
mt, nt = size(Ac)
111+
112+
global_norm = LinearAlgebra.norm2(A)
113+
114+
MP = fill("T", mt, nt)
115+
DMP = view(MP, Blocks(1, 1))
116+
MPc = parent(DMP).chunks
117+
118+
119+
for m in range(1, mt)
120+
for n in range(1, m)
121+
if m == n
122+
MPc[m, n] = Dagger.@spawn tile_precision(
123+
LowerTriangular(Ac[m, n]),
124+
global_norm,
125+
max(mt, nt),
126+
tolerance)
127+
else
128+
MPc[m, n] = Dagger.@spawn tile_precision(
129+
Ac[m, n],
130+
global_norm,
131+
max(mt, nt),
132+
tolerance)
133+
end
134+
135+
end
136+
end
137+
138+
return LowerTriangular(collect(DMP))
139+
end
140+
"""
141+
142+
"""
143+
adapt_precision(A::DArray{T,2}, MP::DArray{String,2}, tolerance::T) where {T}
144+
145+
it iterates over all tiles and calculates the required precision per tile based on formulation from Nicholas J. Higham
146+
147+
### Input
148+
- `A` -- Dagger array of tiles with real values
149+
- `MP` -- Dagger array to associate precision with each tile
150+
- `tolerance` -- User defined tolerance as required aby the application
151+
152+
### Output
153+
The Dagger array shows the required precision of each tile
154+
155+
"""
156+
157+
function adapt_precision(A::DArray{T,2}, tolerance::T) where {T}
158+
159+
Ac = parent(A).chunks
160+
mt, nt = size(Ac)
161+
162+
global_norm = LinearAlgebra.norm2(A)
163+
164+
MP = fill("T", mt, nt)
165+
DMP = view(MP, Blocks(1, 1))
166+
MPc = DMP.chunks
167+
168+
169+
for m in range(1, mt)
170+
for n in range(1, nt)
171+
MPc[m, n] =
172+
Dagger.@spawn tile_precision(
173+
Ac[m, n],
174+
global_norm,
175+
max(mt, nt),
176+
tolerance)
177+
end
178+
end
179+
180+
return collect(DMP)
181+
end

0 commit comments

Comments
 (0)