Scattering length of an attractive square-well potential
                                                                                            7.6.2000 B. Friman

α[k_, v_] := (Sqrt[2*(e[k] + v)])

e[k_] := k^2/2

Compute Wronskians for l=0 phase shift

w1[k_, v_, a_] := 
      k*Cos[k*a]*Sin[α[k, v]*a] - α[k, v]*Sin[k*a]*
          Cos[α[k, v]*a]

w2[k_, v_, a_] := 
      k*Sin[k*a]*Sin[α[k, v]*a] + α[k, v]*Cos[k*a]*
          Cos[α[k, v]*a]

Compute the phase shift naively

δ[k_, v_, a_] := ArcTan[w1[k, v, a]/w2[k, v, a]]

Plot[δ[k, 10, 2], {k, 0, 15}, 
      PlotRange -> {{0, 15}, {-2, 2}}]

The change of branch of ArcTan gives rise to the jumps.
A way to avoid this problem is to compute the phase at large momenta, where it is small, and iterate down to small momenta

kmom[i_] := N[(i - 1)/15]

vv = 10; aa = 2. ;

Phase shift at highest momentum

delta = 
      Arg[(w2[kmom[141], vv, aa] + I*
              w1[kmom[141], vv, aa])]

Change in phase shift as one moves from point i+1 to point i

delincr[i_, vv_, aa_] := 
      Arg[(w2[kmom[i], vv, aa] + I*w1[kmom[i], vv, aa])/(w2[
                kmom[i + 1], vv, aa] + I*
                w1[kmom[i + 1], vv, aa])]

Iterative calculation of the phase shift

phase[i_] := delta = delta  + delincr[i, vv, aa]

phtable = 
        Table[{kmom[i], phase[i]}, {i, 140, 1, -1}];

ListPlot[phtable, PlotRange -> {{0, 10}, {0, 10}}]

The value at k=0 is 3 π. Thus, the potential has 3 bound states