According to Cauchy-Davenport Theorem, if $A,B$ are subsets of a prime field ($F_p$) then we have the following bound on the number of elements within the sumset $A + B = \left\{ {\left. {a + b} \right|a \in A,b \in B} \right\}$:
$$\min \left( {p,\left| A \right| + \left| B \right| - 1} \right) \le \left| {A + B} \right|$$
(There is also an extension of this this theorem, by Károlyi, see here: https://mathoverflow.net/questions/133408/karolyis-theorem-for-finite-groups-and-its-extensions)
My questions are:
Can we say something on average about the number of elements within the sumset $A+B$? that is, what is the probability that $\left| {A + B} \right|=1,2,...,p$? can it be extended to more than two subsets?
In addition to the previous question, can we say that the expectation of the number of the elements is larger/smaller than the expectation when assuming that the $\left| A \right| \cdot \left| B \right|$ sums within the sumset are drawn at uniform? according to some simulations I did it appears to be larger.
Not really a useful answer, just too long to fit into a comment. I was playing a bit with this, and observed the following.
Let us fix $m=|A|$, $n=|B|$ and $p$. Let us fix an element $x\in G=\Bbb{Z}_p$. We calculate the probabilities of $x$ being included in the subsets resulting from the two processes.
Process 1: We draw randomly, with repetition $mn$ elements from $G$. A single attempt will fail to give us $x$ with probability $1-1/p$. Thus the probability that we get $x$ in $mn$ trials is $$ P_1(m,n,p)=1-\left(1-\frac1p\right)^{mn}. $$
Process 2: We select a random subset $A\subset G$ with $m$ elements and another random subset $B\subset G$ with $n$ elements. We record whether $x\in A+B$ or not. Here $x\in A+B$ if and only if the sets $A$ and $x-B=\{x-b\mid b\in B\}$ intersect non-trivially. The set $x-B$ is also a random subset of $G$ with $n$ distinct elements. So we fail to include $x$ in $A+B$, if and only if the elements of $x-B$ were all drawn from the complement of $A$. That complement has $p-m$ elements, so we have ${p-m\choose n}$ ways of drawing the elements of $B$ from the complement. Thus $$ P_2(m,n,p)=1-\frac{p-m\choose n}{p\choose n}=1-\frac{(p-m)(p-m-1)\cdots(p-m-n+1)}{p(p-1)(p-2)\cdots (p-n+1)} $$ is the probability of including $x$ in process 2.
A partial explanation to your observations may be that it is always this the case that $$P_2(m,n,p)\ge P_1(m,n,p).$$
This follows from Bernoulli's inequality: $$ \begin{aligned} \left(1-\frac1p\right)^{mn}&=\prod_{j=0}^{n-1}\left(1-\frac1p\right)^m\\ &>\prod_{j=0}^{n-1}\left(1-\frac1{p-j}\right)^m\\ &>\prod_{j=0}^{n-1}\left(1-\frac{m}{p-j}\right)= \prod_{j=0}^{n-1}\left(\frac{p-m-j}{p-j}\right). \end{aligned} $$ This does not answer your question, because while the probabilities for various elements $x$ being included in process 1 are independent from each other, there is a built-in dependence in process 2. I cannot quantify how this will affect the expected sizes of the two resulting sets.
How pronounced were the differences that you observed? To what extent can you explain them with this?