Showing that surface area is equivalent to $\int_{S}\|\partial_u\phi\times\partial_v\phi\|dudv$, and is there MVT for bijections: $\Bbb R\to\Bbb R^2$? $\newcommand{\d}{\,\mathrm{d}}$It can be shown that arclength, considered as a sum of increasingly fine partitions of the graph, approaches the integral formulation. However, I have only ever seen the surface area integral treated as a definition, or as a truth, and today I wish to prove it - my proof is quite unwieldy, so I'm not sure if it's rigorous or if it can be cleaned up. I also don't study differential geometry, so I'd like to know if my reasoning and terminology is consistent with the field. As soon as I saw the above formula I knew I had to attempt a proof, so here it is:
Definitions:

Let $M$ be a smooth surface if: $M\subset\Bbb R^3$ is bounded and of null Lebesgue measure and there exists a simply connected open $\Omega\subset\Bbb R^2$ together with a $C^2$ bijection $\varphi:\Omega\to M$.
Let $P=\{x_{i,j}\in M:i,j\in\{1,\cdots,N\}\}$ be a partition of a smooth surface $M$ if $\{x_{i,j}:j\in\{1,\cdots,N\}\}$ partitions a contour in $M$ for any fixed $1\le i\le N$ and visa versa. Define the area of $P$ to be the sum of areas of the triangles with vertices $\{x_{i,j},x_{i,j+1},x_{i+1,j}\}$ and $\{x_{i+1,j},x_{i+1,j+1},x_{i,j+1}\}$.

Note: I chose to define it via triangles as triangles are always coplanar, whereas when I attempted this before with quadrilaterals I ran into problems - decomposing the quadrilaterals into two triangles seems to work!

Define the surface area of a smooth surface $M$ to be the limit of the areas of the partitions of $M$ as $\max\{\max_{i,j}\{x_{i+1,j}-x_{i,j}\},\max_{i,j}\{x_{i,j+1}-x_{i,j}\}\}\to0$, if it exists.
If $S\subset\Bbb R^3$ can be decomposed as the disjoint union of finitely many $M_1,M_2,\cdots,M_n$ smooth surfaces with parametristations $(\Omega_1,\varphi_1),\cdots,(\Omega_n,\varphi_n)$, then it has a surface area equivalent to the sum of the surface areas of the $M_k$.

