diff --git a/.github/workflows/build.yml b/.github/workflows/build.yml
index d011a250a..651a31a0c 100644
--- a/.github/workflows/build.yml
+++ b/.github/workflows/build.yml
@@ -57,6 +57,62 @@ jobs:
             exit 1
           fi
         shell: bash
+      - name: Run qobj
+        run: |
+          pip install -U qiskit-terra
+          python tools/generate_qobj.py
+          cd out
+          Release/qasm_simulator ../qobj.json | python ../tools/verify_standalone_results.py
+        shell: bash
+  mpi_standalone:
+    runs-on: ${{ matrix.os }}
+    strategy:
+      matrix:
+        os: ["ubuntu-latest"]
+    steps:
+      - uses: actions/checkout@v2
+      - name: Set up Python 3.7
+        uses: actions/setup-python@v2
+        with:
+          python-version: 3.7
+      - name: Install deps
+        run: pip install "conan>=1.31.2"
+      - name: Install openblas and mpi
+        run: |
+          set -e
+          sudo apt-get update
+          sudo apt-get install -y libopenblas-dev openmpi-bin libopenmpi-dev
+        shell: bash
+      - name: Compile Standalone
+        run: |
+          set -e
+          mkdir out; cd out; cmake .. -DBUILD_TESTS=1 -DAER_MPI=True
+          make
+        shell: bash
+      - name: Run Unit Tests with mpi
+        run: |
+          cd out/bin
+          for test in test*
+          do echo $test
+            if ! /usr/bin/mpirun.openmpi -host localhost:2 -np 2 ./$test
+            then
+              ERR=1
+            fi
+          done
+          if [ ! -z "$ERR" ]
+          then
+            exit 1
+          fi
+        shell: bash
+      - name: Run qobj
+        run: |
+          pip install -U qiskit-terra
+          python tools/generate_qobj.py
+          cd out
+          /usr/bin/mpirun.openmpi -host localhost:2 -np 2 Release/qasm_simulator ../qobj.json | python ../tools/verify_standalone_results.py
+        env:
+          USE_MPI: 1
+        shell: bash
   wheel:
     runs-on: ${{ matrix.os }}
     needs: ["standalone"]
