From 4a0994acebd9e921c552e0430bdd37d77648e0a1 Mon Sep 17 00:00:00 2001 From: Parric007 Date: Mon, 13 May 2024 12:09:59 +0200 Subject: [PATCH] Montecarlo Matlab --- Bounds.m | 27 ++++++++++++++++++++++++++ LinearBounds.m | 25 ++++++++++++++++++++++++ MonteCarloMatlab.m | 47 ++++++++++++++++++++++++++++++++++++++++++++++ Point.m | 12 ++++++++++++ 4 files changed, 111 insertions(+) create mode 100644 Bounds.m create mode 100644 LinearBounds.m create mode 100644 MonteCarloMatlab.m create mode 100644 Point.m diff --git a/Bounds.m b/Bounds.m new file mode 100644 index 0000000..f950e35 --- /dev/null +++ b/Bounds.m @@ -0,0 +1,27 @@ +classdef Bounds + %BOUNDS Summary of this class goes here + % Detailed explanation goes here + + properties + x LinearBounds + y LinearBounds + end + + methods + function obj = Bounds(xIn,yIn) %BOUNDS Construct an instance of this class + % Detailed explanation goes here + obj.x = xIn; + obj.y = yIn; + end + + function outputArg = area(obj) + %METHOD1 Summary of this method goes here + % Detailed explanation goes here + outputArg = obj.x.length * obj.y.length; + end + function outputArg = getRandomPoint(obj) + outputArg = Point(obj.x.getRandomValue, obj.y.getRandomValue); + end + end +end + diff --git a/LinearBounds.m b/LinearBounds.m new file mode 100644 index 0000000..c6b5225 --- /dev/null +++ b/LinearBounds.m @@ -0,0 +1,25 @@ +classdef LinearBounds + %LINEARBOUNDS Summary of this class goes here + % Detailed explanation goes here + + properties + lower double + higher double + end + + methods + function obj = LinearBounds(lowerIn,higherIn) + obj.lower = lowerIn; + obj.higher = higherIn; + end + + function outputArg = length(obj) + outputArg = abs(obj.higher-obj.lower); + end + + function outputArg = getRandomValue(obj) + outputArg = obj.lower + (rand() * obj.length); + end + end +end + diff --git a/MonteCarloMatlab.m b/MonteCarloMatlab.m new file mode 100644 index 0000000..c964a7c --- /dev/null +++ b/MonteCarloMatlab.m @@ -0,0 +1,47 @@ +tic; +samples = 10000000; +xBounds = LinearBounds(0,20); +yBounds = LinearBounds(-10,150000); +bounds = Bounds(xBounds,yBounds); +pointsInside = 0; +realValue = realIntegral(bounds.x); + +for i = 0:samples + toTestPoint = bounds.getRandomPoint; + if getIsInside(toTestPoint) + pointsInside = pointsInside+1; + end +end + + +format longG +integral = (pointsInside/samples)*bounds.area +realValue +error = abs(realValue-integral) +samples +toc + + + +function out = mainFunction(x) + out = x^4-4*x^3+18*x^2-12*x-69; +end + +function out = indefiniteIntegral(x) + out = (1/5)* x^5 - x^4 + 6* x^3 - 6*x^2 - 69*x; +end + +function out = getIsInside(p) + arguments + p Point + end + y_0 = mainFunction(p.x); + out = abs(y_0) > abs(p.y) && y_0 * p.y >= 0; +end + +function out = realIntegral(bounds) + arguments + bounds LinearBounds + end + out = indefiniteIntegral(bounds.higher) - indefiniteIntegral(bounds.lower); +end \ No newline at end of file diff --git a/Point.m b/Point.m new file mode 100644 index 0000000..367c0ba --- /dev/null +++ b/Point.m @@ -0,0 +1,12 @@ +classdef Point + properties + x double + y double + end + methods + function obj = Point(xIn, yIn) + obj.x = xIn; + obj.y = yIn; + end + end +end \ No newline at end of file