Graph.qs

namespace QBitGraph {

    open Microsoft.Quantum.Diagnostics;
    open Microsoft.Quantum.Intrinsic;
    open Microsoft.Quantum.Math;
    open Microsoft.Quantum.Arrays;

    operation QuantumGraphIsomorphism(adjMatrix1 : Bool[][], adjMatrix2 : Bool[][]) : Bool {
        // Check if both matrices have the same size
        let size1 = Length(adjMatrix1);
        let size2 = Length(adjMatrix2);
        if (size1 != size2) {
            fail "Matrices have different sizes.";
        }

        // Calculate the number of vertices and edges in the graph
        let nVertices = size1;
        let nEdges = CountEdges(adjMatrix1);

        // Allocate qubits to store the adjacency matrices and the intermediate results
        use qs = Qubit[2 * nVertices * nVertices + 2 * nVertices + 1];

        // Prepare the initial state
        PrepareInitialState(qs, adjMatrix1, adjMatrix2, nVertices);

        // Apply the isomorphism checking algorithm
        IsomorphismCheck(qs, nVertices, nEdges);

        // Measure the final result
        let result = M(qs[2 * nVertices * nVertices + 2 * nVertices]);

        // Reset all qubits
        ResetAll(qs);

        return result == One;
    }

    operation CountEdges(adjMatrix : Bool[][]) : Int {
        mutable nEdges = 0;
        for i in 0 .. Length(adjMatrix) - 1 {
            for j in i + 1 .. Length(adjMatrix) - 1 {
                if (adjMatrix[i][j]) {
                    set nEdges += 1;
                }
            }
        }
        return nEdges;
    }

    operation PrepareInitialState(qs : Qubit[], adjMatrix1 : Bool[][], adjMatrix2 : Bool[][], nVertices : Int) : Unit {
        // Encode the adjacency matrices as quantum states
        let state1 = EncodeAdjacencyMatrix(adjMatrix1);
        let state2 = EncodeAdjacencyMatrix(adjMatrix2);

        // Apply the Hadamard gate to the first set of qubits
        for i in 0 .. 2 * nVertices * nVertices - 1 {
            if (i < 2 * nVertices * nVertices / 2) {
                H(qs[i]);
            }
        }

        // Apply the controlled-not gate to the second set of qubits
        for i in 0 .. 2 * nVertices * nVertices - 1 {
            if (i >= 2 * nVertices * nVertices / 2) {
                let controlIndex = i - 2 * nVertices * nVertices / 2;
                let targetIndex = controlIndex + 2 * nVertices * nVertices;
                Controlled X([qs[controlIndex]], qs[targetIndex]);
            }
        }

        // Initialize the final qubit in the |1> state
        X(qs[2 * nVertices * nVertices + 2 * nVertices]);
    }

    operation EncodeAdjacencyMatrix(adjMatrix : Bool[][]) : Int[] {
        let nVertices = Length(adjMatrix);
        mutable state = Int[nVertices * nVertices];

        for i in 0 .. nVertices - 1 {
            for j in 0 .. nVertices - 1 {
                if (adjMatrix[i][j]) {
                    set state[(i * nVertices) + j] = 1;
                } else {
                    set state[(i * nVertices) + j] = 0;
                }
            }
        }

        return state;
    }

    operation IsomorphismCheck(qs : Qubit[], nVertices : Int, nEdges : Int) : Unit {
        // Apply the Grover's algorithm
        for i in 0 .. Floor(Sqrt(2.0 * PI() / 4.0 * Sqrt(2.0 * nVertices * nVertices))) - 1 {
            GroverIteration(qs, nVertices, nEdges);
        }
    }

    operation GroverIteration(qs : Qubit[], nVertices : Int, nEdges : Int) : Unit {
        // Apply the oracle
        Oracle(qs, nVertices, nEdges);

        // Apply the diffusion operator
        DiffusionOperator(qs, nVertices);
    }

    operation Oracle(qs : Qubit[], nVertices : Int, nEdges : Int) : Unit {
        // Apply the oracle to the first set of qubits
        for i in 0 .. 2 * nVertices * nVertices - 1 {
            if (i < 2 * nVertices * nVertices / 2) {
                Controlled X([qs[i]], qs[2 * nVertices * nVertices + 2 * nVertices]);
            }
        }

        // Apply the oracle to the second set of qubits
        for i in 0 .. 2 * nVertices * nVertices - 1 {
            if (i >= 2 * nVertices * nVertices / 2) {
                let controlIndex = i - 2 * nVertices * nVertices / 2;
                let targetIndex = controlIndex + 2 * nVertices * nVertices;
                Controlled X([qs[controlIndex]], qs[2 * nVertices * nVertices + 2 * nVertices]);
            }
        }

        // Apply the oracle to the final qubit
        for i in 0 .. 2 * nVertices * nVertices - 1 {
            if (i < 2 * nVertices * nVertices / 2) {
                Controlled X([qs[i]], qs[2 * nVertices * nVertices + 2 * nVertices]);
            }
        }
    }

