Relationship between $(n - 1)$ forms and flux of a vector field across a hypersurface

I am currently studying about differential forms and want to deduce the Divergence Theorem (the one in $\mathbb{R}^n$) from the general Stokes' Theorem, which is obtained by taking $$\omega = \sum_{i=1}^n (-1)^{i+1}F_idx_1\wedge\dots\widehat{dx_i}\dots,\wedge dx_n$$ However, I also want to find the relation between $\omega$ and the flux of $F$ through a surface $\Sigma$. Note we defined flux using the integral $$\text{flux} = \int_\Sigma F\cdot \hat{n} $$ It seems to me that if $g$ is a parameterization of $\Sigma$, then $g^* \omega=F\cdot \hat{n}\ \text{dvol}_\Sigma$, but I cannot prove this. I tried writing the normal explicitly as a cross product

$$ N=\det \begin{bmatrix} e_1 & | & |& |\\ |& \frac{\partial g}{\partial u_1} & \cdots & \frac{\partial g}{\partial u_{n-1}}\\ e_n & | & | & | \end{bmatrix},\ \hat{n}=\frac{N}{||N||}$$

Then if we dot product with $F$, we do get by opening the determinant a sum of the form $$F\cdot \hat{n} = \sum_{i=1}^n (-1)^{i+1}F_i\circ g \cdot \text{det of a weird minor}$$ I couldn't get any further though.

EDIT: I should point out that my knowledge volume forms is basic and stems from the definition

$$\text{vol}_M(x)(v_1,\dots,v_k) = \varepsilon \text{vol}_k (v_1,\dots,v_k) ,\ \forall v_i \in T_xM$$ where $\varepsilon$ is chosen such that $(v_1,\dots,v_k ; \varepsilon)$ is a positivly oriented frame. I also know that if I pull back a volume form I get $\sqrt{\det Dg^T Dg} du_1 \wedge \dots \wedge du_k$

Of course you mean the flux of $\vec F$ across an oriented hypersurface $\Sigma$. You don't need any parametrization to do this.

If $\vec n$ is the unit outward normal vector to $\Sigma$, then the "area" $(n-1)$-form for $\Sigma$ is given by $$\sigma = \sum (-1)^{i-1}n_i dx_1\wedge\dots\wedge\widehat{dx_i}\wedge\dots\wedge dx_n.$$ To see this, note that for any $v_1,\dots,v_{n-1}$ in the tangent plane of $\Sigma$ at a point, we have $\sigma(v_1,\dots,v_{n-1}) = \det(\vec n,v_1,\dots,v_{n-1})$, and this is the $(n-1)$-dimensional (signed) volume of the (oriented) parallelepiped spanned by $v_1,\dots,v_{n-1}$. Next, write $\vec F = (\vec F\cdot\vec n)\vec n + \vec G$, and note that \begin{multline*} (\vec F\cdot\vec n)\sigma(v_1,\dots,v_{n-1}) = \det((\vec F\cdot\vec n)\vec n,v_1,\dots,v_{n-1}) \overset{(*)}= \det(\vec F,v_1,\dots,v_{n-1}) \\ = \big(\sum (-1)^{i-1}F_i dx_1\wedge\dots\wedge\widehat{dx_i}\wedge\dots\wedge dx_n\big)(v_1,\dots,v_{n-1}), \end{multline*} as desired. The equality ($*$) follows because $\det(\vec G,v_1,\dots,v_{n-1})=0$, noting that we have $n$ vectors in the $(n-1)$-dimensional tangent plane of $\Sigma$.

Thus, $\sum (-1)^{i-1}F_i dx_1\wedge\dots\wedge\widehat{dx_i}\wedge\dots\wedge dx_n$ gives the correct $(n-1)$-form whose integral over $\Sigma$ is the flux of $\vec F$ across $\Sigma$.

Hint It's perhaps easier to see this in a coordinate-free way (which thus establishes the Divergence Theorem for any ambient manifold $M$).

Along $\Sigma$ decompose $F$ into parts normal to and tangent to $\Sigma$, respectively, $F\vert_\Sigma = F^\perp + F^\top $, so that $$ {\iota_F \operatorname{vol}_M} \vert_\Sigma = {\iota_{F^\perp} \operatorname{vol}_M} \vert_\Sigma + {\iota_{F^\top} \operatorname{vol}_M} \vert_\Sigma .$$ The second term on the right is zero. To handle the first term on the right, use the definition of orthogonality to write $F^\perp$ in terms of $F\vert_\Sigma$ and $\hat n$.

In our case, $M = \Bbb R^n$ is endowed with the Euclidean metric and the standard coordinates and orientation, giving $\operatorname{vol}_M = dx^1 \wedge \cdots \wedge dx^n$, and so $\omega = \iota_F \operatorname{vol}_M$.