################################ # # TABLE OF CONTENTS # (You can move to any line by using the "lightning bolt" button to the right of the large 'A') # # 65 Distributive property for dot- and cross-product # 112 Miscellaneous examples, done in class on R 8/28/18 # 125 Plot a plane, given three points on it # 170 Plot a plane, given a normal vector and the distance from the plane to the origin # 220 Miscelleaneous computations involving vectors # 240 Interactive 3D grapher: VERY USEFUL!! # 340 Solution to the extra credit for quiz 1 # 391 Plotting parametric curves with normal, tangent, binormal vectors, # osculating plane, acceleration, curvature, etc. # 542 Interactive TANGENT PLANE PLOTTER. Includes partial derivatives. Very useful! # 645 Surface plot compared with two different styles of contour plot. Includes a "fancy" surface plot that is prettier than # the standard one, but you cannot tilt it with the mouse. But is is designed to exactly match the # contour plot. # 680 Chain-rule practice. Great homework helper! # 715 Midpoint approximation and exact computation of simple iterated integrals. Another good homework helper! # 755 Problems 15.6.33-34 # 885 Three views of the element of volume for spherical coordinates # 1022 Images from quiz #5 (examples of cylindrical and spherical plotting) # 1075 Interactive (u,v) --> (x,y) transformation explorer (Section 15.9). Very instructive! Shows Jacobians, too! # 1120 Examples of plotting vector fields # 1158 Examples of parametric equations of surfaces # 1175 INTERACTIVE Parametric Surface EXPLORER #ILLUSTRATION OF DISTRIBUTIVE PROPERTY FOR DOT PRODUCT AND FOR CROSS PRODUCT #In the diagram created, vector A is a unit vector. We draw both parallel and perpendicular projections of vectors B and C onto A. # The main idea is that given any vectors B and A, we can resolve B into two orthogonal components: the PARALLEL component of B onto A is parallel to A, while the # PERPENDICULAR component of B onto A is perpendicular to A and in the same plane as A and B. # The MAGNITUDE of the PARALLEL projection of B onto A is the absolute value of the dot-product of B and A (assume A is a unit vector), # and the PERPENDICULAR projection of B onto A (the actual VECTOR, not the magnitude) is the cross-product of B and A. def perpProj(x,y): #output is the perpendicular component of x wrt y, where #it is understood for convenience that y is a unit vector #and x, y are both vectors return x-y*x.dot_product(y) def parProj(x,y): #output is the parallel component of x wrt y, where y is unit return y*x.dot_product(y) a=vector([0,1,0]) b=vector([-.3,.1,.7]) c=vector([-.7,.5,.3]) t=6 #thickness of some vectors bc=b.cross_product(a) P=plot(a, color='red', thickness=t) P += plot(b, color='blue', thickness=t) P += plot(c, color='black',thickness =t, start=b) #P += plot(c,color='black',start=b) #P += plot(b,color='blue',start=c) P += plot(text3d("B",b*0.7-a*.05,color='blue')) P += plot(text3d("A",a*1.05,color='red')) P += plot(text3d("C",b*1.2+c*0.3,color='black')) P += plot(b+c,color="magenta", thickness=t) P += plot(text3d("B+C", (b+c)*.6, color='magenta')) P += plot(-perpProj(b,a), color='blue',start=b) P += plot(-perpProj(b,a), color='blue',start=b+parProj(c,a)) P += plot(parProj(b,a), color='blue',start=-0.05*perpProj(b,a)) P += plot(parProj(c,a), color='black',start = parProj(b,a)-0.08*perpProj(b,a)) P += plot(parProj(b+c,a), color='magenta',start=-0.12*perpProj(b,a)) P += plot(-perpProj(b+c,a), color='magenta',start=b+c) P += plot(parProj(c,a),color='black',start=b) P += plot(perpProj(c,a), color='black',start=b+parProj(c,a)) #P += plot(parProj(b+c,a), color='magenta',start=perpProj(b+c,a)) P.show(aspect_ratio=1.0,frame=false, spin=5)
3D rendering not yet implemented
244860842756831203896798793038363593764863637715027581243662187665028113084218676375430402110652416
plot(sin,0,pi)
%var x,y,z pic=implicit_plot3d(x^2+y^2+z^2-1,(x,-1,0.5),(y,-1,1),(z,-1,1),color='red',mesh=1) pic.show()
3D rendering not yet implemented
#Plot a plane, given three points on the plane # try changing the points defined in line 27 of this #declare variables %var x, y, z #define a function to draw axes. Very useful for understanding plots def axes(size,xcolor='red',ycolor='green',zcolor='blue'): #returns a plot of axis vectors, labeled, with size xSize=size; ySize= size; zSize=size xAxis = plot(vector((xSize,0,0)),color=xcolor,width=3) yAxis = plot(vector((0,ySize,0)),color=ycolor,width=3) zAxis = plot(vector((0,0,zSize)),color=zcolor,width=3) #Labels for axes #set fontsize (change if you need to) fsize=20 xT=text3d("X",(xSize*1.1,0,0),color=xcolor,fontsize=fsize) yT=text3d("Y",(0,ySize*1.1,0),color=ycolor,fontsize=fsize) zT=text3d("Z",(0,0,zSize*1.1),color=zcolor,fontsize=fsize) return xAxis+xT+yAxis+yT+zAxis+zT asize = 16 #axis size psize = 12 #plane size #define the three points a= (0,0,1); b= (0,1,0); c= (1,0,0) #turn em into vectors A = vector(a); B= vector(b); C=vector(c) AB = B-A; AC = C-A N = AB.cross_product(AC) #note the syntax for cross product f = N[0]*x+N[1]*y + N[2]*z-N.dot_product(A) #set this equal to zero to get equation of plane P = plot(axes(asize)) P += implicit_plot3d(f, (x,-psize, psize), (y,-psize, psize), (z, -psize, psize), color='orange', opacity=0.1) P += plot(points([a,b,c], thickness=40, color='black')) P += plot(N, start = a, width = 3, color='black') show('The equation of the plane is $'+latex(f)+'=0$') show('Normal vector is black.') show(P, spin =5)
The equation of the plane is −x−y−z+1=0
Normal vector is black.
3D rendering not yet implemented
#Plot a plane, given a normal vector and the distance from the origin #declare variables %var x, y, z #define a function to draw axes. Very useful for understanding plots def axes(size,xcolor='red',ycolor='green',zcolor='blue'): #returns a plot of axis vectors, labeled, with size xSize=size; ySize= size; zSize=size xAxis = plot(vector((xSize,0,0)),color=xcolor,width=3) yAxis = plot(vector((0,ySize,0)),color=ycolor,width=3) zAxis = plot(vector((0,0,zSize)),color=zcolor,width=3) #Labels for axes #set fontsize (change if you need to) fsize=20 xT=text3d("X",(xSize*1.1,0,0),color=xcolor,fontsize=fsize) yT=text3d("Y",(0,ySize*1.1,0),color=ycolor,fontsize=fsize) zT=text3d("Z",(0,0,zSize*1.1),color=zcolor,fontsize=fsize) return xAxis+xT+yAxis+yT+zAxis+zT asize = 16 #axis size psize = 12 #plane size #define the normal and distance to origin N=vector((1,1,1)) d=5 #dist to origin #equation of plane is K*<x,y,z>=d ONLY as long as K is a unit normal vector K=N/N.norm() #note syntax for norm of a vector: the norm of A is A.norm() f = K[0]*x+K[1]*y + K[2]*z-d #set this equal to zero to get equation of plane P = plot(axes(asize)) P += implicit_plot3d(f, (x,-psize, psize), (y,-psize, psize), (z, -psize, psize), color='orange', opacity=0.1) P += plot(10*N, width = 3, color='black') #I multiplied N to make it longer. You may want to change this. show('The equation of the plane is $'+latex(f)+'=0$') show('Normal vector is black') show(P, spin =5)
The equation of the plane is 313x+313y+313z−5=0
Normal vector is black
3D rendering not yet implemented
a=(2,-1,5) A=vector(a) b=(1,1,1); c=(0,-2,4) B=vector(b);C=vector(c) A.dot_product(B)
6
A.cross_product(B)
(-6, 3, 3)
B.cross_product(A)
(6, -3, -3)
A.cross_product(B)
(-6, 3, 3)
A.norm()
sqrt(30)
################################################################################ # # INTERACTIVE 3D GRAPHER # # This allows you to graph up to 3 three-dimensional plots at once, choosing whether they # will be parametric, implicit, or z=f(x,y) style. You don't have to do any coding, but instead just fill in # the input boxes. It **should** be self-explanatory and easy to use. # ################################################################################ %var t,x,y,z def axes(size,xcolor='red',ycolor='green',zcolor='blue',asp=(1,1)): #returns a plot of axis vectors, labeled, with size xSize=size; ySize= asp[0]*size; zSize=asp[1]*size xAxis = plot(vector((xSize,0,0)),color=xcolor,width=3) yAxis = plot(vector((0,ySize,0)),color=ycolor,width=3) zAxis = plot(vector((0,0,zSize)),color=zcolor,width=3) #Labels for axes #set fontsize (change if you need to) fsize=20 xT=text3d("X",(xSize*1.1,0,0),color=xcolor,fontsize=fsize) yT=text3d("Y",(0,ySize*1.1,0),color=ycolor,fontsize=fsize) zT=text3d("Z",(0,0,zSize*1.1),color=zcolor,fontsize=fsize) return xAxis+xT+yAxis+yT+zAxis+zT #return xAxis+yAxis+zAxis ############################################################################### def plotPlane(N,P,s,clr='grey',op=0.2,meshvalue=0): x,y,z=var('x y z') #N=normal vec, P = point (vector), s=size #This returns a plot of plane with the given normal vector and point, with x,y,z going plus or minus s units from the point. #optional arguments for color and opacity (default is nearly transparent grey) return implicit_plot3d(N.dot_product(vector((x,y,z)))-N.dot_product(P),(x,P[0]-s,P[0]+s),(y,P[1]-s,P[1]+s),(z,P[2]-s,P[2]+s),color=clr,opacity=op,mesh=meshvalue) #default values color1='red';color2='lightgreen';color3='lightblue' paraf=(sin(t),cos(t),t) implf=x^2+y^2-z^2-4 zfxy =x^2-y^2 tlims =(0,4*pi);xlims = (-3,3);ylims=(-3,3);zlims=(-3,3) html('''<font color='red'>f1 is RED</font> <font color='green'>f2 is GREEN</font> <font color='blue'>f3 is BLUE</font>''') @interact(layout=dict(top=[['f1Show','f1Type','eqn1','tx'],['t1','x1','y1','z1'],['f2Show','f2Type','eqn2','tx2'],['t2','x2','y2','z2'],['f3Show','f3Type','eqn3','tx3'],['t3','x3','y3','z3']])) def _(f1Show=checkbox(false,label= 'f1Show?'), f1Type=selector(['parametric','implicit (mesh)','z=f(x,y) (mesh)'],width=10), eqn1=input_box(paraf,width=50), tx=text_control(' ',label=''), t1=input_box(tlims,width=20,label='t1vals'),x1=input_box(xlims,width=20,label="x1vals"),y1=input_box(ylims,width=18,label="y1vals"),z1=input_box(zlims,width=18,label="z1vals"), f2Show=checkbox(true,label= 'f2Show?'), f2Type=selector(['parametric','implicit (mesh)','z=f(x,y) (mesh)'],default='implicit (mesh)',width=10), eqn2=input_box(implf,width=50), tx2=text_control(' ',label=''), t2=input_box(tlims,width=20,label='t2vals'),x2=input_box(xlims,width=20,label="x2vals"),y2=input_box(ylims,width=18,label="y2vals"),z2=input_box(zlims,width=18,label="z2vals"), f3Show=checkbox(false,label= 'f3Show?'), f3Type=selector(['parametric','implicit (mesh)','z=f(x,y) (mesh)'],default='z=f(x,y) (mesh)',width=10), eqn3=input_box(zfxy,width=50), tx3=text_control(' ',label=''), t3=input_box(tlims,width=20,label='t3vals'),x3=input_box(xlims,width=20,label="x3vals"),y3=input_box(ylims,width=18,label="y3vals"),z3=input_box(zlims,width=18,label="z3vals")): P1=Graphics();P2=Graphics(); P3=Graphics() if f1Show: if f1Type=='parametric': P1=parametric_plot(eqn1,(t,t1[0],t1[1]),color=color1,plot_points=1000,thickness=4) elif f1Type=='implicit (mesh)': P1=implicit_plot3d(eqn1,(x,x1[0],x1[1]),(y,y1[0],y1[1]),(z,z1[0],z1[1]),mesh=1,color=color1) elif f1Type=='z=f(x,y) (mesh)': P1=plot3d(eqn1,(x,x1[0],x1[1]),(y,y1[0],y1[1]),mesh=1,color=color1) if f2Show: if f2Type=='parametric': P2=parametric_plot(eqn2,(t,t2[0],t2[1]),color=color2,plot_points=1000,thickness=4) elif f2Type=='implicit (mesh)': P2=implicit_plot3d(eqn2,(x,x2[0],x2[1]),(y,y2[0],y2[1]),(z,z2[0],z2[1]),mesh=1,color=color2) elif f2Type=='z=f(x,y) (mesh)': P2=plot3d(eqn2,(x,x2[0],x2[1]),(y,y2[0],y2[1]),mesh=1,color=color2) if f3Show: if f3Type=='parametric': P3=parametric_plot(eqn3,(t,t3[0],t3[1]),color=color3,plot_points=1000,thickness=4) elif f3Type=='implicit (mesh)': P3=implicit_plot3d(eqn3,(x,x3[0],x3[1]),(y,y3[0],y3[1]),(z,z3[0],z2[1]),mesh=1,color=color3) elif f3Type=='z=f(x,y) (mesh)': P3=plot3d(eqn3,(x,x3[0],x3[1]),(y,y3[0],y3[1]),mesh=1,color=color3) #axes set up. need size, but don't use limits if not shown xsize=max(f1Show*x1[1],f2Show*x2[1],f3Show*x3[1]) ysize=max(f1Show*y1[1],f2Show*y2[1],f3Show*y3[1]) zsize=max(f1Show*z1[1],f2Show*z2[1],f3Show*z3[1]) ax = axes(xsize,asp=[ysize/xsize, zsize/xsize]) show(P1+P2+P3+ax, spin =2)
f1 is RED f2 is GREEN f3 is BLUE
Interact: please open in CoCalc
#Solution to the extra credit for quiz 1 %var x, y, z def axes(size,xcolor='red',ycolor='green',zcolor='blue'): #returns a plot of axis vectors, labeled, with size xSize=size; ySize= size; zSize=size xAxis = plot(vector((xSize,0,0)),color=xcolor,width=3) yAxis = plot(vector((0,ySize,0)),color=ycolor,width=3) zAxis = plot(vector((0,0,zSize)),color=zcolor,width=3) #Labels for axes #set fontsize (change if you need to) fsize=20 xT=text3d("X",(xSize*1.1,0,0),color=xcolor,fontsize=fsize) yT=text3d("Y",(0,ySize*1.1,0),color=ycolor,fontsize=fsize) zT=text3d("Z",(0,0,zSize*1.1),color=zcolor,fontsize=fsize) return xAxis+xT+yAxis+yT+zAxis+zT asize = 16 psize = 12 a= (1,0,-2); b= (3,1,-2); c= (-5,-1,0); d= (1,1,1) A = vector(a); B= vector(b); C=vector(c); D= vector(d) AB = B-A; AC = C-A N = AB.cross_product(AC) #equation of plane f = N[0]*x+N[1]*y + N[2]*z-N.dot_product(A) #the displacement vector which takes D to its reflection about the plane displacement = N*2*(N.dot_product(A-D))/(N.dot_product(N)) P = plot(axes(asize)) P += implicit_plot3d(f, (x,-psize, psize), (y,-psize, psize), (z, -psize, psize), color='orange', opacity=0.1) P += plot(points([a,b,c,d], thickness=40, color='black')) P += plot(displacement, start = D, width = 3, color='black') #display the coordinates of E and the amount of distance from D to E. print D+displacement, displacement.norm() show(P, spin =5)
(1/9, 25/9, -7/9) 8/3
3D rendering not yet implemented
############################################################################################ # # An example of a parametric curve of a vector function, with tangent vectors, osculating plane, etc. # # ON LINES 69--72, You define a parametric curve r(t), and set the starting and ending t-values. You also include a # value of t in the middle (the constant "c"), a display size for the osculating plane, and you will be rewarded with a parametric plot that shows the # velocity, acceleration, unit tangent, and unit normal vectors, along with the osculating plane, and computes values for # these vectors as well as the curvature. # # The acceleration vector (computed and displayed in lines 80--83, 112--116) is often too large to display nicely. Try experimenting with keeping it commented out # (the default) and instead displaying a vector with magnitude 3 with the ACTUAL magnitude displayed, or commenting out this magnitude display (lines 113--116) and # displaying the actual acceleration vector (line 112). # # THAT'S NOT ALL! You also get a graph showing speed, scalar acceleration, and curvature vs. t, which will help you to understand the important # a = v'T + (kv^2)N equation (the decomposition of the acceleration vector into unit tangent and unit normal components) # ############################################################################################# # NOTE: This example uses functions that produce axes and planes functions defined above. So we will include them in our code here, starting on the next line ######################### # # AXES function # ######################### def axes(size,xcolor='red',ycolor='green',zcolor='blue',asp=(1,1)): #returns a plot of axis vectors, labeled, with length equal to the 'size' argument #Note the OPTIONAL color arguments. If you leave them out, the xyz axes are red, green, blue, respectively #But, for example, the command axes(12,xcolor='black') would produce axes of length 12 with x in black and the others #still in green and blue. #also notice the OPTIONAL axp argument, which is two numbers to adjust the size of the y- and z- axis. For #example, axes(12, asp=(2,.5)) will produce axes with the x-axis of length 12, the y-axis has length 25, and the z-axis is length 6. xSize=size; ySize= asp[0]*size; zSize=asp[1]*size xAxis = plot(vector((xSize,0,0)),color=xcolor,width=3) yAxis = plot(vector((0,ySize,0)),color=ycolor,width=3) zAxis = plot(vector((0,0,zSize)),color=zcolor,width=3) #Labels for axes #set fontsize (change if you need to; not sure if this works) fsize=20 xT=text3d("X",(xSize*1.1,0,0),color=xcolor,fontsize=fsize) yT=text3d("Y",(0,ySize*1.1,0),color=ycolor,fontsize=fsize) zT=text3d("Z",(0,0,zSize*1.1),color=zcolor,fontsize=fsize) return xAxis+xT+yAxis+yT+zAxis+zT #return xAxis+yAxis+zAxis (if we didn't want labels) ############################################################################### # # PLANE drawing function # ###################################### def plotPlane(N,P,s,clr='grey',op=0.2,meshvalue=0): x,y,z=var('x y z') #N=normal vec, P = point (vector), s=size #This returns a plot of plane with the given normal vector and point, with x,y,z going plus or minus s units from the point. #optional arguments for color and opacity and meshvalue (default is nearly transparent grey) return implicit_plot3d(N.dot_product(vector((x,y,z)))-N.dot_product(P),(x,P[0]-s,P[0]+s),(y,P[1]-s,P[1]+s),(z,P[2]-s,P[2]+s),color=clr,opacity=op,mesh=meshvalue) # # end of function definitions # ##################################### #you need to define symbolic variables other than x, to do useful symbolic math like plotting and calculus %var t ################# DEFINE YOUR r(t) PARAMETRIZATION BELOW, along with start and end t and intermediate t-value #define some parametric functions for our curve #this example uses a wiggly curve with variable curvature helix x(t) = (4*cos(t)+0.5)*cos(t); y(t) = (4*cos(t)+0.5)*sin(t); z(t) = t t0=0; t1= 4*pi #starting and ending t values for our curve c= 0.3*pi #a t value in between to draw tangent and normal etc size=1 #the size for the osculating plane. Experiment with different values ######################################### # # FROM THIS POINT ON, EVERYTHING IS AUTOMATED. READ THE CODE TO SEE HOW IT IS DONE # ######################################## r(t) = (x(t),y(t),z(t)) #the position r(t) at time t v=r.diff(t) #velocity! The "diff" suffix does symbolic differentiation; it is a function of t a=r.diff(t,2) #acceleration (the "2" means second derivative); also a function of t #Our convention is to use upper-case for vectors R = vector(r); V= vector(v); A = vector(a) #We need to turn them into vectors in order to use vector operations like cross-product, etc. #Now compute unit tangent and unit normal vectors. #We will compute T, N, B, and T' vectors from textbook #recall that norm() computes magnitude of a vector T = V/norm(V) #unit tangent vector (function of t) Tprime = T.diff(t) #derivative of unit tangent; this is the T' vector in your textbook #turn this into a unit normal N = Tprime/norm(Tprime) # Finally, the little-used binormal vector B B=T.cross_product(N) # So we have produced the following vectors, all functions of t: # R, V, A, T, Tprime, N, B # # Let's also compute the curvature, using the formula k = |T'|/|ds/dt|=|T'|/|V|. #curvature (two ways; should be equal) K = norm(Tprime)/norm(V) # This is a scalar function of t. #create the parametric plot of curve r(t) P = Graphics() # start with a blank slate # Next, include the parametric plot P += parametric_plot3d(r(t),(t,t0,t1),plot_points=400,thickness=6) P += plot(T(t=c),start=R(t=c),color='red') #include the unit tangent at t=c; note starting point P += plot(N(t=c),start=R(t=c),color='black') #add unit normal; note starting point P += plot(B(t=c),start=R(t=c),color='green') #the binormal vector is normal to the osculating plane # the next line allows you to see the acceleration vector, which, like T and N, lies on the osculating plane # However, it is often too long, and you may want to comment out the next line or replace it with the line after it which displays a shortened vector but with the value of the magnitude labeled. #P += plot(A(c),start=R(c),color='magenta') #the ACTUAL acceleration vector P += plot(3.0*A(t=c)/norm(A(t=c)),start=R(t=c),color='magenta') #same direction as A, but with magnitude 3 accel_label=1.1*(3.0*A(c)/norm(A(c))+R(c)) #a location near the endpoint of this magnitude-3 vector to label the actual magnitude P += text3d(str(float(floor(100*norm(A(t=c)))/100)),accel_label,color='magenta') #this messy code eliminates all but a few decimals #plot the osculating plane and choose color and size and opacity P += plotPlane(B(t=c),R(t=c),size,'grey',0.2) #plot the point as well. Note the lower-case r, since it is a point, not a vector. we just want a dot, not an arrow. P += point(r(t=c), thickness=30, color= 'green') #Now compute some stuff about how fast the object is traveling # recall that the speed is the magnitude of velocity vector speed = V(t=c).norm() accel = A(t=c).norm() k = K(t=c) #add axes, sized to be a little taller than the helix reaches P += plot(axes(1.1*t1)) show(P,spin=1) #Now create plots of K, A, V with color-coded legend pK=plot(K,(t,t0,t1),legend_label='curvature',legend_color = 'green', color='green') pA=plot(norm(A(t)),(t,t0,t1),legend_label='acceleration',legend_color = 'magenta', color='magenta') pV=plot(norm(V),(t,t0,t1),legend_label='velocity',legend_color = 'blue', color='blue') #print values at our t=c point. print('At t = %f:'%c) print('speed = %f'%speed) print('scalar acceleration = %f'%accel) print('curvature = %f' %k) #finally, display the K-V-A plot show(pA+pK+pV)
3D rendering not yet implemented
At t = 0.942478:
speed = 4.427318
scalar acceleration = 8.303751
curvature = 0.423226
######################################################### # # interactive TANGENT PLANE PLOTTER. Shows a z=f(x,y) surface # and lets you choose x=a, y=b and draws tangent vectors and tangent plane at (a,b,f(a,b)) # and as a bonus, shows you the partial derivatives and their values. # ########################################################### %var t,x,y,z def axes(size,xcolor='red',ycolor='green',zcolor='blue',asp=(1,1)): #returns a plot of axis vectors, labeled, with size xSize=size; ySize= asp[0]*size; zSize=asp[1]*size xAxis = plot(vector((xSize,0,0)),color=xcolor,width=3) yAxis = plot(vector((0,ySize,0)),color=ycolor,width=3) zAxis = plot(vector((0,0,zSize)),color=zcolor,width=3) #Labels for axes #set fontsize (change if you need to) fsize=20 xT=text3d("X",(xSize*1.1,0,0),color=xcolor,fontsize=fsize) yT=text3d("Y",(0,ySize*1.1,0),color=ycolor,fontsize=fsize) zT=text3d("Z",(0,0,zSize*1.1),color=zcolor,fontsize=fsize) return xAxis+xT+yAxis+yT+zAxis+zT #return xAxis+yAxis+zAxis def plotPlane(N,P,s=1,clr='grey',op=0.2,meshvalue=0): x,y,z=var('x y z') #N=normal vec, P = point (vector), s=size #This returns a plot of plane with the given normal vector and point, with x,y,z going plus or minus s units from the point. #optional arguments for color and opacity (default is nearly transparent grey) return implicit_plot3d(N.dot_product(vector((x,y,z)))-N.dot_product(P),(x,P[0]-s,P[0]+s),(y,P[1]-s,P[1]+s),(z,P[2]-s,P[2]+s),color=clr,opacity=op,mesh=meshvalue) #global values for interactive function defined below xvals = (x,-2,2); yvals = (y,-2,2); ptvals=(1,n(1.5,12)) @interact(layout=dict(top=[['f','xrange','yrange'],['pnt','showTP','showDerivs','showAxes']])) def _(f=input_box(sqrt(9-x^2-y^2),width=40),xrange=input_box(xvals, width=20),yrange=input_box(yvals, width = 20),pnt=input_box(ptvals,label="(a,b)",width=20),showTP=checkbox(false,label="Show tangent plane"),showDerivs=checkbox(false,label="Show derivatives"), showAxes=checkbox(true,label="Show axes")): (a,b) = (pnt[0],pnt[1]) # we will draw tangent plane on surface at (a,b,f(a,b)) #(a,b) = (n(a,12),n(b,12)) #set options for surface plot etc pointColor ='black' surfaceColor = 'grey'; meshvalue = 0; opvalue= 0.4 #make the surface plot P = plot3d(f,xrange, yrange, mesh = meshvalue, color = surfaceColor, opacity = opvalue) # include axes if checked if showAxes: P += axes(1.2*max(xrange[2],yrange[2]),xcolor='grey', ycolor='grey',zcolor='grey') #compute partial derivatives fx(x,y) = f.diff(x) fy(x,y) = f.diff(y) # define point and its vector form pt = (a,b,f(a,b)) Pt = vector(pt) P += plot(point(pt,thickness = 10, color=pointColor)) #compute normal to tangent plane N = vector([fx(a,b),fy(a,b),-1]) #add in the new plane if option checked if showTP: P += plotPlane(N,Pt, s=0.75, clr='yellow', op=0.4) # draw normal vector #P += plot(N/norm(N), start=Pt, color='black') #include traces for y=b, x=a on xy-plane and on surface #the xtrace is parallel to the x-axis, ytrace is parallel to y xtraceColor='magenta'; ytraceColor='green' P += parametric_plot3d((a,t,0),(t,yrange[1],yrange[2]),color=ytraceColor) P += parametric_plot3d((t,b,0),(t,xrange[1],xrange[2]),color=xtraceColor ) P += parametric_plot3d((a,t,f(a,t)),(t,yrange[1],yrange[2]),color=ytraceColor) P += parametric_plot3d((t,b,f(t,b)),(t,xrange[1],xrange[2]),color=xtraceColor ) #next, create the two unit tangent vectors #for the ytrace, we take derivative of the parametric curve (a,t,f(a,t)) ytraceTangent = vector((0,1,fy(a,b))) # for xtrace, differentiate (t, b, f(t,b)) xtraceTangent = vector((1, 0, fx(a,b))) unitYTraceTan=ytraceTangent/norm(ytraceTangent) unitXTraceTan=xtraceTangent/norm(xtraceTangent) P += plot(unitYTraceTan, start = Pt, color='darkgreen') P += plot(unitXTraceTan, start = Pt, color='darkred') if showDerivs: show('$f(x,y) ='+ latex(f)+ ',\quad\quad{\partial f\over\partial x} =' + latex(fx(x,y))+',\quad\quad {\partial f\over\partial y}=' + latex(fy(x,y))+'$') show('At $(a,b)='+latex((a,b))+', \partial f/\partial x ' +'\doteq' + latex(n(fx(a,b),12))+',\quad \partial f/\partial y ' +'\doteq' + latex(n(fy(a,b),12))+'$') #print('At (a,b) = (%s,%s), f_x=%f'%(a,b,fx(a,b))) show(P)
Interact: please open in CoCalc
################################ # # surface plot compared with contour plot # ################################## %var x,y,z,t import matplotlib as mp # pick a function and range to plot it #f=y^2-x^2 f=6*x*y*exp(-(x^2+y^2)/2) #f= x^3 - 3*x*y^2 #monkey saddle!! s=3 xrange=(x,-s,s); yrange=(y,-s,s) #set a colormap for the contour and the matplotlib plot cmap ='coolwarm' #cmap = mp.colors.ListedColormap([0,0,0]) #the first one is a standard 3d plot that we can rotate,etc. plot3d(f, xrange, yrange,aspect_ratio=[1.,1.,1.], mesh=1,color='lightblue') # next plot is fixed, but fancy and prettier plot3d_using_matplotlib(f, xrange, yrange,azim=300,cmap=cmap) #next plot is contour plot which matches second plot in color map cp =contour_plot(f, xrange, yrange,colorbar=True,contours=20,labels=False,fill=true,figsize=6,linestyles='-',cmap=cmap) #let's add in the gradient vectors gv = plot_vector_field(f.gradient(),xrange, yrange) show(cp+gv) #here is another contour style that you may like: black lines with labels import matplotlib as mp blackCmap=mp.colors.ListedColormap([0,0,0]) contour_plot(f, xrange, yrange,colorbar=False,contours=20,labels=True,fill=false,figsize=6,linestyles='-',cmap=blackCmap)
3D rendering not yet implemented
# chain rule practice 1 #DONT FORGET TO DECLARE VARIABLES!! %var x,y,z,t #First define F(x,y,z) F(x,y,z)= x+y^2+y*z #next, express x,y,z as functions of another variable t #in other words, G(t) := F(f(t), g(t), h(t)) f=t^2; g = cos(t); h=3*exp(t^2) G = F(x=f, y=g, z = h) show('$F(x,y,z)='+latex(F(x,y,z))+'$') show('$G(t)=F('+latex(f)+','+latex(g)+','+latex(h)+')='+latex(G(t=t))+'$') show('$\partial G/\partial t='+latex(G.diff(t))+'$')
F(x,y,z)=y2+yz+x
G(t)=F(t2,cos(t),3e(t2))=t2+cos(t)2+3cos(t)e(t2)
∂G/∂t=6tcos(t)e(t2)−2cos(t)sin(t)−3e(t2)sin(t)+2t
# chain rule practice 2 #DONT FORGET TO DECLARE VARIABLES!! %var x,y,z,t,s #First define F(x,y,z) F(x,y,z)= x+y^2+y*z #next, express x,y,z as functions of two other variables t and s #in other words, G(t) := F(f(t,s), g(t,s), h(t,s)) f=t^2+s; g = cos(t*s); h=3*exp(t^2-s^3) G = F(x=f, y=g, z = h) show('$F(x,y,z)='+latex(F(x,y,z))+'$') show('$G(t,s)=F('+latex(f)+','+latex(g)+','+latex(h)+')='+latex(G(t=t,s=s))+'$') show('$\partial G/\partial t='+latex(G.diff(t))+'$') show('$\partial G/\partial s='+latex(G.diff(s))+'$')
F(x,y,z)=y2+yz+x
G(t,s)=F(t2+s,cos(st),3e(−s3+t2))=t2+cos(st)2+3cos(st)e(−s3+t2)+s
∂G/∂t=6tcos(st)e(−s3+t2)−2scos(st)sin(st)−3se(−s3+t2)sin(st)+2t
∂G/∂s=−9s2cos(st)e(−s3+t2)−2tcos(st)sin(st)−3te(−s3+t2)sin(st)+1
#MIDPOINT APPROXIMATION OF SIMPLE DOUBLE INTEGRALS %var x,y # the code for drawing boxes needs to be imported from sage.plot.plot3d.shapes import Box #define the function f=(sin(x+y))^2 #f=13-x^2-y^2 #the function will be integrated over a rectangle centered at (cx, cy) with size sx, sy sx,sy=2,3 cx,cy=0,0 #how much to subdivide in the x direction (m) and y direction (n) m=30 n= 90 boxes=[] dx=2*sx/m dy=2*sy/n vol=0 for i in range(m): xi=(i+0.5)*dx-sx for j in range(n): yj=(j+0.5)*dy-sy h=f(x=xi, y=yj) boxes.append(Box([dx/2, dy/2,abs(h)/2],color=(1+(-1)^(i+j),-1-(-1)^(i+j),.5),opacity=.8,specular=0,shininess=0).translate(xi,yj,h/2)) vol += h*dx*dy #Display the boxes show(sum(boxes),spin=1) #This is the code to actually do the symbolic integral ii=integrate(integrate(f,(x,cx-sx, cx+sx)),(y,cy-sy, cy+sy)) #Display the integral value and the approximation value show('$\int_{'+latex(cy-sy)+'}^{'+latex(cy+sy)+'}\int_{'+latex(cx-sx)+'}^{'+latex(cx+sx)+'}(' +latex(f)+')dx dy='+latex(ii)+'='+latex(ii.n())+'$') show('Midpoint approximation '+'$ = '+latex(vol)+'$')
3D rendering not yet implemented
∫−33∫−22(sin(x+y)2)dxdy=41cos(10)−41cos(2)+12=11.8942688268677
Midpoint approximation =11.8938763056596
#problem 15.6.33 %var x,y,z #reuse axes function from earlier def axes(size,xcolor='red',ycolor='green',zcolor='blue',asp=(1,1)): #returns a plot of axis vectors, labeled, with size xSize=size; ySize= asp[0]*size; zSize=asp[1]*size xAxis = plot(vector((xSize,0,0)),color=xcolor,width=3) yAxis = plot(vector((0,ySize,0)),color=ycolor,width=3) zAxis = plot(vector((0,0,zSize)),color=zcolor,width=3) #Labels for axes #set fontsize (change if you need to) fsize=20 xT=text3d("X",(xSize*1.1,0,0),color=xcolor,fontsize=fsize) yT=text3d("Y",(0,ySize*1.1,0),color=ycolor,fontsize=fsize) zT=text3d("Z",(0,0,zSize*1.1),color=zcolor,fontsize=fsize) return xAxis+xT+yAxis+yT+zAxis+zT #return xAxis+yAxis+zAxis ############################################################################### def plotPlane(N,P,s,clr='grey',op=0.2,meshvalue=0): x,y,z=var('x y z') #N=normal vec, P = point (vector), s=size #This returns a plot of plane with the given normal vector and point, with x,y,z going plus or minus s units from the point. #optional arguments for color and opacity (default is nearly transparent grey) return implicit_plot3d(N.dot_product(vector((x,y,z)))-N.dot_product(P),(x,P[0]-s,P[0]+s),(y,P[1]-s,P[1]+s),(z,P[2]-s,P[2]+s),color=clr,opacity=op,mesh=meshvalue) #default values color1='red';color2='lightgreen';color3='lightblue' eqn1=1-y-z eqn2=y-sqrt(x) xlims = (x,0,1);ylims=(y,0,1);zlims=(z,0,1) P1=Graphics();P2=Graphics(); P3=Graphics() P1+= implicit_plot3d(eqn1,xlims, ylims, zlims,mesh=1,color=color1) P1+= implicit_plot3d(eqn2,xlims, ylims, zlims,mesh=1,color=color2) #axes set up. need size, but don't use limits if not shown ax = axes(1.25) #make some translucent gray planes to bound the figure yzPlane=plotPlane(vector((1,0,0)),vector((0,.5,.5)),.5) xyPlane=plotPlane(vector((0,0,1)),vector((.5,.5,0)),.5) print("Problem 15.6.33") show(P1+P2+P3+ax+yzPlane+xyPlane, spin = 0)
Problem 15.6.33
3D rendering not yet implemented
# 15.6.34 # we will compute some integrals using the function f(x,y,z)=1 and hence find the volume. %var x,y,z #reuse axes function from earlier def axes(size,xcolor='red',ycolor='green',zcolor='blue',asp=(1,1)): #returns a plot of axis vectors, labeled, with size xSize=size; ySize= asp[0]*size; zSize=asp[1]*size xAxis = plot(vector((xSize,0,0)),color=xcolor,width=3) yAxis = plot(vector((0,ySize,0)),color=ycolor,width=3) zAxis = plot(vector((0,0,zSize)),color=zcolor,width=3) #Labels for axes #set fontsize (change if you need to) fsize=20 xT=text3d("X",(xSize*1.1,0,0),color=xcolor,fontsize=fsize) yT=text3d("Y",(0,ySize*1.1,0),color=ycolor,fontsize=fsize) zT=text3d("Z",(0,0,zSize*1.1),color=zcolor,fontsize=fsize) return xAxis+xT+yAxis+yT+zAxis+zT #return xAxis+yAxis+zAxis ############################################################################### def plotPlane(N,P,s,clr='grey',op=0.2,meshvalue=0): x,y,z=var('x y z') #N=normal vec, P = point (vector), s=size #This returns a plot of plane with the given normal vector and point, with x,y,z going plus or minus s units from the point. #optional arguments for color and opacity (default is nearly transparent grey) return implicit_plot3d(N.dot_product(vector((x,y,z)))-N.dot_product(P),(x,P[0]-s,P[0]+s),(y,P[1]-s,P[1]+s),(z,P[2]-s,P[2]+s),color=clr,opacity=op,mesh=meshvalue) #default values color1='red';color2='lightgreen';color3='lightblue' eqn1=y-1+x eqn2=z+x^2-1 xlims = (x,0,1);ylims=(y,0,1);zlims=(z,0,1) P1=Graphics();P2=Graphics(); P3=Graphics() P1+= implicit_plot3d(eqn1,xlims, ylims, zlims,mesh=1,color=color1) P1+= implicit_plot3d(eqn2,xlims, ylims, zlims,mesh=1,color=color2) #axes set up. need size, but don't use limits if not shown ax = axes(1.25) #make some translucent gray planes to bound the figure yzPlane=plotPlane(vector((1,0,0)),vector((0,.5,.5)),.5) xyPlane=plotPlane(vector((0,0,1)),vector((.5,.5,0)),.5) #also create planes parallel to the three coordinate planes to visualize different ways to set up integrals #we make these planes yellow and meshed with high opacity so they are very visible xyPlane2= plotPlane(vector((0,0,1)),vector((.5,.5,.5)),.5,clr='yellow',meshvalue=1, op=.7) yzPlane2= plotPlane(vector((1,0,0)),vector((.5,.5,.5)),.5,clr='yellow',meshvalue=1, op=.7) xzPlane2= plotPlane(vector((0,1,0)),vector((.5,.5,.5)),.5,clr='yellow',meshvalue=1, op=.7) print("Problem 15.6.34") # this is the "dy dz dx" version in the problem statement Iyzx =integrate(integrate(integrate(1,y,0,1-x),z,0,1-x^2),x,0,1) # this is the "dy dx dz" version in the problem statement Iyxz =integrate(integrate(integrate(1,y,0,1-x),x,0,sqrt(1-z)),z,0,1) # this is the "dx dy dz" version Ixyz=integrate(integrate(integrate(1,x,0,1-y),y,1-sqrt(1-z),1),z,0,1)+integrate(integrate(integrate(1,x,0,sqrt(1-z)),y,0,1-sqrt(1-z)),z,0,1) #they should be the same value! print("Iyzx = %s"%Iyzx) print("Ixyz = %s"%Ixyz) print("Iyxz = %s"%Iyxz) show(P1+P2+P3+ax+xyPlane+yzPlane+xyPlane2, spin = 0)
Problem 15.6.34
Iyzx = 5/12
Ixyz = 5/12
Iyxz = 5/12
3D rendering not yet implemented
########################### # These three cells illustrate ways to see how the element of volume depends on dRho, dTheta, and dPhi # # the first two use the spherical_plot3d function which automatically draws a surface of the form rho=f(theta, phi). # the third cell just uses the regular plot3d function, but explicitly specifies the transformation that takes the spherical variables and converts them # into cartesian coordinates (i.e., takes (rho, theta, phi)--->(x,y,z)). This method is very flexible but also pretty sophisticated. # ############################################################################# # this first picture reminds you that longitudes get close to one another at the poles! %var r,theta, phi f(theta,phi)=1 s=0.1 thetaRange=(theta,0,2*pi) phiRange=(phi,0,pi/2) thetaRange1=(theta,0,pi*s) phiRange1=(phi,0,pi/2) S0 = spherical_plot3d(f,thetaRange, phiRange,mesh=0,opacity=0.2, color='yellow') S1=spherical_plot3d(f,thetaRange1, phiRange1,mesh=0, color='lightblue',opacity=0.9) S2 =spherical_plot3d(f*1.1,thetaRange1, phiRange1,color='red',mes0=1,opacity=0.5) show(S0+S1+S2)
3D rendering not yet implemented
# In this second picture, try altering the values for theta and dtheta, etc. to see how the element of volume changes %var theta, phi ########## try altering these values th =0 #starting value of theta deltaTh =0.1*pi ph=pi/4 #starting phi value deltaPh = 0.1*pi r=1 #starting rho deltaR=0.1 ############################## def onSphere(rho,theta, phi): #coords on a sphere return((rho*sin(phi)*cos(theta),rho*sin(phi)*sin(theta), rho*cos(phi))) f(theta,phi)=r #equation of sphere rho=1 thetaRange=(theta,th,th+ deltaTh) phiRange=(phi,ph,ph+ deltaPh) #for a hemisphere thetaRange0=(theta,0,pi*2) phiRange0=(phi,0,pi/2) L1= line([(0,0,1),(0,0,0),(cos(th),sin(th),0),(0,0,0),(cos(th+deltaTh),sin(th+deltaTh),0)]) L2=line([(0,0,0),onSphere(r,th,ph),(0,0,0),onSphere(r, th, ph+deltaPh),(0,0,0),onSphere(r,th+deltaTh,ph),(0,0,0),onSphere(r, th+deltaTh, ph+deltaPh)]) L3 = line([(sin(ph)*cos(th),sin(ph)*sin(th),0),onSphere(r,th,ph)],color='red',linestyle=':') L4 = line([(sin(ph)*cos(th+deltaTh),sin(ph)*sin(th+deltaTh),0),onSphere(r,th+deltaTh,ph)],color='red',linestyle=':') Arc1= parametric_plot3d((sin(ph)*cos(theta),sin(ph)*sin(theta),0),thetaRange) Arc2= parametric_plot3d((sin(ph+ deltaPh)*cos(theta),sin(ph+ deltaPh)*sin(theta),0),thetaRange) S0 = spherical_plot3d(r,thetaRange0, phiRange0,mesh=0,opacity=0.2, color='yellow') #yellow hemisphere S1=spherical_plot3d(r,thetaRange, phiRange, color='lightblue',opacity=0.7) #inner surface S2 =spherical_plot3d(r+deltaR,thetaRange, phiRange,color='red',mesh=0,opacity=0.7) #outer surface show(S0+S1+S2+L1+L2+L3+L4+Arc1+Arc2,spin=5)
3D rendering not yet implemented
%var rho,theta, phi,x,y,z ############ # # This final picture shows two different elements of volume, sitting on a "taco" of "theta wedges" # It also illustrates very sophisticated ways to use the spherical-to-cartesian transformation to draw interesting shapes, where we # declare which variables are independent (the undecleared one, that is not in the final brackets, is the dependent variable) # ###################### ################################## #define three types of spherical trans, depending on the dependent and indep variable #rho is dep, theta and phi are indep SphereRho = (rho*cos(theta)*sin(phi), rho*sin(theta)*sin(phi), rho*cos(phi), [theta,phi]) #phi is dep, theta and rho are indep SpherePhi = (rho*cos(theta)*sin(phi), rho*sin(theta)*sin(phi), rho*cos(phi), [theta,rho]) #theta is dep, phi and rho are indep SphereTheta = (rho*cos(theta)*sin(phi), rho*sin(theta)*sin(phi), rho*cos(phi), [phi,rho]) ################################## #specify starting theta and phi and delta theta and delta phi. th =0 #starting value of theta deltaTh =0.1*pi ph=.1*pi #starting phi value ph2=.4*pi #second starting phi value deltaPh = 0.05*pi #start rho and delta rho rh=1 deltaRh=0.1 thetaRange=(theta,th,th+ deltaTh) phiRange=(phi,ph,ph+ deltaPh) phiRange2=(phi,ph2,ph2+ deltaPh) rhoRange=(rho, rh,rh+deltaRh) ################################### #yellow hemisphere, radius = 1 P0=implicit_plot3d(x^2+y^2+z^2-1,(x,-1,1),(y,-1,1),(z,0,1),color='yellow',opacity=0.1) #"theta wedge" P1=plot3d(th, (phi,0,pi/2),(rho,0,1),transformation=SphereTheta, color='grey', opacity =0.5) P2=plot3d(th+deltaTh, (phi,0,pi/2),(rho,0,1),transformation=SphereTheta, color='brown', opacity =0.3) #first "phi wedge" P3=plot3d(rh, thetaRange,phiRange,transformation=SphereRho, color='green',opacity=0.3) P6=plot3d(rh+deltaRh, thetaRange,phiRange,transformation=SphereRho, color='green',opacity=0.3) P4=plot3d(th,phiRange,rhoRange, transformation=SphereTheta, color='red', opacity =0.3) P5=plot3d(th+deltaTh,phiRange,rhoRange, transformation=SphereTheta, color='red', opacity =0.3) PhiWedge=P3+P6+P4+P5 #second "phi wedge" P7=plot3d(rh, thetaRange,phiRange2,transformation=SphereRho, color='green',opacity=0.3) P8=plot3d(rh+deltaRh, thetaRange,phiRange2,transformation=SphereRho, color='green',opacity=0.3) P9=plot3d(th,phiRange2,rhoRange, transformation=SphereTheta, color='red', opacity =0.3) P10=plot3d(th+deltaTh,phiRange2,rhoRange, transformation=SphereTheta, color='red', opacity =0.3) PhiWedge2=P7+P8+P9+P10 show(P1+P2+P0+PhiWedge+PhiWedge2,spin=0)
3D rendering not yet implemented
#pictures from quiz 5 #Problem 1 (cylindrical) %var r,theta,z zRange=(z,0,2) rRange=(r,0,1) thetaRange=(theta,0,2*pi) ##IMPORTANT NOTE: cylindrical_plot3d allows one ONLY to plot a function of the form r=f(theta, z). So #in this problem, z=2r, so we write r=z/2. Cone=cylindrical_plot3d(z/2,thetaRange,zRange,color='blue', mesh=1) Cyl =cylindrical_plot3d(1,thetaRange,zRange,color='red', mesh=0, opacity=0.4) show(Cone+Cyl)
3D rendering not yet implemented
#quiz 5, problem 2 (spherical) %var rho, theta, phi ################################## #define three types of spherical trans, depending on the dependent and indep variable #rho is dep, theta and phi are indep SphereRho = (rho*cos(theta)*sin(phi), rho*sin(theta)*sin(phi), rho*cos(phi), [theta,phi]) #phi is dep, theta and rho are indep SpherePhi = (rho*cos(theta)*sin(phi), rho*sin(theta)*sin(phi), rho*cos(phi), [theta,rho]) #theta is dep, phi and rho are indep SphereTheta = (rho*cos(theta)*sin(phi), rho*sin(theta)*sin(phi), rho*cos(phi), [phi,rho]) ################################## thetaRange=(theta,0,1.9*pi) #really should be 0 to 2pi but this way we can peek inside! rhoRange=(rho, 1,2) phiRange=(phi,0,pi/6) Cone= plot3d(pi/6, thetaRange, rhoRange,transformation=SpherePhi,color='blue', mesh=1,opacity=0.3) S1=plot3d(1,thetaRange,phiRange, transformation=SphereRho, color='red',mesh=1) S2=plot3d(2,thetaRange,phiRange, transformation=SphereRho, color='green',mesh=1) #'bottoms' of objects S1b=plot3d(1,thetaRange,(phi,pi/6,pi/2), transformation=SphereRho, color='red',opacity=0.1) S2b=plot3d(2,thetaRange,(phi,pi/6,pi/2), transformation=SphereRho, color='green',opacity=0.1) Coneb=plot3d(pi/6, thetaRange, (rho,0,1),transformation=SpherePhi,color='blue', opacity=0.1) show(Cone+Coneb+S1+S2+S1b+S2b)
3D rendering not yet implemented
#INTERACTIVE U,V--->X,Y TRANSFORMATION EXPLORER #INPUT FORMULA FOR X=X(U,V), Y=Y(U,V), A STARTING POINT (U0,Y0) AND THE SIZE AND NUMBER OF DIVISIONS #THE TWO GRAPHS ARE DISPLAYED, ALONG WITH THE JACOBIAN WHICH MEASURES THE AREA SCALING FACTOR. %var u,v,s xdef=u^2-v^2 ydef = 2*u*v uvdef= (2,3) @interact(layout=dict(top=[['x','y','D'],['n','pt','ss']])) def _(x=input_box(xdef, label='x(u,v)'), y=input_box(ydef, label='y(u,v)'), D=input_box(1,label='delta',width=20), n=input_box(5,label='nDivisions',width=20),pt=input_box(uvdef,label='(u0,v0)',width=20),ss=checkbox(false,label='same scale')): u0=pt[0]; v0=pt[1] d=D/n UVPlot=Graphics() XYPlot=Graphics() #add dots to each plot UVPlot += point((u0+d/2,v0+d/2)) XYPlot += point((x(u=u0+d/2,v=v0+d/2),y(u=u0+d/2,v=v0+d/2))) for t in srange(u0-D,u0+D+d,d): UVPlot += line([(t, v0-D),(t,v0+D)]) XYPlot += parametric_plot((x(u=t,v=s),y(u=t,v=s)),(s,pt[1]-D, pt[1]+D)) UVPlot += line([(u0,v0-D),(u0,v0+D)],color='red',thickness=3) UVPlot += line([(u0-D,v0),(u0+D,v0)],color='green',thickness=3) for t in srange(v0-D,v0+D+d,d): UVPlot += line([(u0-D,t),(u0+D,t)]) XYPlot += parametric_plot((x(u=s,v=t),y(u=s,v=t)),(s,pt[0]-D, pt[0]+D)) XYPlot += parametric_plot((x(u=u0,v=s),y(u=u0,v=s)),(s,pt[1]-D, pt[1]+D),color='red',thickness=3) XYPlot += parametric_plot((x(u=s,v=v0),y(u=s,v=v0)),(s,pt[0]-D, pt[0]+D),color='green',thickness=3) L = line([(-3,2), (3,2)],thickness=0) L.axes(false) if ss:#if same scale give uniform xmins, etc xmn=min(XYPlot.xmin(),UVPlot.xmin()) XYPlot.xmin(xmn) UVPlot.xmin(xmn) xmx=max(XYPlot.xmax(),UVPlot.xmax()) XYPlot.xmax(xmx) UVPlot.xmax(xmx) graphics_array([UVPlot,L,XYPlot]).show(aspect_ratio=1) jacobian=x.diff(u)*y.diff(v)-x.diff(v)*y.diff(u) print('UV is left graph, XY is right graph. Jacobian scale factor is %f' %abs(jacobian(u=u0,v=v0))) show(' $\partial x/\partial u='+latex(x.diff(u))+',\quad\partial x/\partial v ='+latex(x.diff(v))+'$') show(' $\partial y/\partial u='+latex(y.diff(u))+',\quad\partial y/\partial v ='+latex(y.diff(v))+'$') show(' $\partial(x,y)/\partial (u,v)='+latex(jacobian)+'$')
Interact: please open in CoCalc
gradient field =
(2x,6y)
3D rendering not yet implemented
vector field =
(x, y, 0)
# PLOTTING VECTOR FIELDS # %var x,y,z,u,v,t f = x^2+3*y^2 #ellipse g = f.gradient() #duh a = 3 # scale cp = contour_plot(f,(x,-a,a),(y,-a,a),fill=false) vf = plot_vector_field(g,(x,-a,a),(y,-a,a),color='blue') show(cp+vf) show("gradient field ="); show(g)
gradient field =
(2x,6y)
# a 3-d example %var x,y,z P=x+y Q=y R=z*y a = 3 # scale xrange=(x,-a,a) yrange=(y,-a,a) zrange=(z,-a,a) plot_vector_field3d((P,Q,R),xrange,yrange,zrange) show("vector field ="); show((P,Q,R))
3D rendering not yet implemented
vector field =
(x+y, y, yz)
%var x,y P=1+sin(y); Q=1+x*cos(y) f=y+x+x*sin(y) show(P.diff(y)) show(Q.diff(x)) vp=plot_vector_field((P,Q), (x,0,3), (y, 0,3),aspect_ratio=1) cp = contour_plot(f,(x,0,3), (y, 0,3),fill=false) show(vp+cp)
cos(y)
cos(y)
%var x,y,t P=2*x*y; Q=x^2 f=y+x+x*sin(y) show(P.diff(y)) show(Q.diff(x)) vp=plot_vector_field((P,Q), (x,0,3.5), (y, 0,3.5),aspect_ratio=1) cp = contour_plot(f,(x,0,3), (y, 0,3),fill=false) circ1=parametric_plot((cos(t)+2,sin(t)+2),(t,0,pi),thickness=3) circ2=parametric_plot((cos(t)+2,1.5*sin(t)+2),(t,pi,2*pi),thickness=3) wiggle=parametric_plot((2*t/(2*pi)+1,2+.1*sin(6*t)),(t,0,2*pi),thickness=3) pts=point([(1,2),(3,2)],size=40) ta=text("A",(2,3.1),fontsize=16) tb=text("B",(2,1.6),fontsize=16) tc=text("C",(1,1),fontsize=16) show(vp+circ1+wiggle+pts+circ2+ta+tb+tc)
2x
2x
#Examples of parametric equations of surfaces %var u,v r=1;R=4; h=.5 nu=10; nv =10 u0=0;u1=4*pi v0=.5; v1=1 urange=(u,u0,u1); vrange=(v,v0,v1) r(u,v) = (R*v*cos(u),R*v*sin(u),u) parametric_plot3d(r(u,v),urange,vrange,mesh=0,opacity=0.2,plotpoints=10)
3D rendering not yet implemented
#INTERACTIVE Parametric Surface EXPLORER #INPUT FORMULA FOR X=X(U,V), Y=Y(U,V), Z(U,V); also input a "central" point (U0,V0) AND THE number of gridlines #TWO GRAPHS ARE DISPLAYED: the rectangular UV "domain space" ALONG WITH THE 3D surface in xyz space. %var u,v,s #defaults : change them if you wish! xdef=cos(u)*v ydef = sin(u)*v zdef = u urangedef=(u,0,2*pi) vrangedef=(v,1,2) uvdef= (pi/4,1.5) #define a function to draw axes. Very useful for understanding plots def axes(size,xcolor='red',ycolor='green',zcolor='blue'): #returns a plot of axis vectors, labeled, with size xSize=size; ySize= size; zSize=size xAxis = plot(vector((xSize,0,0)),color=xcolor,width=3) yAxis = plot(vector((0,ySize,0)),color=ycolor,width=3) zAxis = plot(vector((0,0,zSize)),color=zcolor,width=3) #Labels for axes #set fontsize (change if you need to) fsize=20 xT=text3d("X",(xSize*1.1,0,0),color=xcolor,fontsize=fsize) yT=text3d("Y",(0,ySize*1.1,0),color=ycolor,fontsize=fsize) zT=text3d("Z",(0,0,zSize*1.1),color=zcolor,fontsize=fsize) return xAxis+xT+yAxis+yT+zAxis+zT @interact(layout=dict(top=[['x','y','z','n'],['pt','urange','vrange','showax']])) def _(x=input_box(xdef, label='x(u,v)'), y=input_box(ydef, label='y(u,v)'),z=input_box(zdef, label='z(u,v)'), n=input_box(6,label='nDivisions',width=20),pt=input_box(uvdef,label='(u0,v0)',width=40),urange=input_box(urangedef,label='u range',width=36),vrange=input_box(vrangedef,label='v range',width=36),showax=checkbox(false,label='show axes')): u0=pt[0]; v0=pt[1] Du=urange[2]-urange[1] du=Du/n Dv=vrange[2]-vrange[1] dv=Dv/n UVPlot=Graphics() SPlot= Graphics() #add dots to each plot UVPlot += point((u0+du/2,v0+dv/2),color='black') SPlot += point((x(u=u0+du/2,v=v0+dv/2),y(u=u0+du/2,v=v0+dv/2),z(u=u0+du/2,v=v0+dv/2)),color='black',thickness=5) for t in srange(urange[1],urange[2]+du,du): UVPlot += line([(t, vrange[1]),(t,vrange[2])],color='red') SPlot += parametric_plot((x(u=t,v=s),y(u=t,v=s),z(u=t,v=s)),(s,vrange[1], vrange[2]),color='red') UVPlot += line([(u0,vrange[1]),(u0,vrange[2])],color='red',thickness=3) SPlot += parametric_plot((x(u=u0,v=s),y(u=u0,v=s),z(u=u0,v=s)),(s,vrange[1], vrange[2]),color='red',thickness=5) UVPlot += line([(urange[1],v0),(urange[2],v0)],color='green',thickness=3) SPlot += parametric_plot((x(u=s,v=v0),y(u=s,v=v0),z(u=s,v=v0)),(s,urange[1], urange[2]),color='green',thickness=5) for t in srange(vrange[1],vrange[2]+dv,dv): UVPlot += line([(urange[1],t),(urange[2],t)],color='green') SPlot += parametric_plot((x(u=s,v=t),y(u=s,v=t),z(u=s,v=t)),(s,urange[1], urange[2]),color='green') show('U axis is green, V axis is red. Axes cross at (u0,v0)') SPlot += parametric_plot3d((x,y,z),urange, vrange,color='grey',opacity=0.1) UVPlot.show(aspect_ratio=1) if showax: SPlot += axes(4) show(SPlot)
Interact: please open in CoCalc
ic=icosahedron(size=4, frame=false,frame_size=3, frame_color='red',opacity=0.5) show(ic)
3D rendering not yet implemented
ic
3D rendering not yet implemented