Green's alternative formula for integration by parts I'm currently implementing a method to solve usual elliptic problems where the classical form is the following:
$$-\text{div}(k\nabla u) + \vec \beta \cdot \nabla u + \gamma u = f$$
Due to an implementation task, I need to integrate by parts the following integral:
$$\int_E k \nabla \varphi \cdot \nabla m$$
where $E$ is a poligonal element, $k$ is the diffusion term and $\varphi, m$ are respectively one element of the basis of $V_h \subseteq H^1_0(E)$ a finite dimensional space and $\mathbb{P}_1(E)$.
Since we are working with polynomials of degree one, in the case $k \equiv 1$ I can use integration by parts, leave the laplacian on $m$ such that the integral on $E$ is $=0$ and I'm left with an integral on the boundary of $E$.

What can I do in the case $k \not \equiv 1$? Are there any manipulations of Green's identities to do such a job?

 A: $
\newcommand\grade[1]{\left\langle#1\right\rangle}
\newcommand\d{\mathrm d}
\newcommand\vec\mathbf
\newcommand\R{\mathbb R}
\newcommand\PD[2]{\frac{\partial#1}{\partial#2}}
$
I don't know if this will be helpful or not, but here I give a possible "integration by parts", and give a proof further below.
The boundary $\partial E = \partial E_1 \cup \partial E_2 \cup \cdots \cup \partial E_n$, where each $\partial E_i$ is an edge of the polygon $E$. Each $\partial E_i$ has it's own (constant) unit vector $\vec E_i$ along the edge and unit vector $\vec O_i$ orthogonal to it but still in the plane. If $\vec N$ is a unit normal vector of $V$, then $\vec O_i = \vec E_i\times \vec N$.
In coordinates $\vec x = x_i\vec E_i + y_i\vec O_i + z\vec N$, your integral becomes
$$
  \int_E\d S\,k\,(\nabla\phi){\cdot}(\nabla m)
= -\sum_i\int_{\partial E_i}\d x_i\,k\phi\PD m{y_i}
    + \int_E\d S\Bigl[
      k\PD\phi z\PD mz
      - \phi\,(\partial k){\cdot}(\partial m)
      - k\phi(\partial z^2)
    \Bigr],
$$
where $\d S = \d x_i\d y_i$ is the scalar surface measure and $\partial$ is to be interpreted as the 2D gradient within the plane of $E$; explicitly
$$
  \nabla = \partial + \vec N\PD{}z,\quad
  \nabla^2 = \partial^2 + \PD{^2}{z^2}.
$$

Below I use geometric calculus, which is calculus done with geometric (Clifford) algebra. I know that likely makes it difficult to understand, which is unfortunate, but I don't know how to derive this with just vector calculus. Anyone is free to translate this into vector calculus if they are able to.
The integration by parts is accomplished using the Fundamental Theorem of Geometric Calculus (closely related to the generalized Stokes' Theorem), which says
$$
  \int_V\d^n\vec x\,\partial f(\vec x)
= \int_{\partial V}\d^{n-1}\vec x\,f(\vec x),
$$
where $V$ is an $n$-dimensional manifold, $\d^n\vec x$ is the oriented $n$-vector valued measure on $V$, $\partial V$ is the $(n-1)$-dimensional boundary of $V$ with $(n-1)$-vector valued measure $\d^{n-1}\vec x$, and $\partial$ is the tangential vector derivative. By convention, this derivative does not differentiate any position dependence of $\d^n\vec x$. This derivative can be thought of as the gradient projected onto $V$, and if $\nabla$ is the vector derivative of the ambient space than $\nabla = \partial + \partial_\perp$, where $\partial_\perp$ is the vector derivative orthogonal to $V$.
We work in $\R^3$ with pseudoscalar $I = e_1e_2e_3$, where $\{e_1,e_2,e_3\}$ is the standard basis.
First we rewrite the integral as
$$
  J := \int_E|\d^2\vec x|\,k\,(\nabla\phi){\cdot}(\nabla m) = \grade{\int_E|\d^2\vec x|\dot\nabla\hat\nabla k\dot\phi\hat m}.
$$
To integrate by parts, we make use of the product rule
$$
  \nabla\hat\nabla k\phi\hat m
= \dot\nabla\hat\nabla\dot k\phi\hat m
    + \dot\nabla\hat\nabla k\dot\phi\hat m
    + \dot\nabla\hat\nabla k\phi\hat{\dot m},
$$
rearranging to get
$$
  \dot\nabla\hat\nabla k\dot\phi\hat m
= \nabla\hat\nabla k\phi\hat m
  - \dot\nabla\hat\nabla\dot k\phi\hat m
  - \dot\nabla\hat\nabla k\phi\hat{\dot m}.
$$
The dots and hats indicate what each derivative is differentiating.
Now $J$ becomes
$$
  J
= \grade{\int_E|\d^2\vec x|\,\nabla\hat\nabla k\phi\hat m}
    - \grade{
      \int_E|\d^2\vec x|\Bigl[\dot\nabla\hat\nabla\dot k\phi\hat m
      + \dot\nabla\hat\nabla k\phi\hat{\dot m}\Bigr]
    }.
$$
We rearrange the last two terms and find
$$
  \grade{\dot\nabla\hat\nabla\dot k\phi\hat m}
= \phi(\nabla k)\cdot(\nabla m),\quad
  \grade{\dot\nabla\hat\nabla k\dot\phi\hat{\dot m}}
= k\phi(\nabla^2 m),
$$
so now $J$ is
$$
  J
= \grade{\int_E|\d^2\vec x|\,\nabla\hat\nabla k\phi\hat m}
    - \int_E|\d^2\vec x|\Bigl[
      \phi\,(\nabla k){\cdot}(\nabla m) + k\phi(\nabla^2 m)\Bigr
    ].
$$
We now focus on the first term. We split $\nabla$ orthogonally into
$$
  \nabla = \partial + \partial_\perp \qquad\text{where}\quad
  \partial = \vec N\,\vec N{\wedge}\nabla,\quad
  \partial_\perp = \vec N\,\vec N{\cdot}\nabla
$$
and get
$$
  \grade{\int_E|\d^2\vec x|\,\nabla\hat\nabla k\phi\hat m}
= \grade{\int_E|\d^2\vec x|\,\partial\hat\nabla k\phi\hat m}
    + \grade{\int_E|\d^2\vec x|\,\partial_\perp\hat\nabla k\phi\hat m}.
\tag{$*$}
$$
Since the ambient space is $\R^3$, we have $\d^2\vec x\,\vec N = |\d^2\vec x|I$, so the first integral on the RHS becomes
$$
  \grade{\int_E|\d^2\vec x|\,\partial\hat\nabla k\phi\hat m}
= \grade{I^{-1}\int_E\d^2\vec x\,\vec N\partial\hat\nabla k\phi\hat m}
= -\grade{I^{-1}\int_E\d^2\vec x\,\partial\vec N\hat\nabla k\phi\hat m},
$$
where in the last line $\vec N\partial = -\partial\vec N$ since $\vec N$ is orthogonal to $E$ while $\partial$ is tangential. Since $E$ is flat, $\vec N$ is a constant and doesn't need to be differentiated; so the Fundamental Theorem applies as is and we get
$$
  -\grade{I^{-1}\int_E\d^2\vec x\,\partial\vec N\hat\nabla k\phi\hat m}
= -\grade{I^{-1}\int_{\partial E}\d^1\vec x\,\vec N\hat\nabla k\phi\hat m}
= -\grade{I^{-1}\int_{\partial E}|\d^1\vec x|\,\vec E\vec N\hat\nabla k\phi\hat m},
$$
where $\vec E$ is the unit vector along the edge of $\partial E$ being integrated. Since $\vec E$ and $\vec N$ are orthogonal, we can factor $I = \vec E\vec N\vec O$ for $\vec O$ the unit vector orthogonal to $\vec E,\vec N$ within the plane of $E$ (i.e. $\vec O = \vec E\times\vec N$, the cross-product). Then $I^{-1} = \vec O\vec N\vec E$ and
$$
  -\grade{I^{-1}\int_{\partial E}|\d^1\vec x|\,\vec E\vec N\hat\nabla k\phi\hat m}
= -\grade{\int_{\partial E}|\d^1\vec x|\,\vec O\hat\nabla k\phi\hat m}
= -\int_{\partial E}|\d^1\vec x|\,k\phi(\vec O\cdot\nabla m).
$$
Now working on the RHS of ($*$), applying the product rule for $\partial_\perp$ gives
$$
  \grade{\partial_\perp\hat\nabla k\phi\hat m}
= \grade{\partial_\perp k\phi(\nabla m)}
= \grade{\phi(\partial_\perp k)(\nabla m)
    + k(\partial_\perp\phi)(\nabla m)
    + k\phi(\partial_\perp\nabla m)
  }
$$
$$
= \phi\,(\partial_\perp k){\cdot}(\partial_\perp m)
    + k\,(\partial_\perp\phi){\cdot}(\partial_\perp m)
    + k\phi(\partial_\perp^2 m).
$$
Returning to $J$, altogether we now have
$$
  J
= -\int_{\partial E}|\d^1\vec x|\,k\phi(\vec O\cdot\nabla m)
    + \int_E|\d^2\vec x|\Bigl[
      \phi\,(\partial_\perp k){\cdot}(\partial_\perp m)
      + k\,(\partial_\perp\phi){\cdot}(\partial_\perp m)
      + k\phi(\partial_\perp^2 m)
    \Bigr]
    - \int_E|\d^2\vec x|\Bigl[
      \phi\,(\nabla k){\cdot}(\nabla m) + k\phi(\nabla^2 m)
    \Bigr].
$$
In the last term, we have
$$
  (\nabla k)\cdot(\nabla m) = (\partial k)\cdot(\partial m) + (\partial_\perp k)\cdot(\partial_\perp m),\quad
  \nabla^2 m = \partial^2 m + \partial^2_\perp m,
$$
and these cancel with the other terms, giving
$$
  J
= -\int_{\partial E}|\d^1\vec x|\,k\phi(\vec O\cdot\nabla m)
    + \int_E|\d^2\vec x|\Bigl[
      k\,(\partial_\perp\phi){\cdot}(\partial_\perp m)
      - \phi\,(\partial k){\cdot}(\partial m)
      - k\phi(\partial^2 m)
    \Bigr].
$$
The boundary $\partial E = \partial E_1 \cup \partial E_2 \cup \cdots \cup \partial E_n$, where each $\partial E_i$ is an edge of the polygon $E$. Each $\partial E_i$ has it's own (constant) $\vec E_i$ along the edge and $\vec O_i$ orthogonal to it but still in the plane. In coordinates $\vec x = x_i\vec E_i + y_i\vec O_i + z\vec N$, we may write
$$
  J
= -\sum_i\int_{\partial E_i}\d x_i\,k\phi\PD m{y_i}
    + \int_E\d S\Bigl[
      k\PD\phi z\PD mz
      - \phi(\partial k){\cdot}(\partial m)
      - k\phi(\partial^2z)
    \Bigr],
$$
where $\d S$ is the scalar surface measure and $\partial$ is to be interpreted as the 2D gradient within the plane of $E$.
