Preliminary remark. I am not certain whether Bendixon-Dulac grants the global attractiveness of the equilibria. However, via sheer leveraging on the (strict) monotonicity of $g$ and symmetry of the ODE, we can prove that the equilibria is indeed the global attractor which leaves no room for limit cycles or other nontrivial attractors.
The right region in the phase space to be studied is the invariant set $\mathcal{I}=\left\{\left(x,y\right)\in \left[0,1\right]^2\,:\,g(x)\leq y \leq g^{-1}(x)\mbox{ or }g^{-1}(x)\leq y \leq g(x)\right\}$, i.e., the region between the graphs $g$ and $g^{-1}$. This is the right region in the sense that trajectories necessarily accumulate onto it (exponentially fast). Further, observe that $\mathcal{I}$ contains the equilibria given by $\mathcal{E}=\left\{(x,y)\in\left[0,1\right]^{2}\,:\,y=g(x)=g^{-1}(x)\right\}$.
Let $\left(x(t),y(t)\right)_{t\geq 0}$ be the solution to the ODE with initial condition $(x(0),y(0))\in\left[0,1\right]^2$. Define $T$ as the hitting time to hit the set $\mathcal{I}$, i.e., $T\overset{\Delta}=\inf\left\{T\geq 0\,:\, (x(T),y(T))\in \mathcal{I}\right\}$. Let ${\sf d}(w,z)\overset{\Delta}=\|w-z\|_2$.
Claim $1$. If $T=\infty$, then ${\sf d}((x(t),y(t)),\mathcal{E})\overset{t\rightarrow\infty}\longrightarrow 0$ (exponentially fast).
Proof. Define $f_1(t)=\frac{1}{2}\left(x(t)-g(y(t))\right)^2$ and $f_2(t)=\frac{1}{2}\left(x(t)-g^{-1}(y(t))\right)^2$. Assume that $(x(0),y(0))$ is in the lower triangular part of the phase space, i.e., $y(0)<g(x(0))$ and $y(0)<g^{-1}(x(0))$ -- everything that follows applies similarly if we assume an initial condition at the upper triangular part. If $t<T$, then we have that
$$f'_1(t)=(x(t)-g(y(t)))(\dot{x}(t)-g'(y(t))\dot{y}(t))=(x(t)-g(y(t)))(g(y(t))-x(t)-g'(y(t))(g(x(t))-y(t)))\leq -(x(t)-g(y(t)))^2,$$
where the inequality follows since $-g'(y(t))(x(t)-g(y(t)))(g(x(t))-y(t))<0$ as $(g(x(t))-y(t))>0$, $-g'(y(t))>0$ and $(x(t)-g(y(t)))<0$ since $y(t)<g^{-1}(x(t))$.
Therefore, from Grönwall's inequality, $f_1(t)\leq f_1(0) e^{-2t}$. If $T=\infty$, then $f_1(t)\overset{t\rightarrow \infty}\longrightarrow 0$ exponentially fast. Similarly, if $T=\infty$, we conclude that $f_2(t)\overset{t\rightarrow \infty}\longrightarrow 0$ exponentially fast. This is equivalent to ${\sf d}((x(t),y(t)),\mathcal{E})\overset{t\rightarrow \infty}\longrightarrow 0$.
Claim $2$. If $T<\infty$, then ${\sf d}((x(t),y(t)),\mathcal{E})\overset{t\rightarrow\infty}\longrightarrow 0$.
Proof. Now, $(x(T),y(T))\in \mathcal{I}$. Let $g^{-1}(x(T))\leq y(T) \leq g(x(T))$ -- the other case $g(x(T))\leq y(T) \leq g^{-1}(x(T))$ can be dealt with similarly. If $g^{-1}(x(T)) = y(T) = g(x(T))$, then $(x(T),y(T))$ is already at equlibrium. Assume $g^{-1}(x(T))\leq y(T) < g(x(T))$. Remark that $g^{-1}(x(t))\leq y(t) \leq g(x(t))$ for all $t\geq T$. Let us refer to this invariant set as $\mathcal{I}_1$. Further, let $(x^{\star},y^{\star})$ be the equilibrium that lies in the left part of $\mathcal{I}_1$, i.e., $x^{\star}=\sup\limits_{(w_1,w_2)\in \mathcal{E}} w_1 < x(T)$. Then, $V((x,y))\overset{\Delta}= \frac{1}{2}\left(\left(x-x^{\star}\right)^2+\left(y-y^{\star}\right)^2\right)$ conforms to a Lyapunov function granting attractiveness to the equilibrium for any $(x(T),y(T))\in\mathcal{I}_1$: $V$ is definite positive and $\dot{V}(x,y)<0$ for all $(x,y)\in \mathcal{I}_1\setminus \left\{{\sf eq}_{{\sf right}} \right\}$, where ${\sf eq}_{{\sf right}}$ is the equilibrium on the right side. In words, $g^{-1}(x(T))\leq y(T) < g(x(T))$ implies that the solution will acumulate onto the left closest equilibrium, whereas $g(x(T))\leq y(T) < g^{-1}(x(T))$ will imply convergence to the right closest equilibrium.
Claims $1$ and $2$ combined yield the global attractiveness of the equilibria.
Theorem $3$. [$\mathcal{E}$ is the global attractor] ${\sf d}((x(t),y(t)),\mathcal{E})\overset{t\rightarrow\infty}\longrightarrow 0$ regardless of the initial condition $\left(x(0),y(0)\right)\in\left[0,1\right]^2$.
