Linear Algebra operations: Jiang et al

src/SecureArithmetic.jl

@@ -32,5 +32,6 @@ include("auxiliary.jl")
 include("openfhe.jl")
 include("unencrypted.jl")
 include("arithmetic.jl")
+include("linear_algebra.jl")
 
 end # module SecureArithmetic

src/linear_algebra.jl

@@ -0,0 +1,56 @@
+"""
+Multiply two d x d matrices
+
+
+See https://eprint.iacr.org/2018/1041.pdf
+
+"""
+function mat_times_mat(A, B)
+
+n = length(A)
+d = Int(sqrt(n))
+shape = A.shape
+ctx = A.context
+
+# Flatten the matrices so that rotations act on flat vectors instead of being 2d rotations
+A = reshape(A, (n))
+B = reshape(B, (n))
+
+# Step 1-1 
+# Compute linear Transformation sigma on A
+
+
+A_0 = PlainArray(zeros(n), ctx)
+
+for k in -d+1:d-1
+    if k >= 0
+        u_k_sigma = PlainArray([0 <= l-d*k && l-d*k < d-k ? 1 : 0 for l in 0:n-1], ctx)
+    else 
+        u_k_sigma = PlainArray([-k <= l-(d+k)*d && l-(d+k)*d < d ? 1 : 0 for l in 0:n-1], ctx)
+    end
+    # Note that Rot(ct; l) in the paper is a leftshift, i.e. circshift(ct, -l)
+    A_0 += circshift(A, -(k)) * u_k_sigma
+end
+
+# Step 1-2
+# Compute linear Transformation tau on B
+
+
+B_0 = PlainArray(zeros(n), ctx)
+
+for k in 0:d-1
+    u_dk_tau = PlainArray([(l - k) / d in 0:d-1 ? 1 : 0 for  l in 0:n-1], ctx)
+    B_0 += circshift(B, -(d*k)) * u_dk_tau
+end
+
+# Step 2 and 3
+AB = A_0 * B_0
+for k in 1:d-1
+    v_k = PlainArray([0 <= (l % d) && (l % d) < (d-k) ? 1 : 0 for  l in 0:n-1], ctx)
+    v_k_d = PlainArray([(d-k) <= (l % d) && (l % d) < d ? 1 : 0 for  l in 0:n-1], ctx)
+    A_k = circshift(A_0, -(k)) * v_k + circshift(A_0, -(k-d)) * v_k_d
+    B_k = circshift(B_0, -(d*k))
+    AB += A_k * B_k
+end
+return reshape(AB, shape) 
+end

src/types.jl

@@ -183,6 +183,40 @@ See also: [`length`](@ref), [`SecureArray`](@ref), [`PlainArray`](@ref)
 """
 capacity(a::Union{PlainArray, SecureArray}) = a.capacity
 
+"""
+    reshape(a::SecureArray, shape)
+
+Return a `SecureArray` with the same data but a different shape.
+The new shape must have the same total number of elements as the original.
+
+See also: [`SecureArray`](@ref)
+"""
+function Base.reshape(a::SecureArray, shape::NTuple{M, Int}) where M
+    if prod(shape) != length(a)
+        throw(DimensionMismatch("new shape $(shape) is incompatible with array of length $(length(a))"))
+    end
+    SecureArray(a.data, shape, a.capacity, a.context)
+end
+
+Base.reshape(a::SecureArray, dims::Int...) = reshape(a, dims)
+
+"""
+    reshape(a::PlainArray, shape)
+
+Return a `PlainArray` with the same data but a different shape.
+The new shape must have the same total number of elements as the original.
+
+See also: [`PlainArray`](@ref)
+"""
+function Base.reshape(a::PlainArray, shape::NTuple{M, Int}) where M
+    if prod(shape) != length(a)
+        throw(DimensionMismatch("new shape $(shape) is incompatible with array of length $(length(a))"))
+    end
+    PlainArray(a.data, shape, a.capacity, a.context)
+end
+
+Base.reshape(a::PlainArray, dims::Int...) = reshape(a, dims)
+
 # Get wrapper name of a potentially parametric type
 # Copied from: https://github.com/ClapeyronThermo/Clapeyron.jl/blob/f40c282e2236ff68d91f37c39b5c1e4230ae9ef0/src/utils/core_utils.jl#L17
 # Original source: https://github.com/JuliaArrays/ArrayInterface.jl/blob/40d9a87be07ba323cca00f9e59e5285c13f7ee72/src/ArrayInterface.jl#L20

test/test_linear_algebra.jl

@@ -0,0 +1,77 @@
+module TestLinearAlgebra
+
+using Test
+using SecureArithmetic
+using OpenFHE
+
+function make_context()
+    level_budget = [4, 4]
+    levels_after_bootstrap = 10
+    depth = levels_after_bootstrap + GetBootstrapDepth(level_budget, UNIFORM_TERNARY)
+
+    parameters = CCParams{CryptoContextCKKSRNS}()
+    SetSecretKeyDist(parameters, UNIFORM_TERNARY)
+    SetSecurityLevel(parameters, HEStd_NotSet)
+    SetRingDim(parameters, 1 << 12)
+    SetScalingModSize(parameters, 59)
+    SetScalingTechnique(parameters, FLEXIBLEAUTO)
+    SetFirstModSize(parameters, 60)
+    SetMultiplicativeDepth(parameters, depth)
+
+    cc = GenCryptoContext(parameters)
+    Enable(cc, PKE)
+    Enable(cc, KEYSWITCH)
+    Enable(cc, LEVELEDSHE)
+    Enable(cc, ADVANCEDSHE)
+    Enable(cc, FHE)
+
+    EvalBootstrapSetup(cc; level_budget)
+
+    context = SecureContext(OpenFHEBackend(cc))
+    public_key, private_key = generate_keys(context)
+
+
+    return (context, public_key, private_key)
+end 
+
+@testset verbose=true showtiming=true "test_linear_algebra.jl" begin
+
+
+@testset verbose=true showtiming=true "square_mat_mat_nul" begin
+    (context, public_key, private_key) = make_context()
+
+    # SecureArithmetic stores matrices in column order, thus we need to transpose before eval mult
+    m1 = collect(transpose([0.25 0.5 0.75;
+         1.0 2.0 3.0;
+         4.0 5.0 6.0]))
+
+    m2 = collect(transpose([6.0 5.0 4.0;
+          3.0 2.0 1.0;
+          0.75 0.5 0.25]))
+
+    pm1 = PlainMatrix(m1, context)
+    pm2 = PlainMatrix(m2, context)
+
+
+    sm1 = encrypt(pm1, public_key)
+    sm2 = encrypt(pm2, public_key)
+
+    init_rotation!(context, private_key, (length(m1)), 1:length(m1)...)
+    init_multiplication!(context, private_key)
+
+
+
+    # the result is in row order. This we need to transpose to get the correct matrix in column order
+    result = transpose(collect(decrypt(SecureArithmetic.mat_times_mat(sm1, sm2), private_key)))
+    @test result ≈   [  3.5625   2.625   1.6875;
+                        14.25    10.5     6.75;
+                        43.5     33.0    22.5]
+
+end
+
+release_context_memory()
+GC.gc()
+
+end # @testset "test_serialization.jl"
+
+end # module