Re: Surface integral
- To: mathgroup at smc.vnet.net
- Subject: [mg9854] Re: [mg9774] Surface integral
- From: Allan Hayes <hay at haystack.demon.co.uk>
- Date: Fri, 28 Nov 1997 05:36:10 -0500
- Sender: owner-wri-mathgroup at wolfram.com
Yves Verbandt wrote:
>
> Deat Mathematica experts,
>
> Can anyone help me with the calculation of a surface in 3D which is
> given by a list of {x,y,z} values? Perhaps there is an adaptation of
> ListIntegrate somewhere out there?
>
Yves:
The area of a triangle with verticies a,b,c in any R^n is
area[{a_,b_,c_}]:=
Sqrt[Abs[Det[#.Transpose[#]]]]/2&[{b-a,c-a}]
Take a matrix of points
data=Table[{x,y,Sin[x-Sin[y]]}//N, {x,1,5},{y,1,6}];
Plot it
ListPlot3D[Map[Last,data,{2}], AxesLabel-> {"x","y","z"}];
Cosider a point, p, in data (not along the bottom or the right of the
matrix), with the neighboring points, e (east of p), s (south p} and
se (south east of p).
p--- e
| |
| |
s--- se
The sum of the areas of the triangles with vertex lists {p,s,e} and {se
e,s} is given by the function
area[{p_,s_,e_,se_}]:=
area[{p,s,e}] + area[{se,s,e}]
We shall replace each such p in the matrix with the quadruple {p,se,s,e}
and then sum the values of area[{p,se,s,e}]
p = Drop[#,-1]&/@Drop[data,-1];
s=Drop[#,-1]&/@Drop[RotateLeft[data],-1];
e=Drop[#,-1]&/@Drop[RotateLeft[data,{0,1}],-1];
se=Drop[#,-1]&/@Drop[RotateLeft[data,{1,1}],-1];
Make the matrix:
data4 = MapThread[List,{p,s,e,se},2];
Make this into a list of quadruples {p,s,e,se}
data4 = Flatten[data4,1]
Now sum these areas
Plus@@(area/@data4)
25.1413
Thes steps could be combined.
We can compare this with the area of the related smooth surface:
Clear[x,y]
NIntegrate[
Evaluate[Sqrt[1+ D[Sin[x-Sin[y]],x]^2 + D[Sin[x-Sin[y]],y]^2]],
{x,1,5},{y,1,6}]
25.6051
Another approximation is given by using the area of the parallelogram
{p,s,e,s+e-p} instead of area[{p,s,e}] + area[{se,s,e}] This is the
approximation underlying theintegral above.
p = Drop[#,-1]&/@Drop[data,-1];
s=Drop[#,-1]&/@Drop[RotateLeft[data],-1];
e=Drop[#,-1]&/@Drop[RotateLeft[data,{0,1}],-1];
data4p=
Flatten[MapThread[List,{p,s,e},2],1];
Plus@@(2 area/@data4p)
25.2283
--
Allan Hayes
Mathematica Training and Consulting
Leicester, UK
hay at haystack.demon.co.uk
http://www.haystack.demon.co.uk
voice: +44 (0)116 271 4198
fax: +44 (0)116 271 4198