From: Daniel V. <dan...@ya...> - 2017-01-05 19:57:37
|
Hi Viktor, I'm trying to calculate the curl using the tensor packages. This is my code: reset()$kill(all)$ if get('itensor,'version)=false then load(itensor)$ if get('ctensor,'version)=false then load(ctensor)$ dim:3$ imetric(g)$ t:contract(expand('levi_civita([i,j,k],[])*g([],[m,j])*g([],[n,k])*covdiff(H([m],[]),n)))$ t:contract(expand(ev(t,'levi_civita)))$ t:contract(expand(ev(t,'kdelta)))$ t:contract(expand(t))$ ishow(t)$ ct_coordsys(spherical,all)$ remcomps(g([],[i,j]))$ components(g([],[i,j]),ug)$ levi_civita[i,j,k]:=levi_civita([i,j,k])/sqrt(gdet)$ V:[V_r,V_theta,V_phi]$ depends(V,ct_coords)$ ldisplay(V:list_matrix_entries(V.sqrt(lg)))$ EQ:ic_convert(B([i],[])=t); The last command complains about H having both co and contravariant component, which is true.What should I do to have only one type of H component in t? Thanks, Daniel Volinski |
From: Valery P. <pi...@is...> - 2017-01-06 04:03:58
|
Hi, Current limitation of ctensor does not allow to calculate components of arbitrary expresions containing levi_civita tensor. We should use covariant components of vector (contravariant levi_civita symbols) The code below example of curl of the axisymmetric magnetic field H=phi B+ curl (a phi) ############### load(facexp); load("itensor.lisp"); if get('ctensor, 'version) = false then load(ctensor); idim(3); remcomps(g); imetric(g); remcomps(emfd); components(emfd([],[i]), 'levi_civita([],[i,n,m])*covdiff(B([m],[]),n) ); eqn :ishow(emf([],[i]) =canform(rename(expand(emfd([],[i])))))$ display2d:false; emfc : ic_convert(eqn); ct_coords:[r,th,ph]; depends([af,br,bt,bf,af_t,af_r],[r,th]); B:[diff(af,th)/r^2/sin(th),-diff(af,r)/sin(th),bf*r*sin(th)]; lg:matrix([1,0,0],[0,r^2,0],[0,0,r^2*sin(th)^2])$ dlt:matrix([1,0,0],[0,1,0],[0,0,1])$ ug:invert(lg)$ sgdet:r^2*sin(th); christof(all); /*Note the replacement 'levi_civita[%1,%2,i] -> lc0([%1,%2,i])/sgdet */ for i thru dim do emf[i]:sum(sum(lc0([%1,%2,i])/sgdet*diff(B[%2],ct_coords[%1]),%1, 1,dim),%2,1,dim) -sum(sum(sum(lc0([%1,%2,i])/sgdet*mcs[%1,%2,%3]*B[%3],%1,1, dim),%2,1,dim),%3,1,dim)$ display2d:true; collecttermsl(ratexpand(trigsimp(ratexpand(emf[1]))), ['diff(bf,th),'diff(bf,r),'diff(af,th),'diff(af,r),'diff(af,th,2),'diff(af,r,2), 'diff('diff(af,th),r)]); collecttermsl(ratexpand(trigsimp(ratexpand(emf[2]*r^2))), ['diff(bf,th),'diff(bf,r),'diff(af,th),'diff(af,r),'diff(af,th,2),'diff(af,r,2), 'diff('diff(af,th),r)]); collecttermsl(ratexpand(trigsimp(ratexpand(emf[3]*r^2*sin(th)^2))), ['diff(bf,th),'diff(bf,r),'diff(af,th),'diff(af,r),'diff(af,th,2),'diff(af,r,2)]); ##################### > Hi Viktor, > > I'm trying to calculate the curl using the tensor packages. > > This is my code: > > reset()$kill(all)$ > if get('itensor,'version)=false then load(itensor)$ > if get('ctensor,'version)=false then load(ctensor)$ > dim:3$ > imetric(g)$ > t:contract(expand('levi_civita([i,j,k],[])*g([],[m,j])*g([],[n,k])*covdiff(H([m],[]),n)))$ > t:contract(expand(ev(t,'levi_civita)))$ > t:contract(expand(ev(t,'kdelta)))$ > t:contract(expand(t))$ > ishow(t)$ > ct_coordsys(spherical,all)$ > remcomps(g([],[i,j]))$ > components(g([],[i,j]),ug)$ > levi_civita[i,j,k]:=levi_civita([i,j,k])/sqrt(gdet)$ > V:[V_r,V_theta,V_phi]$ > depends(V,ct_coords)$ > ldisplay(V:list_matrix_entries(V.sqrt(lg)))$ > EQ:ic_convert(B([i],[])=t); > > The last command complains about H having both co and contravariant > component, which is true. > What should I do to have only one type of H component in t? > > Thanks, > > Daniel Volinski > > > > > > ------------------------------------------------------------------------------ > Check out the vibrant tech community on one of the world's most > engaging tech sites, SlashDot.org! http://sdm.link/slashdot > > > _______________________________________________ > Maxima-discuss mailing list > Max...@li... > https://lists.sourceforge.net/lists/listinfo/maxima-discuss |
From: Daniel V. <dan...@ya...> - 2017-01-06 12:13:57
|
Hi Valery, Thank for you explanations I use you code and then modify to suit my needsI think I got it, at least it is the same result given by the "vect" package. reset()$kill(all)$ derivabbrev:false$ if get('itensor,'version)=false then load(itensor)$ if get('ctensor,'version)=false then load(ctensor)$ dim:3$ imetric(g)$ ishow(B([],[i])='levi_civita([],[i,j,k])*covdiff(V([j]),k))$ Eq:subst(lc(i,j,k),'levi_civita[i,j,k],ic_convert(%)); assume(sin(theta)>=0)$ ct_coordsys(spherical,all)$ lc(i,j,k):=lc0([i,j,k])/sqrt(gdet)$ V:[Vr,Vt,Vp]$ depends(V,ct_coords)$ ldisplay(V:list_matrix_entries(expand(trigsimp(V.sqrt(lg)))))$ B:[0,0,0]$ ev(Eq)$ ldisplay(B:list_matrix_entries(expand(trigsimp(-B.sqrt(lg)))))$ Thanks again, Daniel Volinski El Viernes 6 de enero de 2017 6:06, Valery Pipin <pi...@is...> escribió: Hi, Current limitation of ctensor does not allow to calculate components of arbitrary expresions containing levi_civita tensor. We should use covariant components of vector (contravariant levi_civita symbols) The code below example of curl of the axisymmetric magnetic field H=phi B+ curl (a phi) ############### load(facexp); load("itensor.lisp"); if get('ctensor, 'version) = false then load(ctensor); idim(3); remcomps(g); imetric(g); remcomps(emfd); components(emfd([],[i]), 'levi_civita([],[i,n,m])*covdiff(B([m],[]),n) ); eqn :ishow(emf([],[i]) =canform(rename(expand(emfd([],[i])))))$ display2d:false; emfc : ic_convert(eqn); ct_coords:[r,th,ph]; depends([af,br,bt,bf,af_t,af_r],[r,th]); B:[diff(af,th)/r^2/sin(th),-diff(af,r)/sin(th),bf*r*sin(th)]; lg:matrix([1,0,0],[0,r^2,0],[0,0,r^2*sin(th)^2])$ dlt:matrix([1,0,0],[0,1,0],[0,0,1])$ ug:invert(lg)$ sgdet:r^2*sin(th); christof(all); /*Note the replacement 'levi_civita[%1,%2,i] -> lc0([%1,%2,i])/sgdet */ for i thru dim do emf[i]:sum(sum(lc0([%1,%2,i])/sgdet*diff(B[%2],ct_coords[%1]),%1, 1,dim),%2,1,dim) -sum(sum(sum(lc0([%1,%2,i])/sgdet*mcs[%1,%2,%3]*B[%3],%1,1, dim),%2,1,dim),%3,1,dim)$ display2d:true; collecttermsl(ratexpand(trigsimp(ratexpand(emf[1]))), ['diff(bf,th),'diff(bf,r),'diff(af,th),'diff(af,r),'diff(af,th,2),'diff(af,r,2), 'diff('diff(af,th),r)]); collecttermsl(ratexpand(trigsimp(ratexpand(emf[2]*r^2))), ['diff(bf,th),'diff(bf,r),'diff(af,th),'diff(af,r),'diff(af,th,2),'diff(af,r,2), 'diff('diff(af,th),r)]); collecttermsl(ratexpand(trigsimp(ratexpand(emf[3]*r^2*sin(th)^2))), ['diff(bf,th),'diff(bf,r),'diff(af,th),'diff(af,r),'diff(af,th,2),'diff(af,r,2)]); ##################### > Hi Viktor, > > I'm trying to calculate the curl using the tensor packages. > > This is my code: > > reset()$kill(all)$ > if get('itensor,'version)=false then load(itensor)$ > if get('ctensor,'version)=false then load(ctensor)$ > dim:3$ > imetric(g)$ > t:contract(expand('levi_civita([i,j,k],[])*g([],[m,j])*g([],[n,k])*covdiff(H([m],[]),n)))$ > t:contract(expand(ev(t,'levi_civita)))$ > t:contract(expand(ev(t,'kdelta)))$ > t:contract(expand(t))$ > ishow(t)$ > ct_coordsys(spherical,all)$ > remcomps(g([],[i,j]))$ > components(g([],[i,j]),ug)$ > levi_civita[i,j,k]:=levi_civita([i,j,k])/sqrt(gdet)$ > V:[V_r,V_theta,V_phi]$ > depends(V,ct_coords)$ > ldisplay(V:list_matrix_entries(V.sqrt(lg)))$ > EQ:ic_convert(B([i],[])=t); > > The last command complains about H having both co and contravariant > component, which is true. > What should I do to have only one type of H component in t? > > Thanks, > > Daniel Volinski > > > > > > ------------------------------------------------------------------------------ > Check out the vibrant tech community on one of the world's most > engaging tech sites, SlashDot.org! http://sdm.link/slashdot > > > _______________________________________________ > Maxima-discuss mailing list > Max...@li... > https://lists.sourceforge.net/lists/listinfo/maxima-discuss ------------------------------------------------------------------------------ Check out the vibrant tech community on one of the world's most engaging tech sites, SlashDot.org! http://sdm.link/slashdot _______________________________________________ Maxima-discuss mailing list Max...@li... https://lists.sourceforge.net/lists/listinfo/maxima-discuss |
From: Viktor T. T. <vt...@vt...> - 2017-01-09 05:30:04
|
Daniel, Forgive me for not responding earlier. One thing you could do would be a hack of sorts: t:ev(t,H([m],[],n)=H1([m],[],n)); This replaces H with H1 in some places, eliminating the conflicting indices. But there is a better way: first, simplify your expression further using itensor, using in particular the symmetry properties of the metric. Also, don't forget to declare your H to be a function of the coordinates, and yes, you need to explicitly assign values to kdelta, too. Here it is: reset()$kill(all)$ if get('itensor,'version)=false then load(itensor)$ if get('ctensor,'version)=false then load(ctensor)$ dim:3$ imetric(g)$ t:contract(expand('levi_civita([i,j,k],[])*g([],[m,j])*g([],[n,k])*covdiff(H([m],[]),n)))$ t:contract(expand(ev(t,'levi_civita)))$ t:contract(expand(ev(t,'kdelta)))$ t:contract(expand(t))$ ishow(t)$ decsym(g,2,0,[sym(all)],[]); decsym(g,0,2,[],[sym(all)]); ishow(canform(rename(t)))$ t:ishow(rename(canform(contract(canform(rename(t))))))$ ct_coordsys(spherical,all)$ remcomps(g([],[i,j]))$ components(g([],[i,j]),ug)$ levi_civita[i,j,k]:=levi_civita([i,j,k])/sqrt(gdet)$ V:[V_r,V_theta,V_phi]$ depends(V,ct_coords)$ ldisplay(V:list_matrix_entries(V.sqrt(lg)))$ kdelta:ident(3); depends(H,ct_coords); EQ:ic_convert(B([i],[])=t); B:[x,x,x]; derivabbrev:true; ev(EQ,'kdelta); B; I believe this yields the result that you were expecting to see. Viktor > -----Original Message----- > From: Daniel Volinski [mailto:dan...@ya...] > Sent: Friday, January 6, 2017 6:58 AM > To: Maxima Discuss <max...@li...> > Subject: Re: [Maxima-discuss] Curl in tensor notation > > Hi Valery, > > Thank for you explanations > > I use you code and then modify to suit my needs > I think I got it, at least it is the same result given by the "vect" package. > > reset()$kill(all)$ > derivabbrev:false$ > if get('itensor,'version)=false then load(itensor)$ > if get('ctensor,'version)=false then load(ctensor)$ > dim:3$ > imetric(g)$ > ishow(B([],[i])='levi_civita([],[i,j,k])*covdiff(V([j]),k))$ > Eq:subst(lc(i,j,k),'levi_civita[i,j,k],ic_convert(%)); > > assume(sin(theta)>=0)$ > ct_coordsys(spherical,all)$ > lc(i,j,k):=lc0([i,j,k])/sqrt(gdet)$ > V:[Vr,Vt,Vp]$ > depends(V,ct_coords)$ > ldisplay(V:list_matrix_entries(expand(trigsimp(V.sqrt(lg)))))$ > B:[0,0,0]$ > ev(Eq)$ > ldisplay(B:list_matrix_entries(expand(trigsimp(-B.sqrt(lg)))))$ > > Thanks again, > > Daniel Volinski > > > > > El Viernes 6 de enero de 2017 6:06, Valery Pipin <pi...@is...> escribió: > > > > Hi, > > Current limitation of ctensor does not allow to calculate components of > arbitrary expresions > containing levi_civita tensor. We should use covariant components of > vector (contravariant levi_civita symbols) > The code below example of curl of the axisymmetric magnetic field H=phi > B+ curl (a phi) > > ############### > load(facexp); > load("itensor.lisp"); > if get('ctensor, 'version) = false then load(ctensor); > > idim(3); > remcomps(g); > imetric(g); > > remcomps(emfd); > components(emfd([],[i]), 'levi_civita([],[i,n,m])*covdiff(B([m],[]),n) ); > > eqn :ishow(emf([],[i]) =canform(rename(expand(emfd([],[i])))))$ > > display2d:false; > emfc : ic_convert(eqn); > > ct_coords:[r,th,ph]; > depends([af,br,bt,bf,af_t,af_r],[r,th]); > B:[diff(af,th)/r^2/sin(th),-diff(af,r)/sin(th),bf*r*sin(th)]; > > > lg:matrix([1,0,0],[0,r^2,0],[0,0,r^2*sin(th)^2])$ > dlt:matrix([1,0,0],[0,1,0],[0,0,1])$ > ug:invert(lg)$ > sgdet:r^2*sin(th); > christof(all); > > /*Note the replacement 'levi_civita[%1,%2,i] -> lc0([%1,%2,i])/sgdet */ > for i thru dim do > emf[i]:sum(sum(lc0([%1,%2,i])/sgdet*diff(B[%2],ct_coords[%1]),%1, > 1,dim),%2,1,dim) > -sum(sum(sum(lc0([%1,%2,i])/sgdet*mcs[%1,%2,%3]*B[%3],%1,1, > dim),%2,1,dim),%3,1,dim)$ > display2d:true; > > collecttermsl(ratexpand(trigsimp(ratexpand(emf[1]))), > ['diff(bf,th),'diff(bf,r),'diff(af,th),'diff(af,r),'diff(af,th,2),'diff(af,r,2), > 'diff('diff(af,th),r)]); > > collecttermsl(ratexpand(trigsimp(ratexpand(emf[2]*r^2))), > ['diff(bf,th),'diff(bf,r),'diff(af,th),'diff(af,r),'diff(af,th,2),'diff(af,r,2), > 'diff('diff(af,th),r)]); > > collecttermsl(ratexpand(trigsimp(ratexpand(emf[3]*r^2*sin(th)^2))), > ['diff(bf,th),'diff(bf,r),'diff(af,th),'diff(af,r),'diff(af,th,2),'diff(af,r,2)]); > > ##################### > > Hi Viktor, > > > > I'm trying to calculate the curl using the tensor packages. > > > > This is my code: > > > > reset()$kill(all)$ > > if get('itensor,'version)=false then load(itensor)$ > > if get('ctensor,'version)=false then load(ctensor)$ > > dim:3$ > > imetric(g)$ > > t:contract(expand('levi_civita([i,j,k],[])*g([],[m,j])*g([],[n,k])*covdiff(H([m],[]),n)))$ > > t:contract(expand(ev(t,'levi_civita)))$ > > t:contract(expand(ev(t,'kdelta)))$ > > t:contract(expand(t))$ > > ishow(t)$ > > ct_coordsys(spherical,all)$ > > remcomps(g([],[i,j]))$ > > components(g([],[i,j]),ug)$ > > levi_civita[i,j,k]:=levi_civita([i,j,k])/sqrt(gdet)$ > > V:[V_r,V_theta,V_phi]$ > > depends(V,ct_coords)$ > > ldisplay(V:list_matrix_entries(V.sqrt(lg)))$ > > EQ:ic_convert(B([i],[])=t); > > > > The last command complains about H having both co and contravariant > > component, which is true. > > What should I do to have only one type of H component in t? > > > > Thanks, > > > > Daniel Volinski > > > > > > > > > > > > > ------------------------------------------------------------------------------ > > Check out the vibrant tech community on one of the world's most > > engaging tech sites, SlashDot.org! http://sdm.link/slashdot > > > > > > _______________________________________________ > > Maxima-discuss mailing list > > Max...@li... <mailto:Max...@li...> > > https://lists.sourceforge.net/lists/listinfo/maxima-discuss > > > > ------------------------------------------------------------------------------ > Check out the vibrant tech community on one of the world's most > engaging tech sites, SlashDot.org! http://sdm.link/slashdot > _______________________________________________ > Maxima-discuss mailing list > Max...@li... <mailto:Max...@li...> > https://lists.sourceforge.net/lists/listinfo/maxima-discuss > > |
From: Daniel V. <dan...@ya...> - 2017-01-10 14:31:44
|
Hi Viktor, Thank you for your response. 1) I found a bug in Maxima for Windows regarding the first lines of your code which produces the following error: Maxima encountered a Lisp error: 3 is not a string designator. The bug report is #3272. 2) Valery Pipin helped me with another approach that is to snatch the 'levi_civita symbol right at the beginning and substitute it by a function that is defined later on when the metric is set within the ctensor package. 3) I was able to test your code using CESGA - Maxima on line | | | CESGA - Maxima on line | | | and it works fine except for the fact that the last line is not the curl of the vector H.I think that we had this conversation before when someone was dealing with hodge. Thanks, Daniel Volinski |
From: Viktor T. T. <vt...@vt...> - 2017-01-10 22:02:22
|
Daniel, That bug is really "by design", so to speak: a limitation of levi_civita, which does something the itensor package isn't designed to handle, replacing a symbolic (dummy) index with an explicit numerical value. That said, curiously the code actually works in GCL, though it fails in all other lisps to which I have access. Viktor > -----Original Message----- > From: Daniel Volinski [mailto:dan...@ya...] > Sent: Tuesday, January 10, 2017 9:28 AM > To: 'Maxima Discuss' <max...@li...> > Subject: Re: [Maxima-discuss] Curl in tensor notation > > > > Hi Viktor, > > Thank you for your response. > > 1) I found a bug in Maxima for Windows regarding the first lines of your code which produces the following error: > > Maxima encountered a Lisp error: > 3 is not a string designator. > > The bug report is #3272. > > 2) Valery Pipin helped me with another approach that is to snatch the 'levi_civita symbol right at the beginning and > substitute it by a function that is defined later on when the metric is set within the ctensor package. > > 3) I was able to test your code using CESGA - Maxima on line > <http://maxima.cesga.es/index.php?c=pl5yjoo4dklztijxhavuv&n=24> > > <http://maxima.cesga.es/index.php?c=pl5yjoo4dklztijxhavuv&n=24> > > CESGA - Maxima on line > > > > and it works fine except for the fact that the last line is not the curl of the vector H. > I think that we had this conversation before when someone was dealing with hodge. > > Thanks, > > Daniel Volinski > > > > > > |
From: Daniel V. <dan...@ya...> - 2017-01-11 17:02:14
|
Hi Viktor, Would you say that the best strategy in this case is to substitute the levi_civita symbol right at the beginning by some undefined function and later on define that function when the metric has been defined with ctensor? Thanks, Daniel Volinski El Miércoles 11 de enero de 2017 0:02, Viktor T. Toth <vt...@vt...> escribió: Daniel, That bug is really "by design", so to speak: a limitation of levi_civita, which does something the itensor package isn't designed to handle, replacing a symbolic (dummy) index with an explicit numerical value. That said, curiously the code actually works in GCL, though it fails in all other lisps to which I have access. Viktor > -----Original Message----- > From: Daniel Volinski [mailto:dan...@ya...] > Sent: Tuesday, January 10, 2017 9:28 AM > To: 'Maxima Discuss' <max...@li...> > Subject: Re: [Maxima-discuss] Curl in tensor notation > > > > Hi Viktor, > > Thank you for your response. > > 1) I found a bug in Maxima for Windows regarding the first lines of your code which produces the following error: > > Maxima encountered a Lisp error: > 3 is not a string designator. > > The bug report is #3272. > > 2) Valery Pipin helped me with another approach that is to snatch the 'levi_civita symbol right at the beginning and > substitute it by a function that is defined later on when the metric is set within the ctensor package. > > 3) I was able to test your code using CESGA - Maxima on line > <http://maxima.cesga.es/index.php?c=pl5yjoo4dklztijxhavuv&n=24> > > <http://maxima.cesga.es/index.php?c=pl5yjoo4dklztijxhavuv&n=24> > > CESGA - Maxima on line > > > > and it works fine except for the fact that the last line is not the curl of the vector H. > I think that we had this conversation before when someone was dealing with hodge. > > Thanks, > > Daniel Volinski > > > > > > |
From: Viktor T. T. <vt...@vt...> - 2017-01-12 03:51:37
|
Yes, Daniel, that is indeed a very viable strategy, one that I employed myself on a few occasions. Viktor > -----Original Message----- > From: Daniel Volinski [mailto:dan...@ya...] > Sent: Wednesday, January 11, 2017 11:58 AM > To: 'Maxima Discuss' <max...@li...> > Subject: Re: [Maxima-discuss] Curl in tensor notation > > Hi Viktor, > > Would you say that the best strategy in this case is to substitute the levi_civita symbol right at the beginning by > some undefined function and later on define that function when the metric has been defined with ctensor? > > Thanks, > > Daniel Volinski > > > > El Miércoles 11 de enero de 2017 0:02, Viktor T. Toth <vt...@vt...> escribió: > > > > Daniel, > > That bug is really "by design", so to speak: a limitation of levi_civita, which does something the itensor package > isn't designed to handle, replacing a symbolic (dummy) index with an explicit numerical value. > > That said, curiously the code actually works in GCL, though it fails in all other lisps to which I have access. > > > Viktor > > > > > > -----Original Message----- > > From: Daniel Volinski [mailto:dan...@ya... <mailto:dan...@ya...> ] > > Sent: Tuesday, January 10, 2017 9:28 AM > > To: 'Maxima Discuss' <max...@li... <mailto:max...@li...> > > > Subject: Re: [Maxima-discuss] Curl in tensor notation > > > > > > > > Hi Viktor, > > > > Thank you for your response. > > > > 1) I found a bug in Maxima for Windows regarding the first lines of your code which produces the following error: > > > > Maxima encountered a Lisp error: > > 3 is not a string designator. > > > > The bug report is #3272. > > > > 2) Valery Pipin helped me with another approach that is to snatch the 'levi_civita symbol right at the beginning > and > > substitute it by a function that is defined later on when the metric is set within the ctensor package. > > > > 3) I was able to test your code using CESGA - Maxima on line > > <http://maxima.cesga.es/index.php?c=pl5yjoo4dklztijxhavuv&n=24> > > > > <http://maxima.cesga.es/index.php?c=pl5yjoo4dklztijxhavuv&n=24> > > > > > CESGA - Maxima on line > > > > > > > > and it works fine except for the fact that the last line is not the curl of the vector H. > > I think that we had this conversation before when someone was dealing with hodge. > > > > Thanks, > > > > Daniel Volinski > > > > > > > > > > > > > > > |