Let function $u\in BV(\Omega)$ be a function of bounded variation and $\Omega\subset \mathbb R^2$ be a smooth domain. I know it is possible to approximate function $u$ with polynomials, i.e., $$ u = \sum_{n=1}^{+\infty} P_n, $$ where $P_n$ is of polynomial of order $n$.
My question: would it be possible to design polynomial $P_n$ so that $$ TV(u) = \sum_{n=1}^{\infty}TV(P_n), $$ where by $TV(u)$ we mean the total variation of $u$.
The answer is no. Orthogonality is with respect to the $L^2$ norm and the $TV$ norm (of a smooth function) is the $L^1$ norm of the derivative and for this norm the $L^2$ orthogonality does not mean much.
Let $\Vert\cdot\Vert_1$ denote the $L^1$ norm on $\Omega$. Smooth functions are dense in BV (Theorem 2, Section 5.2.2 in [1]) so I think you can approximate u by polynomials $Q_n$ and then you can write $P_n=Q_n-Q_{n-1}$. However, for smooth functions, and in particular for polynomials, $TV(P)=|DP|(\Omega)=\Vert \nabla P\Vert_1$.
Let $P$ and $Q$ be polynomials. If $\nabla P$ is not parallel to $\nabla Q$ at every point of $\Omega$, then $\Vert \nabla P + \nabla Q\Vert_1<\Vert \nabla P\Vert_1 + \Vert\nabla Q\Vert_1$. Note that the gradients are parallel iff $\nabla P = \lambda\nabla Q$, $\nabla (P-\lambda Q)=0$, $P=\lambda Q+c$ which is a very restrictive condition. Therefore, unless you have this linear relation between all polynomials you have: $$ TV(u)=|Du|(\Omega)=\lim_{k\to\infty}\Vert\sum_{n=1}^k \nabla P_n\Vert_1<\sum_{n=1}^\infty \Vert\nabla P_n\Vert_1=\sum_{n=1}^\infty TV(P_n). $$ Therefore you cannot get the equality that you want, but instead you get a sharp inequality.
A good introduction to BV functions is Chapter 5 in:
[1] L. Evans, R. F. Gariepy, Measure theory and fine properties of functions. Revised edition. Textbooks in Mathematics. CRC Press, Boca Raton, FL, 2015.