I am a touch worried that there may be pathological surfaces for which this doesn't work but I think this is all ok so far.
Let $Q$ be an arbitrary rectangular partition into points $(u_i,v_i)_{i=1}^N$ of a smooth surface $\Omega$, and put $\Delta Q=\max\{\max_i\{u_{i+1}-u_i\},\max_i\{v_{i+1}-v_i\}\}$. As $\varphi$ is a continuous bijection, $\varphi(Q)$ will be a valid partition of $M=\varphi(\Omega)$. Consider for any fixed $i$ the two triangles described by the points $(u_i,v_j),(u_{i+1},v_j),(u_i,v_{j+1}),(u_{i+1},v_{j+1})$. The corresponding area of the two triangles in the partition of $M$, $P=\varphi(Q)$, will be: $$\frac{1}{2}\|(\varphi(u_{i+1},v_j)-\varphi(u_i,v_j))\times(\varphi(u_i,v_{j+1})-\varphi(u_i,v_j))\|$$And: $$\frac{1}{2}\|(\varphi(u_{i+1},v_{j+1})-\varphi(u_i,v_{j+1}))\times(\varphi(u_{i+1},v_{j+1})-\varphi(u_{i+1},v_j))\|$$
EDIT: As per Matteo's comment, I've now realised that the MVT fails not only for functions with multidimensional domain, but also for functions $\Bbb R\to\Bbb R^n$. I think the bijectivity of $\varphi$ saves the proof, but it will require more work that I cannot straight away think of.
CONT:
As $\varphi$ is $C^2$, one may employ the mean value theorem twice by holding $u$ or $v$ fixed and letting the other vary, making the function defined $\Bbb R\to\Bbb R^3$, so defining $\Delta u_i=u_{i+1}-u_i$ and similarly for $\Delta v_j$ gives: $$\begin{align}\varphi(u_{i+1},v_{j+1})-\varphi(u_i,v_{j+1})&=\varphi(u_{i+1},v_j)-\varphi(u_i,v_j)+\varphi(u_{i+1},v_{j+1})\\&\hspace{10pt}-\varphi(u_{i+1},v_j)-(\varphi(u_i,v_{j+1})-\varphi(u_i,v_j))\\&=\varphi(u_{i+1},v_j)-\varphi(u_i,v_j)+\Delta v_j\cdot(\partial_v\varphi(u_{i+1},\theta_v)-\partial_v\varphi(u_i,\theta_v))\\&=\varphi(u_{i+1},v_j)-\varphi(u_i,v_j)+\Delta v_j\Delta u_i\cdot\partial^2_{uv}(\theta_u,\theta_v)\\\varphi(u_{i+1},v_{j+1})-\varphi(u_{i+1},v_j)&=\cdots\\&=\varphi(u_i,v_{j+1})-\varphi(u_i,v_j)+\Delta u_i\Delta v_j\cdot\partial^2_{vu}\varphi(\vartheta_u,\vartheta_v)\end{align}$$
For mean values $\theta_u,\vartheta_u\in(u_i,u_{i+1}),\,\theta_v,\vartheta_v\in(v_j,v_{j+1})$.
Then the second area becomes: $$\frac{1}{2}\|(\varphi(u_{i+1},v_j)-\varphi(u_i,v_j))\times(\varphi(u_i,v_{j+1})-\varphi(u_i,v_j))\|\\\hspace{80pt}\\\hspace{20pt}+\frac{1}{2}\Delta u_i\Delta v_j\left\{\|\partial^2_{uv}\varphi(\theta_u,\theta_v)\times(\varphi(u_i,v_{j+1})-\varphi(u_i,v_j))\|\\\hspace{40pt}+\|(\varphi(u_{i+1},v_j)-\varphi(u_i,v_j))\times\partial^2_{vu}\varphi(\vartheta_u,\vartheta_v)\|\\\hspace{30pt}+\Delta u_i\Delta v_j\cdot\|\partial^2_{uv}\varphi(\theta_u,\theta_v)\times\partial^2_{vu}\varphi(\vartheta_u,\vartheta_v)\|\right\}$$
Employ the mean value theorem once more with feeling to the expressions of the form $\varphi(u_i,v_{j+1})-\varphi(u_i,v_j)$, using mean values $\lambda_u,\lambda_v$ in the appropriate ranges, and the two areas become:
$$\frac{1}{2}\Delta u_i\Delta v_j\|\partial_u\varphi(\lambda_u,v_j)\times\partial_v\varphi(u_i,\lambda_v)\|$$
And:
$$\frac{1}{2}\Delta u_i\Delta v_j\|\partial_u\varphi(\lambda_u,v_j)\times\partial_v\varphi(u_i,\lambda_v)\|+\\\hspace{20pt}\frac{1}{2}\Delta u_i\Delta v_j\left\{\Delta v_j\cdot\|\partial^2_{uv}\varphi(\theta_u,\theta_v)\times\partial_v\varphi(u_i,\lambda_v)\|+\Delta u_i\cdot\|\partial_u\varphi(\lambda_u,v_j)\times\partial^2_{vu}\varphi(\vartheta_u,\vartheta_v)\|\\\hspace{30pt}+\Delta u_i\Delta v_j\cdot\|\partial^2_{uv}\varphi(\theta_u,\theta_v)\times\partial^2_{vu}\varphi(\vartheta_u,\vartheta_v)\|\right\}$$
The area of the partition $P=\varphi(Q)$ is then expressible as (by summing the two triangles per summand):
$$\sum_{i,j=1}^{N-1}\|\partial_u\varphi(\lambda_u,v_j)\times\partial_v(u_i,\lambda_v)\|\cdot\Delta u_i\Delta v_j\\+\frac{1}{2}\sum_{i,j=1}^{N-1}\left\{\Delta v_j\cdot\|\partial^2_{uv}\varphi(\theta_u,\theta_v)\times\partial_v\varphi(u_i,\lambda_v)\|+\Delta u_i\cdot\|\partial_u\varphi(\lambda_u,v_j)\times\partial^2_{vu}\varphi(\vartheta_u,\vartheta_v)\|\\\hspace{30pt}+\Delta u_i\Delta v_j\cdot\|\partial^2_{uv}\varphi(\theta_u,\theta_v)\times\partial^2_{vu}\varphi(\vartheta_u,\vartheta_v)\|\right\}\cdot\Delta u_i\Delta v_j$$
As $\Delta Q\to0$, the latter sum clearly goes to $0$ as it features quadratic $\Delta$ terms, and by considering the first sum as a Riemann sum over a double partition tagged with the values $\lambda$, in the passage to the limit (if it exists) we find that the surface area is:
$$A(M)=\iint_{\Omega}\left\|\frac{\partial\varphi}{\partial u}\times\frac{\partial\varphi}{\partial v}\right\|\d u\d v$$
And:
$$A(S)=\sum_{k=1}^n\iint_{\Omega_k}\left\|\frac{\partial\varphi_k}{\partial u}\times\frac{\partial\varphi_k}{\partial v}\right\|\d u\d v$$
Is this a correct proof? And is the preliminary definition a sensible one? Furthermore, I expect that the $C^2$ condition can be dropped but I don't see how. It would have certainly been easier to consider quadrilaterals rather than triangles, but my attempts at that failed on the realisation that four points need not be coplanar, and thus need not be a good approximation to the surface - finding the area and creating good partitions was tricky. Hopefully my triangle approach works!
 A: Here is the famous Schwarz example of a triangular mesh inscribed in a finite cylinder and the total surface area goes to infinity as we move the planes of the hexagons closer and closer. I am not sure if this gives a counterexample to your procedure or not, but I'll let you ponder it.

