How to find the inverse metrics? I know one can calculate the inverse of metric tensor $g$ in coordinates as the inverse of it's matrix $g_{ij}$. However what I really liked about differential geometry is how one can actually avoid writing matrices, one just use expressions like
$$g = dx \otimes dx + dy \otimes dy + \frac{1}{2} \left( dx \otimes dy + dy \otimes dx \right)$$
Honestly, in the case of matrices I was never able to remember whether $i$ in $g_{ij}$ is a row or a column, though for a metrics it doesn't matter since it is symmetric. The next best thing is that all the skills needed to do many of the transformations are in essence polynomial handling and algebraic differentiation.
So I ask if there is a technique to reverse the tensor above without dealing with the matrix explicitly? Just some fancy rewritings, substitutions, applications to certain vectors etc.

P.S. I've just realized that the way I wrote here $g$ can be viewed as a sparse matrix representation, whereas the ordinary notation as a table is more like a dense representation.
 A: The inverse metric is defined by the property that:
$$ g\circ g^{-1} = \delta $$
where the Kronecker delta is defined to be the section of $T^*M\otimes TM$ that is acts as the identity on $T^*M$:
$$ \delta := \sum_{i = 1}^{n} \mathrm{d}x^i \otimes \partial_{x^i} $$
In other words, if there were a general way to compute the inverse metric without computing the inverse matrix, then the method can be applied to the case of a constant metric on a linear space (that is, a simple inner product space) to compute the inverse matrix of a matrix without, erh, computing the inverse matrix. (And I think you see the problem with this...)

Now, in your special case your metric can be written as the following:
$$ g = \frac12 \left( \mathrm{d}x+\mathrm{d}y\right)\otimes \left(\mathrm{d}x + \mathrm{d}y\right) + \frac12 \left( \mathrm{d}x\otimes\mathrm{d}x + \mathrm{d}y \otimes \mathrm{d}y\right) $$
the second half is the "identity matrix". So we see that your metric can be easily put into diagonal form
$$ g = \frac34 \mathrm{d}(x + y)^2 + \frac14\mathrm{d}(x - y)^2 $$
for which inversion is easy. 

The OP asks:

But every symmetric matrix is diagonalizable, right? Then one can try to rewrite the metrics in such a simple form.

Yes, every symmetric matrix is diagonalisable. But that is besides the point. When we speak of symmetric matrices being diagonalisable, we are usually speaking in the context of spectral theorem of diagonalising an operator $A:V\to V$ from a vector space $V$ to itself. Diagonalisation is the process of choosing a basis of $V$ consisting of only eigenvectors of $A$. The statement "$A$ is symmetric" is interpreted as: "$A$ is symmetric with respect to some positive definite inner product on $V$". 
In this case diagonalisation of $A$ means writing $A$ as
$$ A = \sum_{i = 1}^K \lambda_i e_i \otimes (e_i)^*~. $$
In the context of the metric tensor you are dealing with a "symmetric bilinear form" on the vector space $V$. It is a simple algebraic fact that 

If $B$ is a symmetric bilinear form on a vector space $V$, there exists a basis of $V$ such that $B$ is diagonal. 

In this context $B$ is written as
$$ B = \sum_{i = 1}^K \nu_i (e_i)^* \otimes (e_i)^*~.$$
The algebraic content of the two statements are very, very different. 
With that said: given an arbitrary (Riemannian) metric tensor field
$$ g = \sum g_{ij}(x) \mathrm{d}x^i \otimes \mathrm{d}x^j $$
by the above fact it is true that you can find a family of vector fields $\vec{E_\mu} = \sum_i E_\mu^i \partial_i$ such that
$$ g(E_\mu,E_\nu) = \delta_{\mu\nu} $$
However, the vector fields are constructed at each point $x$ by finding an orthogonal basis for $g$ at the point $x$, and there are two major problems:


*

*The vector fields need not be continuous (this, at least, in practice can be solved locally). 

*The vector fields need not be integrable: it is a necessary condition that a family of coordinate vector fields $\partial_i$ be commutative: that is $[\partial_i,\partial_j] = 0$ where $[,]$ is the Lie bracket. (This is just the statement that coordinate partial derivatives must commute.) In general there is absolutely no reason that your chosen vector fields $E_\mu$, or some $x$-dependent multiples thereof, will satisfy this condition. 


Given a (pseudo)Riemannian metric, whether one can find a local coordinate system in which the metric can be diagonalised is a very classical problem. In two dimensions, smooth Riemannian and Lorentzian metrics can always be diagonalised: the former by isothermal (conformal) coordinates, the latter by null coordinates. It is also know that in three dimensions any Riemannian metric can be diagonalised locally. In four and higher dimensions, metrics which admit a coordinate system in which they are diagonal must force certain algebraic conditions on the Weyl curvature tensor and its derivatives. For references on these facts and more, please see this paper of Paul Tod's. 
Even in the two dimensional case the question can be difficult: finding isothermal coordinates given arbitrary smooth metric is a hard problem in general. You case admits a simple solution because your metric is "constant" in the coordinates you gave, which automatically implies that it is flat, and hence admits a linear transformation of the coordinates making it isometric to the standard Euclidean space. 
