# Solve equation of motion with Dirac-fermions?

Posted 2 months ago
579 Views
|
12 Replies
|
1 Total Likes
|
 Dear Wolfram team:I am a beginner of Mathematica.My Problem is that I want solve a System with n equation of Motion in first order. In this equation of motion are creation and annihilation operators of Dirac-fermions. I don't know and don't find a it, how I can describes the creation and annihilation operators of Dirac-fermions in Mathematica. The equation of motion have the form: $$\dot{c}_i^\dagger [t]=f*c_i^\dagger[t]+g[t]*c_{i+1}[t]-g[t]*c_{i+1}^\dagger[t]+h[t]*c_{i-1}[t]-h[t]*c_{i+1}^\dagger[t]\\ \dot{c}_i [t]=f*c_i[t]+g[t]*c_{i+1}^\dagger[t]-g[t]*c_{i+1}[t]-h[t]*c_{i-1}^\dagger[t]+h[t]*c_{i+1}[t],$$where $c_i^\dagger,c_i$ are creation and annihilation operators and f,g,h are functions.Then I want use DSolve or NDSolve to solve the equation of motion. Thanks, for your help.
12 Replies
Sort By:
Posted 2 months ago
 The FeynCalc package may help.https://feyncalc.github.io/FeynCalcBook/guide/FeynCalc.html
Posted 2 months ago
 Thanks for your answer. I found a lot of nice commands but Mathematica don't work. My code is for the equation are eoMkitaevchaindiracnormal = { -I c1'[t] == -\[Omega] (c2[t] + ConjugateTranspose[c2[t]]) + Exp[-((x1 - v*t)^2/(2 \[Sigma]^2))]/ 2 c1[t], -I ConjugateTranspose[c1]'[ t] == \[Omega] (c2[t] + ConjugateTranspose[c2[t]]) + Exp[-((x1 - v*t)^2/(2 \[Sigma]^2))]/ 2 ConjugateTranspose[c1[t]], -I c2'[ t] == -\[Omega] (c1[t] - ConjugateTranspose[c1[t]]) + Exp[-((x2 - v*t)^2/(2 \[Sigma]^2))]/ 2 c2[t], -I ConjugateTranspose[c2]'[ t] == \[Omega] (-c1[t] + ConjugateTranspose[c1[t]]) + Exp[-((x2 - v*t)^2/(2 \[Sigma]^2))]/2 ConjugateTranspose[c2[t]], , c1[-100] == 1, c2[-100] == 1, ConjugateTranspose[c1][-100] == 1, ConjugateTranspose[c2][-100] == 1}; Now, I want use NDSolve to solve the equation and use this code solutionkitaevchaindiracnormal = NDSolve[eoMkitaevchaindiracnormal /. {v0 -> 1, v -> 1, \[Sigma] -> 0.75, x1 -> 1, x2 -> 2, \[Omega] -> 1}, {c1[ t], c2[t], ConjugateTranspose[c1][t], ConjugateTranspose[c2][t]}, {t, -100, 100}]; After the second step Mathematics give me this Output: NDSolve::deqn: Equation or list of equations expected instead of Null in the first argument {-I c1^\[Prime](t)==1/2 E^(-0.888889 Plus[<<2>>]^2) c1(t)-c2(t)-c2(t)^\[ConjugateTranspose],-I ((c1^\[ConjugateTranspose])^\[Prime])\[InvisibleApplication](t)==c2(t)+1/2 E^(-0.888889 Plus[<<2>>]^2) c1(t)^\[ConjugateTranspose]+c2(t)^\[ConjugateTranspose],-I c2^\[Prime](t)==-c1(t)+1/2 E^(-0.888889 Plus[<<2>>]^2) c2(t)+c1(t)^\[ConjugateTranspose],-I ((c2^\[ConjugateTranspose])^\[Prime])\[InvisibleApplication](t)==-c1(t)+c1(t)^\[ConjugateTranspose]+1/2 E^(-0.888889 Plus[<<2>>]^2) c2(t)^\[ConjugateTranspose],Null,c1(-100)==1,c2(-100)==1,(c1^\[ConjugateTranspose])\[InvisibleApplication](-100)==1,(c2^\[ConjugateTranspose])\[InvisibleApplication](-100)==1}. >> What is my mistake? Can I use ConjugateTranspose too describe a creation operator or I have use an other command ?
