dart/src/wAssign.h

@@ -83,4 +83,4 @@ public:
 
 } // namespace dart
 
-#endif //WASSIGN_H
+#endif // WASSIGN_H

dart/src/wBlowing.cpp

@@ -59,14 +59,6 @@ void Blowing::create()
     }
 }
 
-/**
- * @brief Set blowing velocity
- */
-void Blowing::setU(size_t i, double ue)
-{
-    uE[i].second->update(ue);
-}
-
 void Blowing::write(std::ostream &out) const
 {
     out << " Blowing boundary condition on " << *tag << std::endl;

dart/src/wBlowing.h

@@ -41,7 +41,8 @@ public:
     Blowing(std::shared_ptr<tbox::MshData> _msh, std::string const &names);
     virtual ~Blowing();
 
-    void setU(size_t i, double ue);
+    void setU(size_t i, double ue) { uE[i].second->update(ue); }
+    double getU(size_t i) const { return uE[i].second->eval(*uE[i].first, std::vector<double>(), 0); }
 
 #ifndef SWIG
     virtual void write(std::ostream &out) const override;
@@ -53,4 +54,4 @@ private:
 
 } // namespace dart
 
-#endif //WBLOWING_H
+#endif // WBLOWING_H

dart/src/wBlowingResidual.cpp

@@ -27,7 +27,7 @@ using namespace dart;
 /**
  * @brief Build the residual vector, on one boundary element
  *
- * b = \int psi * vn dV
+ * b = \int psi * vn dS
  */
 Eigen::VectorXd BlowingResidual::build(Element const &e, std::vector<double> const &phi, F0El const &f)
 {
@@ -44,3 +44,44 @@ 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
+ *
+ * b = \int psi dS
+ */
+Eigen::VectorXd BlowingResidual::buildGradientBlowing(tbox::Element const &e)
+{
+    // Get pre-computed values
+    Cache &cache = e.getVCache();
+    Gauss &gauss = cache.getVGauss();
+
+    // Build velocity gradient
+    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
+ *
+ * A = \int psi * vn ddS
+ */
+Eigen::MatrixXd BlowingResidual::buildGradientMesh(tbox::Element const &e, std::vector<double> const &phi, F0El const &f, int const nDim)
+{
+    // Get pre-computed values
+    Cache &cache = e.getVCache();
+    Gauss &gauss = cache.getVGauss();
+
+    // Build velocity gradient
+    double vn = f.eval(e, phi, 0);
+    Eigen::MatrixXd A = Eigen::MatrixXd::Zero(e.nodes.size(), nDim * e.nodes.size());
+    for (size_t k = 0; k < gauss.getN(); ++k)
+    {
+        std::vector<double> dDetJ = e.getGradDetJ(k);
+        for (size_t i = 0; i < dDetJ.size(); i++)
+            A.col(i) += cache.getSf(k) * vn * gauss.getW(k) * dDetJ[i];
+    }
+    return A;
+}

dart/src/wBlowingResidual.h

@@ -33,7 +33,9 @@ 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
-#endif //WBLOWINGRESIDUAL_H
+#endif // WBLOWINGRESIDUAL_H

dart/src/wBody.cpp

@@ -22,8 +22,6 @@
 #include "wElement.h"
 #include "wNode.h"
 #include <unordered_set>
-#include <iomanip>
-#include <fstream>
 using namespace tbox;
 using namespace dart;
 
@@ -111,84 +109,7 @@ void Body::create()
     nLoads.resize(nodes.size(), Eigen::Vector3d::Zero());
 }
 
-/**
- * @brief Save nodal data for further post-processing
- */
-void Body::save(std::string const &name, Results const &res)
-{
-    // Write to file
-    std::cout << "writing file: " << name + ".dat"
-              << "... " << std::flush;
-    std::ofstream outfile;
-    outfile.open(name + ".dat");
-    // Header
-    outfile << "$Body data" << std::endl;
-    // Aerodynamic coefficients
-    outfile << "$Aerodynamic coefficients" << std::endl;
-    outfile << "!" << std::fixed
-            << std::setw(14) << "Cl"
-            << std::setw(15) << "Cd"
-            << std::setw(15) << "Cs"
-            << std::setw(15) << "Cm"
-            << std::endl;
-    outfile << std::fixed
-            << std::setw(15) << Cl
-            << std::setw(15) << Cd
-            << std::setw(15) << Cs
-            << std::setw(15) << Cm
-            << std::endl;
-    // Elements (connectvity)
-    outfile << "$Elements" << std::endl;
-    outfile << groups[0]->tag->elems.size() << std::endl;
-    for (auto e : groups[0]->tag->elems)
-    {
-        outfile << std::fixed
-                << std::setw(10) << e->no;
-        for (auto n : e->nodes)
-            outfile << std::setw(10) << n->no;
-        outfile << std::endl;
-    }
-    //Nodes (data)
-    outfile << "$Nodal data" << std::endl;
-    outfile << nodes.size() << std::endl;
-    outfile << "!" << std::fixed
-            << std::setw(9) << "no"
-            << std::setw(15) << "x"
-            << std::setw(15) << "y"
-            << std::setw(15) << "z";
-    for (auto &p : res.scalars_at_nodes)
-        outfile << std::setw(15) << p.first;
-    for (auto &p : res.vectors_at_nodes)
-    {
-        outfile << std::setw(14) << p.first << "x"
-                << std::setw(14) << p.first << "y"
-                << std::setw(14) << p.first << "z";
-    }
-    outfile << std::endl;
-    for (auto n : nodes)
-    {
-        outfile << std::fixed
-                << std::setw(10) << n->no
-                << std::scientific << std::setprecision(6)
-                << std::setw(15) << n->pos(0)
-                << std::setw(15) << n->pos(1)
-                << std::setw(15) << n->pos(2);
-        for (auto &p : res.scalars_at_nodes)
-            outfile << std::setw(15) << (*(p.second))[n->row];
-        for (auto &p : res.vectors_at_nodes)
-            outfile << std::setw(15) << (*(p.second))[n->row](0)
-                    << std::setw(15) << (*(p.second))[n->row](1)
-                    << std::setw(15) << (*(p.second))[n->row](2);
-        outfile << std::endl;
-    }
-    // Footer
-    outfile << std::endl;
-    // Close file
-    outfile.close();
-    std::cout << "done" << std::endl;
-}
-
 void Body::write(std::ostream &out) const
 {
     out << *groups[0]->tag << " is a Body immersed in " << *groups[1]->tag << std::endl;
-}
+}
\ No newline at end of file

