Calculating exact volume using maple

369 Views Asked by At

For an assignment I have to use numerical integration technique to calculate volume with cylindrical surface

$\Omega = \{(x−0.5)^2 +(y−0.5)^2 \leq 0.25\}$ and height is $\ln(x+y)$

I have used Monte Carlo technique to calculate the volume. But to be sure the answer is correct I want to check the exact volume using Maple. I have been searching online on how to do it but couldn't find it.

So the question is, is there a way to calculate exact volume using Maple for that object or integral like $$ \iint_\Omega \lvert \ln(x+y)\rvert \, \mathrm dx \, \mathrm dy $$ with $$ \Omega = \{(x-0.5) ^2 + (y-0.5)^2 \leq 0.25\} $$

3

There are 3 best solutions below

0
On BEST ANSWER

Note that this answer does not contain an exact solution

For a general solution, you can use indicator functions, implemented in Maple with $1_A$ = piecewise(A,1) to integrate over your set $\Omega$.

f:=abs(log(x+y));                            ## function
omega := (x-1/2)^2 + (y-1/2)^2 <= 1/4;       ## set
one__omega := piecewise(omega,1);            ## indicator fn of set

a:=-1;b:=1;                                  ## integral bounds
int(int(f*one__omega,x = a..b),y = a..b);  
value(%);                                    ## exact value
evalf(%);                                    ## approximation

The integral bounds form a cartesian square $(x,y)\in [a,b]^2$. Note that the bounds must cover $\Omega$ entirely. As $\omega$ in your case is a disk with center $(1/2,1/2)$ and radius $1/4$ then $a = 0$ and $b = 1$ will suffice. Since Maple was not able to find an exact solution I used evalf to get a numerical approximation. Maple approximates it to be $\approx 0.2550197391$ in agreement with your MC method.

1
On

I was able to corroborate your result by two independent methods, using Matlab. In the first instance, I developed my own Monte Caro simulation. Using a million points, the fraction $\pi/4$ land in the circular cylinder, I find that $V\approx 0.25513$. (The calculation takes less than 1 second.)

The second method was to create matrix of $(x,y)$ points on the plane and calculate the height of the surface $z=|\ln(x+y)|$. Then I sum the incremental volume elements, $dV=z~dx~dy$. In this way I found the volume to be $V\approx 0.25501$ using an $(x,y)$ array of $2000\times2000.$

0
On

Using Maple, you can use a formula (in y say) for the range of the inner integral in x. Thus you can integrate from one side of the circular domain to the other (for any given y). That makes it easier for a numeric quadrature method than using a piecewise function which is not smooth at the circular boundary. It is easy enough to isolate the equation for the circular edge, to obtain a formula (in y) for the x-value at either side of the circular domain.

The same kind of thing can be done for making a 3D plot.

BTW, why do you take the absolute value of the integrand given in the first line of your question? It makes a difference to the answer (See below, where I do it both with and without the abs.)

restart;
f := abs(log(x+y)):

sols:=[ solve((x-1/2)^2 + (y-1/2)^2 = 1/4, x ) ];

                  [             (1/2)               (1/2)]
                  [1   /  2    \       1   /  2    \     ]
          sols := [- + \-y  + y/     , - - \-y  + y/     ]
                  [2                   2                 ]

plot3d(f, x=sols[2]..sols[1], y=0..1,
       filled, style=surface, grid=[100,100],
       lightmodel=Light1, orientation=[-50,70,0] );

enter image description here

evalf( Int( Int( f, x=sols[2]..sols[1] ), y=0..1 ) );

                            0.2550197391

restart;
f := log(x+y):

sols:=[ solve((x-1/2)^2 + (y-1/2)^2 = 1/4, x ) ]:

plot3d(f, x=sols[2]..sols[1], y=0..1,
       filled, style=surface, grid=[100,100],
       lightmodel=Light1, orientation=[-50,70,0] );

enter image description here

evalf( Int( Int( f, x=sols[2]..sols[1] ), y=0..1 ) );

                           -0.05698907679