diff --git a/.gitignore b/.gitignore
index 33565c14f8..7c9decc636 100644
--- a/.gitignore
+++ b/.gitignore
@@ -86,6 +86,7 @@ _version.py
 
 # third party integrations
 simulators/
+apps/
 
 # macOS
 .DS_Store
diff --git a/docs/sphinx/api/languages/python_api.rst b/docs/sphinx/api/languages/python_api.rst
index 160d2f2e9f..fb9b2b715c 100644
--- a/docs/sphinx/api/languages/python_api.rst
+++ b/docs/sphinx/api/languages/python_api.rst
@@ -33,6 +33,7 @@ Program Construction
     .. automethod:: rz
     .. automethod:: r1
     .. automethod:: swap
+    .. automethod:: exp_pauli
     .. automethod:: mx
     .. automethod:: my
     .. automethod:: mz
diff --git a/docs/sphinx/conf.py b/docs/sphinx/conf.py
index f2298889f0..08d7d83d3e 100644
--- a/docs/sphinx/conf.py
+++ b/docs/sphinx/conf.py
@@ -169,6 +169,7 @@ def setup(app):
     ('cpp:identifier', 'mlir::ImplicitLocOpBuilder'),
     ('cpp:identifier', 'BinarySymplecticForm'),
     ('cpp:identifier', 'CountsDictionary'),
+    ('cpp:identifier', 'QuakeValueOrNumericType'),
     ('py:class', 'function'),
     ('py:class', 'type'),
     ('py:class', 'cudaq::spin_op'),
diff --git a/include/cudaq/Frontend/nvqpp/ASTBridge.h b/include/cudaq/Frontend/nvqpp/ASTBridge.h
index 822ab2a47e..62ca38c724 100644
--- a/include/cudaq/Frontend/nvqpp/ASTBridge.h
+++ b/include/cudaq/Frontend/nvqpp/ASTBridge.h
@@ -704,4 +704,16 @@ inline bool isCallOperator(clang::OverloadedOperatorKind kindValue) {
   return kindValue == clang::OverloadedOperatorKind::OO_Call;
 }
 
+// Is \p t of type `char *`?
+inline bool isCharPointerType(mlir::Type t) {
+  if (auto ptrTy = dyn_cast<cc::PointerType>(t)) {
+    mlir::Type eleTy = ptrTy.getElementType();
+    if (auto arrTy = dyn_cast<cc::ArrayType>(eleTy))
+      eleTy = arrTy.getElementType();
+    if (auto intTy = dyn_cast<mlir::IntegerType>(eleTy))
+      return intTy.getWidth() == 8;
+  }
+  return false;
+}
+
 } // namespace cudaq
diff --git a/include/cudaq/Optimizer/Dialect/CC/CCOps.td b/include/cudaq/Optimizer/Dialect/CC/CCOps.td
index ef5514065f..98dc2215c4 100644
--- a/include/cudaq/Optimizer/Dialect/CC/CCOps.td
+++ b/include/cudaq/Optimizer/Dialect/CC/CCOps.td
@@ -1369,4 +1369,22 @@ def cc_CallableClosureOp : CCOp<"callable_closure", [Pure]> {
   }];
 }
 
+def cc_CreateStringLiteralOp : CCOp<"string_literal"> {
+  let summary = "Create a constant string literal.";
+  let description = [{
+    This operation creates a ASCIIZ string literal value. It's argument is a
+    constant MLIR String Attribute. The literal will have a null character
+    appended automatically.
+
+    ```mlir
+      %0 = cc.string_literal "Quantum Computing" : !cc.ptr<!cc.array<i8 x 18>>
+    ```
+  }];
+
+  let arguments = (ins StrAttr:$stringLiteral);
+  let results = (outs cc_PointerType:$result);
+  let assemblyFormat = [{
+     $stringLiteral `:` qualified(type(results)) attr-dict
+  }];
+}
 #endif // CUDAQ_OPTIMIZER_DIALECT_CC_OPS
diff --git a/include/cudaq/Optimizer/Dialect/Quake/QuakeOps.td b/include/cudaq/Optimizer/Dialect/Quake/QuakeOps.td
index 8aa52e9024..882a993ec5 100644
--- a/include/cudaq/Optimizer/Dialect/Quake/QuakeOps.td
+++ b/include/cudaq/Optimizer/Dialect/Quake/QuakeOps.td
@@ -786,6 +786,23 @@ class TwoTargetOp<string mnemonic, list<Trait> traits = []> :
 // Quantum operators (gates)
 //===----------------------------------------------------------------------===//
 
