I think I figured this out.
The bottom line is that if you transform the problem into the language of vector fields according to the correspondance described in example 1.2, then this is all well known classical stuff.
I found it a little easier to do the case R² - 0 first.
So we are identifying 1-forms and vector fields in the following way: \alpha=Mdx+Ldy \leftrightarrow \mathbf{F}=M\hat{x}+L\hat{y}. Under this identification, integrating a one form over a curve is the same as computing the work of the vector field F over that curve:
\int_{\gamma}\alpha = \int_{\gamma}\mathbf{F}\cdot d\mathbf{r}
Also, if \alpha is closed, then it means that
\frac{\partial L}{\partial x}-\frac{\partial M}{\partial y}=0
and in particular, by Green's thm, if C is a contractible closed curve (i.e. here, any curve not enclosing the origin), then
\int_{C}\mathbf{F}\cdot d\mathbf{r}=0
And \alpha is exact iff \mathbf{F}=\nabla\phi for some scalar function \phi.
Note also that there is in this setting a particular case of DeRham's theorem that is well known, and that is that a vector field \mathbf{F} is the gradient of some function iff
\int_C\mathbf{F}\cdot d\mathbf{r}=0
for all closed curves C. (That is, iff F is conservative.)
This said, notice that there is at least one closed 1-form on R² - 0 that is not exact, and that is the 1-form \alpha represented by the vector field of length 1/r curling around 0 in the anti-clockwise direction:
\mathbf{F}=\frac{y\hat{x}+x\hat{y}}{x^2+y^2}
It is straightforward to check that it is closed, and it is not exact because obviously, for C_0 the unit circle around the origin oriented in the anti-clockwise direction,
\int_{C_0}\mathbf{F}\cdot d\mathbf{r}=2\pi \neq 0
Very well, so H^1_{DR}(R^2 - 0) is not trivial and obviously, for each real number r, [r\alpha] is a different cohomology class.
Now, suppose \omega is a closed but inexact 1-form. Since the integral of 1-form around a loop depends only on the winding number of the loop around 0, then it must be that
\int_{C_0}\omega \neq 0
In particular, there exists a real number r with
\int_{C_0}\omega =\int_{C_0}r\alpha
And so for C any loop of winding number n about 0, we have
\int_{C}\omega =n\int_{C_0}\omega=n\int_{C_0}r\alpha=\int_{C}r\alpha
So
\int_C(\omega-r\alpha) = 0
which is to say, \omega -r\alpha is an exact 1-form. Thus, the cohomology classes r\alpha are the only things in the first cohomology group: H^1_{DR}(R^2 - 0) = \mathbb{R}\left<[\alpha]\right>\approx \mathbb{R}.
For the other groups, it's easy:
-H^k(R² - 0) is 0 for k>2
-H^0(R² - 0) are the functions whose derivatives vanish everywhere, i.e. the constant functions. So H^0(R² - 0)=R.
-H²(R² - 0) is the space of all 2-forms fdxdy on R² - 0 modulo the exact ones, which can be identified with the space of all smooth maps f on R² - 0 modulo those of the form
f = \frac{\partial L}{\partial x}-\frac{\partial M}{\partial y}
for some smooth maps M and L. But all smooth maps are of this form since one can always take M(x,y)=0 and
L(x,y):=\int_{a}^x f(x',y)dx'
(a =\=0). Thus H²(R² - 0)=0.
Now, for R² - p - q (or for R² - {any finite set of points}), it is the same thing, except in the computation of the first cohomology group, the argument is that H^1 now has at least two distinct non zero elements. Namely, those 1-forms \alpha,\beta corresponding to the vector fields rotating anti-clockwise around p and q respectively.
And now a closed non exact 1-form \omega and C_p, C_q loops around p and q of winding number 1, there are real numbers r and s such that
\int_{C_p}\omega =\int_{C_p}r\alpha
\int_{C_q}\omega =\int_{C_q}s\beta
And so for C any loop of winding number n about p and m about q, we have
\int_{C}\omega =n\int_{C_p}\omega+m\int_{C_q}\omega =n\int_{C_p}r\alpha + m\int_{C_q}s\beta = \int_{C}r\alpha+ \int_{C}s\beta,
so \omega is cohomologous to r\alpha + s\beta and thus H^1 is the free vector space generated by [\alpha] and [\beta] : H^1(\mathbb{R}^2-p-q)\approx \mathbb{R}^2.