1
$\begingroup$

I am in the process of trying to use Mathematica to find the coefficients of transmission and reflection for the wavefunction of a beam of particles tunnelling across any number of identical potential barriers.

I have got to the stage where I am able to generate the equations for n numbers of barriers, but I am not able to equate the appropriate equations to each other to solve them for the coefficients. The module I have written to generate the relavant equations is:

GeneratePhiEqs[N_] :=
    Module[{inceqs, diffinceqs, bareqs, diffbareqs, finaleq, difffinaleq,
    IncidentPhi, BarrierPhi, FinalTransmittedPhi, phieqs}, 
    IncidentPhi[x_, n_] := T[n]*E^(I*k*x) + R[n]*E^(-I*k*x);
    BarrierPhi[x_, n_] := A[n]*E^(r*x) + B[n]*E^(-r*x);
    FinalTransmittedPhi[x_, n_] := T[n]*E^(I*k*x); T[0] = 1;

    inceqs = Table[Phi[j] = {IncidentPhi[x, j]}, {j, 0, N}];
    diffinceqs = 
    Table[Phi[N + j] = {D[IncidentPhi[x, j], x]}, {j, 0, N - 1}];
    bareqs = Table[Phi[2 N + j - 1] = {BarrierPhi[x, j]}, {j, N}];
    diffbareqs = 
    Table[Phi[3 N + j - 1] = {D[BarrierPhi[x, j], x]}, {j, N}];
    finaleq = Phi[4 N] = FinalTransmittedPhi[ x, N];
    difffinaleq = Phi[4 N + 1] = D[FinalTransmittedPhi[x, N], x];
    phieqs = {inceqs, diffinceqs, bareqs, diffbareqs, finaleq, 
    difffinaleq}; phieqs]

The module I have written (or what it looks like so far) to equate the relevant equations is:

SolvePhiEqs[number_] := 
    Module[{eqs, simeqs}, eqs = GeneratePhiEqs[number]; 
    simeqs = Table[{Phi[j] /. phieqs == Phi[2*number + j] /. 
    phieqs}, {j, 2*number - 1}]; simeqs] 

Basically, the equations that I generated with the first module will not substitute in correctly, they just remain isolated from each other surrounded by curly brackets.

Any suggestions would be much appreciated!

Thanks.

$\endgroup$
6
  • $\begingroup$ I get a collection of stuff by running GeneratePhiEqs[1]. But these are not equations, yet, they have no equality. Can you explain what terms are to be equated with what? $\endgroup$
    – bill s
    Commented Jul 10, 2013 at 13:38
  • $\begingroup$ Can you write what explicitly A[n], B[n], R[n], T[n] are? Also, the usage Table[Phi[j] = {IncidentPhi[x, j]}, {j, 0, N}] looks so weird to me. I would rather go with something much simpler like Phi[j_]:={IncidentPhi[x, j]}. $\endgroup$
    – Ali
    Commented Jul 10, 2013 at 14:03
  • $\begingroup$ @bills GeneratePhiEqs effectively gives me a set of functions and their derivatives, some of which I know are equal to each other at a certain value of x. I need to equate particular functions to each other and substitute in a value for x in order to solve for the coefficients. The SolvePhiEqs module should take care of equating them to each other and solving them. $\endgroup$
    – CJO
    Commented Jul 10, 2013 at 14:23
  • $\begingroup$ In answer to your question @Ali the A[n], B[n] are coefficients. For example, if I run GeneratePhiEqs[1], I get a certain number of equations, and it numbers the coefficients for me. I am trying to solve for the coefficients. $\endgroup$
    – CJO
    Commented Jul 10, 2013 at 14:24
  • $\begingroup$ I get a function like $e^{ikx}$ at one barrier, but then I get the same thing at another barrier. I want to differentiate between them whilst still being able to see that they are related. $\endgroup$
    – CJO
    Commented Jul 10, 2013 at 14:29

1 Answer 1

3
$\begingroup$

Possibly this?

Clear[Phi];
GeneratePhiEqs[N_] := (* <- Don't use N *)
 Module[{inceqs, diffinceqs, bareqs, diffbareqs, finaleq, difffinaleq,
    IncidentPhi, BarrierPhi, FinalTransmittedPhi}, 
  IncidentPhi[x_, n_] := T[n]*E^(I*k*x) + R[n]*E^(-I*k*x);
  BarrierPhi[x_, n_] := A[n]*E^(r*x) + B[n]*E^(-r*x);
  FinalTransmittedPhi[x_, n_] := T[n]*E^(I*k*x); T[0] = 1;
  inceqs      = Table[Phi[j] -> IncidentPhi[x, j], {j, 0, N}];
  diffinceqs  = Table[Phi[N + j] -> D[IncidentPhi[x, j], x], {j, 0, N - 1}];
  bareqs      = Table[Phi[2 N + j - 1] -> BarrierPhi[x, j], {j, N}];
  diffbareqs  = Table[Phi[3 N + j - 1] -> D[BarrierPhi[x, j], x], {j, N}];
  finaleq     = Phi[4 N] -> FinalTransmittedPhi[x, N];
  difffinaleq = Phi[4 N + 1] -> D[FinalTransmittedPhi[x, N], x];
  Flatten @ {inceqs, diffinceqs, bareqs, diffbareqs, finaleq, difffinaleq}]

SolvePhiEqs[number_] := 
 Module[{eqs},
  eqs = GeneratePhiEqs[number];
  Table[Phi[j] == Phi[2*number + j] /. eqs, {j, 2*number - 1}]
 ]

If it's right, then you should

Example:

SolvePhiEqs[2]
(* {E^(-I k x) R[1] + E^(I k x) T[1] == E^(r x) A[2] + E^(-r x) B[2], 
    E^(-I k x) R[2] + E^(I k x) T[2] == E^(r x) r A[1] - E^(-r x) r B[1],
    -I E^(-I k x) k R[1] + I E^(I k x) k T[1] == E^(r x) r A[2] - E^(-r x) r B[2]} *)
$\endgroup$
1
  • $\begingroup$ It works - Thank you so much! Please excuse my cobbled together syntax, I'm a bit of a Mathematica beginner. $\endgroup$
    – CJO
    Commented Jul 10, 2013 at 16:01

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.