Path: news.mathworks.com!not-for-mail
From: "zedong 
" <zdongwu@gmail.com>
Newsgroups: comp.soft-sys.matlab
Subject: Re: How to integration a polynomial on a polygon?
Date: Mon, 22 Dec 2008 01:20:05 +0000 (UTC)
Organization: The MathWorks, Inc.
Lines: 54
Message-ID: <gimq05$2gv$1@fred.mathworks.com>
References: <gikgkj$4qp$1@fred.mathworks.com> <gil3rj$ss6$1@fred.mathworks.com> <gimie6$gr7$1@fred.mathworks.com>
Reply-To: "zedong 
" <zdongwu@gmail.com>
NNTP-Posting-Host: webapp-05-blr.mathworks.com
Content-Type: text/plain; charset="ISO-8859-1"
Content-Transfer-Encoding: 8bit
X-Trace: fred.mathworks.com 1229908805 2591 172.30.248.35 (22 Dec 2008 01:20:05 GMT)
X-Complaints-To: news@mathworks.com
NNTP-Posting-Date: Mon, 22 Dec 2008 01:20:05 +0000 (UTC)
X-Newsreader: MATLAB Central Newsreader 1639198
Xref: news.mathworks.com comp.soft-sys.matlab:508266

"Roger Stafford" <ellieandrogerxyzzy@mindspring.com.invalid> wrote in message <gimie6$gr7$1@fred.mathworks.com>...
> "Roger Stafford" <ellieandrogerxyzzy@mindspring.com.invalid> wrote in message <gil3rj$ss6$1@fred.mathworks.com>...
> >   For linear functions such as in your f(x,y), the integral over a polygonal domain can be evaluated by a comparatively simple summation.  I'll just give the method for f(x,y) = x.
> > 
> >   Let x and y be vectors of the polygon vertices' coordinates.  Then do:
> > 
> >  x1 = x([2:end,1]);
> >  y1 = y([2:end,1]);
> >  I = sum((x.^2+x.*x1+x1.^2).*(y1-y))/6;
> > 
> > I is the desired integral of f(x,y) = x over the polygonal area.
> > 
> >   There are similar formulas for f(x,y) = y and f(x,y) = 1 (area), and you could combine them in a linear combination to evaluate the function above.
> > 
> >   Note that it is assumed that the polygon vertices in the vectors are sequenced in counterclockwise order around the polygon.  Also it is assumed that the polygon boundary never crosses itself, as for example in a figure eight shape.  The polygon does not have to be convex.
> > 
> >   Such formulas become considerably more complicated for functions like x^2, x*y, x^2*y, etc. or linear combinations of these, but they are possible.
> > 
> > Roger Stafford
> 
>   The following allows you to generalize that formula I gave you, Zedong.  Unfortunately you will have to convert its format to that of up-to-date symbolic toolboxes.  Mine is an older version.  Suppose you want to integrate the function f(x,y) = x*y.  Then do this:
> 
>  a = 1; b = 1;
>  f = subs(subs('x^a*y^b',a,'a'),b,'b');
>  I1 = int(int(f,'y',0,'y2/x2*x'),'x',0,'x2');
>  I2 = int(int(f,'y',0,'y1/x1*x'),'x',0,'x1');
>  I3 = int(int(f,'y',0,'(y1*x2-y2*x1+(y2-y1)*x)/(x2-x1)'),'x','x2','x1');
>  I = factor(simplify(symop(I1,'-',I2,'+',I3)));
>  pretty(I)
>  
>        1/24 (- y1 x2 + y2 x1) (2 x1 y1 + y2 x1 + y1 x2 + 2 x2 y2)
> 
> Corresponding expression for I:
> 
>  I = '(y2*x1-y1*x2)*(y1*(2*x1+x2)+y2*(2*x2+x1))/24';
> 
>   Now substitute x for x1, x1 for x2, y for y1, and y1 for y2 in this expression and use the method I described earlier with polygon vectors x and y.
> 
>  x1 = x([2:end,1]);
>  y1 = y([2:end,1]);
>  K = sum((y1.*x-y.*x1).*(y.*(2.*x+x1)+y1.*(2.*x1+x)))/24;
> 
> K is your integral. 
> 
>   Note 1: "symop(I1,'-',I2,'+',I3)" is just my symbolic toolbox's crude way of combining the expressions: I1-I2+I3.
> 
>   Note 2: The above expression for I is the integral of f(x,y) = x*y over the area of the triangle (0,0), (x1,y1), (x2,y2).  In the summation all areas of integration lying outside the polygon eventually cancel out.
> 
>   Note 3: At about degree 5, these formulas begin to be uncomfortably large and cumbersome so you might prefer numerical methods for polynomials of that high a degree.
> 
> Roger Stafford

I am sorry.I think your method must be very beautiful.But I can't follow your algorithm.Could you tell me why this kind of combination summation is right?Thank you very much.If it's hard to explain, Could you tell me where can I find the theory?Thank you.