Normal vector on a surface
Matthi...@oppenheim.de
I have a function in the variables x1 and y1:
Out[27]=
19.74211746962547 - 61.78321746073334*
x1 + 70.84823523445556*x1^2 -
34.64681309362152*x1^3 +
5.822595947190386*x1^4 -
61.783217460733795*y1 +
188.56171712734522*x1*y1 -
208.7457484391798*x1^2*y1 +
99.21114279328117*x1^3*y1 -
16.223098505477388*x1^4*y1 +
70.8482352344551*y1^2 -
208.7457484391805*x1*y1^2 +
225.08774661852397*x1^2*y1^2 -
103.5151716236312*x1^3*y1^2 +
16.351931921608763*x1^4*y1^2 -
34.64681309362163*y1^3 +
99.21114279328117*x1*y1^3 -
103.51517162363109*x1^2*y1^3 +
45.654124756950296*x1^3*y1^3 -
6.928857192755963*x1^4*y1^3 +
5.822595947190411*y1^4 -
16.22309850547743*x1*y1^4 +
16.351931921608763*x1^2*y1^4 -
6.928857192755952*x1^3*y1^4 +
1.0137658500940734*x1^4*y1^4
This function yields a surface very similar to Sin[x1*y1] for 1<x1<3 and
1<y1<3.
Now I want to calculate (how?) and draw (how?) several "Normalenvektors"
(sorry, I do not know the English termini technici) which should sit smugly
- like palisades - on the plane tangential to the surface.
The "Normalenvektor" N in point P - according to Bronstein-Semendjajew - is
a unity vector perpendicular to the tangential plane; its accompanying
vectors e1 and e2 on the plane form a "right-handed system". N, e1 and e2
are referred to as the "accompanying tripod". - I understand the words but
not their meaning.
My attempts with Calculus`VectorAnalysis` and PlotVectorField3D &c. failed
dismally.
Thank you for your assistance,
Matthias Bode
Sal. Oppenheim jr. & Cie. KGaA
Koenigsberger Strasse 29
D-60487 Frankfurt am Main
GERMANY
Tel.: +49(0)69 71 34 53 80
Mobile: +49(0)172 6 74 95 77
Fax: +49(0)69 71 34 6380
E-mail: matthi...@oppenheim.de
Internet: http://www.oppenheim.de
David Park
Here is one method for obtaining a unit vector field normal to your surface.
Parameterize your surface as what is called a patch in differential
geometry.
X[u_,v_]:= {u,v,(your expression with u and v for x1 and y1)}
Then you can obtain two independent tangent vectors to the surface by taking
the derivative of X, first with respect to u and then with respect to v.
eu[u_,v_] = D[X[u,v],u]
ev[u_,v_] = D[X[u,v],v]
These two vectors are not orthogonal, but for the moment that doesn't
matter. You can then calculate the normal vector field by using a cross
product and making a unit vector of the result.
Cross[eu[u_,v_],ev[u_,v_]
eN[u_,v_]=%/Sqrt[%.%]
Now, for the plotting. This is a little difficult with regular Mathematica
because if has no Arrow3D routine. But I have the routines for making a nice
plot and these are available at my web site (except for one which I will
send to you if you wish). The packages are DrawingCube, ParametricDrawing3D
and DrawingArrows. DrawingArrows has an Arrow3D routine in it, which is
similar to the regular Arrow routine. DrawingCube makes it easy to combine
various elements into one plot, such as your surface and the arrows.
If you download the packages, they should be placed in the
AddOns/ExtraPackages/Graphics directory. You may have to create the Graphics
directory there. Then they behave just like the standard packages.
I made up a notebook which draws the surface with the normal vector field on
it using these packages. I will send it to you if you are interested in
trying the packages. I can also show you how to draw a tangent plane at one
point along with the vector triad that goes with that point.
One other package that I used was one that was posted by Allan Hayes &
Hartmut Wolf a month or so ago on MathGroup. It is called Smooth3D. It is
useful in this way. I wanted to put a normal vector at each intersection of
the "mesh" lines. However we don't want too many vectors because it will
take too long and the plot will be cluttered. On the other hand, we want to
use enough points in plotting the surface to make it smooth. Smooth3D allows
you to suppress the regular "mesh" lines and then plot new, more widely
spaced lines. This makes a much better looking plot. I can send you the
Smooth3D package along with the notebook. I slightly modified it so that it
goes in the ExtraPackages/Graphics directory, just like the other graphics
packages.
So, if you are interested in pursuing this approach I will send you the
material.
David Park
dj...@earthlink.net
http://home.earthlink.net/~djmp/
jmt
Needs["Calculus`VectorAnalysis`"]
n[u_, v_] := Module[{dx, dy, p},
dx[x_, y_] = {1, 0, D[f[x, y], {x, 1}]};
dy[x_, y_] = {0, 1, D[f[x, y], {y, 1}]};
p[x_, y_] = CrossProduct[dx[x, y], dy[x, y]];
With[{m = p[u, v]}, m/Sqrt[m.m]]
]
Albert Retey
> Now I want to calculate (how?) and draw (how?) several "Normalenvektors"
> (sorry, I do not know the English termini technici) which should sit smugly
> - like palisades - on the plane tangential to the surface.
I think it is as simple as "normal vector" :-)
> The "Normalenvektor" N in point P - according to Bronstein-Semendjajew - is
> a unity vector perpendicular to the tangential plane; its accompanying
> vectors e1 and e2 on the plane form a "right-handed system". N, e1 and e2
> are referred to as the "accompanying tripod". - I understand the words but
> not their meaning.
One simple way to determine the normal vector of a surface is to find
two (linear independent) tangential vectors to the surface. Then find a
vector that is perpendicular to those two and has length one. If you
have expressed the surface as z=z[x,y], as in your case, the most simple
choice for the two tangential vectors is the ones in "x-" and
"y-direction", that is e1 and e2 (If you choose other coordinates, e1
and e2 might change, but not the normal vector). The fact that these
three vectors should form a right-handed-system is then just a
convention, because there are actually two choices of sign: for your
surface the normal vectors could either point "up" or "down".
Here is a simple way to calculate and visualize the normal vectors of
your surface (basically using the definitions in Bronstein):
In[200]:=
surf[x1_,y1_]:=
19.74211746962547-61.78321746073334*x1+70.84823523445556*
x1^2-34.64681309362152*x1^3+5.822595947190386*x1^4-61.783217460733795*
y1+188.56171712734522*x1*y1-208.7457484391798*x1^2*
y1+99.21114279328117*x1^3*y1-16.223098505477388*x1^4*
y1+70.8482352344551*y1^2-208.7457484391805*x1*y1^2+225.08774661852397*
x1^2*y1^2-103.5151716236312*x1^3*y1^2+16.351931921608763*x1^4*
y1^2-34.64681309362163*y1^3+99.21114279328117*x1*
y1^3-103.51517162363109*x1^2*y1^3+45.654124756950296*x1^3*
y1^3-6.928857192755963*x1^4*y1^3+5.822595947190411*
y1^4-16.22309850547743*x1*y1^4+16.351931921608763*x1^2*
y1^4-6.928857192755952*x1^3*y1^4+1.0137658500940734*x1^4*y1^4;
In[202]:=
surfplot=
Plot3D[surf[x,y],{x,1,3},{y,1,3},PlotRange->All,PlotPoints->100,
Mesh->False];
(* This are the tangential vectors r_u and r_v in Bronstein, for the
special case u=x and v=y *)
In[203]:=
rx[x_,y_]:={1,0,D[surf[x,y],x]}
ry[x_,y_]:={0,1,D[surf[x,y],y]}
(* the cross product (will be perpendicular to rx and ry) *)
In[205]:=
rxcrossry[x_,y_]:=Cross[rx[x,y],ry[x,y]]
(* divide by length for normalization *)
In[206]:=
normalvector[x_,y_]:=
Evaluate[rxcrossry[x,y]/Sqrt[rxcrossry[x,y].rxcrossry[x,y]]]
(* make a table of normal vectors in the form appropriate for
ListPlotVectorField *)
In[207]:=
normvectable=
Flatten[Table[{{xx,yy,surface[xx,yy]},normalvector[xx,yy]},{xx,1,3,
0.3},{yy,1,3,0.3}],1];
In[208]:=
<<Graphics`PlotField3D`
In[209]:=
?ListPlotVectorField3D
"ListPlotVectorField3D[{{pt,vec},...},(options)] plots a list \
of vectors in three dimensions, each vector based at a corresponding
point \
pt."
(* I tried to make the plot look somewhat better by scaling the length
and plot the vectors in red = Hue[1] *)
In[210]:=
normvecplot=ListPlotVectorField3D[normvectable,VectorHeads->True,
ScaleFactor->0.3,ColorFunction->(Hue[1]&)];
In[211]:=
Show[surfplot,normvecplot,AspectRatio->1];
You will see that the vectors plotted do not exactly meet what you
probably expect concerning perpendicularity, I guess it should be
possible to find a ViewPoint that makes the plot look somewhat better...
Hope that helps,
Albert
--
Visual Analysis AG Internet: www.visualanalysis.com
Neumarkter Str. 87 Telefon: 089 / 431 981 36
D-81673 Muenchen Telefax: 089 / 431 981 1
Albert Retey
There was a misprint in my first post, the line
In[207]:=
normvectable=
Flatten[Table[{{xx,yy,surface[xx,yy]},normalvector[xx,yy]},{xx,1,3,
0.3},{yy,1,3,0.3}],1];
should be:
In[207]:=
normvectable=
Flatten[Table[{{xx,yy,surf[xx,yy]},normalvector[xx,yy]},{xx,1,3,
0.3},{yy,1,3,0.3}],1];
Daniel Lichtblau
> (sorry, I do not know the English termini technici) which should sit smugly
> - like palisades - on the plane tangential to the surface.
>
> a unity vector perpendicular to the tangential plane; its accompanying
> vectors e1 and e2 on the plane form a "right-handed system". N, e1 and e2
> are referred to as the "accompanying tripod". - I understand the words but
> not their meaning.
>
> dismally.
>
> Thank you for your assistance,
>
> Matthias Bode
> Sal. Oppenheim jr. & Cie. KGaA
> Koenigsberger Strasse 29
> D-60487 Frankfurt am Main
> GERMANY
> Tel.: +49(0)69 71 34 53 80
> Mobile: +49(0)172 6 74 95 77
> Fax: +49(0)69 71 34 6380
> E-mail: matthi...@oppenheim.de
> Internet: http://www.oppenheim.de
Your surface is given in the form z = f[x,y] so you can take
-{D[f[x,y],x], Df[x,y],y], -1}
as an upward-pointing normal vector, and then normalize to get unit
length. The code below will do this.
normalvector[x_,y_]:= With[
{vec = -{D[f[x1,y1],x1],D[f[x1,y1],y1],-1} /. {x1->x,y1->y}},
vec / Sqrt[vec.vec]]
Daniel Lichtblau
Wolfram Research