procedure 2 to investigate velocity fields Normal form of Greens theorem
This procedure takes a vector field F1 (a list of two expressions in x and y), a radius R, and a center P and draws the circle C of radius R centered at P. Then it draws the the magenta graph of the normal component of F over the curve and shades vertically the magenta fence trapped between the graph and the circle C. This enables you to visually estimate the net outflow of fluid thru C. Then it also draws the graph of the div of F over a square containg D and draw the navy projection of C up onto that graph together with a navy fence of the solid. This enables you to visually estimate the net divergencel of F over D. Green's theorem says that these two quantities are equal numerically.
> | with(plots): |
> | greenpic2 := proc(R,P,F1) #R=radius of D, P=center of D, F1=[y,-x] velocity field local xrng,yrng,trng,rp,F,Fdn,PQ,r,prpwk,Rd,dv; use plottools in F:=unapply([F1[1],F1[2],0],x,y): xrng:=(P[1]-1.1*R)..(P[1]+1.1*R): yrng:=(P[2]-1.1*R)..(P[2]+1.1*R):trng:=0..2*Pi: r:=[P[1]+R*cos(t),P[2]+R*sin(t),0]: rp:= [diff(r[1],t),diff(r[2],t),0]: Fdn := F(r[1],r[2])[1]*rp[2]-F(r[1],r[2])[2]*rp[1]: PQ := diff(F(x,y)[1],x)+diff(F(x,y)[2],y): prpwk:=evalf(Int(Fdn,t=trng),6): dv:=evalf(Int(Int(subs({x=P[1]+Rd*cos(theta),y=P[2]+Rd*sin(theta)},PQ*Rd),Rd=0..R),theta=0..2*Pi),6); |
> | display( plot3d(PQ,x=xrng,y=yrng,style=wireframe,color=navy), fieldplot3d(F(x,y),x=xrng,y=yrng,z=0..(.01),arrows=SLIM), seq(line([P[1]+R*cos(i*op(trng)[2]/36),P[2]+R*sin(i*op(trng)[2]/36),0],[P[1]+R*cos(i*op(trng)[2]/36),P[2]+R*sin(i*op(trng)[2]/36),subs(t=i*op(trng)[2]/36,Fdn)],color=magenta),i=0..36), rotate(display(seq(line([P[1]+R*cos(i*op(trng)[2]/36),P[2]+R*sin(i*op(trng)[2]/36),0],[P[1]+R*cos(i*op(trng)[2]/36),P[2]+R*sin(i*op(trng)[2]/36),subs(t=i*op(trng)[2]/36,subs({x=r[1],y=r[2]},PQ))],color=maroon),i=0..36)),Pi/36,[[P[1],P[2],0],[P[1],P[2],1]]), spacecurve([r[1],r[2],Fdn],t=trng,color=magenta), spacecurve([r[1],r[2],subs({x=r[1],y=r[2]},PQ)],t=trng,thickness=1,color=navy), spacecurve(r,t=trng),color=black,axes=boxed,labels=[x,y,z],title=cat("net outflow (magenta area) = ",convert(prpwk,string)," and net divergence(navy volume) = ",convert(dv,string))); end use; |
> | end: |
> |
> |