Need help using Matlab for an approximation of Pi via the Mandelbrot set.

Hello, for my final project in my computational physics class, I wanted to compare the Mandelbrot approximation of pi to Matlab's built in pi. I got the idea from a NumberPhile video: https://www.youtube.com/watch?v=d0vY0CKYhPY and here is my code so far:

%%%%%%%%%%%%%

close all;clear all; clc %Prepared by Lauren Stearns

% Project: Using Matlab to evaluate pie via the Mandlebrot set and comparing the approximation to % Matlab's built in pi.

prompt = 'Input a small, real, positive value for epsilon. Input:';

c = 1/4 + input(prompt) ; % c needs to be greater than 1/4 such as 1/4 + epsilon. % Where epsilon is a very small, real, positive, number. 1/4 represents the

% cusp of the mandlebrot set.

z = 0; n = 0;

while z < 2 % The number of iterations required to surpass 2 is the approximation of pi.

z = z.^2 + c; n = n+1 ;

end

n format long

Pi = pi

%%%%%%%%%%%%%

Table of Values obtained from running this code:

1) I wanted to know if there was a way to insert a period in front of the first digit in my n output. For example, I would like the n that corresponds to epsilon 1e-6 to be displayed as 3.1410 instead of 31410

2) My next question, is why won't anything to a higher power than epsilon = 1e-16 generate an output? 1e-16 gives me an out put (as seen in the table), but I let 1e-17 run for over 20 minutes and got nothing.

3) I also hoped to do a comparison at the end of my code between the output of n and Matlab's built in value for pi. I tried adding the following lines at the end of my code, but I couldn't figure out how to account for the changing magnitude of the iterations.

%%%%%%%%

Comparison = abs(((pi-n)/pi*10^(-4)))*100;

display(Comparison)

%%%%%%%

4) Lastly, I noticed that only the even powers are converging to pi and the odd powers appear to be converging to 1. So why do only even powers converge to pi?

Any and all help would be greatly appreciated. Thank you for taking the time to read my post :)

Answer 1

[David Goodmanson]

Hello Lauren,

For all but your second question, the reasons are because if c = 1/4 + eps and there are n steps to get past 2, then what the video didn't say is

pi = n*sqrt(eps), approximately.

With that expression you will get successive approximations to pi without the decimal point problem.

The situation in your second question is predictable due to machine precision. If you compare (1/4 + 1e-17) to 1/4 by subtraction you will find no difference, because Matlab has run out of significant figures at that point. For 1e-16 you do get a difference (although it is not exactly 1e-16). For the 1e-17 case, c = 1/4 and the process gets never gets past a value that you can go check on.

[Lauren Stearns]

Thank you so very much!!

[David Goodmanson]

You're welcome. Here is a reference on how it works: http://www.pi314.net/eng/mandelbrot.php

[Lauren Stearns]

Your answer to my post and the link you provided answered all of my previous questions! However I have a new question regarding the code. Here is my updated script:

%%%%%%

% Third Version

close all;clear all; clc

%Prepared by Lauren Stearns

prompt = 'Input a small, real, positive value for epsilon. Input:';

c = 1/4 + input(prompt) ;

n = 0;

while z < 2

z = z.^2 + c;

n = n+1 ;

end

Approximation = n*sqrt(input(prompt));

display(Approximation)

format long

Pi = pi

%%%%%%%%%%%

The code works well with exception to being required to type epsilon in twice. I'm assuming it is because I have used input(prompt) a second time in the Approximation variable, but I don't know how to get around it. Im sure this is a simple fix, but this is my first term working with Matlab and i'm still learning. =^-^=

Best Regards Lauren

[Stephen23]

@Lauren Stearns: store the value in a variable and then you do not need to call input twice:

epsilon = str2double(input(prompt,'s')); 
c = 1/4 + epsilon;
...
Approximation = n*sqrt(epsilon);

[Lauren Stearns]

Thank you for your help. :)