(* :Context: versora` *) (* :Title: Unitary versors *) (* :Author: Renan Cabrera *) (* :Summary: This package has definitions and properties of unitary vectors in cartesian, cylindrical and spherical coordinates systems. There are some examples of use at the end of this package!!!!!!!!! *) (* Restrictions: In order to use this package the default coordinates cannot be changed to others, but I expect to solve this problem in the next version. *) (* :Package Version: (1.1 a) September 2000 This is one of a set of packages to be published with a book at the end of the year 2000. This package was done by the Author Lic. Renan Cabrera in collaboration with the "Academia Nacional de Ciencias de Bolivia" and can be used and distributed free always refering to the author if it is used. *) (* :Mathematica Version: 3.0 *) (* :Warning: This package was not fully tested, so it may have many bugs. If you want to help me and join the developing group of this package or you have any suggestion please contact me at "renanbo@hotmail.com" or "rencabla@ceibo.entelnet.bo " *) BeginPackage["versora`"] e::usage = " e[ _] unitary versor " Spherical::usage = " Spherical Symbol" Cartesian::usage = " Cartesian Symbol" Cylindrical::usage = " Cylindrical Symbol" SetVariables::usage = " SetVariables[] Set default variables " TransformVersors::usage = " TransformVersor[ initial , final ]" TransformCoordinates::usage = " TransformCoordinates[ initial , final]" PE::usage = " PE[a,b] , Escalar Product between a y b " PV::usage = " PV[a,b] , Vectorial Product between a y b " Vec::usage = " Vec[f] , Define Vec[f] as a vector " Rot::usage = " Rotational DifferentiaL Operator " Div::usage = " Divergence Differential Operator " Grad::usage = " Gradient Differential Operator " Laplacian::usage = "Laplacian Differential Operator" CollectVersors::usage = " CollectVersors[ exp] Collects exp in term like {e[x],e[z]....} " SetVariables[] := Block[{}, Print["The default variables are: {x,y,z} ,{r,th,ph} , {rhc,thc,z} for cartesian, spherical and cylindrical coordinates systems"]; x::usage = " defalut cartesian coordinate x"; y::usage = " default cartesian coordinate y"; z::usage = " default cartesian coordinate z"; r::usage = " default radial spherical coordinate r "; th::usage = " default spherical coordinate th"; ph::usage = " default spherical coordiante ph"; rhc::usage = " default radial cylindrical coordinate rho"; thc::usage = " default angular cylindrical coordinate "; ] Begin["`Private`"] VersorsCar = {e[x],e[y],e[z]}; VersorsCil = {e[rhc],e[thc],e[z]}; VersorsEsf = {e[r],e[th],e[ph]}; Format[e[x]]= Subscript[e,x]; Format[e[y]]=Subscript[e,y]; Format[e[z]]=Subscript[e,z]; Format[e[rhc]]=Subscript[e,rhc]; Format[e[thc]]=Subscript[e,thc]; Format[e[r]]=Subscript[e,r]; Format[e[th]]=Subscript[e,th]; Format[e[ph]]=Subscript[e,ph]; CollectVersors[w_]:= Collect[ w , {e[_]...} ]; (*....Transformaciones.........*) MEsf2Car={ {Sin[th]Cos[ph],Sin[th]Sin[ph],Cos[th] }, {Cos[th]Cos[ph],Cos[th]Sin[ph],-Sin[th] }, {-Sin[ph],Cos[ph],0} }; MCil2Car = { {Cos[thc],Sin[thc],0}, {-Sin[thc],Cos[thc],0}, {0,0,1} } ; (*...*) MCar2Esf=Inverse[MEsf2Car]//Simplify; MCar2Cil=Inverse[MCil2Car]//Simplify; TransGeneral[V1_,M_,V2_]:= MapThread[Rule,{ V1 ,M.V2}]; TransformVersors[Spherical,Cartesian]= TransGeneral[VersorsEsf,MEsf2Car,VersorsCar]; TransformVersors[Cartesian,Spherical] = TransGeneral[VersorsCar,MCar2Esf,VersorsEsf]; TransformVersors[Cylindrical,Cartesian] = TransGeneral[VersorsCil,MCil2Car,VersorsCar]; TransformVersors[Cartesian,Cylindrical] = TransGeneral[VersorsCar,MCar2Cil,VersorsCil]; (*..........*) TransformCoordinates[Cartesian,Spherical] = MapThread[Rule,{ { x, y, z} , {r Sin[th]Cos[ph], r Sin[th]Sin[ph] , r Cos[th] } }]; TransformCoordinates[Spherical,Cartesian]= MapThread[Rule,{ {r ,th,ph }, {Sqrt[x^2+y^2+z^2], ArcCos[z/Sqrt[x^2+y^2+z^2]], ArcCos[x/Sqrt[x^2+y^2]] } }]; (*......*) TransformCoordinates[Cartesian,Cylindrical] = MapThread[Rule,{ { x, y} , {rhc Cos[thc], rhc Sin[thc] } }]; TransformCoordinates[Cylindrical,Cartesian] = MapThread[Rule,{ {rhc,thc}, { Sqrt[x^2+y^2], ArcCos[x/Sqrt[x^2+y^2]]} }]; (*...........*) (*...Diferenciación de versores......*) Unprotect[Dt]; Unprotect[D]; Dt[x,x]=1; Dt[x,y]=0; Dt[x,z]=0; Dt[y,x]=0; Dt[y,y]=1; Dt[y,z]=0; Dt[z,x]=0; Dt[z,y]=0; Dt[z,z]=1; Dt[r,r]=1; Dt[r,th]=0; Dt[r,ph]=0; Dt[th,r]=0; Dt[th,th]=1; Dt[th,ph]=0; Dt[ph,r]=0; Dt[ph,th]=0; Dt[ph,ph]=1; Dt[rhc,rhc]=1; Dt[rhc,thc]=0; Dt[rhc,z]=0; Dt[thc,rhc]=0; Dt[thc,thc]=1; Dt[thc,z]=0; Dt[z,rhc]=0; Dt[z,thc]=0; Dt[z,z]=1; (*...Cartersianas....*) Dt/:Dt[ e[x],_]=0; Dt/:Dt[ e[y],_]=0; Dt/:Dt[ e[z],_]=0; Derivative[_][e][x] = 0; Derivative[_][e][y] = 0; Derivative[_][e][z] = 0; (*...Cilinfricas...*) Dt/:Dt[ e[rhc],thc]= e[thc]; Dt/:Dt[ e[thc],thc]= -e[rhc]; Dt/:Dt[ e[rhc],t_]= Dt[thc,t] e[thc]; Dt/:Dt[ e[thc],t_]= - Dt[thc,t]e[rhc]; Dt/:Dt[rhc,t_] Derivative[1][e][rhc] = Dt[ e[rhc],t]; Dt/:Dt[thc,t_] Derivative[1][e][thc] = Dt[ e[thc],t]; Dt/:Dt[rhc] Derivative[1][e][rhc] = Dt[ e[rhc],thc]Dt[thc]; Dt/:Dt[thc] Derivative[1][e][thc] = Dt[ e[thc] ,thc] Dt[thc]; (*...Esfericas...*) Dt/:Dt[ e[r],th]= e[th]; Dt/:Dt[ e[r],ph]= Sin[th] e[ph]; Dt/:Dt[ e[th],th]= - e[r]; Dt/:Dt[ e[th],ph]= Cos[th] e[ph]; Dt/:Dt[ e[ph],ph]= -Sin[th] e[r]-Cos[th] e[th]; Dt/:Dt[ e[r],t_]=Dt[ e[r],th] Dt[th,t]+ Dt[ e[r],ph] Dt[ph,t]; Dt/:Dt[ e[th],t_]=Dt[ e[th],th]Dt[th,t]+ Dt[ e[th],ph] Dt[ph,t]; Dt/: Dt[e[ph],t_]=Dt[ e[ph],ph] Dt[ph,t]; Dt/:Dt[r,t_] Derivative[1][e][r]= Dt[ e[r],t]; Dt/:Dt[th,t_] Derivative[1][e][th]= Dt[ e[th],t]; Dt/:Dt[ph,t_] Derivative[1][e][ph]= Dt[ e[ph],t]; Dt/:Dt[r] Derivative[1][e][r]= Dt[ e[r],th] Dt[th]+ Dt[ e[r],ph] Dt[ph]; Dt/:Dt[th] Derivative[1][e][th]= Dt[ e[th],th] Dt[th]+ Dt[ e[th],ph] Dt[ph]; Dt/:Dt[ph] Derivative[1][e][ph]= Dt[ e[ph],ph] Dt[ph]; D/: D[ e[rhc],rhc] = 0; D/: D[ e[rhc],thc] = e[thc]; D/: D[ e[rhc],z] = 0; D/: D[ e[thc],rhc] = 0; D/: D[ e[thc],thc] = -e[rhc]; D/: D[ e[thc],z] = 0; D/: D[ e[z],rch] = 0; D/: D[ e[z],thc] = 0; D/: D[ e[z],z] = 0; D/: D[ e[r],r] = 0; D/: D[ e[r],th] = e[th]; D/: D[ e[r],ph] = e[ph] Sin[th]; D/: D[ e[th],r] = 0; D/: D[ e[th],th] = - e[r] ; D/: D[ e[th],ph] = e[ph] Cos[th] ; D/: D[ e[ph],r] = 0; D/: D[ e[ph],th] = 0; D/: D[ e[ph],ph] = - e[th] Cos[th] - e[r] Sin[th]; D/: D[a_,{_,0}...] := a; D/: D[ a_ e[w_] , v__ ] := D[a,v]e[w] + a D[ e[w] , v ]; D/: D[a_.( b_. e[w_] + K_ ),v__] := a b D[e[w],v] + D[ a b ,v] e[w] + a D[K,v] + D[a,v]K; D/: D[ e[r_] ,v_ ,w__] := D[ D[e[r],v] , w]; D/: D[ e[r_],{v_,n_} ]:= Nest[D[#,v]& , e[r] , n ]; Protect[D]; Protect[Dt]; (*.........Producto Escalar........*) SetAttributes[PE,{Orderless}] Format[PE[x_,y_]]:=SmallCircle[x,y] PE[a_.e[x_],b_.e[x_]]=a*b; PE[a_.e[x],b_.e[y]]=0; PE[a_.e[x],b_.e[z]]=0; PE[a_.e[y],b_.e[z]]=0; PE[a_.e[r],b_.e[th]]=0; PE[a_.e[th],b_.e[ph]]=0; PE[a_.e[ph],b_.e[r]]=0; PE[a_.e[rhc],b_.e[thc]]=0; PE[a_.e[thc],b_.e[z]]=0; PE[a_.e[z],b_.e[rhc]]=0; PE[a_ e[w_],b_ e[v_] ]:= a*b PE[e[w],e[v]]; PE[a_ e[w_],e[v_] ]:= a PE[e[w],e[v]]; PE[0,_]=0; (*......*) PE[a_.*w_Plus,u_]:= Distribute[PEX[w,u//Expand]]/. PEX[g__]->a*PE[g] (*.......Producto Vectorial...*) Format[PV[x_,y_]]:=Cross[x,y] PV[a_.e[x_],b_.e[x_]]=0; PV[a_.e[x],b_.e[y]]=a*b e[z]; PV[a_.e[y],b_.e[z]]=a*b e[x]; PV[a_.e[z],b_.e[x]]=a*b e[y]; PV[a_.e[y],b_.e[x]]=-a*b e[z]; PV[a_.e[z],b_.e[y]]=-a*b e[x]; PV[a_.e[x],b_.e[z]]=-a*b e[y]; (*...*) PV[a_.e[r],b_.e[th]]=a*b e[ph]; PV[a_.e[th],b_.e[ph]]=a*b e[r]; PV[a_.e[ph],b_.e[r]]=a*b e[th]; PV[a_.e[th],b_.e[r]]=-a*b e[ph]; PV[a_.e[ph],b_.e[th]]=-a*b e[r]; PV[a_.e[r],b_.e[ph]]=-a*b e[th]; (*...*) PV[a_.e[rhc],b_.e[thc]]=a*b e[z]; PV[a_.e[z],b_.e[rhc]]=a*b e[thc]; PV[a_.e[thc],b_.e[z]]=a*b e[rhc]; PV[a_.e[thc],b_.e[rhc]]=-a*b e[z]; PV[a_.e[rhc],b_.e[z]]=-a*b e[thc]; PV[a_.e[z],b_.e[thc]]=-a*b e[rhc]; (*...*) PV[a_.*w_Plus,u_]:= Distribute[PVX[w,u//Expand]]/. PVX[g__]->a*PV[g] PV[u_,a_.*w_Plus]:= Distribute[PVX[u//Expand,w]]/. PVX[g__]->a*PV[g] PV[0,_]=0; PV[_,0]=0; (*...............................................................*) (*...............................................................*) Format[Grad[f_]]= Del[f]; Format[ Vec[f_] ]= Overscript[f,"\[Rule]"]; Format[ Div[ f_ ]]= SmallCircle[ "\[Del]" , f ]; Format[Rot[ f_ ] ] = Cross[ "\[Del]" ,f ]; Format[Laplacian[f_]]= Power["\[Del]",2] f; PEDel[A_,B_] = SequenceForm["(", SmallCircle[A,"\[Del]"] ,")" B] Grad[f_ + g_]:= Grad[f]+Grad[g]; Grad[ f_*g_]:=f Grad[g]+Grad[f] g; Div[ b_.(a_.Vec[F_] + K_) ]:= Div[ b a Vec[F]]+ Div[b K ] Div[f_ Vec[F_]] := PE[ Grad[f] , Vec[F] ]+ f Div[ Vec[F] ] Div[ Grad[f_] ]:= Laplacian[f] Rot[ b_.(a_.Vec[F] + K_) ]:= Rot[b a Vec[F]]+Rot[b K ] Rot[f_ Vec[F_]] := PV[ Grad[f] , Vec[F] ]+ f Rot[ Vec[F] ] Div[ a_. Rot[ Vec[A_] , b_.Vec[B_] ]] := PE[b Vec[B], Rot[a Vec[A]]]-PE[a A,Rot[ b Vec[B] ]] Rot[Grad[f_]]= 0; Rot[PV[ Vec[A_] ,Vec[B_] ]]:= PEDel[ Vec[B] , Vec[A] ]-B Div[Vec[A]]-PEDel[ Vec[A] , Vec[B] ]+ A Div[Vec[B]] Grad[ PE[Vec[A_] ,Vec[B_]] ]:= PEDel[Vec[B],Vec[A]]+ PEDel[Vec[A],Vec[B]]+PV[Vec[B],Rot[Vec[A]]]+ PV[Vec[A],Rot[Vec[B]]] Div[ Rot[A_] ]=0; Grad[Cartesian]= Grad[f_]:> D[f,x] e[x]+ D[f,y] e[y]+ D[f,z] e[z]; Grad[Spherical]= Grad[f_]:> D[f,r] e[r]+ D[f,th]/r e[th]+ D[f,ph]/(r Sin[th]) e[ph]; Grad[Cylindrical]= Grad[f_]:> D[f,rhc] e[rhc]+ D[f,thc]/rhc e[thc]+ D[f,z] e[z]; Div[Cartesian] :={ Div[ a_.e[x] ]:>D[a,x], Div[a_. e[y] ]:>D[a,y], Div[ a_. e[z] ]:>D[a,z]} Div[Spherical]:={ Div[ a_.e[r] ]:>D[ r^2 Sin[th] a,r]/(r^2 Sin[th]), Div[a_. e[th] ]:> D[r Sin[th] a,th]/(r^2 Sin[th]), Div[ a_. e[ph] ]:> D[ r a,ph]/(r^2 Sin[th]) } Div[Cylindrical]:={ Div[ a_.e[rhc] ]:>D[ rhc a,rhc]/rhc , Div[a_. e[thc] ]:>D[a,th]/rhc, Div[ a_. e[z] ]:>D[a,z] } Div[ b_.(a_.e[x_] + K_)]:= Div[b a e[x] ]+ Div[ b K] Rot[Cartesian]:= { Rot[a_. e[x]]:> e[y]D[a, z] - e[z]D[a,y], Rot[a_. e[y]]:> e[z]D[a, x] - e[x]D[a,z], Rot[a_. e[z]]:> e[x]D[a, y] - e[y]D[a,x] } Rot[Spherical]:={ Rot[a_. e[r]] :> e[th]D[a,ph]/(r Sin[th]) - e[ph]D[ a,th]/r , Rot[a_. e[th]] :> e[ph]D[ r a,r]/r - e[r]D[ r a,ph]/(r^2 Sin[th] ), Rot[a_. e[ph]] :> e[r]D[r Sin[th] a,th]/(r^2 Sin[th] ) - e[th]D[r Sin[th] a,r]/(r Sin[th]), } Rot[Cylindrical]:={ Rot[ a_. e[rhc] ] :> e[thc] D[ a,z] - e[z]/rhc D[ a ,thc], Rot[ a_. e[thc] ] :> e[z]/rhc D[ rhc a,r] - e[r]/rhc D[rhc a ,z], Rot[ a_. e[z] ] :> e[rhc]/rhc D[ a,thc] - e[thc] D[ a ,rhc] } Rot[ b_.(a_.e[x_] + K_)]:= Rot[b a e[x] ]+ Rot[ b K] Laplacian[Cartesian]:= {Laplacian[f_]:> D[ f , {x,2} ] + D[ f , {y,2} ] + D[ f , {z,2} ] } (*xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx*) (*........Examples......... Needs["versora`","f:\\renan\\fiscom\\book\\versora.m"] SetVariables[] Format[ph]=\[CurlyPhi]; Format[th]=\[Theta]; Format[rhc]=\[Rho]; Format[thc]= \[Theta]c; x e[x]+y e[y]+z e[z] /. TransformVersors[Cartesian,Spherical] /. TransformCoordinates[Cartesian,Spherical] //Simplify PowerExpand[Expand[ r e[r]/. TransformVersors[Spherical,Cartesian]/.TransformCoordinates[ Spherical,Cartesian]]] x e[ x]+y e[y]+z e[ z] /. TransformVersors[Cartesian,Cylindrical] /. TransformCoordinates[Cartesian,Cylindrical] //Simplify PowerExpand[Expand[ rhc e[rhc]+z e[z]/. TransformVersors[Cylindrical,Cartesian]/.TransformCoordinates[ Cylindrical,Cartesian]]] Dt[ rhc e[rhc]+z e[z],{t,2}] //CollectVersors Dt[ r e[r],{t,2}] //CollectVersors D[b e[th],r] D[c(b e[th]+e[r]),th] D[b e[th],ph] PE[ f(k e[x]+ c e[z]),b e[y]+d e[z]+2e[x]] PE[ e[r], e[x]+e[ph] +2e[r]] PV[ e[r],2( e[th]+e[z])] PV[ a e[rhc] , e[z]+e[thc] ] PV[e[x],e[r]] Grad[F*G*H]//Expand Grad[x*y*z+F] /. Grad[Cartesian] Grad[ 1/r ]/.Grad[Spherical] Div[ b Vec[F] +a(Vec[G]+Vec[H]) ] Div[Grad[f]] Rot[a Vec[F]+Vec[G]] Div[ a Rot[ Vec[A] ,Vec[B] ] ] Div[ x ^2 e[x] +q(y^2 e[y]+z e[z]) ] % /.Div[Cartesian] Div[ 1/r e[r] ] %/.Div[Spherical] Rot[ PV[ Wx e[x] + Wy e[y]+ Wz e[z] , x e[x]+y e[y]+z e[z]] ] % /.Rot[Cartesian] BSun = e[r]/r^2 - a Sin[th]/r e[ph] PV[Dt[r e[r]] , BSun ]//CollectVersors PV[Vec[A],Vec[B]] Rot[%] Div[Grad[T[x,y,z]]] %/.Laplacian[Cartesian] PE[Grad[F[x,y,z]],Grad[G[x,y,z]]] %/.Grad[Cartesian] .....................................*) Protect[PE,PV,Grad,Div,Rot,Laplacian]; Protect[e]; End[] EndPackage[]