+def ExpPauliOp : QuakeOp<"exp_pauli", []> {
+  let summary = "General Pauli tensor product rotation";
+  let description = [{
+    This operation affects a general Pauli tensor product rotation on 
+    the input qubits. The number of Pauli characters in the input Pauli word 
+    string must equal the number of qubits in the veq. Mathematically, this operation 
+    applies exp(i theta P) where P is a general Pauli tensor product. 
+  }];
+  
+  let arguments = (ins AnyFloat:$parameter, VeqType:$qubits, cc_PointerType:$pauli);
+  let results = (outs );
+
+  let assemblyFormat = [{
+    `(` $parameter `)` $qubits `,` $pauli `:` functional-type(operands, results) attr-dict
+  }];
+}
+
 def HOp : OneTargetOp<"h", [Hermitian]> {
   let summary = "Hadamard operation";
   let description = [{
@@ -815,7 +832,7 @@ def PhasedRxOp : QuakeOperator<"phased_rx",
     Matrix representation:
     ```
     PhasedRx(θ,φ) = |        cos(θ/2)        -iexp(-iφ) * sin(θ/2) |
-                    | -iexp(iφ)) * sin(θ/2)         cos(θ/2)       |
+                    |  -iexp(iφ) * sin(θ/2)         cos(θ/2)       |
     ```
 
     Circuit symbol:
diff --git a/lib/Frontend/nvqpp/ConvertDecl.cpp b/lib/Frontend/nvqpp/ConvertDecl.cpp
index b116383498..a1cfae0eb8 100644
--- a/lib/Frontend/nvqpp/ConvertDecl.cpp
+++ b/lib/Frontend/nvqpp/ConvertDecl.cpp
@@ -90,9 +90,10 @@ void QuakeBridgeVisitor::addArgumentSymbols(
       // Transform pass-by-value arguments to stack slots.
       auto loc = toLocation(argVal);
       auto parmTy = entryBlock->getArgument(index).getType();
-      if (isa<cc::CallableType, cc::StdvecType, cc::ArrayType, cc::StructType,
-              LLVM::LLVMStructType, FunctionType, quake::RefType,
-              quake::VeqType>(parmTy)) {
+      if (isa<FunctionType, cc::ArrayType, cc::CallableType, cc::PointerType,
+              cc::StdvecType, cc::StructType, LLVM::LLVMStructType,
+              quake::ControlType, quake::RefType, quake::VeqType,
+              quake::WireType>(parmTy)) {
         symbolTable.insert(name, entryBlock->getArgument(index));
       } else {
         auto stackSlot = builder.create<cc::AllocaOp>(loc, parmTy);
diff --git a/lib/Frontend/nvqpp/ConvertExpr.cpp b/lib/Frontend/nvqpp/ConvertExpr.cpp
index 19ae358392..d7e377af3a 100644
--- a/lib/Frontend/nvqpp/ConvertExpr.cpp
+++ b/lib/Frontend/nvqpp/ConvertExpr.cpp
@@ -1254,6 +1254,47 @@ bool QuakeBridgeVisitor::VisitCallExpr(clang::CallExpr *x) {
               isAdjoint = structTypeAsRecord->getName() == "adj";
       }
 
+    if (funcName.equals("exp_pauli")) {
+      assert(args.size() > 2);
+      SmallVector<Value> processedArgs;
+      auto addTheString = [&](Value v) {
+        // The C-string argument (char*) may be loaded by an lvalue to rvalue
+        // cast. Here, we must pass the pointer and not the first character's
+        // value.
+        if (isCharPointerType(v.getType())) {
+          processedArgs.push_back(v);
+        } else if (auto load = v.getDefiningOp<cudaq::cc::LoadOp>()) {
+          processedArgs.push_back(load.getPtrvalue());
+        } else {
+          reportClangError(x, mangler, "could not determine string argument");
+        }
+      };
+      if (args.size() == 3 && isa<quake::VeqType>(args[1].getType())) {
+        // Have f64, veq, string
+        processedArgs.push_back(args[0]);
+        processedArgs.push_back(args[1]);
+        addTheString(args[2]);
+      } else {
+        // should have f64, string, qubits...
+        // need f64, veq, string, so process here
+
+        // add f64 value
+        processedArgs.push_back(args[0]);
+
+        // concat the qubits to a veq
+        SmallVector<Value> quantumArgs;
+        for (std::size_t i = 2; i < args.size(); i++)
+          quantumArgs.push_back(args[i]);
+        processedArgs.push_back(builder.create<quake::ConcatOp>(
+            loc, quake::VeqType::get(builder.getContext(), quantumArgs.size()),
+            quantumArgs));
+        addTheString(args[1]);
+      }
+
+      builder.create<quake::ExpPauliOp>(loc, TypeRange{}, processedArgs);
+      return true;
+    }
+
     if (funcName.equals("mx") || funcName.equals("my") ||
         funcName.equals("mz")) {
       // Measurements always return a bool or a std::vector<bool>.
@@ -2140,7 +2181,7 @@ bool QuakeBridgeVisitor::VisitCXXConstructExpr(clang::CXXConstructExpr *x) {
 
     // TODO: remove this when we can handle ctors more generally.
     if (!ctor->isDefaultConstructor()) {
-      LLVM_DEBUG(llvm::dbgs() << "unhandled ctor:\n"; x->dump());
+      LLVM_DEBUG(llvm::dbgs() << ctorName << " - unhandled ctor:\n"; x->dump());
       TODO_loc(loc, "C++ ctor (not-default)");
     }
 
@@ -2206,8 +2247,10 @@ bool QuakeBridgeVisitor::VisitDeclRefExpr(clang::DeclRefExpr *x) {
 }
 
 bool QuakeBridgeVisitor::VisitStringLiteral(clang::StringLiteral *x) {
-  TODO_x(toLocation(x->getSourceRange()), x, mangler, "string literal");
-  return false;
+  auto strLitTy = cc::PointerType::get(cc::ArrayType::get(
+      builder.getContext(), builder.getI8Type(), x->getString().size() + 1));
+  return pushValue(builder.create<cc::CreateStringLiteralOp>(
+      toLocation(x), strLitTy, builder.getStringAttr(x->getString())));
 }
 
 } // namespace cudaq::details
diff --git a/lib/Frontend/nvqpp/ConvertType.cpp b/lib/Frontend/nvqpp/ConvertType.cpp
index 5b260aaaab..b9e1f14932 100644
--- a/lib/Frontend/nvqpp/ConvertType.cpp
+++ b/lib/Frontend/nvqpp/ConvertType.cpp
@@ -81,6 +81,7 @@ QuakeBridgeVisitor::findCallOperator(const clang::CXXRecordDecl *decl) {
 
 bool QuakeBridgeVisitor::TraverseRecordType(clang::RecordType *t) {
   auto *recDecl = t->getDecl();
+
   if (ignoredClass(recDecl))
     return true;
   auto reci = records.find(t);
@@ -311,7 +312,7 @@ bool QuakeBridgeVisitor::doSyntaxChecks(const clang::FunctionDecl *x) {
     // device kernels may take veq and/or ref arguments.
     if (isArithmeticType(t) || isArithmeticSequenceType(t) ||
         isQuantumType(t) || isKernelCallable(t) || isFunctionCallable(t) ||
-        isReferenceToCallableRecord(t, p))
+        isCharPointerType(t) || isReferenceToCallableRecord(t, p))
       continue;
     reportClangError(p, mangler, "kernel argument type not supported");
     return false;
diff --git a/lib/Optimizer/CodeGen/LowerToQIR.cpp b/lib/Optimizer/CodeGen/LowerToQIR.cpp
index 555bf58071..568c2bc545 100644
--- a/lib/Optimizer/CodeGen/LowerToQIR.cpp
+++ b/lib/Optimizer/CodeGen/LowerToQIR.cpp
@@ -316,16 +316,14 @@ class SubveqOpRewrite : public ConvertOpToLLVMPattern<quake::SubVeqOp> {
 };
 
 /// Lower the quake.reset op to QIR
-template <typename ResetOpType>
-class ResetRewrite : public ConvertOpToLLVMPattern<ResetOpType> {
+class ResetRewrite : public ConvertOpToLLVMPattern<quake::ResetOp> {
 public:
-  using Base = ConvertOpToLLVMPattern<ResetOpType>;
-  using Base::Base;
+  using ConvertOpToLLVMPattern::ConvertOpToLLVMPattern;
 
   LogicalResult
-  matchAndRewrite(ResetOpType instOp, typename Base::OpAdaptor adaptor,
+  matchAndRewrite(quake::ResetOp instOp, OpAdaptor adaptor,
                   ConversionPatternRewriter &rewriter) const override {
-    auto parentModule = instOp->template getParentOfType<ModuleOp>();
+    auto parentModule = instOp->getParentOfType<ModuleOp>();
     auto context = parentModule->getContext();
     std::string qirQisPrefix(cudaq::opt::QIRQISPrefix);
     std::string instName = instOp->getName().stripDialect().str();
@@ -348,6 +346,37 @@ class ResetRewrite : public ConvertOpToLLVMPattern<ResetOpType> {
   }
 };
 
+/// Lower exp_pauli(f64, veq, cc.string) to __quantum__qis__exp_pauli
+class ExpPauliRewrite : public ConvertOpToLLVMPattern<quake::ExpPauliOp> {
+public:
+  using ConvertOpToLLVMPattern::ConvertOpToLLVMPattern;
+
+  LogicalResult
+  matchAndRewrite(quake::ExpPauliOp instOp, OpAdaptor adaptor,
+                  ConversionPatternRewriter &rewriter) const override {
+    auto loc = instOp->getLoc();
+    auto parentModule = instOp->getParentOfType<ModuleOp>();
+    auto *context = rewriter.getContext();
+    std::string qirQisPrefix(cudaq::opt::QIRQISPrefix);
+    auto qirFunctionName = qirQisPrefix + "exp_pauli";
+    FlatSymbolRefAttr symbolRef = cudaq::opt::factory::createLLVMFunctionSymbol(
+        qirFunctionName, /*return type=*/LLVM::LLVMVoidType::get(context),
+        {rewriter.getF64Type(), cudaq::opt::getArrayType(context),
+         cudaq::opt::factory::getPointerType(context)},
+        parentModule);
+    SmallVector<Value> operands = adaptor.getOperands();
+    // Make sure to drop any length information from the type of the Pauli word.
+    auto pauliWord = operands.back();
+    operands.pop_back();
+    auto castedPauli = rewriter.create<LLVM::BitcastOp>(
+        loc, cudaq::opt::factory::getPointerType(context), pauliWord);
+    operands.push_back(castedPauli);
+    rewriter.replaceOpWithNewOp<LLVM::CallOp>(instOp, TypeRange{}, symbolRef,
+                                              operands);
+    return success();
+  }
+};
+
 /// Lower single target Quantum ops with no parameter to QIR:
 /// h, x, y, z, s, t
 template <typename OP>
@@ -1310,6 +1339,42 @@ class StdvecSizeOpPattern
   }
 };
 
+class CreateStringLiteralOpPattern
+    : public ConvertOpToLLVMPattern<cudaq::cc::CreateStringLiteralOp> {
+public:
+  using Base = ConvertOpToLLVMPattern<cudaq::cc::CreateStringLiteralOp>;
+  using Base::Base;
+
+  LogicalResult
+  matchAndRewrite(cudaq::cc::CreateStringLiteralOp stringLiteralOp,
+                  OpAdaptor adaptor,
+                  ConversionPatternRewriter &rewriter) const override {
+    auto loc = stringLiteralOp.getLoc();
+    auto parentModule = stringLiteralOp->getParentOfType<ModuleOp>();
+    StringRef stringLiteral = stringLiteralOp.getStringLiteral();
+
+    // Write to the module body
+    auto insertPoint = rewriter.saveInsertionPoint();
+    rewriter.setInsertionPointToStart(parentModule.getBody());
+
+    // Create the register name global
+    auto builder = cudaq::IRBuilder::atBlockEnd(parentModule.getBody());
+    auto slGlobal =
+        builder.genCStringLiteralAppendNul(loc, parentModule, stringLiteral);
+
+    // Shift back to the function
+    rewriter.restoreInsertionPoint(insertPoint);
+
+    // Get the string address
+    rewriter.replaceOpWithNewOp<LLVM::AddressOfOp>(
+        stringLiteralOp,
+        cudaq::opt::factory::getPointerType(slGlobal.getType()),
+        slGlobal.getSymName());
+
+    return success();
+  }
+};
+
 class StoreOpPattern : public ConvertOpToLLVMPattern<cudaq::cc::StoreOp> {
 public:
   using Base = ConvertOpToLLVMPattern<cudaq::cc::StoreOp>;
@@ -1420,25 +1485,26 @@ class QuakeToQIRRewrite : public cudaq::opt::QuakeToQIRBase<QuakeToQIRRewrite> {
 
     patterns.insert<GetVeqSizeOpRewrite, MxToMz, MyToMz, ReturnBitRewrite>(
         context);
-    patterns.insert<
-        AllocaOpRewrite, AllocaOpPattern, CallableClosureOpPattern,
-        CallableFuncOpPattern, CallCallableOpPattern, CastOpPattern,
-        ComputePtrOpPattern, ConcatOpRewrite, DeallocOpRewrite,
-        ExtractQubitOpRewrite, ExtractValueOpPattern, FuncToPtrOpPattern,
-        InsertValueOpPattern, InstantiateCallableOpPattern, LoadOpPattern,
-        OneTargetRewrite<quake::HOp>, OneTargetRewrite<quake::XOp>,
-        OneTargetRewrite<quake::YOp>, OneTargetRewrite<quake::ZOp>,
-        OneTargetRewrite<quake::SOp>, OneTargetRewrite<quake::TOp>,
-        OneTargetOneParamRewrite<quake::R1Op>,
-        OneTargetTwoParamRewrite<quake::PhasedRxOp>,
-        OneTargetOneParamRewrite<quake::RxOp>,
-        OneTargetOneParamRewrite<quake::RyOp>,
-        OneTargetOneParamRewrite<quake::RzOp>,
-        OneTargetTwoParamRewrite<quake::U2Op>,
-        OneTargetTwoParamRewrite<quake::U3Op>, ResetRewrite<quake::ResetOp>,
-        StdvecDataOpPattern, StdvecInitOpPattern, StdvecSizeOpPattern,
-        StoreOpPattern, SubveqOpRewrite, TwoTargetRewrite<quake::SwapOp>,
-        UndefOpPattern>(typeConverter);
+    patterns
+        .insert<AllocaOpRewrite, AllocaOpPattern, CallableClosureOpPattern,
+                CallableFuncOpPattern, CallCallableOpPattern, CastOpPattern,
+                ComputePtrOpPattern, ConcatOpRewrite, DeallocOpRewrite,
+                CreateStringLiteralOpPattern, ExtractQubitOpRewrite,
+                ExtractValueOpPattern, FuncToPtrOpPattern, InsertValueOpPattern,
+                InstantiateCallableOpPattern, LoadOpPattern, ExpPauliRewrite,
+                OneTargetRewrite<quake::HOp>, OneTargetRewrite<quake::XOp>,
+                OneTargetRewrite<quake::YOp>, OneTargetRewrite<quake::ZOp>,
+                OneTargetRewrite<quake::SOp>, OneTargetRewrite<quake::TOp>,
+                OneTargetOneParamRewrite<quake::R1Op>,
+                OneTargetTwoParamRewrite<quake::PhasedRxOp>,
+                OneTargetOneParamRewrite<quake::RxOp>,
+                OneTargetOneParamRewrite<quake::RyOp>,
+                OneTargetOneParamRewrite<quake::RzOp>,
+                OneTargetTwoParamRewrite<quake::U2Op>,
+                OneTargetTwoParamRewrite<quake::U3Op>, ResetRewrite,
+                StdvecDataOpPattern, StdvecInitOpPattern, StdvecSizeOpPattern,
+                StoreOpPattern, SubveqOpRewrite,
+                TwoTargetRewrite<quake::SwapOp>, UndefOpPattern>(typeConverter);
     patterns.insert<MeasureRewrite<quake::MzOp>>(typeConverter, measureCounter);
 
     target.addLegalDialect<LLVM::LLVMDialect>();
diff --git a/lib/Optimizer/Dialect/CC/CCTypes.cpp b/lib/Optimizer/Dialect/CC/CCTypes.cpp
index 3a2ea1e014..95fe5e1fd4 100644
--- a/lib/Optimizer/Dialect/CC/CCTypes.cpp
+++ b/lib/Optimizer/Dialect/CC/CCTypes.cpp
@@ -133,14 +133,14 @@ void cc::ArrayType::print(AsmPrinter &printer) const {
 #define GET_TYPEDEF_CLASSES
 #include "cudaq/Optimizer/Dialect/CC/CCTypes.cpp.inc"
 
+//===----------------------------------------------------------------------===//
+
 namespace cudaq {
 
 cc::CallableType cc::CallableType::getNoSignature(MLIRContext *ctx) {
   return CallableType::get(ctx, FunctionType::get(ctx, {}, {}));
 }
 
-//===----------------------------------------------------------------------===//
-
 void cc::CCDialect::registerTypes() {
   addTypes<ArrayType, CallableType, PointerType, StdvecType, StructType>();
 }
diff --git a/python/runtime/cudaq/builder/py_kernel_builder.cpp b/python/runtime/cudaq/builder/py_kernel_builder.cpp
index c5984beefe..8884f47e39 100644
--- a/python/runtime/cudaq/builder/py_kernel_builder.cpp
+++ b/python/runtime/cudaq/builder/py_kernel_builder.cpp
@@ -893,7 +893,24 @@ provided `function` will be applied within `self` at each iteration.
       .def("to_quake", &kernel_builder<>::to_quake, "See :func:`__str__`.")
       .def("__str__", &kernel_builder<>::to_quake,
            "Return the :class:`Kernel` as a string in its MLIR representation "
-           "using the Quake dialect.\n");
+           "using the Quake dialect.\n")
+      .def(
+          "exp_pauli",
+          [](kernel_builder<> &self, py::object theta, const QuakeValue &qubits,
+             const std::string &pauliWord) {
+            if (py::isinstance<py::float_>(theta))
+              self.exp_pauli(theta.cast<double>(), qubits, pauliWord);
+            else if (py::isinstance<QuakeValue>(theta))
+              self.exp_pauli(theta.cast<QuakeValue &>(), qubits, pauliWord);
+            else
+              throw std::runtime_error(
+                  "Invalid `theta` argument type. Must be a "
+                  "`float` or a `QuakeValue`.");
+          },
+          "Apply a general Pauli tensor product rotation, `exp(i theta P)`, on "
+          "the specified qubit register. The Pauli tensor product is provided "
+          "as a string, e.g. `XXYX` for a 4-qubit term. The angle parameter "
+          "can be provided as a concrete float or a `QuakeValue`.");
 }
 
 void bindBuilder(py::module &mod) {
diff --git a/python/runtime/cudaq/spin/py_spin_op.cpp b/python/runtime/cudaq/spin/py_spin_op.cpp
index 8154d85b38..d85f58dcff 100644
--- a/python/runtime/cudaq/spin/py_spin_op.cpp
+++ b/python/runtime/cudaq/spin/py_spin_op.cpp
@@ -79,6 +79,14 @@ void bindSpinOperator(py::module &mod) {
       .def(py::init(
                [](py::object o) { return fromOpenFermionQubitOperator(o); }),
            "Create from OpenFermion QubitOperator.")
+      .def(py::init<std::vector<double> &, std::size_t>(), py::arg("data"),
+           py::arg("num_qubits"),
+           "Construct a :class:`SpinOperator` from a list of numeric values. "
+           "The encoding is as follows: for each term, a list of doubles where "
+           "element `i` is a 3.0 for a `Y`, a 1.0 for a `X`, and a 2.0 for a "
+           "`Z` on qubit `i`, followed by the real and imaginary part of the "
+           "coefficient. Each set of term elements is appended to the one "
+           "list. The list is ended with the total number of terms.")
 
       /// @brief Bind the member functions.
       .def("get_raw_data", &cudaq::spin_op::get_raw_data,
diff --git a/python/tests/unittests/test_kernel_builder.py b/python/tests/unittests/test_kernel_builder.py
index f2e67cabbd..d36d16a156 100644
--- a/python/tests/unittests/test_kernel_builder.py
+++ b/python/tests/unittests/test_kernel_builder.py
@@ -281,3 +281,34 @@ def test_from_state():
     ss = cudaq.get_state(kernel)
     for i in range(4):
         assert np.isclose(ss[i], state[i], 1e-3)
+
+def test_exp_pauli():
+    cudaq.reset_target()
+    kernel = cudaq.make_kernel()
+    qubits = kernel.qalloc(4)
+    kernel.x(qubits[0])
+    kernel.x(qubits[1])
+    print(type(-.22))
+    kernel.exp_pauli(-.22, qubits, "XXXY")
+    print(kernel)
+    h2_data = [
+      3, 1, 1, 3, 0.0454063,  0,  2,  0, 0, 0, 0.17028,    0,
+      0, 0, 2, 0, -0.220041,  -0, 1,  3, 3, 1, 0.0454063,  0,
+      0, 0, 0, 0, -0.106477,  0,  0,  2, 0, 0, 0.17028,    0,
+      0, 0, 0, 2, -0.220041,  -0, 3,  3, 1, 1, -0.0454063, -0,
+      2, 2, 0, 0, 0.168336,   0,  2,  0, 2, 0, 0.1202,     0,
+      0, 2, 0, 2, 0.1202,     0,  2,  0, 0, 2, 0.165607,   0,
+      0, 2, 2, 0, 0.165607,   0,  0,  0, 2, 2, 0.174073,   0,
+      1, 1, 3, 3, -0.0454063, -0, 15
+    ]
+    h = cudaq.SpinOperator(h2_data, 4)
+    want_exp = cudaq.observe(kernel, h).expectation_z()
+    assert np.isclose(want_exp, -1.13, atol=1e-2)
+
+    kernel, theta = cudaq.make_kernel(float)
+    qubits = kernel.qalloc(4)
+    kernel.x(qubits[0])
+    kernel.x(qubits[1])
+    kernel.exp_pauli(theta, qubits, "XXXY")
+    want_exp = cudaq.observe(kernel, h, -.22).expectation_z()
+    assert np.isclose(want_exp, -1.13, atol=1e-2)
\ No newline at end of file
diff --git a/runtime/cudaq/builder/kernel_builder.cpp b/runtime/cudaq/builder/kernel_builder.cpp
index 3df248d359..676d9bc0b0 100644
--- a/runtime/cudaq/builder/kernel_builder.cpp
+++ b/runtime/cudaq/builder/kernel_builder.cpp
@@ -163,6 +163,39 @@ bool isArgStdVec(std::vector<QuakeValue> &args, std::size_t idx) {
   return args[idx].isStdVec();
 }
 
+void exp_pauli(ImplicitLocOpBuilder &builder, const QuakeValue &theta,
+               const std::vector<QuakeValue> &qubits,
+               const std::string &pauliWord) {
+  Value qubitsVal;
+  if (qubits.size() == 1)
+    qubitsVal = qubits.front().getValue();
+  else {
+    // we have a vector of quake value qubits, need to concat them
+    SmallVector<Value> values;
+    for (auto &v : qubits)
+      values.push_back(v.getValue());
+
+    qubitsVal = builder.create<quake::ConcatOp>(
+        quake::VeqType::get(builder.getContext(), qubits.size()), values);
+  }
+
+  auto thetaVal = theta.getValue();
+  if (!isa<quake::VeqType>(qubitsVal.getType()))
+    throw std::runtime_error(
+        "exp_pauli must take a QuakeValue of veq type as second argument.");
+  if (!thetaVal.getType().isIntOrFloat())
+    throw std::runtime_error("exp_pauli must take a QuakeValue of float/int "
+                             "type as first argument.");
+  cudaq::info("kernel_builder apply exp_pauli {}", pauliWord);
+
+  auto strLitTy = cc::PointerType::get(cc::ArrayType::get(
+      builder.getContext(), builder.getI8Type(), pauliWord.size() + 1));
+  Value stringLiteral = builder.create<cc::CreateStringLiteralOp>(
+      strLitTy, builder.getStringAttr(pauliWord));
+  SmallVector<Value> args{thetaVal, qubitsVal, stringLiteral};
+  builder.create<quake::ExpPauliOp>(TypeRange{}, args);
+}
+
 /// @brief Search the given `FuncOp` for all `CallOps` recursively.
 /// If found, see if the called function is in the current `ModuleOp`
 /// for this `kernel_builder`, if so do nothing. If it is not found,
diff --git a/runtime/cudaq/builder/kernel_builder.h b/runtime/cudaq/builder/kernel_builder.h
index 4314f790d8..4e181a8ac1 100644
--- a/runtime/cudaq/builder/kernel_builder.h
+++ b/runtime/cudaq/builder/kernel_builder.h
@@ -191,6 +191,10 @@ CUDAQ_DETAILS_MEASURE_DECLARATION(mx)
 CUDAQ_DETAILS_MEASURE_DECLARATION(my)
 CUDAQ_DETAILS_MEASURE_DECLARATION(mz)
 
+void exp_pauli(ImplicitLocOpBuilder &builder, const QuakeValue &theta,
+               const std::vector<QuakeValue> &qubits,
+               const std::string &pauliWord);
+
 void swap(ImplicitLocOpBuilder &builder, const std::vector<QuakeValue> &ctrls,
           const std::vector<QuakeValue> &targets, bool adjoint = false);
 
@@ -231,17 +235,17 @@ void control(ImplicitLocOpBuilder &builder, std::string &name,
 void adjoint(ImplicitLocOpBuilder &builder, std::string &name,
              std::string &quakeCode, std::vector<QuakeValue> &values);
 
-/// @brief Add a for loop that starts from the given `start` integer index, ends
-/// at the given `end` integer index, and applies the given `body` as a callable
-/// function. This callable function must take as input an index variable that
-/// can be used within the body.
+/// @brief Add a for loop that starts from the given `start` integer index,
+/// ends at the given `end` integer index, and applies the given `body` as a
+/// callable function. This callable function must take as input an index
+/// variable that can be used within the body.
 void forLoop(ImplicitLocOpBuilder &builder, std::size_t start, std::size_t end,
              std::function<void(QuakeValue &)> &body);
 
-/// @brief Add a for loop that starts from the given `start` integer index, ends
-/// at the given `end` index, and applies the given `body` as a
-/// callable function. This callable function must take as input an index
-/// variable that can be used within the body.
+/// @brief Add a for loop that starts from the given `start` integer index,
+/// ends at the given `end` index, and applies the given `body` as a callable
+/// function. This callable function must take as input an index variable that
+/// can be used within the body.
 void forLoop(ImplicitLocOpBuilder &builder, std::size_t start, QuakeValue &end,
              std::function<void(QuakeValue &)> &body);
 
@@ -541,6 +545,32 @@ class kernel_builder : public details::kernel_builder_base {
     details::c_if(*opBuilder, result, thenFunctor);
   }
 
+  /// @brief Apply a general pauli rotation, exp(i theta P),
+  /// takes a QuakeValue representing a register of qubits.
+  template <QuakeValueOrNumericType ParamT>
+  void exp_pauli(const ParamT &theta, const QuakeValue &qubits,
+                 const std::string &pauliWord) {
+    std::vector<QuakeValue> qubitValues{qubits};
+    if constexpr (std::is_floating_point_v<ParamT>)
+      details::exp_pauli(*opBuilder, QuakeValue(*opBuilder, theta), qubitValues,
+                         pauliWord);
+    else
+      details::exp_pauli(*opBuilder, theta, qubitValues, pauliWord);
+  }
+
+  /// @brief Apply a general pauli rotation, exp(i theta P),
+  /// takes a variadic list of QuakeValues representing a individual qubits.
+  template <QuakeValueOrNumericType ParamT, typename... QubitArgs>
+  void exp_pauli(const ParamT &theta, const std::string &pauliWord,
+                 QubitArgs &&...qubits) {
+    std::vector<QuakeValue> qubitValues{qubits...};
+    if constexpr (std::is_floating_point_v<ParamT>)
+      details::exp_pauli(*opBuilder, QuakeValue(*opBuilder, theta), qubitValues,
+                         pauliWord);
+    else
+      details::exp_pauli(*opBuilder, theta, qubitValues, pauliWord);
+  }
+
   /// @brief Apply the given `otherKernel` with the provided `QuakeValue`
   /// arguments.
   template <typename OtherKernelBuilder>
diff --git a/runtime/cudaq/domains/chemistry/uccsd.h b/runtime/cudaq/domains/chemistry/uccsd.h
index 29ce66338c..d3dc14729e 100644
--- a/runtime/cudaq/domains/chemistry/uccsd.h
+++ b/runtime/cudaq/domains/chemistry/uccsd.h
@@ -46,84 +46,16 @@ void singletExcitation(KernelBuilder &&kernel, QuakeValue &qubits,
                        QuakeValue &theta, const SingleIndices &indices) {
   auto r = indices.front();
   auto p = indices.back();
-
-  kernel.rx(-M_PI_2, qubits[r]);
-  kernel.h(qubits[p]);
-
-  std::vector<std::pair<std::size_t, std::size_t>> cnots;
-
-  for (std::size_t i = r; i < p; i++)
-    cnots.emplace_back(i, i + 1);
-
-  auto reversed = cnots;
-  std::reverse(reversed.begin(), reversed.end());
-
-  for (auto &[i, j] : cnots)
-    kernel.template x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  kernel.rz(0.5 * theta, qubits[p]);
-
-  for (auto &[i, j] : reversed)
-    kernel.template x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  kernel.rx(M_PI_2, qubits[r]);
-  kernel.h(qubits[p]);
-
-  kernel.h(qubits[r]);
-  kernel.rx(-M_PI_2, qubits[p]);
-
-  for (auto &[i, j] : cnots)
-    kernel.template x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  kernel.rz(-0.5 * theta, qubits[p]);
-
-  for (auto &[i, j] : reversed)
-    kernel.template x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  kernel.h(qubits[r]);
-  kernel.rx(M_PI_2, qubits[p]);
+  kernel.exp_pauli(0.5 * theta, "YX", qubits[r], qubits[p]);
+  kernel.exp_pauli(-0.5 * theta, "XY", qubits[r], qubits[p]);
 }
 
 __qpu__ void singletExcitation(cudaq::qspan<> qubits, double theta,
                                const SingleIndices &indices) {
   auto r = indices.front();
   auto p = indices.back();
-
-  rx(-M_PI_2, qubits[r]);
-  h(qubits[p]);
-
-  std::vector<std::pair<std::size_t, std::size_t>> cnots;
-
-  for (std::size_t i = r; i < p; i++)
-    cnots.emplace_back(i, i + 1);
-
-  auto reversed = cnots;
-  std::reverse(reversed.begin(), reversed.end());
-
-  for (auto &[i, j] : cnots)
-    x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  rz(theta / 2., qubits[p]);
-
-  for (auto &[i, j] : reversed)
-    x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  rx(M_PI_2, qubits[r]);
-  h(qubits[p]);
-
-  h(qubits[r]);
-  rx(-M_PI_2, qubits[p]);
-
-  for (auto &[i, j] : cnots)
-    x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  rz(-theta / 2., qubits[p]);
-
-  for (auto &[i, j] : reversed)
-    x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  h(qubits[r]);
-  rx(M_PI_2, qubits[p]);
+  exp_pauli(theta / 2., "YX", qubits[r], qubits[p]);
+  exp_pauli(-theta / 2., "XY", qubits[r], qubits[p]);
 }
 
 template <typename KernelBuilder>
@@ -135,168 +67,30 @@ void doubletExcitation(KernelBuilder &kernel, QuakeValue &qubits,
   auto q = d2.front();
   auto p = d2.back();
 
-  std::vector<std::pair<std::size_t, std::size_t>> cnots;
-  for (auto &i : cudaq::range(d1.size() - 1))
-    cnots.emplace_back(d1[i], d1[i + 1]);
-
-  cnots.emplace_back(r, q);
-
-  for (auto &i : cudaq::range(d2.size() - 1))
-    cnots.emplace_back(d2[i], d2[i + 1]);
-
-  auto reversed_cnots = cnots;
-  std::reverse(reversed_cnots.begin(), reversed_cnots.end());
-
-  kernel.h(qubits[s]);
-  kernel.h(qubits[r]);
-  kernel.rx(-M_PI_2, qubits[q]);
-  kernel.h(qubits[p]);
-
-  for (auto &[i, j] : cnots)
-    kernel.template x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  kernel.rz((1. / 8.0) * theta, qubits[p]);
-
-  for (auto &[i, j] : reversed_cnots)
-    kernel.template x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  kernel.h(qubits[s]);
-  kernel.h(qubits[r]);
-  kernel.rx(M_PI_2, qubits[q]);
-  kernel.h(qubits[p]);
-
+  // layer 1
+  kernel.exp_pauli((1. / 8.) * theta, "XXYX", qubits[s], qubits[r], qubits[q],
+                   qubits[p]);
   // layer 2
-  kernel.rx(-M_PI_2, qubits[s]);
-  kernel.h(qubits[r]);
-  kernel.rx(-M_PI_2, qubits[q]);
-  kernel.rx(-M_PI_2, qubits[p]);
-
-  for (auto &[i, j] : cnots)
-    kernel.template x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  kernel.rz((1. / 8.0) * theta, qubits[p]);
-
-  for (auto &[i, j] : reversed_cnots)
-    kernel.template x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  kernel.rx(M_PI_2, qubits[s]);
-  kernel.h(qubits[r]);
-  kernel.rx(M_PI_2, qubits[q]);
-  kernel.rx(M_PI_2, qubits[p]);
-
+  kernel.exp_pauli((1. / 8.) * theta, "YXYY", qubits[s], qubits[r], qubits[q],
+                   qubits[p]);
   // layer 3
-  kernel.h(qubits[s]);
-  kernel.rx(-M_PI_2, qubits[r]);
-  kernel.rx(-M_PI_2, qubits[q]);
-  kernel.rx(-M_PI_2, qubits[p]);
-
-  for (auto &[i, j] : cnots)
-    kernel.template x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  kernel.rz((1. / 8.0) * theta, qubits[p]);
-
-  for (auto &[i, j] : reversed_cnots)
-    kernel.template x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  kernel.h(qubits[s]);
-  kernel.rx(M_PI_2, qubits[r]);
-  kernel.rx(M_PI_2, qubits[q]);
-  kernel.rx(M_PI_2, qubits[p]);
-
+  kernel.exp_pauli((1. / 8.) * theta, "XYYY", qubits[s], qubits[r], qubits[q],
+                   qubits[p]);
   // layer 4
-  kernel.h(qubits[s]);
-  kernel.h(qubits[r]);
-  kernel.h(qubits[q]);
-  kernel.rx(-M_PI_2, qubits[p]);
-
-  for (auto &[i, j] : cnots)
-    kernel.template x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  kernel.rz((1. / 8.0) * theta, qubits[p]);
-
-  for (auto &[i, j] : reversed_cnots)
-    kernel.template x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  kernel.h(qubits[s]);
-  kernel.h(qubits[r]);
-  kernel.h(qubits[q]);
-  kernel.rx(M_PI_2, qubits[p]);
-
+  kernel.exp_pauli((1. / 8.) * theta, "XXXY", qubits[s], qubits[r], qubits[q],
+                   qubits[p]);
   // layer 5
-  kernel.rx(-M_PI_2, qubits[s]);
-  kernel.h(qubits[r]);
-  kernel.h(qubits[q]);
-  kernel.h(qubits[p]);
-
-  for (auto &[i, j] : cnots)
-    kernel.template x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  kernel.rz((-1. / 8.0) * theta, qubits[p]);
-
-  for (auto &[i, j] : reversed_cnots)
-    kernel.template x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  kernel.rx(M_PI_2, qubits[s]);
-  kernel.h(qubits[r]);
-  kernel.h(qubits[q]);
-  kernel.h(qubits[p]);
-
+  kernel.exp_pauli((-1. / 8.) * theta, "YXXX", qubits[s], qubits[r], qubits[q],
+                   qubits[p]);
   // layer 6
-  kernel.h(qubits[s]);
-  kernel.rx(-M_PI_2, qubits[r]);
-  kernel.h(qubits[q]);
-  kernel.h(qubits[p]);
-
-  for (auto &[i, j] : cnots)
-    kernel.template x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  kernel.rz((-1. / 8.0) * theta, qubits[p]);
-
-  for (auto &[i, j] : reversed_cnots)
-    kernel.template x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  kernel.h(qubits[s]);
-  kernel.rx(M_PI_2, qubits[r]);
-  kernel.h(qubits[q]);
-  kernel.h(qubits[p]);
-
+  kernel.exp_pauli((-1. / 8.) * theta, "XYXX", qubits[s], qubits[r], qubits[q],
+                   qubits[p]);
   // layer 7
-  kernel.rx(-M_PI_2, qubits[s]);
-  kernel.rx(-M_PI_2, qubits[r]);
-  kernel.rx(-M_PI_2, qubits[q]);
-  kernel.h(qubits[p]);
-
-  for (auto &[i, j] : cnots)
-    kernel.template x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  kernel.rz((-1. / 8.0) * theta, qubits[p]);
-
-  for (auto &[i, j] : reversed_cnots)
-    kernel.template x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  kernel.rx(M_PI_2, qubits[s]);
-  kernel.rx(M_PI_2, qubits[r]);
-  kernel.rx(M_PI_2, qubits[q]);
-  kernel.h(qubits[p]);
-
+  kernel.exp_pauli((-1. / 8.) * theta, "YYYX", qubits[s], qubits[r], qubits[q],
+                   qubits[p]);
   // layer 8
-  kernel.rx(-M_PI_2, qubits[s]);
-  kernel.rx(-M_PI_2, qubits[r]);
-  kernel.h(qubits[q]);
-  kernel.rx(-M_PI_2, qubits[p]);
-
-  for (auto &[i, j] : cnots)
-    kernel.template x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  kernel.rz((-1. / 8.0) * theta, qubits[p]);
-
-  for (auto &[i, j] : reversed_cnots)
-    kernel.template x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  kernel.rx(M_PI_2, qubits[s]);
-  kernel.rx(M_PI_2, qubits[r]);
-  kernel.h(qubits[q]);
-  kernel.rx(M_PI_2, qubits[p]);
+  kernel.exp_pauli((-1. / 8.) * theta, "YYXY", qubits[s], qubits[r], qubits[q],
+                   qubits[p]);
 }
 
 __qpu__ void doubletExcitation(cudaq::qspan<> qubits, double theta,
@@ -307,168 +101,30 @@ __qpu__ void doubletExcitation(cudaq::qspan<> qubits, double theta,
   auto q = d2.front();
   auto p = d2.back();
 
-  std::vector<std::pair<std::size_t, std::size_t>> cnots;
-  for (auto &i : cudaq::range(d1.size() - 1))
-    cnots.emplace_back(d1[i], d1[i + 1]);
-
-  cnots.emplace_back(r, q);
-
-  for (auto &i : cudaq::range(d2.size() - 1))
-    cnots.emplace_back(d2[i], d2[i + 1]);
-
-  auto reversed_cnots = cnots;
-  std::reverse(reversed_cnots.begin(), reversed_cnots.end());
-
-  h(qubits[s]);
-  h(qubits[r]);
-  rx(-M_PI_2, qubits[q]);
-  h(qubits[p]);
-
-  for (auto &[i, j] : cnots)
-    x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  rz(theta / 8., qubits[p]);
-
-  for (auto &[i, j] : reversed_cnots)
-    x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  h(qubits[s]);
-  h(qubits[r]);
-  rx(M_PI_2, qubits[q]);
-  h(qubits[p]);
-
+  // layer 1
+  exp_pauli((1. / 8.) * theta, "XXYX", qubits[s], qubits[r], qubits[q],
+            qubits[p]);
   // layer 2
-  rx(-M_PI_2, qubits[s]);
-  h(qubits[r]);
-  rx(-M_PI_2, qubits[q]);
-  rx(-M_PI_2, qubits[p]);
-
-  for (auto &[i, j] : cnots)
-    x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  rz(theta / 8., qubits[p]);
-
-  for (auto &[i, j] : reversed_cnots)
-    x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  rx(M_PI_2, qubits[s]);
-  h(qubits[r]);
-  rx(M_PI_2, qubits[q]);
-  rx(M_PI_2, qubits[p]);
-
+  exp_pauli((1. / 8.) * theta, "YXYY", qubits[s], qubits[r], qubits[q],
+            qubits[p]);
   // layer 3
-  h(qubits[s]);
-  rx(-M_PI_2, qubits[r]);
-  rx(-M_PI_2, qubits[q]);
-  rx(-M_PI_2, qubits[p]);
-
-  for (auto &[i, j] : cnots)
-    x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  rz(theta / 8., qubits[p]);
-
-  for (auto &[i, j] : reversed_cnots)
-    x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  h(qubits[s]);
-  rx(M_PI_2, qubits[r]);
-  rx(M_PI_2, qubits[q]);
-  rx(M_PI_2, qubits[p]);
-
+  exp_pauli((1. / 8.) * theta, "XYYY", qubits[s], qubits[r], qubits[q],
+            qubits[p]);
   // layer 4
-  h(qubits[s]);
-  h(qubits[r]);
-  h(qubits[q]);
-  rx(-M_PI_2, qubits[p]);
-
-  for (auto &[i, j] : cnots)
-    x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  rz(theta / 8., qubits[p]);
-
-  for (auto &[i, j] : reversed_cnots)
-    x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  h(qubits[s]);
-  h(qubits[r]);
-  h(qubits[q]);
-  rx(M_PI_2, qubits[p]);
-
+  exp_pauli((1. / 8.) * theta, "XXXY", qubits[s], qubits[r], qubits[q],
+            qubits[p]);
   // layer 5
-  rx(-M_PI_2, qubits[s]);
-  h(qubits[r]);
-  h(qubits[q]);
-  h(qubits[p]);
-
-  for (auto &[i, j] : cnots)
-    x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  rz(-theta / 8., qubits[p]);
-
-  for (auto &[i, j] : reversed_cnots)
-    x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  rx(M_PI_2, qubits[s]);
-  h(qubits[r]);
-  h(qubits[q]);
-  h(qubits[p]);
-
+  exp_pauli((-1. / 8.) * theta, "YXXX", qubits[s], qubits[r], qubits[q],
+            qubits[p]);
   // layer 6
-  h(qubits[s]);
-  rx(-M_PI_2, qubits[r]);
-  h(qubits[q]);
-  h(qubits[p]);
-
-  for (auto &[i, j] : cnots)
-    x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  rz(-theta / 8., qubits[p]);
-
-  for (auto &[i, j] : reversed_cnots)
-    x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  h(qubits[s]);
-  rx(M_PI_2, qubits[r]);
-  h(qubits[q]);
-  h(qubits[p]);
-
+  exp_pauli((-1. / 8.) * theta, "XYXX", qubits[s], qubits[r], qubits[q],
+            qubits[p]);
   // layer 7
-  rx(-M_PI_2, qubits[s]);
-  rx(-M_PI_2, qubits[r]);
-  rx(-M_PI_2, qubits[q]);
-  h(qubits[p]);
-
-  for (auto &[i, j] : cnots)
-    x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  rz(-theta / 8., qubits[p]);
-
-  for (auto &[i, j] : reversed_cnots)
-    x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  rx(M_PI_2, qubits[s]);
-  rx(M_PI_2, qubits[r]);
-  rx(M_PI_2, qubits[q]);
-  h(qubits[p]);
-
+  exp_pauli((-1. / 8.) * theta, "YYYX", qubits[s], qubits[r], qubits[q],
+            qubits[p]);
   // layer 8
-  rx(-M_PI_2, qubits[s]);
-  rx(-M_PI_2, qubits[r]);
-  h(qubits[q]);
-  rx(-M_PI_2, qubits[p]);
-
-  for (auto &[i, j] : cnots)
-    x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  rz(-theta / 8., qubits[p]);
-
-  for (auto &[i, j] : reversed_cnots)
-    x<cudaq::ctrl>(qubits[i], qubits[j]);
-
-  rx(M_PI_2, qubits[s]);
-  rx(M_PI_2, qubits[r]);
-  h(qubits[q]);
-  rx(M_PI_2, qubits[p]);
+  exp_pauli((-1. / 8.) * theta, "YYXY", qubits[s], qubits[r], qubits[q],
+            qubits[p]);
 }
 
 std::size_t uccsd_num_parameters(std::size_t nElectrons, std::size_t nQubits) {
diff --git a/runtime/cudaq/qis/execution_manager.h b/runtime/cudaq/qis/execution_manager.h
index daf223f728..79f6f56119 100644
--- a/runtime/cudaq/qis/execution_manager.h
+++ b/runtime/cudaq/qis/execution_manager.h
@@ -108,11 +108,13 @@ class ExecutionManager {
 
   /// Apply the quantum instruction with the given name, on the provided
   /// target qudits. Supports input of control qudits and rotational parameters.
+  /// Can also optionally take a spin_op as input to affect a general
+  /// Pauli rotation.
   virtual void apply(const std::string_view gateName,
                      const std::vector<double> &params,
                      const std::vector<QuditInfo> &controls,
                      const std::vector<QuditInfo> &targets,
-                     bool isAdjoint = false) = 0;
+                     bool isAdjoint = false, const spin_op op = spin_op()) = 0;
 
   /// Reset the qubit to the |0> state
   virtual void reset(const QuditInfo &target) = 0;
diff --git a/runtime/cudaq/qis/managers/BasicExecutionManager.h b/runtime/cudaq/qis/managers/BasicExecutionManager.h
index 8d79537ad9..89d36a4c35 100644
--- a/runtime/cudaq/qis/managers/BasicExecutionManager.h
+++ b/runtime/cudaq/qis/managers/BasicExecutionManager.h
@@ -34,11 +34,11 @@ class BasicExecutionManager : public cudaq::ExecutionManager {
 
 protected:
   /// @brief An instruction is composed of a operation name,
-  /// a optional set of rotation parameters, control qudits, and
-  /// target qudits.
-  using Instruction =
-      std::tuple<std::string, std::vector<double>,
-                 std::vector<cudaq::QuditInfo>, std::vector<cudaq::QuditInfo>>;
+  /// a optional set of rotation parameters, control qudits,
+  /// target qudits, and an optional spin_op.
+  using Instruction = std::tuple<std::string, std::vector<double>,
+                                 std::vector<cudaq::QuditInfo>,
+                                 std::vector<cudaq::QuditInfo>, spin_op>;
 
   /// @brief `typedef` for a queue of instructions
   using InstructionQueue = std::queue<Instruction>;
@@ -207,7 +207,7 @@ class BasicExecutionManager : public cudaq::ExecutionManager {
   void apply(const std::string_view gateName, const std::vector<double> &params,
              const std::vector<cudaq::QuditInfo> &controls,
              const std::vector<cudaq::QuditInfo> &targets,
-             bool isAdjoint = false) override {
+             bool isAdjoint = false, spin_op op = spin_op()) override {
 
     // Make a copy of the name that we can mutate if necessary
     std::string mutable_name(gateName);
@@ -240,21 +240,21 @@ class BasicExecutionManager : public cudaq::ExecutionManager {
     if (!adjointQueueStack.empty()) {
       // Add to the adjoint instruction queue
       adjointQueueStack.top().emplace(std::make_tuple(
-          mutable_name, mutable_params, mutable_controls, mutable_targets));
+          mutable_name, mutable_params, mutable_controls, mutable_targets, op));
       return;
     }
 
     // Add to the instruction queue
     instructionQueue.emplace(std::make_tuple(std::move(mutable_name),
                                              mutable_params, mutable_controls,
-                                             mutable_targets));
+                                             mutable_targets, op));
   }
 
   void synchronize() override {
     while (!instructionQueue.empty()) {
       auto instruction = instructionQueue.front();
       if (isInTracerMode()) {
-        auto [gateName, params, controls, targets] = instruction;
+        auto [gateName, params, controls, targets, op] = instruction;
         std::vector<std::size_t> controlIds;
         std::transform(controls.begin(), controls.end(),
                        std::back_inserter(controlIds),
diff --git a/runtime/cudaq/qis/managers/default/DefaultExecutionManager.cpp b/runtime/cudaq/qis/managers/default/DefaultExecutionManager.cpp
index 27f2d317cd..c27f6f09d2 100644
--- a/runtime/cudaq/qis/managers/default/DefaultExecutionManager.cpp
+++ b/runtime/cudaq/qis/managers/default/DefaultExecutionManager.cpp
@@ -102,7 +102,7 @@ class DefaultExecutionManager : public cudaq::BasicExecutionManager {
     flushRequestedAllocations();
 
     // Get the data, create the Qubit* targets
-    auto [gateName, parameters, controls, targets] = instruction;
+    auto [gateName, parameters, controls, targets, op] = instruction;
 
     // Map the Qudits to Qubits
     std::vector<std::size_t> localT;
@@ -139,6 +139,8 @@ class DefaultExecutionManager : public cudaq::BasicExecutionManager {
               })
         .Case("swap",
               [&]() { simulator()->swap(localC, localT[0], localT[1]); })
+        .Case("exp_pauli",
+              [&]() { simulator()->applyExpPauli(parameters[0], localT, op); })
         .Default([&]() {
           throw std::runtime_error("[DefaultExecutionManager] invalid gate "
                                    "application requested " +
diff --git a/runtime/cudaq/qis/qubit_qis.h b/runtime/cudaq/qis/qubit_qis.h
index e4625af83b..54d2e3feae 100644
--- a/runtime/cudaq/qis/qubit_qis.h
+++ b/runtime/cudaq/qis/qubit_qis.h
@@ -14,6 +14,7 @@
 #include "qarray.h"
 #include "qreg.h"
 #include "qvector.h"
+#include <cstring>
 #include <functional>
 
 #define __qpu__ __attribute__((annotate("quantum")))
@@ -276,6 +277,35 @@ inline void cs(qubit &q, qubit &r) { s<cudaq::ctrl>(q, r); }
 inline void ct(qubit &q, qubit &r) { t<cudaq::ctrl>(q, r); }
 inline void ccx(qubit &q, qubit &r, qubit &s) { x<cudaq::ctrl>(q, r, s); }
 
+/// @brief Apply a general Pauli rotation, takes a qubit register and the
+/// size must be equal to the pauli word length.
+template <typename QubitRange>
+  requires(std::ranges::range<QubitRange>)
+inline void exp_pauli(double theta, QubitRange &&qubits,
+                      const char *pauliWord) {
+  std::vector<QuditInfo> quditInfos;
+  std::transform(qubits.begin(), qubits.end(), std::back_inserter(quditInfos),
+                 [](auto &q) { return cudaq::qubitToQuditInfo(q); });
+  getExecutionManager()->apply("exp_pauli", {theta}, {}, quditInfos, false,
+                               spin_op::from_word(pauliWord));
+}
+
+/// @brief Apply a general Pauli rotation, takes a variadic set of
+/// qubits, and the number of qubits must be equal to the pauli word length.
+template <typename... QubitArgs>
+inline void exp_pauli(double theta, const char *pauliWord,
+                      QubitArgs &...qubits) {
+
+  if (sizeof...(QubitArgs) != std::strlen(pauliWord))
+    throw std::runtime_error(
+        "Invalid exp_pauli call, number of qubits != size of pauliWord.");
+
+  // Map the qubits to their unique ids and pack them into a std::array
+  std::vector<QuditInfo> quditInfos{qubitToQuditInfo(qubits)...};
+  getExecutionManager()->apply("exp_pauli", {theta}, {}, quditInfos, false,
+                               spin_op::from_word(pauliWord));
+}
+
 /// @brief Measure an individual qubit, return 0,1 as `bool`
 inline measure_result mz(qubit &q) {
   return getExecutionManager()->measure({q.n_levels(), q.id()});
diff --git a/runtime/nvqir/CircuitSimulator.h b/runtime/nvqir/CircuitSimulator.h
index c37ce74aa0..9c952e69e2 100644
--- a/runtime/nvqir/CircuitSimulator.h
+++ b/runtime/nvqir/CircuitSimulator.h
@@ -62,6 +62,15 @@ class CircuitSimulator {
     // do nothing
   }
 
+  /// @brief Apply exp(-i theta PauliTensorProd) to the underlying state.
+  /// This must be provided by subclasses.
+  virtual void applyExpPauli(double theta,
+                             const std::vector<std::size_t> &qubitIds,
+                             const cudaq::spin_op &op) {
+    throw std::runtime_error("CircuitSimulator::applyExpPauli not implemented, "
+                             "must be implemented by subclasses.");
+  }
+
   /// @brief Compute the expected value of the given spin op
   /// with respect to the current state, <psi | H | psi>.
   virtual cudaq::ExecutionResult observe(const cudaq::spin_op &term) = 0;
diff --git a/runtime/nvqir/NVQIR.cpp b/runtime/nvqir/NVQIR.cpp
index 311e648ee8..57f37790ce 100644
--- a/runtime/nvqir/NVQIR.cpp
+++ b/runtime/nvqir/NVQIR.cpp
@@ -363,6 +363,14 @@ Result *__quantum__qis__mz__to__register(Qubit *q, const char *name) {
   return b ? ResultOne : ResultZero;
 }
 
+void __quantum__qis__exp_pauli(double theta, Array *qubits, char *pauliWord) {
+  std::string pauliWordStr(pauliWord);
+  auto qubitsVec = arrayToVectorSizeT(qubits);
+  nvqir::getCircuitSimulatorInternal()->applyExpPauli(
+      theta, qubitsVec, cudaq::spin_op::from_word(pauliWordStr));
+  return;
+}
+
 void __quantum__rt__result_record_output(Result *, int8_t *) {}
 
 /// @brief Map an Array pointer containing Paulis to a vector of Paulis.
diff --git a/runtime/nvqir/custatevec/CuStateVecCircuitSimulator.cu b/runtime/nvqir/custatevec/CuStateVecCircuitSimulator.cu
index 16e963da24..d18f690a05 100644
--- a/runtime/nvqir/custatevec/CuStateVecCircuitSimulator.cu
+++ b/runtime/nvqir/custatevec/CuStateVecCircuitSimulator.cu
@@ -363,6 +363,34 @@ public:
     }
   }
 
+  /// @brief Override base class functionality for a general Pauli
+  /// rotation to delegate to the performant custatevecApplyPauliRotation.
+  void applyExpPauli(double theta, const std::vector<std::size_t> &qubits,
+                     const cudaq::spin_op &op) override {
+    flushGateQueue();
+    cudaq::info(" [cusv decomposing] exp_pauli({}, {})", theta,
+                op.to_string(false));
+    std::vector<int> controls, targets;
+    std::vector<custatevecPauli_t> paulis;
+    op.for_each_pauli([&](cudaq::pauli p, std::size_t i) {
+      if (p == cudaq::pauli::I)
+        paulis.push_back(custatevecPauli_t::CUSTATEVEC_PAULI_I);
+      else if (p == cudaq::pauli::X)
+        paulis.push_back(custatevecPauli_t::CUSTATEVEC_PAULI_X);
+      else if (p == cudaq::pauli::Y)
+        paulis.push_back(custatevecPauli_t::CUSTATEVEC_PAULI_Y);
+      else
+        paulis.push_back(custatevecPauli_t::CUSTATEVEC_PAULI_Z);
+
+      targets.push_back(qubits[i]);
+    });
+
+    HANDLE_ERROR(custatevecApplyPauliRotation(
+        handle, deviceStateVector, cuStateVecCudaDataType, nQubitsAllocated,
+        theta, paulis.data(), targets.data(), targets.size(), controls.data(),
+        nullptr, controls.size()));
+  }
+
   /// @brief Compute the operator expectation value, with respect to
   /// the current state vector, directly on GPU with the
   /// given the operator matrix and target qubit indices.
diff --git a/runtime/nvqir/qpp/QppCircuitSimulator.cpp b/runtime/nvqir/qpp/QppCircuitSimulator.cpp
index bf3b3e225f..1ed0669aad 100644
--- a/runtime/nvqir/qpp/QppCircuitSimulator.cpp
+++ b/runtime/nvqir/qpp/QppCircuitSimulator.cpp
@@ -181,6 +181,49 @@ class QppCircuitSimulator : public nvqir::CircuitSimulatorBase<double> {
     qpp::RandomDevices::get_instance().get_prng().seed(seed);
   }
 
+  void applyExpPauli(double theta, const std::vector<std::size_t> &qubitIds,
+                     const cudaq::spin_op &op) override {
+    flushGateQueue();
+    cudaq::info(" [qpp decomposing] exp_pauli({}, {})", theta,
+                op.to_string(false));
+    std::vector<std::size_t> qubitSupport;
+    std::vector<std::function<void(bool)>> basisChange;
+    op.for_each_pauli([&](cudaq::pauli type, std::size_t qubitIdx) {
+      if (type != cudaq::pauli::I)
+        qubitSupport.push_back(qubitIds[qubitIdx]);
+
+      if (type == cudaq::pauli::Y)
+        basisChange.emplace_back([&, qubitIdx](bool reverse) {
+          rx(!reverse ? M_PI_2 : -M_PI_2, qubitIds[qubitIdx]);
+        });
+      else if (type == cudaq::pauli::X)
+        basisChange.emplace_back(
+            [&, qubitIdx](bool) { h(qubitIds[qubitIdx]); });
+    });
+
+    if (!basisChange.empty())
+      for (auto &basis : basisChange)
+        basis(false);
+
+    std::vector<std::pair<std::size_t, std::size_t>> toReverse;
+    for (std::size_t i = 0; i < qubitSupport.size() - 1; i++) {
+      x({qubitSupport[i]}, qubitSupport[i + 1]);
+      toReverse.emplace_back(qubitSupport[i], qubitSupport[i + 1]);
+    }
+
+    rz(theta, qubitSupport.back());
+
+    std::reverse(toReverse.begin(), toReverse.end());
+    for (auto &[i, j] : toReverse)
+      x({i}, j);
+
+    if (!basisChange.empty()) {
+      std::reverse(basisChange.begin(), basisChange.end());
+      for (auto &basis : basisChange)
+        basis(true);
+    }
+  }
+
   bool canHandleObserve() override {
     // Do not compute <H> from matrix if shots based sampling requested
     if (executionContext &&
diff --git a/test/AST-Quake/exp-pauli-1.cpp b/test/AST-Quake/exp-pauli-1.cpp
new file mode 100644
index 0000000000..b1e17236b1
--- /dev/null
+++ b/test/AST-Quake/exp-pauli-1.cpp
@@ -0,0 +1,35 @@
+/*******************************************************************************
+ * Copyright (c) 2022 - 2023 NVIDIA Corporation & Affiliates.                  *
+ * All rights reserved.                                                        *
+ *                                                                             *
+ * This source code and the accompanying materials are made available under    *
+ * the terms of the Apache License 2.0 which accompanies this distribution.    *
+ ******************************************************************************/
+
+// RUN: cudaq-quake %s | cudaq-opt | FileCheck %s
+
+#include <cudaq.h>
+
+int main() {
+  auto kernel = [](double theta) __qpu__ {
+    cudaq::qreg q(4);
+    x(q[0]);
+    x(q[1]);
+    exp_pauli(theta, q, "XXXY");
+  };
+}
+
+// CHECK-LABEL:   func.func @__nvqpp__mlirgen__Z4mainE3$_0(
+// CHECK-SAME:                                             %[[VAL_0:.*]]: f64) attributes {"cudaq-entrypoint", "cudaq-kernel"} {
+// CHECK:           %[[VAL_1:.*]] = cc.alloca f64
+// CHECK:           cc.store %[[VAL_0]], %[[VAL_1]] : !cc.ptr<f64>
+// CHECK:           %[[VAL_2:.*]] = quake.alloca !quake.veq<4>
+// CHECK:           %[[VAL_3:.*]] = quake.extract_ref %[[VAL_2]][0] : (!quake.veq<4>) -> !quake.ref
+// CHECK:           quake.x %[[VAL_3]] : (!quake.ref) -> ()
+// CHECK:           %[[VAL_4:.*]] = quake.extract_ref %[[VAL_2]][1] : (!quake.veq<4>) -> !quake.ref
+// CHECK:           quake.x %[[VAL_4]] : (!quake.ref) -> ()
+// CHECK:           %[[VAL_5:.*]] = cc.load %[[VAL_1]] : !cc.ptr<f64>
+// CHECK:           %[[VAL_6:.*]] = cc.string_literal "XXXY" : !cc.ptr<!cc.array<i8 x 5>>
+// CHECK:           quake.exp_pauli(%[[VAL_5]]) %[[VAL_2]], %[[VAL_6]] : (f64, !quake.veq<4>, !cc.ptr<!cc.array<i8 x 5>>) -> ()
+// CHECK:           return
+// CHECK:         }
diff --git a/test/AST-Quake/exp-pauli-2.cpp b/test/AST-Quake/exp-pauli-2.cpp
new file mode 100644
index 0000000000..346e8d3664
--- /dev/null
+++ b/test/AST-Quake/exp-pauli-2.cpp
@@ -0,0 +1,36 @@
+/*******************************************************************************
+ * Copyright (c) 2022 - 2023 NVIDIA Corporation & Affiliates.                  *
+ * All rights reserved.                                                        *
+ *                                                                             *
+ * This source code and the accompanying materials are made available under    *
+ * the terms of the Apache License 2.0 which accompanies this distribution.    *
+ ******************************************************************************/
+
+// RUN: cudaq-quake %s | cudaq-opt | FileCheck %s
+
+#include <cudaq.h>
+
+int main() {
+  auto kernel2 = [](double theta, const char *pauli) __qpu__ {
+    cudaq::qreg q(4);
+    x(q[0]);
+    x(q[1]);
+    exp_pauli(theta, q, pauli);
+  };
+}
+
+// CHECK-LABEL:   func.func @__nvqpp__mlirgen__Z4mainE3$_0(
+// CHECK-SAME:                                             %[[VAL_0:.*]]: f64,
+// CHECK-SAME:                                             %[[VAL_1:.*]]: !cc.ptr<i8>) attributes {"cudaq-entrypoint", "cudaq-kernel"} {
+// CHECK:           %[[VAL_2:.*]] = cc.alloca f64
+// CHECK:           cc.store %[[VAL_0]], %[[VAL_2]] : !cc.ptr<f64>
+// CHECK:           %[[VAL_3:.*]] = quake.alloca !quake.veq<4>
+// CHECK:           %[[VAL_4:.*]] = quake.extract_ref %[[VAL_3]][0] : (!quake.veq<4>) -> !quake.ref
+// CHECK:           quake.x %[[VAL_4]] : (!quake.ref) -> ()
+// CHECK:           %[[VAL_5:.*]] = quake.extract_ref %[[VAL_3]][1] : (!quake.veq<4>) -> !quake.ref
+// CHECK:           quake.x %[[VAL_5]] : (!quake.ref) -> ()
+// CHECK:           %[[VAL_6:.*]] = cc.load %[[VAL_2]] : !cc.ptr<f64>
+// CHECK:           quake.exp_pauli(%[[VAL_6]]) %[[VAL_3]], %[[VAL_1]] : (f64, !quake.veq<4>, !cc.ptr<i8>) -> ()
+// CHECK:           return
+// CHECK:         }
+
diff --git a/test/AST-Quake/exp-pauli-3.cpp b/test/AST-Quake/exp-pauli-3.cpp
new file mode 100644
index 0000000000..aa5aee93d0
--- /dev/null
+++ b/test/AST-Quake/exp-pauli-3.cpp
@@ -0,0 +1,40 @@
+/*******************************************************************************
+ * Copyright (c) 2022 - 2023 NVIDIA Corporation & Affiliates.                  *
+ * All rights reserved.                                                        *
+ *                                                                             *
+ * This source code and the accompanying materials are made available under    *
+ * the terms of the Apache License 2.0 which accompanies this distribution.    *
+ ******************************************************************************/
+
+// RUN: cudaq-quake %s | cudaq-opt | FileCheck %s
+
+#include <cudaq.h>
+
+int main() {
+  auto kernel = [](double theta) __qpu__ {
+    cudaq::qreg q(4);
+    x(q[0]);
+    x(q[1]);
+    exp_pauli(theta, "XXXY", q[0], q[1], q[2], q[3]);
+  };
+}
+
+// CHECK-LABEL:   func.func @__nvqpp__mlirgen__Z4mainE3$_0(
+// CHECK-SAME:         %[[VAL_0:.*]]: f64) attributes
+// CHECK:           %[[VAL_1:.*]] = cc.alloca f64
+// CHECK:           cc.store %[[VAL_0]], %[[VAL_1]] : !cc.ptr<f64>
+// CHECK:           %[[VAL_2:.*]] = quake.alloca !quake.veq<4>
+// CHECK:           %[[VAL_3:.*]] = quake.extract_ref %[[VAL_2]][0] : (!quake.veq<4>) -> !quake.ref
+// CHECK:           quake.x %[[VAL_3]] : (!quake.ref) -> ()
+// CHECK:           %[[VAL_4:.*]] = quake.extract_ref %[[VAL_2]][1] : (!quake.veq<4>) -> !quake.ref
+// CHECK:           quake.x %[[VAL_4]] : (!quake.ref) -> ()
+// CHECK:           %[[VAL_5:.*]] = cc.load %[[VAL_1]] : !cc.ptr<f64>
+// CHECK:           %[[VAL_6:.*]] = cc.string_literal "XXXY" : !cc.ptr<!cc.array<i8 x 5>>
+// CHECK:           %[[VAL_7:.*]] = quake.extract_ref %[[VAL_2]][0] : (!quake.veq<4>) -> !quake.ref
+// CHECK:           %[[VAL_8:.*]] = quake.extract_ref %[[VAL_2]][1] : (!quake.veq<4>) -> !quake.ref
+// CHECK:           %[[VAL_9:.*]] = quake.extract_ref %[[VAL_2]][2] : (!quake.veq<4>) -> !quake.ref
+// CHECK:           %[[VAL_10:.*]] = quake.extract_ref %[[VAL_2]][3] : (!quake.veq<4>) -> !quake.ref
+// CHECK:           %[[VAL_11:.*]] = quake.concat %[[VAL_7]], %[[VAL_8]], %[[VAL_9]], %[[VAL_10]] : (!quake.ref, !quake.ref, !quake.ref, !quake.ref) -> !quake.veq<4>
+// CHECK:           quake.exp_pauli(%[[VAL_5]]) %[[VAL_11]], %[[VAL_6]] : (f64, !quake.veq<4>, !cc.ptr<!cc.array<i8 x 5>>) -> ()
+// CHECK:           return
+// CHECK:         }
diff --git a/test/Quake-QIR/exp-pauli-qir-1.qke b/test/Quake-QIR/exp-pauli-qir-1.qke
new file mode 100644
index 0000000000..c646cbf6f0
--- /dev/null
+++ b/test/Quake-QIR/exp-pauli-qir-1.qke
@@ -0,0 +1,39 @@
+// ========================================================================== //
+// Copyright (c) 2022 - 2023 NVIDIA Corporation & Affiliates.                 //
+// All rights reserved.                                                       //
+//                                                                            //
+// This source code and the accompanying materials are made available under   //
+// the terms of the Apache License 2.0 which accompanies this distribution.   //
+// ========================================================================== //
+
+// RUN: cudaq-opt %s --canonicalize --quake-add-deallocs | cudaq-translate --convert-to=qir | FileCheck %s
+
+module attributes {quake.mangled_name_map = {__nvqpp__mlirgen__Z4mainE3$_0 = "_ZZ4mainENK3$_0clEd"}} {
+  func.func @__nvqpp__mlirgen__Z4mainE3$_0(%arg0: f64) attributes {"cudaq-entrypoint", "cudaq-kernel"} {
+    %0 = cc.alloca f64
+    cc.store %arg0, %0 : !cc.ptr<f64>
+    %1 = quake.alloca !quake.veq<4>
+    %2 = quake.extract_ref %1[0] : (!quake.veq<4>) -> !quake.ref
+    quake.x %2 : (!quake.ref) -> ()
+    %3 = quake.extract_ref %1[1] : (!quake.veq<4>) -> !quake.ref
+    quake.x %3 : (!quake.ref) -> ()
+    %4 = cc.load %0 : !cc.ptr<f64>
+    %5 = cc.string_literal "XXXY" : !cc.ptr<!cc.array<i8 x 5>>
+    quake.exp_pauli(%4) %1, %5 : (f64, !quake.veq<4>, !cc.ptr<!cc.array<i8 x 5>>) -> ()
+    return
+  }
+}
+// CHECK:    %[[VAL_0:.*]] = tail call
+// CHECK:    %[[VAL_1:.*]]* @__quantum__rt__qubit_allocate_array(i64 4)
+// CHECK:         %[[VAL_2:.*]] = tail call %[[VAL_1]]* @__quantum__rt__array_slice(%[[VAL_1]]* %[[VAL_0]], i32 1, i64 0, i64 1, i64 3)
+// CHECK:         %[[VAL_3:.*]] = tail call i8* @__quantum__rt__array_get_element_ptr_1d(%[[VAL_1]]* %[[VAL_0]], i64 0)
+// CHECK:         %[[VAL_4:.*]] = bitcast i8* %[[VAL_3]] to %[[VAL_5:.*]]**
+// CHECK:         %[[VAL_6:.*]] = load %[[VAL_5]]*, %[[VAL_5]]** %[[VAL_4]], align 8
+// CHECK:         tail call void @__quantum__qis__x(%[[VAL_5]]* %[[VAL_6]])
+// CHECK:         %[[VAL_7:.*]] = tail call i8* @__quantum__rt__array_get_element_ptr_1d(%[[VAL_1]]* %[[VAL_0]], i64 1)
+// CHECK:         %[[VAL_8:.*]] = bitcast i8* %[[VAL_7]] to %[[VAL_5]]**
+// CHECK:         %[[VAL_9:.*]] = load %[[VAL_5]]*, %[[VAL_5]]** %[[VAL_8]], align 8
+// CHECK:         tail call void @__quantum__qis__x(%[[VAL_5]]* %[[VAL_9]])
+// CHECK:         tail call void @__quantum__qis__exp_pauli(double %[[VAL_10:.*]], %[[VAL_1]]* %[[VAL_2]], i8* nonnull getelementptr inbounds ([5 x i8], [5 x i8]* @cstr.5858585900, i64 0, i64 0))
+// CHECK:         tail call void @__quantum__rt__qubit_release_array(%[[VAL_1]]* %[[VAL_0]])
+// CHECK:         ret void
diff --git a/test/Quake-QIR/exp-pauli-qir-2.qke b/test/Quake-QIR/exp-pauli-qir-2.qke
new file mode 100644
index 0000000000..8a8583f656
--- /dev/null
+++ b/test/Quake-QIR/exp-pauli-qir-2.qke
@@ -0,0 +1,40 @@
+// ========================================================================== //
+// Copyright (c) 2022 - 2023 NVIDIA Corporation & Affiliates.                 //
+// All rights reserved.                                                       //
+//                                                                            //
+// This source code and the accompanying materials are made available under   //
+// the terms of the Apache License 2.0 which accompanies this distribution.   //
+// ========================================================================== //
+
+// RUN: cudaq-opt %s --canonicalize --quake-add-deallocs | cudaq-translate --convert-to=qir | FileCheck %s
+
+
+module attributes {quake.mangled_name_map = {__nvqpp__mlirgen__Z4mainE3$_0 = "_ZZ4mainENK3$_0clEdNSt7__cxx1112basic_stringIcSt11char_traitsIcESaIcEEE"}} {
+  func.func @__nvqpp__mlirgen__Z4mainE3$_0(%arg0: f64, %arg1: !cc.ptr<i8>) attributes {"cudaq-entrypoint", "cudaq-kernel"} {
+    %0 = cc.alloca f64
+    cc.store %arg0, %0 : !cc.ptr<f64>
+    %1 = quake.alloca !quake.veq<4>
+    %2 = quake.extract_ref %1[0] : (!quake.veq<4>) -> !quake.ref
+    quake.x %2 : (!quake.ref) -> ()
+    %3 = quake.extract_ref %1[1] : (!quake.veq<4>) -> !quake.ref
+    quake.x %3 : (!quake.ref) -> ()
+    %4 = cc.load %0 : !cc.ptr<f64>
+    quake.exp_pauli(%4) %1, %arg1 : (f64, !quake.veq<4>, !cc.ptr<i8>) -> ()
+    return
+  }
+}
+
+// CHECK:    %[[VAL_0:.*]] = tail call
+// CHECK:    %[[VAL_1:.*]]* @__quantum__rt__qubit_allocate_array(i64 4)
+// CHECK:         %[[VAL_2:.*]] = tail call %[[VAL_1]]* @__quantum__rt__array_slice(%[[VAL_1]]* %[[VAL_0]], i32 1, i64 0, i64 1, i64 3)
+// CHECK:         %[[VAL_3:.*]] = tail call i8* @__quantum__rt__array_get_element_ptr_1d(%[[VAL_1]]* %[[VAL_0]], i64 0)
+// CHECK:         %[[VAL_4:.*]] = bitcast i8* %[[VAL_3]] to %[[VAL_5:.*]]**
+// CHECK:         %[[VAL_6:.*]] = load %[[VAL_5]]*, %[[VAL_5]]** %[[VAL_4]], align 8
+// CHECK:         tail call void @__quantum__qis__x(%[[VAL_5]]* %[[VAL_6]])
+// CHECK:         %[[VAL_7:.*]] = tail call i8* @__quantum__rt__array_get_element_ptr_1d(%[[VAL_1]]* %[[VAL_0]], i64 1)
+// CHECK:         %[[VAL_8:.*]] = bitcast i8* %[[VAL_7]] to %[[VAL_5]]**
+// CHECK:         %[[VAL_9:.*]] = load %[[VAL_5]]*, %[[VAL_5]]** %[[VAL_8]], align 8
+// CHECK:         tail call void @__quantum__qis__x(%[[VAL_5]]* %[[VAL_9]])
+// CHECK:         tail call void @__quantum__qis__exp_pauli(double %[[VAL_13:.*]], %[[VAL_1]]* %[[VAL_2]], i8* %[[VAL_12:.*]])
+// CHECK:         tail call void @__quantum__rt__qubit_release_array(%[[VAL_1]]* %[[VAL_0]])
+// CHECK:         ret void
diff --git a/test/Quake-QIR/exp-pauli-qir-3.qke b/test/Quake-QIR/exp-pauli-qir-3.qke
new file mode 100644
index 0000000000..5270316b4e
--- /dev/null
+++ b/test/Quake-QIR/exp-pauli-qir-3.qke
@@ -0,0 +1,70 @@
+// ========================================================================== //
+// Copyright (c) 2022 - 2023 NVIDIA Corporation & Affiliates.                 //
+// All rights reserved.                                                       //
+//                                                                            //
+// This source code and the accompanying materials are made available under   //
+// the terms of the Apache License 2.0 which accompanies this distribution.   //
+// ========================================================================== //
+
+// RUN: cudaq-opt %s --canonicalize --quake-add-deallocs | cudaq-translate --convert-to=qir | FileCheck %s
+
+
+module attributes {quake.mangled_name_map = {__nvqpp__mlirgen__Z4mainE3$_0 = "_ZZ4mainENK3$_0clEd"}} {
+  func.func @__nvqpp__mlirgen__Z4mainE3$_0(%arg0: f64) attributes {"cudaq-entrypoint", "cudaq-kernel"} {
+    %0 = cc.alloca f64
+    cc.store %arg0, %0 : !cc.ptr<f64>
+    %1 = quake.alloca !quake.veq<4>
+    %2 = quake.extract_ref %1[0] : (!quake.veq<4>) -> !quake.ref
+    quake.x %2 : (!quake.ref) -> ()
+    %3 = quake.extract_ref %1[1] : (!quake.veq<4>) -> !quake.ref
+    quake.x %3 : (!quake.ref) -> ()
+    %4 = cc.load %0 : !cc.ptr<f64>
+    %5 = cc.string_literal "XXXY" : !cc.ptr<!cc.array<i8 x 5>>
+    %6 = quake.extract_ref %1[0] : (!quake.veq<4>) -> !quake.ref
+    %7 = quake.extract_ref %1[1] : (!quake.veq<4>) -> !quake.ref
+    %8 = quake.extract_ref %1[2] : (!quake.veq<4>) -> !quake.ref
+    %9 = quake.extract_ref %1[3] : (!quake.veq<4>) -> !quake.ref
+    %10 = quake.concat %6, %7, %8, %9 : (!quake.ref, !quake.ref, !quake.ref, !quake.ref) -> !quake.veq<4>
+    quake.exp_pauli(%4) %10, %5 : (f64, !quake.veq<4>, !cc.ptr<!cc.array<i8 x 5>>) -> ()
+    return
+  }
+}
+
+// CHECK:    %[[VAL_0:.*]] = tail call
+// CHECK:    %[[VAL_1:.*]]* @__quantum__rt__qubit_allocate_array(i64 4)
+// CHECK:         %[[VAL_2:.*]] = tail call i8* @__quantum__rt__array_get_element_ptr_1d(%[[VAL_1]]* %[[VAL_0]], i64 0)
+// CHECK:         %[[VAL_3:.*]] = bitcast i8* %[[VAL_2]] to %[[VAL_4:.*]]**
+// CHECK:         %[[VAL_5:.*]] = load %[[VAL_4]]*, %[[VAL_4]]** %[[VAL_3]], align 8
+// CHECK:         tail call void @__quantum__qis__x(%[[VAL_4]]* %[[VAL_5]])
+// CHECK:         %[[VAL_6:.*]] = tail call i8* @__quantum__rt__array_get_element_ptr_1d(%[[VAL_1]]* %[[VAL_0]], i64 1)
+// CHECK:         %[[VAL_7:.*]] = bitcast i8* %[[VAL_6]] to %[[VAL_4]]**
+// CHECK:         %[[VAL_8:.*]] = load %[[VAL_4]]*, %[[VAL_4]]** %[[VAL_7]], align 8
+// CHECK:         tail call void @__quantum__qis__x(%[[VAL_4]]* %[[VAL_8]])
+// CHECK:         %[[VAL_9:.*]] = tail call i8* @__quantum__rt__array_get_element_ptr_1d(%[[VAL_1]]* %[[VAL_0]], i64 2)
+// CHECK:         %[[VAL_10:.*]] = bitcast i8* %[[VAL_9]] to i8**
+// CHECK:         %[[VAL_11:.*]] = load i8*, i8** %[[VAL_10]], align 8
+// CHECK:         %[[VAL_12:.*]] = tail call i8* @__quantum__rt__array_get_element_ptr_1d(%[[VAL_1]]* %[[VAL_0]], i64 3)
+// CHECK:         %[[VAL_13:.*]] = bitcast i8* %[[VAL_12]] to i8**
+// CHECK:         %[[VAL_14:.*]] = load i8*, i8** %[[VAL_13]], align 8
+// CHECK:         %[[VAL_15:.*]] = tail call %[[VAL_1]]* @__quantum__rt__array_create_1d(i32 8, i64 1)
+// CHECK:         %[[VAL_16:.*]] = tail call i8* @__quantum__rt__array_get_element_ptr_1d(%[[VAL_1]]* %[[VAL_15]], i64 0)
+// CHECK:         %[[VAL_17:.*]] = bitcast i8* %[[VAL_16]] to %[[VAL_4]]**
+// CHECK:         store %[[VAL_4]]* %[[VAL_5]], %[[VAL_4]]** %[[VAL_17]], align 8
+// CHECK:         %[[VAL_18:.*]] = tail call %[[VAL_1]]* @__quantum__rt__array_create_1d(i32 8, i64 1)
+// CHECK:         %[[VAL_19:.*]] = tail call i8* @__quantum__rt__array_get_element_ptr_1d(%[[VAL_1]]* %[[VAL_18]], i64 0)
+// CHECK:         %[[VAL_20:.*]] = bitcast i8* %[[VAL_19]] to %[[VAL_4]]**
+// CHECK:         store %[[VAL_4]]* %[[VAL_8]], %[[VAL_4]]** %[[VAL_20]], align 8
+// CHECK:         %[[VAL_21:.*]] = tail call %[[VAL_1]]* @__quantum__rt__array_concatenate(%[[VAL_1]]* %[[VAL_15]], %[[VAL_1]]* %[[VAL_18]])
+// CHECK:         %[[VAL_22:.*]] = tail call %[[VAL_1]]* @__quantum__rt__array_create_1d(i32 8, i64 1)
+// CHECK:         %[[VAL_23:.*]] = tail call i8* @__quantum__rt__array_get_element_ptr_1d(%[[VAL_1]]* %[[VAL_22]], i64 0)
+// CHECK:         %[[VAL_24:.*]] = bitcast i8* %[[VAL_23]] to i8**
+// CHECK:         store i8* %[[VAL_11]], i8** %[[VAL_24]], align 8
+// CHECK:         %[[VAL_25:.*]] = tail call %[[VAL_1]]* @__quantum__rt__array_concatenate(%[[VAL_1]]* %[[VAL_21]], %[[VAL_1]]* %[[VAL_22]])
+// CHECK:         %[[VAL_26:.*]] = tail call %[[VAL_1]]* @__quantum__rt__array_create_1d(i32 8, i64 1)
+// CHECK:         %[[VAL_27:.*]] = tail call i8* @__quantum__rt__array_get_element_ptr_1d(%[[VAL_1]]* %[[VAL_26]], i64 0)
+// CHECK:         %[[VAL_28:.*]] = bitcast i8* %[[VAL_27]] to i8**
+// CHECK:         store i8* %[[VAL_14]], i8** %[[VAL_28]], align 8
+// CHECK:         %[[VAL_29:.*]] = tail call %[[VAL_1]]* @__quantum__rt__array_concatenate(%[[VAL_1]]* %[[VAL_25]], %[[VAL_1]]* %[[VAL_26]])
+// CHECK:         tail call void @__quantum__qis__exp_pauli(double %[[VAL_30:.*]], %[[VAL_1]]* %[[VAL_29]], i8* nonnull getelementptr inbounds ([5 x i8], [5 x i8]* @cstr.5858585900, i64 0, i64 0))
+// CHECK:         tail call void @__quantum__rt__qubit_release_array(%[[VAL_1]]* %[[VAL_0]])
+// CHECK:         ret void
diff --git a/test/Quake/roundtrip-ops.qke b/test/Quake/roundtrip-ops.qke
index 818534eb4b..432c8327d8 100644
--- a/test/Quake/roundtrip-ops.qke
+++ b/test/Quake/roundtrip-ops.qke
@@ -143,6 +143,11 @@ func.func @quantum_ops() {
   quake.u2 (%f, %g) [%46] %7 : (f32, f32, !quake.veq<3>, !quake.ref) -> ()
   quake.u3 (%f, %g, %47) [%46] %7 : (f32, f32, f64, !quake.veq<3>, !quake.ref) -> ()
 
+  // Exp Pauli and StringLiteral
+  %f2 = arith.constant 12.0 : f64
+  %48 = cc.string_literal "XXY" : !cc.ptr<!cc.array<i8 x 4>>
+  quake.exp_pauli (%f2) %46, %48 : (f64, !quake.veq<3>, !cc.ptr<!cc.array<i8 x 4>>) -> ()
+  
   return
 }
 
@@ -250,6 +255,9 @@ func.func @quantum_ops() {
 // CHECK:           quake.rz (%[[VAL_20]]) {{\[}}%[[VAL_54]]] %[[VAL_10]] : (f32, !quake.veq<3>, !quake.ref) -> ()
 // CHECK:           quake.u2 (%[[VAL_20]], %[[VAL_21]]) {{\[}}%[[VAL_54]]] %[[VAL_10]] : (f32, f32, !quake.veq<3>, !quake.ref) -> ()
 // CHECK:           quake.u3 (%[[VAL_20]], %[[VAL_21]], %[[VAL_53]]) {{\[}}%[[VAL_54]]] %[[VAL_10]] : (f32, f32, f64, !quake.veq<3>, !quake.ref) -> ()
+// CHECK:           %[[VAL_59:.*]] = arith.constant 1.200000e+01 : f64
+// CHECK:           %[[VAL_60:.*]] = cc.string_literal "XXY" : !cc.ptr<!cc.array<i8 x 4>>
+// CHECK:           quake.exp_pauli(%[[VAL_59]]) %[[VAL_54]], %[[VAL_60]] : (f64, !quake.veq<3>, !cc.ptr<!cc.array<i8 x 4>>) -> ()
 // CHECK:           return
 // CHECK:         }
 
diff --git a/unittests/domains/ChemistryTester.cpp b/unittests/domains/ChemistryTester.cpp
index 11f5f3b31d..d0e8383e79 100644
--- a/unittests/domains/ChemistryTester.cpp
+++ b/unittests/domains/ChemistryTester.cpp
@@ -53,6 +53,29 @@ CUDAQ_TEST(H2MoleculeTester, checkHamiltonian) {
   }
 }
 
+CUDAQ_TEST(H2MoleculeTester, checkExpPauli) {
+  auto kernel = [](double theta) __qpu__ {
+    cudaq::qreg q(4);
+    x(q[0]);
+    x(q[1]);
+    exp_pauli(theta, q, "XXXY");
+  };
+
+  cudaq::molecular_geometry geometry{{"H", {0., 0., 0.}},
+                                     {"H", {0., 0., .7474}}};
+  auto molecule = cudaq::create_molecule(geometry, "sto-3g", 1, 0);
+  cudaq::observe(kernel, molecule.hamiltonian, 1.1);
+
+  cudaq::optimizers::cobyla optimizer;
+  auto [e, opt] = optimizer.optimize(1, [&](std::vector<double> x) -> double {
+    double e = cudaq::observe(kernel, molecule.hamiltonian, x[0]);
+    printf("E = %lf\n", e);
+    return e;
+  });
+
+  EXPECT_NEAR(-1.137, e, 1e-3);
+}
+
 CUDAQ_TEST(H2MoleculeTester, checkUCCSD) {
   {
     cudaq::molecular_geometry geometry{{"H", {0., 0., 0.}},
diff --git a/unittests/integration/builder_tester.cpp b/unittests/integration/builder_tester.cpp
index dd8ae63f04..4ed4d6800a 100644
--- a/unittests/integration/builder_tester.cpp
+++ b/unittests/integration/builder_tester.cpp
@@ -654,6 +654,61 @@ CUDAQ_TEST(BuilderTester, checkEntryPointAttribute) {
   EXPECT_TRUE(std::regex_search(quake, functionDecleration));
 }
 
+CUDAQ_TEST(BuilderTester, checkExpPauli) {
+  {
+    auto [kernel, theta] = cudaq::make_kernel<double>();
+    auto qubits = kernel.qalloc(4);
+    kernel.x(qubits[0]);
+    kernel.x(qubits[1]);
+    kernel.exp_pauli(theta, qubits, "XXXY");
+    std::cout << kernel << "\n";
+    std::vector<double> h2_data{
+        3, 1, 1, 3, 0.0454063,  0,  2, 0, 0, 0, 0.17028,    0,
+        0, 0, 2, 0, -0.220041,  -0, 1, 3, 3, 1, 0.0454063,  0,
+        0, 0, 0, 0, -0.106477,  0,  0, 2, 0, 0, 0.17028,    0,
+        0, 0, 0, 2, -0.220041,  -0, 3, 3, 1, 1, -0.0454063, -0,
+        2, 2, 0, 0, 0.168336,   0,  2, 0, 2, 0, 0.1202,     0,
+        0, 2, 0, 2, 0.1202,     0,  2, 0, 0, 2, 0.165607,   0,
+        0, 2, 2, 0, 0.165607,   0,  0, 0, 2, 2, 0.174073,   0,
+        1, 1, 3, 3, -0.0454063, -0, 15};
+    cudaq::spin_op h(h2_data, 4);
+    cudaq::optimizers::cobyla optimizer;
+    optimizer.max_eval = 30;
+    auto [e, opt] = optimizer.optimize(1, [&](std::vector<double> x) -> double {
+      double e = cudaq::observe(kernel, h, x[0]);
+      printf("E = %lf, %lf\n", e, x[0]);
+      return e;
+    });
+    EXPECT_NEAR(e, -1.13, 1e-2);
+  }
+  {
+    auto [kernel, theta] = cudaq::make_kernel<double>();
+    auto qubits = kernel.qalloc(4);
+    kernel.x(qubits[0]);
+    kernel.x(qubits[1]);
+    kernel.exp_pauli(theta, "XXXY", qubits[0], qubits[1], qubits[2], qubits[3]);
+    std::cout << kernel << "\n";
+    std::vector<double> h2_data{
+        3, 1, 1, 3, 0.0454063,  0,  2, 0, 0, 0, 0.17028,    0,
+        0, 0, 2, 0, -0.220041,  -0, 1, 3, 3, 1, 0.0454063,  0,
+        0, 0, 0, 0, -0.106477,  0,  0, 2, 0, 0, 0.17028,    0,
+        0, 0, 0, 2, -0.220041,  -0, 3, 3, 1, 1, -0.0454063, -0,
+        2, 2, 0, 0, 0.168336,   0,  2, 0, 2, 0, 0.1202,     0,
+        0, 2, 0, 2, 0.1202,     0,  2, 0, 0, 2, 0.165607,   0,
+        0, 2, 2, 0, 0.165607,   0,  0, 0, 2, 2, 0.174073,   0,
+        1, 1, 3, 3, -0.0454063, -0, 15};
+    cudaq::spin_op h(h2_data, 4);
+    cudaq::optimizers::cobyla optimizer;
+    optimizer.max_eval = 30;
+    auto [e, opt] = optimizer.optimize(1, [&](std::vector<double> x) -> double {
+      double e = cudaq::observe(kernel, h, x[0]);
+      printf("E = %lf, %lf\n", e, x[0]);
+      return e;
+    });
+    EXPECT_NEAR(e, -1.13, 1e-2);
+  }
+}
+
 #ifndef CUDAQ_BACKEND_DM
 
 CUDAQ_TEST(BuilderTester, checkCanProgressivelyBuild) {
diff --git a/unittests/qudit/simple_qudit/SimpleQuditExecutionManager.cpp b/unittests/qudit/simple_qudit/SimpleQuditExecutionManager.cpp
index 5ea3d1b48a..34640dd793 100644
--- a/unittests/qudit/simple_qudit/SimpleQuditExecutionManager.cpp
+++ b/unittests/qudit/simple_qudit/SimpleQuditExecutionManager.cpp
@@ -105,7 +105,7 @@ class SimpleQuditExecutionManager : public cudaq::BasicExecutionManager {
     instructions.emplace("plusGate", [&](const Instruction &inst) {
       qpp::cmat u(3, 3);
       u << 0, 0, 1, 1, 0, 0, 0, 1, 0;
-      auto &[gateName, params, controls, qudits] = inst;
+      auto &[gateName, params, controls, qudits, op] = inst;
       auto target = qudits[0];
       cudaq::info("Applying plusGate on {}<{}>", target.id, target.levels);
       state = qpp::applyCTRL(state, u, {}, {target.id}, target.levels);