P.S. The bijectivity of the mapping won't rescue the Mean Value Equality. Consider $\phi(t) = (\cos t,\sin t)$. Instead of the domain $[0,2\pi]$, take $[0,2\pi-\epsilon]$. Now $\|\phi(2\pi-\epsilon)-\phi(0)\|<\epsilon$, but $\|(2\pi-\epsilon)\phi'(\xi)\| = 2\pi-\epsilon$ for any $\xi\in (0,2\pi-\epsilon)$, so we cannot have $\phi(2\pi-\epsilon)-\phi(0) = (2\pi-\epsilon)\phi'(\xi)$ for any $\xi$.
A: I think you are bumping into the Schwarz paradox, which is a counter-example demonstrating the failure of  the naive approach to defining surface area by inscribed polyhedra.
The counter-example consists of a circular cylinder that is formed from a rectangular sheet of paper that is partitioned into triangular cells. The polyhedral approximations fail to converge to the expected area unless very specific conditions are imposed on how the partitioning is performed in each of the two directions on the rectangular sheet.
The following is a quote from "On the circumradius condition for piecewise linear triangular elements" by
Kenta Kobayashi · Takuya Tsuchiya
Japan J. Indust. Appl. Math. (2015) 32:65–76
N.B. The mathematical formulas were re-typeset by hand, and they differ slightly in layout from that of the original paper. I hope that no errors were introduced in the process.

As is well understood, the length of a curve is defined as the limit
of the length of the inscribed polygonal edges. Hence, one might think
that the area of a surface could be defined in a similar manner.
Actually, mathematicians in the 19th century believed that the area of
surface is the limit of the areas of inscribed polygonal surfaces. In
the 1880s, Schwarz and Peano independently showed, however, that this
definition does not work.

Let $R$ be a rectangle of height $H$ and width $2πr$.

Let $m, n$ be positive integers. Suppose that this rectangle is
divided into $m$ equal strips, each of height $H/m$. Each strip is
then divided into isosceles triangles whose base length is $2πr/n$.
Then, the cylinder is approximated by the inscribed polygonal surface
which consists of $2mn$ congruent isosceles triangles. Because the
height of ­ each triangle is $\sqrt{(H/m)^2 + r^2(1 − \cos(\pi/n))^2}$
and the base length is $2r \sin(\pi/n)$,
the area  of the inscribed polygonal surface is  ­
$A= 2mnr\sin(\pi/n) \sqrt{ (H/m)^2 +r^2 (1−\cos(\pi/n))^2}$
$= 2πr  \frac{ \sin(\pi/n)}{\pi/n} \sqrt{H^2+ \frac{\pi^4 r^2 }{4} \frac{ m^2}{n^4} (\frac{ \sin(\pi/(2n))}{(\pi/2n)})^4 }$

which tends to the traditional surface area definition $A= 2\pi RH$ if and only if $(m/n^2)\to 0$.
­