Posted 2 months ago
 Umm, did you ever actually look at eoMkitaevchaindiracnormal? Because there is a Null in it, and the error message is telling you that, and you should find and remove it before appealing to an entire forum.Next there is the business of double naming. Some places use ConjugateTranspose[c1[t]], others use ConjugateTranspose[c1][t]. I will propose getting rid of both and just using ctc1[t].NDSolve will handle the variant below. eoMkitaevchaindiracnormal = {-I c1'[ t] == -\[Omega] (c2[t] + ctc2[t]) + Exp[-((x1 - v*t)^2/(2 \[Sigma]^2))]/2 c1[t], -I ctc1'[ t] == \[Omega] (c2[t] + ctc2[t]) + Exp[-((x1 - v*t)^2/(2 \[Sigma]^2))]/2 ctc1[t], -I c2'[ t] == -\[Omega] (c1[t] - ctc1[t]) + Exp[-((x2 - v*t)^2/(2 \[Sigma]^2))]/2 c2[t], -I ctc2'[ t] == \[Omega] (-c1[t] + ctc1[t]) + Exp[-((x2 - v*t)^2/(2 \[Sigma]^2))]/2 ctc2[t], c1[-100] == 1, c2[-100] == 1, ctc1[-100] == 1, ctc2[-100] == 1} Bottom line: this could and should have been debugged off-line. If the result from NDSolve is not satisfactory, that's the sort of issue to bring to the forum.
