Rosetta Code

Total circles area: Difference between revisions

Page Discussion

Total circles area (view source)

Revision as of 02:52, 23 September 2012

1,001 bytes added ,  11 years ago
Improved C entry (partially from Ledrug code)
(reflow floating thumbnails)
(Improved C entry (partially from Ledrug code))
Line 49:
<lang c>#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <time.h>
#include <stdbool.h>
 
typedef struct { double x, y, r; } CircleFp;
typedef struct { Fp x, y, r; } Circle;
 
const Circle circles[] = {
{ 1.6417233788, 1.6121789534, 0.0848270516},
{-1.4944608174, 1.2077959613, 1.1039549836},
Line 79 ⟶ 83:
{ 0.0152957411, 0.0638919221, 0.9771215985}};
 
const size_t n_circles = sizeof(circles) / sizeof(Circle);▼
double min(const double a, const double b) { return a <= b ? a : b; }▼
 
doublestatic maxinline Fp min(const doubleFp a, const doubleFp b) { return a ><= b ? a : b; }
 
doublestatic uniforminline Fp max(const doubleFp a, const doubleFp b) { return a >= b ? a : b; }
 
static inline Fp sq(const Fp a) { return a * a; }
 
// Return an uniform random value in [a, b).
▲static inline double minuniform(const double a, const double b) { return a <= b ? a : b; }
const double r01 = rand() / (double)RAND_MAX;
return a + (b - a) * r01;
}
 
static inline bool is_inside_circles(const Fp x, const Fp y) {
for (size_t i = 0; i < n_circles; i++)
if (sq(dxx *- dxcircles[i].x) + dysq(y *- dy) <= (circles[ji].ry) *< circles[ji].r)) {▼
return hits++true;▼
return false;
}
 
int main() {
// Initialize the bounding box (bbox) of the circles.▼
▲ const size_t n_circles = sizeof(circles) / sizeof(Circle);
doubleFp x_min = 1/0.0, x_max = -1e1001/0.0;▼
 
Fp y_min = x_min, y_max = x_max;
▲ // Initialize the bounding box of the circles.
double x_min = 1e100;
▲ double x_max = -1e100;
double y_min = 1e100;
double y_max = -1e100;
 
// Compute the bounding box of the circles.
for (size_t i = 0; i < n_circles; i++) {
const Circle *c = &circles[i];
x_min = min(x_min, c.->x - c.->r);
x_max = max(x_max, c.->x + c.->r);
y_min = min(y_min, c.->y - c.->r);
y_max = max(y_max, c.->y + c.->r);
 
c->r *= c->r; // Square the radii to speed up testing.
}
 
const Fp bbox_area = (x_max - x_min) * (y_max - y_min);
 
// Montecarlo sampling.
srand(1time(0));
const size_t n_samplesto_try = 100 * 10001U *<< 100016;
size_t n_tries = 0;
size_t n_hits = 0;
 
size_twhile hits(true) = 0;{
n_hits += is_inside_circles(uniform(x_min, x_max),
uniform(y_min, y_max));
n_tries++;
 
for (size_t i = 0;if i(n_tries < n_samples;== i++to_try) {
const doubleFp xarea = uniform(x_min,bbox_area * n_hits / x_max)n_tries;
const doubleFp yr = uniform(y_min, y_maxFp)n_hits / n_tries;
for (size_t j const Fp s = 0;area j* <sqrt(r n_circles;* j++(1 - r) {/ n_tries);
constprintf("%.4f double+/- dx%.4f =(%zd xsamples)\n", -area, circles[j].xs, n_tries);
constif double(s dy* 3 <= y 1e-3) // Stop at 3 circles[j]sigmas.y;
▲ if ((dx * dx + dy * dy) <= (circles[j].r * circles[j].r)) {
▲ hits++;
break;
}to_try *= 2;
}
}
 
printf("Approximated area: %.8f\n",
(double)(x_max - x_min) *
(double)(y_max - y_min) *
((double)hits / n_samples));
 
return 0;
}</lang>
{{out}}
<pre>21.4498 +/- 0.0370 (65536 samples)
<pre>Approximated area: 21.56262288</pre>
21.5031 +/- 0.0262 (131072 samples)
21.5170 +/- 0.0185 (262144 samples)
21.5442 +/- 0.0131 (524288 samples)
21.5477 +/- 0.0093 (1048576 samples)
21.5531 +/- 0.0065 (2097152 samples)
21.5624 +/- 0.0046 (4194304 samples)
21.5631 +/- 0.0033 (8388608 samples)
21.5602 +/- 0.0023 (16777216 samples)
21.5632 +/- 0.0016 (33554432 samples)
21.5617 +/- 0.0012 (67108864 samples)
21.5628 +/- 0.0008 (134217728 samples)
21.5639 +/- 0.0006 (268435456 samples)
21.5637 +/- 0.0004 (536870912 samples)
21.5637 +/- 0.0003 (1073741824 samples)</pre>
 
=={{header|D}}==
Retrieved from "https://rosettacode.org/wiki/Special:MobileDiff/219421"