added: Advection to ResidualOps

Apps/Common/EqualOrderOperators.C

@@ -13,6 +13,7 @@
 #include "EqualOrderOperators.h"
 #include "FiniteElement.h"
 #include "Vec3.h"
+#include "Vec3Oper.h"
 
 //! \brief Helper for adding an element matrix to several components.
 //! \param[out] EM The element matrix to add to.
@@ -211,6 +212,19 @@ void EqualOrderOperators::Weak::Source(Vector& EV, const FiniteElement& fe,
 }
 
 
+void EqualOrderOperators::Residual::Advection(Vector& EV, const FiniteElement& fe,
+                                              const Vec3& AC, const Tensor& g,
+                                              double scale, int basis)
+{
+  size_t nsd = fe.grad(basis).cols();
+  for (size_t k = 1; k <= nsd; ++k) {
+    double ag = g[k-1]*AC;
+    for (size_t i = 1; i <= fe.basis(basis).size(); ++i)
+        EV((i-1)*nsd+k) = ag*scale*fe.basis(basis)(i)*fe.detJxW;
+  }
+}
+
+
 void EqualOrderOperators::Residual::Convection(Vector& EV, const FiniteElement& fe,
                                                 const Vec3& U, const Tensor& dUdX,
                                                const Vec3& UC, double scale,

Apps/Common/EqualOrderOperators.h

@@ -142,6 +142,16 @@ public:
   //! \brief Common weak residual operators using equal-ordered discretizations.
   class Residual {
   public:
+    //! \brief Compute an advection term.
+    //! \param[out] EM The element vector to add contribution to
+    //! \param[in] fe The finite element to evaluate for
+    //! \param[in] AC Advecting field
+    //! \param[in] scale Scaling factor for contribution
+    //! \param[in] basis Basis to use
+    static void Advection(Vector& EV, const FiniteElement& fe,
+                          const Vec3& AC, const Tensor& g,
+                          double scale = 1.0, int basis=1);
+
     //! \brief Compute a convection term in a residual vector.
     //! \param EV The element vector to add contribution to
     //! \param[in] fe The finite element to evaluate for