Given a manifold $M$, a differential $k$-form $\omega$ assigns to each point $p$ $\epsilon$ $M$ a $k$-covector $\omega_p$ $\epsilon$ $\bigwedge^k \left(T_p^*M\right)$, where $\bigwedge^k \left(T_p^*M\right)$ is the space of alternating $k$-tensors on the tangent space $T_pM$.
If $\omega = \sum f_i dx^i$ denotes a 1-form of $M$, then the exterior derivative $d:\bigwedge^k \left(T_p^*M\right) \rightarrow \bigwedge^{k+1} \left(T_p^*M\right)$ acts on $\omega$ to give a 2-form:
$dw = \sum \frac{\partial f_i}{\partial x^j} dx^j \wedge dx^i$.
I wonder if there is an analogous construction of "symmetric" $k$-forms. By a symmetric $k$-form I mean a map $\sigma$ that assigns to each point $p$ $\epsilon$ $M$ a $k$-covector $\sigma_p$ that is a symmetric $k$-tensor of the tangent space $T_pM$. If $\sigma = \sum g_i dx^i$, maybe the exterior derivative of this symmetric 1-form could be
$d\sigma = \frac{\partial g_i}{\partial x^j} dx^j \vee dx^i$,
where $dx^j \vee dx^i$ could be analogous to the wedge product with the property $dx^j \vee dx^i = dx^i \vee dx^j$.
Do you have any idea if such a construction exists? Also, could there be an analogous to the de Rham complex?
Not without more structure, no. Computing directly the transformation of the map under a change of coordinates shows that the proposed map in fact depends on the choice of coordinates.
If a manifold $M$ (of general dimension) is equipped with a linear connection $\nabla$, however, one gets for free a map $$\Gamma(\operatorname{Sym}^k T^*M) \to \Gamma(\operatorname{Sym}^k T^*M \otimes T^*M), \qquad \Phi \mapsto \nabla \Phi ,$$ and symmetrizing gives a manifestly coordinate-independent map $$D : \Gamma(\operatorname{Sym}^k T^*M) \to \Gamma(\operatorname{Sym}^{k + 1} T^* M), \qquad \Phi \mapsto \operatorname{Sym} \nabla \Phi .$$ In abstract index notation, $D$ is $\Phi_{i_1 \cdots i_k} \mapsto \Phi_{(i_1 \cdots i_k, i_{k + 1})}$; in index-free notation, $D$ is $$D(\Phi)(X_1, \ldots, X_{k + 1}) = \frac1{k + 1} [(\nabla_{X_1} \Phi)(X_2, \ldots, X_{k + 1}) + \cdots + (\nabla_{X_{k + 1}} \Phi)(X_1, \ldots, X_k)] ;$$ in a coordinate frame, $D$ is $$D(\Phi)_{i_1 \ldots i_{k + 1}} = \partial_{(i_1} \beta_{i_2 \cdots i_{k + 1})} - k \Gamma_{(i_1 i_2}^j \Phi_{i_3 \cdots i_{k + 1})j}^{\phantom{j}},$$ where the $\Gamma_{ij}^l$ are the Christoffel symbols of $\nabla$. In particular, in contrast to the exterior derivative, $d$:
The first two operators are of special interest:
Provided that $\nabla$ is torsion-free, for any function $f$ on $M$ the Hessian $\nabla^2 f$ is symmetric, so it coincides with $D^2 f$. If $\nabla$ is moreover the Levi-Civita connection of a Riemannian metric $g$, $\operatorname{tr}_g D^2 : C^\infty(M) \to C^{\infty}(M)$ is the Laplace–Beltrami operator, $\Delta_g$.