square root of pi

I was banging my head against the wall looking for the error in my code, that kept returning "false" when it should have returned true. I finally realised that my error was because in MATLAB (or at least in my copy) pi^0.5^2==pi returns "false". The same goes for sqrt(pi)^2==pi. I assume this is because the machine precision of pi is more exact then the machine precison of the square root function, even though they are identical up to 16 digits.

I have some unelegant solutions to this. e.g writing down pi=3.141592 at the begining of my code. is there a better solution? a more percise calculaion of the square root?

many thanks

[dpb]

Where would one have this in code, specifically? Show the use case.

NB: reliance on exact equality is always subject to machine precision/rounding for floating point and code should not be written to rely on such.

One possible solution with MATLAB is ismembertol instead of "=="

[Nathan Bblanc]

thanks for your answer dpb. i am attaching my code. the relevant one is TA_Round_Duct it is a class with several dependent properties, representing a round duct or tube. never mind most of what is written there, the important thing is that the area and diameter are dependent via A=pi*d^2/4 or D=(A/pi*4)^0.5. to prevent recursive calling between the set functions of both properties i set "abortset=true" on both properties. however, in some cases this does not work and recursion calling occurs due to the issue i described earlier. I think that this happens only when D is a round number but i am not sure.

I think your suggestion is not relevant to this case because the == is tested for inside the "preset" event and not in an accessible part of my code.

if you want to test the code, you can instantiate the objects with "null" as the "system" input parameter and any string as the "name"

[Adam Danz]

pi - (pi^0.5^2)
ans =
      4.44089209850063e-16

NASA only uses around 15 digits of pi in its calculations for sending rockets into space (source).

[dpb]

Your initialization code has warnings of one property making reference to another and explicit instructions to not leave code that has such an issue.

I don't do objects so not positive of the manner in which this should be resolved, but redesign appears in your future... :)

[Nathan Bblanc]

dpb i actualy couldn't see a way around this that allows setting both properties. in any case this is a seperate issue from my question

Answer 1

[Steven Lord]

I would have both the Area and Diameter properties depend upon a third (maybe public, maybe internal) property. Setting either Area or Diameter computes what the value of that internal property should be and getting them computes the value using that internal property.

classdef Circle715788
    properties(Dependent)
        Area
        Diameter
    end
    
    properties
        Radius
    end
    
    methods % Constructor
        function obj = Circle715788(r)
            obj.Radius = r;
        end
    end
    
    methods % Property accessors
        function obj = set.Area(obj, newArea)
            R = sqrt(newArea./pi);
            obj.Radius = R;
        end
        function obj = set.Diameter(obj, newDiameter)
            obj.Radius = newDiameter./2;
        end
        function A = get.Area(obj)
            A = pi*obj.Radius.^2;
        end
        function D = get.Diameter(obj)
            D = 2*obj.Radius;
        end
    end
end 

Use as:

>> Q = Circle715788(1)
Q = 
  Circle715788 with properties:

        Area: 3.1416
    Diameter: 2
      Radius: 1

>> Q = Circle715788(4)
Q = 
  Circle715788 with properties:

        Area: 50.2655
    Diameter: 8
      Radius: 4

>> Q.Area = pi
Q = 
  Circle715788 with properties:

        Area: 3.1416
    Diameter: 2
      Radius: 1       

[Nathan Bblanc]

Thank you Steven and sorry for the late reply. I thought of doing this but I am worried that this might slow down my code. I am calling these properties about 10^6 times during a run. What do you think?

[Steven Lord]

On my Windows 10 machine running release R2020b:

>> C = Circle715788(7);
>> tic; for k = 1:1e6, A = C.Area; end, toc
Elapsed time is 0.132099 seconds.
>> tic; for k = 1:1e6, A = C.Area; end, toc
Elapsed time is 0.112353 seconds.
>> tic; for k = 1:1e6, A = C.Area; end, toc
Elapsed time is 0.108369 seconds.
>> tic; for k = 1:1e6, A = C.Area; end, toc
Elapsed time is 0.105188 seconds.
>> tic; for k = 1:1e6, A = C.Area; end, toc
Elapsed time is 0.102445 seconds.

So about a tenth of a second for 1e6 queries.

>> tic; for k = 1:1e6, C.Radius = k; A = C.Area; end, toc
Elapsed time is 0.231171 seconds.
>> tic; for k = 1:1e6, C.Radius = k; A = C.Area; end, toc
Elapsed time is 0.210652 seconds.
>> tic; for k = 1:1e6, C.Radius = k; A = C.Area; end, toc
Elapsed time is 0.208293 seconds.

Two-tenths for 1e6 property sets and gets.

You might be able to speed this up, particularly if the calculations of Area / Diameter are expensive, using a memoize object inside your Area / Diameter calculations to cache the values for repeated Radius values.