This is really a question about numerical integration because it reduces immediately to the problem of finding the mean distance between two line segments in the plane.
Monte-Carlo integration would be wasteful and inefficient compared to other methods that are readily available, including
Direct analytic integration. This is messy but can be done.
Numeric quadrature. Adaptive methods will be accurate.
Analytic integration of the mean distance between a point and a line segment and numeric quadrature of that mean distance.
Discrete approximations, such as Riemann sums, the Trapezoidal Rule, etc.
Because the mean distance function in (3) is differentiable and varies relatively slowly, simple numerical methods such as Simpson's Rule work well. It's best to do the numeric integration over the shorter of the two line segments. In experiments I find that using the endpoints and midpoint in Simpson's Rule typically gives better than 1% accuracy and using five equally spaced points is far better than that. (The toughest cases occur where two sides are equally long, parallel, and close together.) If you need to, you can do the analysis to estimate the Simpson's Rule error and adaptively decrease the point spacing.
For the original problem, the approximations sometimes overestimate and sometimes underestimate the mean distance, so we could expect some cancellation of errors.
At any rate, you will definitely achieve 5% accuracy using method (3). For an $n$ sided polygon it requires $n$ computations (each of constant cost) per vertex and $n-1$ computations per midpoint for a cost of $O(n^2)$. The individual computations of mean point-segment distances require some basic arithmetic, four square roots, and two logarithms, so they're quick and easy. (Mathematica computes about 7,500 mean segment-segment distances per core per second.)
If you want to avoid all but the simplest arithmetic, you can approximate each segment-segment integral with a double Simpson's Rule pass, one for each segment (method 4). With three division points (requiring the computation of nine distances) you will achieve 5% accuracy. With five division points (25 distances) you will be better than 1% accurate (usually better than 0.05%) except for extremely close parallel segments. In Mathematica, with its interpreted overhead, this technique actually is only half as fast as method (3).
This plot shows (negative) log (base 10) relative errors for 1000 segments randomly distributed near the endpoint of a unit segment (which is the larger of the two). E.g., a value of 2 is a 1% error; a value of 6 is a 0.0001% error. This is a severe test, because the worst performances occur when the segments are crossing or close to doing so, situations that are impossible or unlikely in convex polygons. (Note that the smallest possible mean distance is 1/4 (for an infinitesimally short segment near the middle of the unit segment).)
Where the approximation is high, the error is plotted on the right, and where the approximation is low, the error is plotted on the left. Blue dots are the three-point by three-point Simpson's Rule calculation (method 4) and red dots are the three-point Simpson's Rule calculation (method 3).
Typically method 3 (red) has much better than 2% accuracy and is better than method 4 (blue), which in a few cases has almost 10% error. Method 3 has a tendency to overestimate slightly (more of its values are at the right). Both methods almost always obtain better than 0.1% accuracy (a height of 3 or better on the plot), especially at mean distances greater than the longer side. This implies that if you want to improve the accuracy, focus first on nearby edges of a convex polygon.
In this plot of the 1000 segments, the brighter, thicker, redder ones have the greatest relative error for Method 3. The black segment is the reference (unit) segment.