    operation DiffusionOperator(qs : Qubit[], nVertices : Int) : Unit {
        // Apply the Hadamard gate to the first set of qubits
        for i in 0 .. 2 * nVertices * nVertices - 1 {
            if (i < 2 * nVertices * nVertices / 2) {
                H(qs[i]);
            }
        }

        // Apply the controlled-not gate to the second set of qubits
        for i in 0 .. 2 * nVertices * nVertices - 1 {
            if (i >= 2 * nVertices * nVertices / 2) {
                let controlIndex = i - 2 * nVertices * nVertices / 2;
                let targetIndex = controlIndex + 2 * nVertices * nVertices;
                Controlled X([qs[controlIndex]], qs[targetIndex]);
            }
        }

        // Apply the Hadamard gate to the first set of qubits
        for i in 0 .. 2 * nVertices * nVertices - 1 {
            if (i < 2 * nVertices * nVertices / 2) {
                H(qs[i]);
            }
        }
    }

    operation ResetAll(qs : Qubit[]) : Unit {
        for i in 0 .. Length(qs) - 1 {
            Reset(qs[i]);
        }
    }

    operation Length(adjMatrix : Bool[][]) : Int {
        return Length(adjMatrix[0]);
    }

    operation Size(adjMatrix1 : Bool[][], adjMatrix2 : Bool[][]) : Int {
        return Length(adjMatrix1);
    }

    operation PrintAdjacencyMatrix(adjMatrix : Bool[][]) : Unit {
        for i in 0 .. Length(adjMatrix) - 1 {
            for j in 0 .. Length(adjMatrix) - 1 {
                if (adjMatrix[i][j]) {
                    Message($"1 ");
                } else {
                    Message($"0 ");
                }
            }
            MessageLine("");
        }
    }

    operation PrintAdjacencyMatrix(adjMatrix : Int[]) : Unit {
        let nVertices = Int(Sqrt(Length(adjMatrix)));
        for i in 0 .. nVertices - 1 {
            for j in 0 .. nVertices - 1 {
                Message($"{adjMatrix[(i * nVertices) + j]} ");
            }
            MessageLine("");
        }
    }

    operation PrintAdjacencyMatrix(adjMatrix : Double[]) : Unit {
        let nVertices = Int(Sqrt(Length(adjMatrix)));
        for i in 0 .. nVertices - 1 {
            for j in 0 .. nVertices - 1 {
                Message($"{adjMatrix[(i * nVertices) + j]} ");
            }
            MessageLine("");
        }
    }

    operation PrintAdjacencyMatrix(adjMatrix : Complex[]) : Unit {
        let nVertices = Int(Sqrt(Length(adjMatrix)));
        for i in 0 .. nVertices - 1 {
            for j in 0 .. nVertices - 1 {
                Message($"{adjMatrix[(i * nVertices) + j]} ");
            }
            MessageLine("");
        }
    }

    operation PrintAdjacencyMatrix(adjMatrix : Qubit[]) : Unit {
        let nVertices = Int(Sqrt(Length(adjMatrix)));
        for i in 0 .. nVertices - 1 {
            for j in 0 .. nVertices - 1 {
                Message($"{adjMatrix[(i * nVertices) + j]} ");
            }
            MessageLine("");
        }
    }

    operation PrintAdjacencyMatrix(adjMatrix : Result[]) : Unit {
        let nVertices = Int(Sqrt(Length(adjMatrix)));
        for i in 0 .. nVertices - 1 {
            for j in 0 .. nVertices - 1 {
                Message($"{adjMatrix[(i * nVertices) + j]} ");
            }
            MessageLine("");
        }
    }

    operation PrintAdjacencyMatrix(adjMatrix : String[]) : Unit {
        let nVertices = Int(Sqrt(Length(adjMatrix)));
        for i in 0 .. nVertices - 1 {
            for j in 0 .. nVertices - 1 {
                Message($"{adjMatrix[(i * nVertices) + j]} ");
            }
            MessageLine("");
        }
    }

    operation PrintAdjacencyMatrix(adjMatrix : Int[][]) : Unit {
        for i in 0 .. Length(adjMatrix) - 1 {
            for j in 0 .. Length(adjMatrix) - 1 {
                Message($"{adjMatrix[i][j]} ");
            }
            MessageLine("");
        }
    }

    operation PrintAdjacencyMatrix(adjMatrix : Double[][]) : Unit {
        for i in 0 .. Length(adjMatrix) - 1 {
            for j in 0 .. Length(adjMatrix) - 1 {
                Message($"{adjMatrix[i][j]} ");
            }
            MessageLine("");
        }
    }

    operation PrintAdjacencyMatrix(adjMatrix : Complex[][]) : Unit {
        for i in 0 .. Length(adjMatrix) - 1 {
            for j in 0 .. Length(adjMatrix) - 1 {
                Message($"{adjMatrix[i][j]} ");
            }
            MessageLine("");
        }
    }
}