diff --git a/tools/generate_qobj.py b/tools/generate_qobj.py
new file mode 100755
index 000000000..8fa856854
--- /dev/null
+++ b/tools/generate_qobj.py
@@ -0,0 +1,77 @@
+# -*- coding: utf-8 -*-
+
+# Copyright 2021, IBM.
+#
+# This source code is licensed under the Apache License, Version 2.0 found in
+# the LICENSE.txt file in the root directory of this source tree.
+
+import json
+import os
+
+from qiskit import ClassicalRegister
+from qiskit.compiler import assemble, transpile
+from qiskit import QuantumCircuit
+from qiskit import QuantumRegister
+
+
+def grovers_circuit(final_measure=True, allow_sampling=True):
+    """Testing a circuit originated in the Grover algorithm"""
+
+    circuits = []
+
+    # 6-qubit grovers
+    qr = QuantumRegister(6)
+    if final_measure:
+        cr = ClassicalRegister(2)
+        regs = (qr, cr)
+    else:
+        regs = (qr, )
+    circuit = QuantumCircuit(*regs)
+
+    circuit.h(qr[0])
+    circuit.h(qr[1])
+    circuit.x(qr[2])
+    circuit.x(qr[3])
+    circuit.x(qr[0])
+    circuit.cx(qr[0], qr[2])
+    circuit.x(qr[0])
+    circuit.cx(qr[1], qr[3])
+    circuit.ccx(qr[2], qr[3], qr[4])
+    circuit.cx(qr[1], qr[3])
+    circuit.x(qr[0])
+    circuit.cx(qr[0], qr[2])
+    circuit.x(qr[0])
+    circuit.x(qr[1])
+    circuit.x(qr[4])
+    circuit.h(qr[4])
+    circuit.ccx(qr[0], qr[1], qr[4])
+    circuit.h(qr[4])
+    circuit.x(qr[0])
+    circuit.x(qr[1])
+    circuit.x(qr[4])
+    circuit.h(qr[0])
+    circuit.h(qr[1])
+    circuit.h(qr[4])
+    if final_measure:
+        circuit.barrier(qr)
+        circuit.measure(qr[0], cr[0])
+        circuit.measure(qr[1], cr[1])
+    if not allow_sampling:
+        circuit.barrier(qr)
+        circuit.iden(qr)
+    circuits.append(circuit)
+
+    return circuits
+
+
+if __name__ == '__main__':
+    # Run qasm simulator
+    shots = 4000
+    circuits = grovers_circuit(final_measure=True, allow_sampling=True)
+    if os.getenv('USE_MPI', False):
+        qobj = assemble(transpile(circuits), shots=shots,
+                        blocking_enable=True, blocking_qubits=2)
+    else:
+        qobj = assemble(transpile(circuits), shots=shots)
+    with open('qobj.json', 'wt') as fp:
+        json.dump(qobj.to_dict(), fp)
diff --git a/tools/verify_standalone_results.py b/tools/verify_standalone_results.py
new file mode 100644
index 000000000..4c04c3c16
--- /dev/null
+++ b/tools/verify_standalone_results.py
@@ -0,0 +1,126 @@
+# -*- coding: utf-8 -*-
+
+# Copyright 2021, IBM.
+#
+# This source code is licensed under the Apache License, Version 2.0 found in
+# the LICENSE.txt file in the root directory of this source tree.
+
+import json
+import os
+import sys
+
+from qiskit.result import Result
+
+
+def assertDictAlmostEqual(dict1, dict2, delta=None, msg=None,
+                          places=None, default_value=0):
+    """Assert two dictionaries with numeric values are almost equal.
+
+    Fail if the two dictionaries are unequal as determined by
+    comparing that the difference between values with the same key are
+    not greater than delta (default 1e-8), or that difference rounded
+    to the given number of decimal places is not zero. If a key in one
+    dictionary is not in the other the default_value keyword argument
+    will be used for the missing value (default 0). If the two objects
+    compare equal then they will automatically compare almost equal.
+
+    Args:
+        dict1 (dict): a dictionary.
+        dict2 (dict): a dictionary.
+        delta (number): threshold for comparison (defaults to 1e-8).
+        msg (str): return a custom message on failure.
+        places (int): number of decimal places for comparison.
+        default_value (number): default value for missing keys.
+
+    Raises:
+        TypeError: raises TestCase failureException if the test fails.
+    """
+    if dict1 == dict2:
+        # Shortcut
+        return
+    if delta is not None and places is not None:
+        raise TypeError("specify delta or places not both")
+
+    if places is not None:
+        success = True
+        standard_msg = ''
+        # check value for keys in target
+        keys1 = set(dict1.keys())
+        for key in keys1:
+            val1 = dict1.get(key, default_value)
+            val2 = dict2.get(key, default_value)
+            if round(abs(val1 - val2), places) != 0:
+                success = False
+                standard_msg += '(%s: %s != %s), ' % (key,
+                                                      val1,
+                                                      val2)
+        # check values for keys in counts, not in target
+        keys2 = set(dict2.keys()) - keys1
+        for key in keys2:
+            val1 = dict1.get(key, default_value)
+            val2 = dict2.get(key, default_value)
+            if round(abs(val1 - val2), places) != 0:
+                success = False
+                standard_msg += '(%s: %s != %s), ' % (key,
+                                                      val1,
+                                                      val2)
+        if success is True:
+            return
+        standard_msg = standard_msg[:-2] + ' within %s places' % places
+
+    else:
+        if delta is None:
+            delta = 1e-8  # default delta value
+        success = True
+        standard_msg = ''
+        # check value for keys in target
+        keys1 = set(dict1.keys())
+        for key in keys1:
+            val1 = dict1.get(key, default_value)
+            val2 = dict2.get(key, default_value)
+            if abs(val1 - val2) > delta:
+                success = False
+                standard_msg += '(%s: %s != %s), ' % (key,
+                                                      val1,
+                                                      val2)
+        # check values for keys in counts, not in target
+        keys2 = set(dict2.keys()) - keys1
+        for key in keys2:
+            val1 = dict1.get(key, default_value)
+            val2 = dict2.get(key, default_value)
+            if abs(val1 - val2) > delta:
+                success = False
+                standard_msg += '(%s: %s != %s), ' % (key,
+                                                      val1,
+                                                      val2)
+        if success is True:
+            return
+        standard_msg = standard_msg[:-2] + ' within %s delta' % delta
+
+    raise Exception(standard_msg)
+
+
+def compare_counts(result, target, delta=0):
+    """Compare counts dictionary to targets."""
+    # Don't use get_counts method which converts hex
+    output = result.data(0)["counts"]
+    assertDictAlmostEqual(output, target, delta=delta)
+
+
+if __name__ == '__main__':
+    if len(sys.argv) == 2:
+        with open(sys.argv[1], 'rt') as fp:
+            result_dict = json.load(fp)
+    else:
+        result_dict = json.load(sys.stdin)
+
+    result = Result.from_dict(result_dict)
+    assert result.status == 'COMPLETED'
+    assert result.success is True
+    if os.getenv('USE_MPI', False):
+        assert result.metadata['num_mpi_processes'] > 1
+    shots = result.results[0].shots
+    targets = {'0x0': 5 * shots / 8, '0x1': shots / 8,
+               '0x2': shots / 8, '0x3': shots / 8}
+    compare_counts(result, targets, delta=0.05 * shots)
+    print("Input result JSON is valid!")