diff --git a/CHANGELOG.md b/CHANGELOG.md
new file mode 100644
--- /dev/null
+++ b/CHANGELOG.md
@@ -0,0 +1,5 @@
+# Revision history for hqcsim
+
+## 0.1.0.0 -- YYYY-mm-dd
+
+* First version. Released on an unsuspecting world.
diff --git a/LICENSE b/LICENSE
new file mode 100644
--- /dev/null
+++ b/LICENSE
@@ -0,0 +1,28 @@
+BSD 3-Clause License
+
+Copyright (c) 2024, Sebastian Mihai Ardelean
+
+Redistribution and use in source and binary forms, with or without
+modification, are permitted provided that the following conditions are met:
+
+1. Redistributions of source code must retain the above copyright notice, this
+   list of conditions and the following disclaimer.
+
+2. Redistributions in binary form must reproduce the above copyright notice,
+   this list of conditions and the following disclaimer in the documentation
+   and/or other materials provided with the distribution.
+
+3. Neither the name of the copyright holder nor the names of its
+   contributors may be used to endorse or promote products derived from
+   this software without specific prior written permission.
+
+THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
+AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
+IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
+DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE
+FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
+DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
+SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
+CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,
+OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
+OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
diff --git a/README.md b/README.md
new file mode 100644
--- /dev/null
+++ b/README.md
@@ -0,0 +1,18 @@
+# hqcsim
+
+## Getting started
+
+To get started with the library, clone the repo and then install Haskell Stack.
+
+Next, build the project:
+
+```
+$ cabal build
+
+```
+Documentation can be generated by running the next command:
+
+```
+$ cabal haddock
+
+```
diff --git a/app/Main.hs b/app/Main.hs
new file mode 100644
--- /dev/null
+++ b/app/Main.hs
@@ -0,0 +1,21 @@
+module Main where
+
+import Control.Monad (replicateM)
+import QC
+
+
+main :: IO ()
+main = do
+  let machine = Machine (makeQuantumState 2) 0
+  let qprog = QProgram [
+        (QInstruction hGate [0])
+        , (QInstruction hGate [1])
+        ]
+  machines <- replicateM 10 (runQProg qprog machine)
+  --mapM_ print machines
+  let measurementRegisters = map measurementRegister machines
+  print measurementRegisters
+ -- mapM_ (\a -> putStrLn $ show (measurementRegister a)) machines
+  --putStrLn $ "Measured state: " ++ show (measurementRegister result)
+
+
diff --git a/hqcsim.cabal b/hqcsim.cabal
new file mode 100644
--- /dev/null
+++ b/hqcsim.cabal
@@ -0,0 +1,56 @@
+cabal-version:      2.4
+name:               hqcsim
+version:            0.1.0.0
+
+-- A short (one-line) description of the package.
+synopsis: A library for simulating quantum circuits.
+
+-- A longer description of the package.
+description: A library used for simulating quantum circuits, based on the "A tutorial quantum interpreter in 150 lines of Lisp" from https://www.stylewarning.com/posts/quantum-interpreter/
+homepage: https://github.com/sebastianardelean/hqcsim#readme
+-- A URL where users can report bugs.
+-- bug-reports:
+
+-- The license under which the package is released.
+license: BSD-3-Clause
+license-file: LICENSE
+author:             Sebastian Mihai Ardelean
+maintainer:         ardeleanasm@gmail.com
+
+-- A copyright notice.
+copyright: 2024 Sebastian Mihai Ardelean
+category: Algorithms
+extra-source-files:    README.md
+                     , CHANGELOG.md
+
+library
+  hs-source-dirs:      src
+                       src/Quantum
+                       
+  exposed-modules:     QC
+                     , Quantum.QDataTypes
+                     , Quantum.QProgram
+                     , Quantum.Gates
+  build-depends:       base ^>= 4.17.2.1 
+                     , hmatrix
+                     , random
+  default-language:    Haskell2010
+                    
+executable hqcsim-exe
+    main-is:          Main.hs
+
+    -- Modules included in this executable, other than Main.
+    -- other-modules:
+
+    -- LANGUAGE extensions used by modules in this package.
+    -- other-extensions:
+    build-depends:    base ^>=4.17.2.1
+                    , hmatrix
+                    , random
+                    , hqcsim
+    hs-source-dirs:   app
+    default-language: Haskell2010
+
+source-repository head
+  type:     git
+  location: https://github.com/sebastianardelean/hqcsim
diff --git a/src/QC.hs b/src/QC.hs
new file mode 100644
--- /dev/null
+++ b/src/QC.hs
@@ -0,0 +1,20 @@
+{-|
+ -Module      : QC
+ -Description : Quantum Computing
+ -Copyright   : (c) Mihai Sebastian Ardelean, 2024
+ -License     : BSD3
+ -Maintainer  : ardeleanasm@gmail.com
+ -Portability : POSIX
+ 
+ This module is used to import needed modules for Quantum Computing.
+ 
+ -}
+module QC(
+  module Quantum.QDataTypes
+  , module Quantum.Gates
+  , module Quantum.QProgram
+         ) where
+
+import Quantum.QDataTypes
+import Quantum.Gates
+import Quantum.QProgram
diff --git a/src/Quantum/Gates.hs b/src/Quantum/Gates.hs
new file mode 100644
--- /dev/null
+++ b/src/Quantum/Gates.hs
@@ -0,0 +1,113 @@
+{-|
+ -Module      : Gates
+ -Description : Basic Quantum Gates
+ -Copyright   : (c) Mihai Sebastian Ardelean, 2024
+ -License     : BSD3
+ -Maintainer  : ardeleanasm@gmail.com
+ -Portability : POSIX
+ -}
+module Quantum.Gates
+  (
+    iGate
+  , swapGate
+  , hGate
+  , xGate
+  , yGate
+  , zGate
+  , cNotGate
+  , Gate
+  )where
+
+
+import qualified Numeric.LinearAlgebra as LA
+
+import Quantum.QDataTypes
+
+{-|
+ -  iGate function represent an Identity Matrix
+ 
+ >>>iGate
+ (2><2)
+ [ 1.0 :+ 0.0, 0.0 :+ 0.0
+ , 0.0 :+ 0.0, 1.0 :+ 0.0 ]
+ -}
+iGate :: Gate
+iGate = LA.ident 2 :: Gate
+
+{-|
+ -  swapGate function represent a Swap Gate
+ 
+ >>>swapGate
+ (4><4)
+ [ 1.0 :+ 0.0, 0.0 :+ 0.0, 0.0 :+ 0.0, 0.0 :+ 0.0
+ , 0.0 :+ 0.0, 0.0 :+ 0.0, 1.0 :+ 0.0, 0.0 :+ 0.0
+ , 0.0 :+ 0.0, 1.0 :+ 0.0, 0.0 :+ 0.0, 0.0 :+ 0.0
+ , 0.0 :+ 0.0, 0.0 :+ 0.0, 0.0 :+ 0.0, 1.0 :+ 0.0 ]
+ -}
+swapGate :: Gate
+swapGate = (4LA.><4)
+    [ 1, 0, 0, 0 
+    , 0, 0, 1, 0 
+    , 0, 1, 0, 0 
+    , 0, 0, 0, 1]::Gate
+
+{-|
+ -  hGate function represent a Hadamard Gate
+  
+  >>>hGate
+ (2><2)
+ [ 0.7071067811865475 :+ 0.0,    0.7071067811865475 :+ 0.0
+ , 0.7071067811865475 :+ 0.0, (-0.7071067811865475) :+ 0.0 ]
+ -}
+hGate :: Gate
+hGate = (2LA.><2) [1/sqrt 2,1/sqrt 2,1/sqrt 2,(-1)/sqrt 2] :: Gate
+
+{-|
+ -  yGate function represent a Pauli Y-Gate
+ 
+ >>>yGate
+(2><2)
+ [ 0.0 :+ 0.0, 0.0 :+ (-1.0)
+ , 0.0 :+ 1.0,    0.0 :+ 0.0 ]
+ -}
+yGate :: Gate
+yGate = (2LA.><2) [0.0,0.0LA.:+(-1.0),0.0LA.:+1.0,0.0] :: Gate
+
+{-|
+ -  zGate function represent a Pauli Z-Gate
+ 
+ >>>zGate
+ (2><2)
+ [ 1.0 :+ 0.0,       0.0 :+ 0.0
+ , 0.0 :+ 0.0, (-1.0) :+ (-0.0) ]
+ -}
+zGate :: Gate
+zGate = (2LA.><2) [1,0,0,-1] :: Gate
+
+{-|
+ -  xGate function represent a Pauli X-Gate
+ 
+ >>>xGate
+ (2><2)
+ [ 0.0 :+ 0.0, 1.0 :+ 0.0
+ , 1.0 :+ 0.0, 0.0 :+ 0.0 ]
+ -}
+xGate :: Gate
+xGate = (2LA.><2) [0,1,1,0] :: Gate
+
+{-|
+ -  cNotGate function represent a Controlled-Not Gate
+ 
+ >>>cNotGate
+ (4><4)
+ [ 1.0 :+ 0.0, 0.0 :+ 0.0, 0.0 :+ 0.0, 0.0 :+ 0.0
+ , 0.0 :+ 0.0, 1.0 :+ 0.0, 0.0 :+ 0.0, 0.0 :+ 0.0
+ , 0.0 :+ 0.0, 0.0 :+ 0.0, 0.0 :+ 0.0, 1.0 :+ 0.0
+ , 0.0 :+ 0.0, 0.0 :+ 0.0, 1.0 :+ 0.0, 0.0 :+ 0.0 ]
+ -}
+cNotGate :: Gate
+cNotGate = (4LA.><4)
+  [1,0,0,0,
+   0,1,0,0,
+   0,0,0,1
+  ,0,0,1,0] :: Gate
diff --git a/src/Quantum/QDataTypes.hs b/src/Quantum/QDataTypes.hs
new file mode 100644
--- /dev/null
+++ b/src/Quantum/QDataTypes.hs
@@ -0,0 +1,47 @@
+{-|
+ -Module      : QDataTypes
+ -Description : Definitions of datatypes for quantum state and quantum gates.
+ -Copyright   : (c) Mihai Sebastian Ardelean, 2024
+ -License     : BSD3
+ -Maintainer  : ardeleanasm@gmail.com
+ -Portability : POSIX
+ -}
+module Quantum.QDataTypes
+  (
+    State
+  , Gate
+  ) where
+
+import qualified Numeric.LinearAlgebra as LA
+
+-- | 
+-- The `State` type is an alias for a vector of complex numbers.
+--
+-- In the context of quantum computing, a `State` represents a quantum state
+-- as a column vector where each element is a complex number. The length of
+-- the vector is typically \(2^n\) for a system of `n` qubits.
+--
+-- For example, a `State` could be:
+-- 
+-- @
+--  [1 :+ 0, 0 :+ 0, 0 :+ 0, 0 :+ 0]  -- Represents |0⟩
+-- @
+type State = LA.Vector (LA.Complex Double)
+
+-- |
+-- The `Gate` type is an alias for a matrix of complex numbers.
+--
+-- In quantum computing, a `Gate` is a unitary matrix that represents a quantum operation
+-- applied to qubits. The matrix elements are complex numbers.
+--
+-- The `Gate` type is used to describe quantum gates in algorithms. For example:
+--
+-- @
+-- -- Represents a 2x2 Hadamard Gate
+-- hGate :: Gate
+-- hGate = (2LA.><2)
+--   [1/sqrt 2,1/sqrt 2,1/sqrt 2,(-1)/sqrt 2] :: Gate
+--
+-- @
+type Gate = LA.Matrix (LA.Complex Double)
+
diff --git a/src/Quantum/QProgram.hs b/src/Quantum/QProgram.hs
new file mode 100644
--- /dev/null
+++ b/src/Quantum/QProgram.hs
@@ -0,0 +1,203 @@
+{-|
+ -Module      : QProgram
+ -Description : Definitions of qprogram datatypes and of simulation function.
+ -Copyright   : (c) Mihai Sebastian Ardelean, 2024
+ -License     : BSD3
+ -Maintainer  : ardeleanasm@gmail.com
+ -Portability : POSIX
+ -}
+module Quantum.QProgram
+  (
+    runQProg
+  , makeQuantumState
+  , Machine(..)
+  , QInstruction(..)
+  , QProgram(..)
+  ) where
+
+import Quantum.QDataTypes
+import Quantum.Gates
+
+import Data.List (nub, foldl')
+
+import System.Random (randomRIO)
+
+import qualified Numeric.LinearAlgebra as LA 
+
+import Control.Monad (replicateM)
+
+{-|
+A `Machine` is defined by the quantum state and the measurement register.
+
+It has two fields:
+* `qstate` of type `State`
+* `measurementRegister` of type `Int`
+-}
+data Machine = Machine {
+    qstate :: State            -- ^ Quantum state.
+  , measurementRegister :: Int -- ^ Measurement register.
+  } deriving (Eq, Show)
+
+
+{-|
+A `QInstruction` is defined by the unitary transformation and by the
+qubits' index on which the transformation is applied.
+
+It has two fields:
+* `gateMatrix` of type `Gate` is the unitary matrix that defines the quantum gate.
+* `affectedQubits` of type `[Int]`
+-}
+data QInstruction = QInstruction {
+    gateMatrix ::Gate       -- ^ Quantum gate matrix.
+  , affectedQubits :: [Int] -- ^ List of qubits' index that are affected by the quantum gate.
+    } deriving (Eq,Show)
+
+{-|
+A `QProgram` is defined by the list quantum instructions:
+* `instructions` is of type `[QInstruction]`
+-}
+data QProgram = QProgram {
+  instructions :: [QInstruction] -- ^ List of program instructions.
+                   } deriving (Eq, Show)
+
+{-|
+ - makeQuantumState function initializes a quantum state of `n` qubits
+-}
+makeQuantumState :: Int -> State
+makeQuantumState n = LA.fromList $ 1 : replicate (2 ^ n - 1) 0
+
+
+dimensionQubits :: Int -> Int
+dimensionQubits size = floor $ logBase 2 (fromIntegral size)
+
+apply :: Gate -> State -> State
+apply = (LA.#>)
+
+compose :: Gate -> Gate -> Gate
+compose = (LA.<>)
+
+kroneckerMul :: Gate -> Gate -> Gate
+kroneckerMul a b = LA.kronecker a b
+
+kroneckerExp :: Gate -> Int -> Gate
+kroneckerExp gate n
+  | n < 1 = (1LA.><1) [1]::Gate
+  | n == 1 = gate
+  | otherwise = kroneckerMul (kroneckerExp gate (n - 1)) gate
+
+
+lift :: Gate -> Int -> Int -> Gate
+lift gate i n = liftResult
+  where
+    left = kroneckerExp iGate (n - i - (dimensionQubits $ LA.rows gate))
+    right = kroneckerExp iGate i
+    liftResult = kroneckerMul left $ kroneckerMul gate right
+
+perm2trans :: [Int] -> [(Int, Int)]
+perm2trans permutation = nub (concatMap processIndex [0..length permutation - 1])
+  where
+    -- Process each index to determine necessary swaps
+    updateSrc :: Int -> Int -> [Int] -> Int
+    updateSrc src dest lst
+      |src >= dest = src
+      |otherwise = updateSrc (lst !! src) dest lst
+    processIndex dest
+      | src < dest = [(src, dest)]
+      | src > dest = [(dest,src)]
+      | otherwise = []
+      where
+        originalSrc = permutation !! dest
+        src = updateSrc originalSrc dest permutation
+
+    
+trans2adj :: [(Int, Int)] -> [Int]
+trans2adj transpositions = concatMap expandConsecutive transpositions
+  where 
+    expandConsecutive :: (Int, Int) -> [Int]
+    expandConsecutive (x, y)
+      | y - x == 1 = [x]
+      | otherwise  = trans ++ reverse (init trans)
+      where
+        trans = [x..y-1]
+
+apply1Q :: State -> Gate -> Int -> State
+apply1Q s u qubit = q1State
+  where
+    liftedU = lift u qubit (dimensionQubits $ LA.size s)
+    q1State = apply liftedU s
+
+
+
+applyNQ :: State -> Gate -> [Int] -> State
+applyNQ s u qubits = qubitsNState
+  where
+    swap :: Int -> Int -> Gate
+    swap i n = lift swapGate i n
+
+    trans2op :: [Int] -> Int -> Gate
+    trans2op [] n = LA.ident (2^n)
+    trans2op (t:ts) n = foldl' compose (swap t n) (map (`swap` n) ts)
+
+    n = dimensionQubits $ LA.size s
+    u01 = lift u 0 n
+    fromSpace = reverse qubits ++ [i | i <- [0..n-1], i `notElem` qubits]
+    trans = perm2trans fromSpace
+    adj = trans2adj trans
+    toFrom = trans2op adj n
+    fromTo = trans2op (reverse adj) n
+    upq = compose toFrom (compose u01 fromTo)
+    qubitsNState = apply upq s
+
+applyGate :: State -> Gate -> [Int] -> State
+applyGate s u qubits
+  | qubitsLength == 1 = apply1Q s u (qubits !! 0)
+  | otherwise         = applyNQ s u qubits
+  where
+    qubitsLength = length qubits
+
+
+
+
+sample :: State -> IO Int
+sample s = do
+  r <- randomRIO (0.0, 1.0)
+  return $ sampleIndex (LA.toList s) r 0.0
+  where
+    sampleIndex :: [LA.Complex Double] -> Double -> Double -> Int
+    sampleIndex [] _ _ = error "Invalid state vector"
+    sampleIndex (c:cs) r accProb =
+      let prob = accProb + (LA.magnitude c) ** 2
+      in if r < prob
+         then 0
+         else 1 + sampleIndex cs r prob
+
+collapse :: State -> Int -> State
+collapse st i = LA.fromList collapsedState
+  where
+    stateLength = LA.size st
+    collapsedState = replicate i 0 ++ [1] ++ replicate (stateLength - i - 1) 0    
+
+observe :: Machine -> IO Machine
+observe machine = do
+  let state = qstate machine
+  i <- sample state
+  let newState = collapse state i
+  return machine { qstate = newState, measurementRegister = i}
+
+
+evolveState :: QInstruction -> Machine ->Machine
+evolveState (QInstruction gateMatrix affectedQubits) m = newMachine
+  where
+    newState = applyGate (qstate m) gateMatrix affectedQubits
+    newMachine = Machine newState (measurementRegister m)
+
+{-|
+ - runQProg function executes the program's instruction.
+-}
+runQProg :: QProgram -> Machine -> IO Machine
+runQProg qprog machine = runInstruction (instructions qprog) machine
+  where
+    runInstruction :: [QInstruction] -> Machine -> IO Machine
+    runInstruction [] m = observe m
+    runInstruction (x:xs) m = runInstruction xs (evolveState x m)
+