Posted 2 months ago
 Thank,s for the hint with double naming. That was the problemIt is irrelevant to use ctc1[t] or ConjugateTranspose[c1][t].When I solve and Plot the System with this code : eoMkitaevchaindiracnormal = { -I c1'[t] == -\[Omega] (c2[t] + ConjugateTranspose[c2][t]) + Exp[-((x1 - v*t)^2/(2 \[Sigma]^2))]/ 2 c1[t], -I ConjugateTranspose[c1]'[ t] == \[Omega] (c2[t] + ConjugateTranspose[c2][t]) + Exp[-((x1 - v*t)^2/(2 \[Sigma]^2))]/ 2 ConjugateTranspose[c1][t], -I c2'[ t] == -\[Omega] (c1[t] - ConjugateTranspose[c1][t]) + Exp[-((x2 - v*t)^2/(2 \[Sigma]^2))]/ 2 c2[t], -I ConjugateTranspose[c2]'[ t] == \[Omega] (-c1[t] + ConjugateTranspose[c1][t]) + Exp[-((x2 - v*t)^2/(2 \[Sigma]^2))]/2 ConjugateTranspose[c2][t], c1[-100] == 1, c2[-100] == 1, ConjugateTranspose[c1][-100] == 1, ConjugateTranspose[c2][-100] == 1}; solutionkitaevchaindiracnormal = NDSolve[eoMkitaevchaindiracnormal /. {v0 -> 1, v -> 1, \[Sigma] -> 0.75, x1 -> 1, x2 -> 2, \[Omega] -> 1}, {c1[ t], c2[t], ConjugateTranspose[c1][t], ConjugateTranspose[c2][t]}, {t, -100, 100}] ParametricPlot[{Re[#], Im[#]} &@ Evaluate[{ConjugateTranspose[c1][t]*c1[t]} /. solutionkitaevchaindiracnormal], {t, -10, 10}, PlotRange -> All, PlotStyle -> {Blue,Thick}] I become this solution for $c_1^\dagger c_1$But $c_1^\dagger c_1$ is a density Operator. I think that the imaginary part of $c_1^\dagger c_1$ should be 0 and ConjugateTranspose[c1] and ctc1 is not $c_1^\dagger$. But what is the correct description of $c_1^\dagger$?
Posted 2 months ago
 Thanks you for the hint with the double Name. Now Mathematica solve the System. I use ctc1[t} und ConjugateTranspose[c1[t]] and both have the same solution. But I think that Mathematica don't know that ConjugateTranspose[c1[t]] should be $c_1^\dagger$ and >Mathematica interpret c1 and ConjugateTranspose[c1[t]] as operator because if I plot $c_1^\dagger c_1$ then I getwhere horizontal axis is the real part and vertical axis is the imaginary part. I think that ConjugateTranspose[c1[t]] is not $c_1^\dagger$. I think that the imaginary part are Zero. I use this Code: eoMkitaevchaindiracnormal = { -I c1'[t] == -\[Omega] (c2[t] + ctc2[t]) + v0 Exp[-((x1 - v*t)^2/(2 \[Sigma]^2))]/2 c1[t], -I ctc1'[t] == \[Omega] (c2[t] + ctc2[t]) + v0 Exp[-((x1 - v*t)^2/(2 \[Sigma]^2))]/2 ctc1[t], -I c2'[t] == -\[Omega] (c1[t] - ctc1[t]) + v0 Exp[-((x2 - v*t)^2/(2 \[Sigma]^2))]/2 c2[t], -I ctc2'[t] == \[Omega] (-c1[t] + ctc1[t]) + v0 Exp[-((x2 - v*t)^2/(2 \[Sigma]^2))]/2 ctc2[t], c1[-0] == 1, c2[-0] == 1, ctc1[-0] == 1, ctc2[-0] == 1}; solutionkitaevchaindiracnormal = NDSolve[eoMkitaevchaindiracnormal /. {v0 -> 0.1, v -> 0, \[Sigma] -> 0.75, x1 -> 1, x2 -> 2, \[Omega] -> 1}, {c1[ t], c2[t], ctc1[t], ctc2[t]}, {t, -0, 10}]; ParametricPlot[{Re[#], Im[#]} &@ Evaluate[{ctc1[t]*c1[t]} /. solutionkitaevchaindiracnormal], {t, -0, 10}, PlotRange -> All, PlotStyle -> {{Blue, Thick}, {Black, Dashed, Thick}}] When I include ctc1 = ComplexConjugate[c1]; ctc2 = ComplexConjugate[c2]; then dont solve Mathematica the System.What is my mistake and who can I describe the operators in Mathematica korrect?
Posted 2 months ago
 I think I have it. My solution: eoMkitaevchaindiracnormal = { -I c1'[t] == -\[Omega] (c2[t] + ctc2[t]) - v0 Exp[-((x1 - v*t)^2/(2 \[Sigma]^2))]/2 c1[t], -I ctc1'[t] == \[Omega] (c2[t] + ctc2[t]) + v0 Exp[-((x1 - v*t)^2/(2 \[Sigma]^2))]/2 ctc1[t], -I c2'[t] == -\[Omega] (c1[t] - ctc1[t]) - v0 Exp[-((x2 - v*t)^2/(2 \[Sigma]^2))]/2 c2[t], -I ctc2'[t] == \[Omega] (-c1[t] + ctc1[t]) + v0 Exp[-((x2 - v*t)^2/(2 \[Sigma]^2))]/2 ctc2[t], c1[-10] == 1, c2[-10] == 1, ctc1[-10] == 1, ctc2[-10] == 1}; solutionkitaevchaindiracnormal = NDSolve[eoMkitaevchaindiracnormal /. {v0 -> 1, v -> 1, \[Sigma] -> 0.75, x1 -> 1, x2 -> 2, \[Omega] -> 1}, {c1[ t], c2[t], ctc1[t], ctc2[t]}, {t, -100, 100}]; Plot[Evaluate[{ctc1[t]*c1[t], ctc2[t]*c2[t], ctc1[t]*c1[t] + ctc2[t]*c2[t]} /. solutionkitaevchaindiracnormal], {t, -10, 10}, PlotRange -> All, PlotStyle -> {{Blue, Thick}, {Orange, Thick}, {Black, Dashed, Thick}}] Thanks for your help. Now, I have another problem. I want that the absolute value of $c_1, c_1^\dagger,c_2, c_2^\dagger$ is between 0 and 1. Assumptions and Assuming don't work with NDSolve. I know why I can't use Assumptions and Assuming. But I don't know what is the command for NDSolve?
Posted 2 months ago
 NDSolve is for functions, not operators.
Posted 2 months ago
 Ok, which command can I use for operators?
Posted 2 months ago
 I don't think there is one in Mathematica. That's why I suggested FeynCalc.
Posted 2 months ago
 Ok, I have the package FeynCalc installed but I don't find in the list of FeynCalc an other command which solve my system of equations.