1 /*
2 * Licensed to the Apache Software Foundation (ASF) under one or more
3 * contributor license agreements. See the NOTICE file distributed with
4 * this work for additional information regarding copyright ownership.
5 * The ASF licenses this file to You under the Apache License, Version 2.0
6 * (the "License"); you may not use this file except in compliance with
7 * the License. You may obtain a copy of the License at
8 *
9 * http://www.apache.org/licenses/LICENSE-2.0
10 *
11 * Unless required by applicable law or agreed to in writing, software
12 * distributed under the License is distributed on an "AS IS" BASIS,
13 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14 * See the License for the specific language governing permissions and
15 * limitations under the License.
16 */
17
18 package org.apache.commons.rng.examples.quadrature;
19
20 import org.apache.commons.rng.UniformRandomProvider;
21
22 /**
23 * <a href="https://en.wikipedia.org/wiki/Monte_Carlo_integration">Monte-Carlo method</a>
24 * for approximating an integral on a n-dimensional unit cube.
25 */
26 public abstract class MonteCarloIntegration {
27 /** RNG. */
28 private final UniformRandomProvider rng;
29 /** Integration domain dimension. */
30 private final int dimension;
31
32 /**
33 * Simulation constructor.
34 *
35 * @param rng RNG.
36 * @param dimension Integration domain dimension.
37 */
38 public MonteCarloIntegration(UniformRandomProvider rng,
39 int dimension) {
40 this.rng = rng;
41 this.dimension = dimension;
42 }
43
44 /**
45 * Run the Monte-Carlo integration.
46 *
47 * @param n Number of random points to generate.
48 * @return the integral.
49 */
50 public double integrate(long n) {
51 long inside = 0;
52 for (long i = 0; i < n; i++) {
53 if (isInside(generateU01())) {
54 ++inside;
55 }
56 }
57
58 return inside / (double) n;
59 }
60
61 /**
62 * Indicates whether the given points is inside the region whose
63 * integral is computed.
64 *
65 * @param point Point whose coordinates are random numbers uniformly
66 * distributed in the unit interval.
67 * @return {@code true} if the {@code point} is inside.
68 */
69 protected abstract boolean isInside(double... point);
70
71 /**
72 * @return a value from a random sequence uniformly distributed
73 * in the {@code [0, 1)} interval.
74 */
75 private double[] generateU01() {
76 final double[] rand = new double[dimension];
77
78 for (int i = 0; i < dimension; i++) {
79 rand[i] = rng.nextDouble();
80 }
81
82 return rand;
83 }
84 }