Let the vectors $\vec{w}, \vec{v} \in{\mathbb{R}^n}$ to be given. Further let $\vec{p}$ the vector projection of $\vec{w}$ on $\vec{v}$:
$$ \vec{p} = \vec{w}+\vec{x} = r\vec{v} $$
The unknown $\vec{x}$ and $r$ can be resolved by:
$$ r\vec{v}-\vec{w}=\vec{x} \\ (r\vec{v}-\vec{w})\cdot\vec{v}=\vec{v}\cdot\vec{x} \\ r\vec{v}^2-\vec{v}\cdot\vec{w}=0 \\ r\vec{v}^2=\vec{v}\cdot\vec{w} \\ r=\frac{\vec{v}\cdot\vec{w}}{\vec{v}^2} \implies \vec{x}=\frac{\vec{v}\cdot\vec{w}}{\vec{v}^2}\cdot\vec{v}-\vec{w} $$
But it seems that $\vec{v}\cdot\left(\frac{\vec{v}\cdot\vec{w}}{\vec{v}^2}\cdot\vec{v}-\vec{w}\right)\neq0$ for $\vec{w}=\begin{pmatrix}1\\1\end{pmatrix}$ and $\vec{v}=\begin{pmatrix}2\\1\end{pmatrix}$ using the following code:
import numpy as np
w, v = np.array([1,1]), np.array([2,1])
p = np.dot(v,w)/np.dot(v,v) * v
x = p-w
print(np.dot(v,x))
Where is my misconception?
Running your code, I get $0$ to within a rounding error. The problem is you divide by $v^2=5$, which isn't a power of two. This introduces a very small rounding error because computers use binary. If you got something much larger in size than my result of
-1.1102230246251565e-16, the explanation will still lie in software, not in your understanding of linear algebra.