diff --git a/.github/workflows/openmp.yaml b/.github/workflows/openmp.yaml
new file mode 100644
index 0000000..9c22366
--- /dev/null
+++ b/.github/workflows/openmp.yaml
@@ -0,0 +1,41 @@
+name: OpenMP SIMD
+
+on: [push, pull_request]
+
+env:
+  # The CMake build type
+  BUILD_TYPE: Release
+
+jobs:
+  build:
+    strategy:
+      fail-fast: false
+      matrix:
+        float: [ON, OFF]
+
+    name: OpenMP SIMD Build - ${{ matrix.float == 'ON' && 'single' || 'double' }}
+    runs-on: ubuntu-latest
+
+    steps:
+      - name: Check out repository
+        uses: actions/checkout@v2
+        with:
+          submodules: 'recursive'
+
+      - name: Setup Environment
+        run: cmake -E make_directory ${{ runner.workspace }}/build
+
+      - name: Configure
+        shell: bash
+        working-directory: ${{ runner.workspace }}/build
+        run: cmake --warn-uninitialized -DCMAKE_BUILD_TYPE=$BUILD_TYPE -DQDLDL_UNITTESTS=ON -DQDLDL_FLOAT=${{ matrix.float }} -DCOVERAGE=OFF -DQDLDL_ENABLE_OPENMP_SIMD=ON $GITHUB_WORKSPACE
+
+      - name: Build
+        shell: bash
+        working-directory: ${{ runner.workspace }}/build
+        run: cmake --build . --config $BUILD_TYPE
+
+      - name: Run tests
+        shell: bash
+        working-directory: ${{ runner.workspace }}/build
+        run: ctest -C $BUILD_TYPE
diff --git a/CMakeLists.txt b/CMakeLists.txt
index 26db59c..402b84f 100644
--- a/CMakeLists.txt
+++ b/CMakeLists.txt
@@ -10,6 +10,7 @@ set(QDLDL_VERSION "${QDLDL_VERSION_MAJOR}.${QDLDL_VERSION_MINOR}.${QDLDL_VERSION
 project(qdldl VERSION ${QDLDL_VERSION})
 
 include( CMakeDependentOption )
+include(CheckCCompilerFlag)
 
 option( QDLDL_BUILD_STATIC_LIB "Build the static library" ON )
 option( QDLDL_BUILD_SHARED_LIB "Build the shared library" ON )
@@ -53,6 +54,11 @@ if( NOT (CMAKE_SIZEOF_VOID_P EQUAL 8) )
 endif()
 message(STATUS "Long integers (64bit) are ${QDLDL_LONG}")
 
+option( QDLDL_ENABLE_OPENMP_SIMD "Enable OpenMP SIMD annotations" OFF )
+if( ${QDLDL_ENABLE_OPENMP_SIMD} )
+  message( STATUS "Enabling OpenMP SIMD annotations" )
+endif()
+
 
 # Set Compiler flags
 # ----------------------------------------------
@@ -74,6 +80,30 @@ if (NOT MSVC)
     set(CMAKE_C_FLAGS_DEBUG "${CMAKE_C_FLAGS_DEBUG} -O0 -g")
 endif (NOT MSVC)
 
+# Add flags for OpenMP SIMD optimizations
+if( ${QDLDL_ENABLE_OPENMP_SIMD} )
+  # Flags for Clang/GCC, if supported
+  check_c_compiler_flag( "-fopenmp-simd" HAVE_FOPENMP_SIMD )
+  if( ${HAVE_FOPENMP_SIMD} )
+    set(CMAKE_C_FLAGS "${CMAKE_C_FLAGS} -fopenmp-simd")
+    set(CMAKE_C_FLAGS_DEBUG "${CMAKE_C_FLAGS_DEBUG} -fopenmp-simd")
+  endif()
+
+  # Flags for Intel linux compiler, if supported
+  check_c_compiler_flag( "-qopenmp-simd" HAVE_DASH_QOPENMP_SIMD )
+  if( ${HAVE_DASH_QOPENMP_SIMD} )
+    set(CMAKE_C_FLAGS "${CMAKE_C_FLAGS} -qopenmp-simd")
+    set(CMAKE_C_FLAGS_DEBUG "${CMAKE_C_FLAGS_DEBUG} -qopenmp-simd")
+  endif()
+
+  # Flags for Intel Windows compiler, if supported
+  check_c_compiler_flag( "/Qopenmp-simd" HAVE_SLASH_QOPENMP_SIMD )
+  if( ${HAVE_SLASH_QOPENMP_SIMD} )
+    set(CMAKE_C_FLAGS "${CMAKE_C_FLAGS} /Qopenmp-simd")
+    set(CMAKE_C_FLAGS_DEBUG "${CMAKE_C_FLAGS_DEBUG} /Qopenmp-simd")
+  endif()
+endif()
+
 # Generate header file with the global options
 # ---------------------------------------------
 
diff --git a/src/qdldl.c b/src/qdldl.c
index 6fb17d4..ca60671 100644
--- a/src/qdldl.c
+++ b/src/qdldl.c
@@ -242,6 +242,8 @@ void QDLDL_Lsolve(const QDLDL_int    n,
   QDLDL_int i,j;
   for(i = 0; i < n; i++){
     QDLDL_float val = x[i];
+
+    #pragma omp simd
     for(j = Lp[i]; j < Lp[i+1]; j++){
       x[Li[j]] -= Lx[j]*val;
     }
@@ -258,6 +260,8 @@ void QDLDL_Ltsolve(const QDLDL_int    n,
   QDLDL_int i,j;
   for(i = n-1; i>=0; i--){
     QDLDL_float val = x[i];
+
+    #pragma omp simd reduction(-:val)
     for(j = Lp[i]; j < Lp[i+1]; j++){
       val -= Lx[j]*x[Li[j]];
     }
@@ -276,6 +280,10 @@ void QDLDL_solve(const QDLDL_int       n,
   QDLDL_int i;
 
   QDLDL_Lsolve(n,Lp,Li,Lx,x);
-  for(i = 0; i < n; i++) x[i] *= Dinv[i];
+
+  #pragma omp simd
+  for(i = 0; i < n; i++)
+    x[i] *= Dinv[i];
+
   QDLDL_Ltsolve(n,Lp,Li,Lx,x);
 }