dart/src/wBody.h

@@ -49,8 +49,6 @@ public:
     Body(std::shared_ptr<tbox::MshData> _msh, std::vector<std::string> const &names);
     virtual ~Body() { std::cout << "~Body()\n"; }
 
-    void save(std::string const &name, tbox::Results const &res);
-
 #ifndef SWIG
     virtual void write(std::ostream &out) const override;
 #endif
@@ -61,4 +59,4 @@ private:
 
 } // namespace dart
 
-#endif //WBODY_H
+#endif // WBODY_H

dart/src/wF0Ps.h

@@ -54,4 +54,4 @@ public:
 
 } // namespace dart
 
-#endif //WF0PS_H
+#endif // WF0PS_H

dart/src/wF1Ct.h

@@ -50,7 +50,7 @@ class DART_API F1CtDrag : public F1Ct
 
 public:
     F1CtDrag(int _nDim, double _sRef, double alpha, double _beta = 0);
-    void update(double alpha);
+    virtual void update(double alpha) override;
 
     virtual Eigen::Vector3d eval() const override;
     virtual Eigen::Vector3d evalGrad() const override;
@@ -69,7 +69,7 @@ class DART_API F1CtSide : public F1Ct
 
 public:
     F1CtSide(int _nDim, double _sRef, double alpha, double _beta = 0);
-    void update(double alpha);
+    virtual void update(double alpha) override;
 
     virtual Eigen::Vector3d eval() const override;
     virtual Eigen::Vector3d evalGrad() const override;
@@ -88,7 +88,7 @@ class DART_API F1CtLift : public F1Ct
 
 public:
     F1CtLift(int _nDim, double _sRef, double alpha, double _beta = 0);
-    void update(double alpha);
+    virtual void update(double alpha) override;
 
     virtual Eigen::Vector3d eval() const override;
     virtual Eigen::Vector3d evalGrad() const override;
@@ -96,4 +96,4 @@ public:
 
 } // namespace dart
 
-#endif //WF1CT_H
+#endif // WF1CT_H

dart/src/wF1El.h

@@ -62,4 +62,4 @@ public:
 
 } // namespace dart
 
-#endif //WF1EL_H
+#endif // WF1EL_H

dart/src/wFace.h

@@ -85,4 +85,4 @@ public:
 
 } // namespace dart
 
-#endif //WFACE_H
+#endif // WFACE_H

dart/src/wFluid.cpp

@@ -54,15 +54,15 @@ void Fluid::initFun(double mInf, int dim, double alpha, double beta)
 {
     if (mInf == 0)
     {
-        rho = new F0ElRhoL();
-        mach = new F0ElMachL();
-        cP = new F0ElCpL();
+        rho = new F0ElRhoLin();
+        mach = new F0ElMachLin();
+        cP = new F0ElCpLin();
     }
     else
     {
-        rho = new F0ElRho(mInf);
-        mach = new F0ElMach(mInf);
-        cP = new F0ElCp(mInf);
+        rho = new F0ElRhoClamp(mInf, 3.0);
+        mach = new F0ElMachClamp(mInf, 3.0);
+        cP = new F0ElCpClamp(mInf, 3.0);
     }
     phiInf = new F0PsPhiInf(dim, alpha, beta);
 }