diff --git a/dart/src/wAdjoint.cpp b/dart/src/wAdjoint.cpp
index 4dd690ba30a5bbeee109ff5ffec975555812ef65..367921e19a175e7b83e68d7d2122e6c93fdb26df 100644
--- a/dart/src/wAdjoint.cpp
+++ b/dart/src/wAdjoint.cpp
@@ -27,6 +27,8 @@
 #include "wWakeResidual.h"
 #include "wKuttaResidual.h"
 #include "wLoadFunctional.h"
+#include "wBlowing.h"
+#include "wBlowingResidual.h"
 
 #include "wMshData.h"
 #include "wNode.h"
@@ -561,6 +563,28 @@ void Adjoint::buildGradientResidualMesh(Eigen::SparseMatrix<double, Eigen::RowMa
     });
     // Farfield B.C. are not taken into account because they have no influence on the solution
 
+    // Blowing boundary
+    if (sol->pbl->bBCs.size() > 0)
+    {
+        for (auto b : sol->pbl->bBCs)
+        {
+            for (auto p : b->uE)
+            {
+                Eigen::MatrixXd A = BlowingResidual::buildGradientMesh(*p.first, sol->phi, *p.second, sol->pbl->nDim);
+                for (size_t i = 0; i < p.first->nodes.size(); ++i)
+                {
+                    size_t rowi = p.first->nodes[i]->row;
+                    for (size_t j = 0; j < p.first->nodes.size(); ++j)
+                    {
+                        int rowj = p.first->nodes[j]->row;
+                        for (int n = 0; n < sol->pbl->nDim; ++n)
+                            T.push_back(Eigen::Triplet<double>(sol->rows[rowi], sol->pbl->nDim * rowj + n, A(i, sol->pbl->nDim * j + n)));
+                    }
+                }
+            }
+        }
+    }
+
     // Wake
     for (auto wake : sol->pbl->wBCs)
     {
diff --git a/dart/src/wAdjoint.h b/dart/src/wAdjoint.h
index a8c678eefac9820a60beaa2834cfa83782a0b03b..79f4eee16f311488270943b49b9c4106e6af62d9 100644
--- a/dart/src/wAdjoint.h
+++ b/dart/src/wAdjoint.h
@@ -87,6 +87,11 @@ public:
     std::vector<double> computeGradientCoefficientsAoa();
     std::vector<std::vector<double>> computeGradientCoefficientsMesh();
 
+    Eigen::SparseMatrix<double, Eigen::RowMajor> getRu_U() const { return dRu_U; }
+    Eigen::VectorXd getRu_A() const { return dRu_A; }
+    Eigen::SparseMatrix<double, Eigen::RowMajor> getRu_X() const { return dRu_X; }
+    Eigen::SparseMatrix<double, Eigen::RowMajor> getRx_X() const { return dRx_X; }
+
 #ifndef SWIG
     virtual void write(std::ostream &out) const override;
 #endif
diff --git a/dart/src/wBlowing.cpp b/dart/src/wBlowing.cpp
index a0c2f14b42198d27f16511bd62a8af839104e87c..ef2c86a5d244f0584f337ad74ef91750dd543dd7 100644
--- a/dart/src/wBlowing.cpp
+++ b/dart/src/wBlowing.cpp
@@ -67,6 +67,17 @@ void Blowing::setU(size_t i, double ue)
     uE[i].second->update(ue);
 }
 
+/**
+ * @brief Get blowing velocity
+ */
+double Blowing::getU(size_t i) const
+{
+    // Return a constant value
+    tbox::Element *e = uE[i].first;
+    std::vector<double> phi(0, 0.);
+    return uE[i].second->eval(*e, phi, 0);
+}
+
 void Blowing::write(std::ostream &out) const
 {
     out << " Blowing boundary condition on " << *tag << std::endl;
diff --git a/dart/src/wBlowing.h b/dart/src/wBlowing.h
index 7334194e0c35134323a0da8388f99215f5aff3ed..2a0cd7db1ab63ec4d995fde5b8ccb1bd159d8c20 100644
--- a/dart/src/wBlowing.h
+++ b/dart/src/wBlowing.h
@@ -42,6 +42,7 @@ public:
     virtual ~Blowing();
 
     void setU(size_t i, double ue);
+    double getU(size_t i) const;
 
 #ifndef SWIG
     virtual void write(std::ostream &out) const override;
diff --git a/dart/src/wBlowingResidual.cpp b/dart/src/wBlowingResidual.cpp
index 2a5b35c1fa306d1839dba60e902f6ea634306084..b8b4ea3a1239a29b1e85ca342fbd9cf3f8f3f637 100644
--- a/dart/src/wBlowingResidual.cpp
+++ b/dart/src/wBlowingResidual.cpp
@@ -44,3 +44,35 @@ Eigen::VectorXd BlowingResidual::build(Element const &e, std::vector<double> con
         b += cache.getSf(k) * vn * gauss.getW(k) * e.getDetJ(k);
     return b;
 }
+
+/**
+ * @brief Build the gradient of the residual vector with respect to the blowing velocity, on one blowing boundary element
+ */
+Eigen::VectorXd BlowingResidual::buildGradientBlowing(tbox::Element const &e)
+{
+    Cache &cache = e.getVCache();
+    Gauss &gauss = cache.getVGauss();
+    Eigen::VectorXd b = Eigen::VectorXd::Zero(e.nodes.size());
+    for (size_t k = 0; k < gauss.getN(); ++k)
+        b += cache.getSf(k) * gauss.getW(k) * e.getDetJ(k);
+    return b;
+}
+
+/**
+ * @brief Build the gradient of the residual vector with respect to the mesh, on one blowing boundary element
+ */
+Eigen::MatrixXd BlowingResidual::buildGradientMesh(tbox::Element const &e, std::vector<double> const &phi, F0El const &f, int const nDim)
+{
+    Eigen::MatrixXd A = Eigen::MatrixXd::Zero(e.nodes.size(), nDim * e.nodes.size());
+    Cache &cache = e.getVCache();
+    Gauss &gauss = cache.getVGauss();
+    // Blowing velocity
+    double vn = f.eval(e, phi, 0);
+    for (size_t k = 0; k < gauss.getN(); ++k)
+    {
+        std::vector<double> gradJ = e.getGradDetJ(k);
+        for (size_t i = 0; i < gradJ.size(); i++)
+            A.col(i) += cache.getSf(k) * vn * gauss.getW(k) * gradJ[i];
+    }
+    return A;
+}
diff --git a/dart/src/wBlowingResidual.h b/dart/src/wBlowingResidual.h
index f2857b86e98dc4c0efd8ab244c2c3d519b53cd18..f49b34f2e57d7c866387872f0cd339188fd647de 100644
--- a/dart/src/wBlowingResidual.h
+++ b/dart/src/wBlowingResidual.h
@@ -33,6 +33,8 @@ class DART_API BlowingResidual
 public:
     // Newton
     static Eigen::VectorXd build(tbox::Element const &e, std::vector<double> const &phi, F0El const &f);
+    static Eigen::VectorXd buildGradientBlowing(tbox::Element const &e);
+    static Eigen::MatrixXd buildGradientMesh(tbox::Element const &e, std::vector<double> const &phi, F0El const &f, int const nDim);
 };
 
 } // namespace dart