Calculation of currents around a triangular indentation by the hodograph method
Jonathan I. Avila Benoit Vanderheyden Alejandro V. Silhanek Sorin Melinte






^1 ICTEAM, Université catholique de Louvain, B-1348, Louvain-La-Neuve, Belgium
^2 SUPRATECS and Department of Electrical Engineering and Computer Science, Université de Liège, B-4000 Sart Tilman, Belgium
^3 Experimental Physics of Nanostructured Materials, Q-MAT, CESAM, Université de Liège, B-4000 Sart Tilman, Belgium





Border indentations in non-linear conductors, such as superconducting thin films in the creep regime, alter the distribution of currents and magnetic fields near and far from the indentation. One of such disturbances are the discontinuity lines, or d-lines, a parabolic-like line originating from the indentation where the current density direction changes abruptly. Hodograph series results are obtained for the currents around a triangular indentation and its corresponding d-lines in a conducting stripe of finite width and in an infinite half plane, considering two cases: uniform creep exponent and mixed infinite and ohmic exponents. The mixed creep exponent case presents currents distributions resembling the purely ohmic case, with significant current disturbances only near the indentation. For uniform creep exponent, results similar to a planar indentation are obtained, with far ranged currents features  and  parabolic-like d-lines with shapes depending on the creep exponent. In particular, the same d-line asymptotic behaviour is obtained for the triangle indentation as that of the planar defect in the critical state, a result obtained here just on continuity considerations of the hodograph expansions. This equivalence is due to identical contributions to the Fourier series of the current stream-function in the hodograph space, obtained from an images method expansion.


...


Calculation of currents around a triangular indentation by the hodograph method
    J. I. Avila^1,  B. Venderheyden^2, A. V. Silhanek^3, S. Melinte^1
    December 30, 2023
===============================================================================





§ INTRODUCTION

Defects have significant effects on the electrical properties of type II superconductor. While nanoscopic defects serve as pinning centres that arrest the motion of flux vortices and avoid dissipation <cit.>, larger defects create zones of increased resistivity that can serve as nucleation spots for flux avalanches <cit.>. Border defects such as cracks or voids are accompanied by long range perturbations in the current density and magnetic field similar to the Bean model <cit.> discontinuity lines <cit.>, where the current density J turns sharply. The exact shape of the lines will be directly related to the critical current and constitutive relationship under external applied fields, such as an external electric field, or an external magnetic field ramp creating induction currents. The resistivity in the flux creep regime <cit.> can be modelled as: 

    ρ(J,B,T)=ρ_0 (J/J_c(B,T))^n-1,

where ρ_0 is the resistivity at J=J_c,  B the magnetic field, T the temperature and n is the creep exponent. The steep variations of the resistivity in this regime are also espected to generate thermo-magnetic instabilities <cit.>. By taking no magnetic field or temperature dependence of J_c, expression  (<ref>) has been used by means of the hodograph method <cit.> to obtain analytical expressions for the electric currents around a planar defect, for finite and infinite creep exponent n. A more basic model consists of considering only the case n→∞ and to disregard non-local effects and magnetic field or temperature depence of the critical current. This leads to the Bean model <cit.>, in which the magnetization current density of a sample exposed to a ramping magnetic field is modelled as piece-wise uniform regions, on whose boundaries the current can change of direction discontinuously, tracing the d-lines. From this, parabolic d-lines of the form x = y^2/2h+h/2 are obtained for an indentation of height h on a border located at x=0 <cit.>. Here, the calculations of currents around a planar defect by the hodograph method <cit.> are extended to a triangular defect. In that formalism the calculated currents actually do not change direction in a discontinuous way, however, the d-line concept is extended to a current domain wall that traces a line similar to that of the Bean's model, so the d-line term is also used here. The current stream function is obtained for different creep exponents and geometries and finally compared to the Bean model. 

  


§ HODOGRAPH METHOD

A general stationary solution method was developed in refs <cit.> for currents in films with a non-linear relationship between the electric field E and the current J of the type E=E_0 (J/J_0)^n, where E_0 is a characteristic electric field magnitude and J_0 the corresponding current density and n≥ 1. In this method the Maxwell equations for the current stream function ψ are solved in the space of the current density magnitude J=|∇×ψ| and the angle θ subtended between J and the x axis, i.e. under the transformation ψ (x,y)→ψ (J,θ). It can be equally applied for the electric potential with ϕ (x,y) →ϕ (E,θ), where E=-∇ϕ. In the space (J,θ), namely the hodograph space, the electromagnetic equations become linear and separable when using the power law resistivity, and can be solved by eigenfunction expansions, with general solution:


    ψ(J,θ)/I= A_0+(B_0+D_0 θ)(J/J_0)^(1-n)+C_0 θ
     + ∑_m=1^∞ C_m( J/J_0)^τ_m^+sin (mθ+ϕ_m) 
     + ∑_m=1^∞ D_m( J/J_0)^τ_m^-sin (mθ+ω_m),


with

    τ_m^±=1/2(1-n±√((n-1)^2+4nm^2)).

The boundary conditions for ψ in the hodograph space are simple only in the case of samples with straight borders, in which θ is constant. Once the solution is found in the hodograph plane, it can be transformed back to the x, y space, i.e. (J,θ) → (x,y). By using the complex variable z=x+iy, this transformation is given by:


    z/l=F_0+C_0 e^iθJ_0/J+
    +e^iθ/n(J/J_0)^-n{i B_0(n-1)+D_0 (1+i(n-1)θ}+
    
    -e^iθ∑_m=1^∞D_m/τ_m^- -1( J/J_0)^τ_m^- -1[mcos (mθ+ω_m) .
    . -i τ_m^-sin(mθ+ω_m)]+
    
    - e^iθ∑_m=1^∞C_m/τ_m^+ -1( J/J_0)^τ_m^+ -1[mcos (mθ+ϕ_m) .
    . -i τ_m^+sin(mθ+ϕ_m)],


Additional details are given in the appendix. The technique has been applied to geometries like corners, bridges, current leads and planar defects. In ref. <cit.>  a planar defect of width 2a perpendicular to the current in an infinite superconductor was considered. In that work, the d-line was identified as the core of a current domain wall starting near the defect. In the limit n→∞ it was found there that, asymptotically, the d-line  traces a curve with x≈ y^2 /1.94a, similar but not identical to the d-line in the Bean result that goes as x= y^2 /2a.




§ STRIPE OF FINITE WIDTH

Equations (<ref>) and (<ref>) are solved for a strip carrying a total current I, of width l  and with a triangular constriction of base 2a and height a, as depicted in Fig.<ref>(a). The solution for the stream function ψ(J,θ) is separated in 3 hodograph regions, as shown in Fig.<ref>(b), namely ψ_1 in region (1) for 0≤ J≤ J_0, with J_0 the magnitude of J far from the defect,  ψ_2 in region (2) for J_0≤  J≤ J_m, where J_m is the magnitude of the current density in the border directly opposite to the defect, and ψ_3 in region (3) for J_m≤ J. The solutions fulfilling the boundary conditions described in Fig.<ref>(b) are:


    ψ_1/I=1+∑_m=1^∞ C_4m^(1)( J/J_0)^τ_4m^+sin (4mθ),





    ψ_2/I=-2+4θ/π+∑_m=1^∞[ D_4m^(2)( J/J_0)^τ_4m^- +.
    . +C_4m^(2)( J/J_m)^τ_4m^+] sin (4mθ),





    ψ_3/I=1+∑_m=1^∞  D_4m-2^(3)( J/J_m)^τ_4m-2^-sin ((4m-2)θ+3π/2),


The C_4m^(1), D_4m^(2), C_4m^(2) and D_4m-2^(3) coefficients are obtained by enforcing continuity of ψ and ∂ψ /∂ J at J=J_0 and J=J_m along with Fourier analysis for θ∈[π/2,3π/4]. 
Applying continuity at J_0:


    8/π∫_π/2^3π/4ψ_1(J_0)sin(4mθ)dθ=8/π∫_π/2^3π/4ψ_2(J_0)sin(4mθ)dθ


It gets:

    C_4m^(1)=-2/mπ+D_4m^(2)+r^τ_4m^+C_4m^(2)

Where r=J_0/J_m. From the derivative continuity at J_0, 
it is obtained:

    C_4m^(1)=τ_4m^-/τ_4m^+ D_4m^(2)+r^τ_4m^+C_4m^(2)

From continuity at J_m,
by using:


    8/π∫_π/2^3π/4sin((4k-2)θ+3π/2)sin(4mθ)dθ=8/π4m/(4m)^2-(4k-2)^2


It gets:

    C_4m^(2)=2/mπ-D_4m^(2)r^-τ_4m^-+8/π∑_k=1^∞4mD_4k-2^(3)/(4m)^2-(4k-2)^2

And from continuity of the derivative at J_m:

    τ_4m^- D_4m^(2)r^-τ_4m^-+τ_4m^+ C_4m^(2)=8/π∑_k=1^∞4mD_4k-2^(3)τ_4k-2^-  /(4m)^2-(4k-2)^2

Combining (<ref>) and (<ref>):


    D_4m^(2)=2/mπτ_4m^+/τ_4m^+ - τ_4m^-

And combining (<ref>), (<ref>) and (<ref>) the linear system of equations for D_4m-2^(3)  is obtained:

    τ_4m^+/(4m)^2(r^-τ_4m^--1)=∑_k=1^∞D_4k-2^(3)(τ_4m^+-τ_4k-2^-)/(4m)^2-(4k-2)^2


 The space coordinates (x, y) associated to ψ_i(J, θ) can be obtained from the solution (<ref>) for the complex variable z=x+iy in the different regions:



    z_1/l=-1+ia/l-e^iθ∑_m=1^∞C_4m^(1)/τ_4m^+-1( J/J_0)^τ_4m^+-1[4mcos (4mθ) .
    .-i τ_4m^+sin(4mθ)],




    z_2/l=-i(l-a)/l+4/πJ_0/Je^iθ+
     -e^iθ∑_m=1^∞D_4m^(2)/τ_4m^- -1( J/J_0)^τ_4m^- -1[4mcos (4mθ) .
    .-i τ_4m^-sin(4mθ)]
    -r e^iθ∑_m=1^∞C_4m^(2)/τ_4m^+ -1( J/J_m)^τ_4m^+ -1[4mcos (4mθ).
    . -i τ_4m^+sin(4mθ)],




    z_3/l=-(l-a)/l+
    -r e^iθ∑_m=1^∞D_4m-2^(3)/τ_4m-2^- -1( J/J_m)^τ_4m-2^- -1[ .
    .(4m-2)cos ((4m-2)θ+3π/2).
    . -i τ_4m-2^-sin((4m-2)θ+3π/2) ..],


Where the integration constants have been chosen to guarantee space continuity of ψ. 
The parameters are all dependent on D_4k-2^(3), and those depend at the same time on r, that is found numerically by enforcing space continuity of the solutions.
In the limit n→∞ the coefficients τ_m^+ → m^2 and τ_m^- → -∞. This greatly simplify the stream function calculation. In particular the hodograph region (3) disappears and r=(l-a)/l due to current conservation. Then C_4m^(2)=2/mπ from (<ref>) and C_4m^(1)=2(r^(4m)^2-1)/mπ from (<ref>). 

The current streamlines, or the same the contour or equipotential lines of the stream function, can be seen for different n values in FIG. <ref>. Those noticably change direction around the triangle constriction much far from the constriction than in the ohmic case. The non-local magnetic field H_d given by Biot-Savart's law has a well behaved local minimum around the same area of current direction change, and as in <cit.>, this minimum was used here as the d-lines, that are  displayed in FIG. <ref>. This approach is practical as in magneto-optical experiments it is the magnetic field that is measured rather than the current distribution. The field H_d was calculated by a Fourier transform method <cit.>, namely, H_d=ℱ^-1[k ℱ(ψ)/2], where ℱ is the 2D Fourier transform, and k the wave-vectors associated to the domain discretization.  As it can be seen in FIG. <ref>, the d-lines are always narrower than the one corresponding to the Bean case. This is expected as in a stripe it is impossible to keep constant current density with constrictions present. Furthermore,  if J> J_c, as expected for hodograph region (3), a flux flow regime must be considered near the indentation, where for finite n the current density J diverges. 




 §.§ Flux Flow

The previous hodograph calculation assumes that the power law resistivity is valid even for currents above the current J_0, as in refs. <cit.>. If we consider the flux flow resistance <cit.>, then a better representation of the resistivity is given by <cit.>:

    ρ={[ ρ_0 (J/J_c)^n-1,           J≤ J_c;            ρ_0 ,          J> J_c,;                  ].
 
Assuming that J_0=J_c for this case, here it is combined the creep regime for regions with J<J_0, and ohmic regions for J>J_0 using n=1. Lets obtain the relationship for the series expansions coefficients for region (1) with arbitrary creep exponent ñ with corresponding τ_m^+=1/2(1-ñ+ √((ñ-1)^2+4ñm^2)), along with region (2) and region (3) both with creep exponent n and τ_m^±=1/2(1-n±√((n-1)^2+4nm^2)).
The continuity relationships  (<ref>), (<ref>) and (<ref>) stay unmodified. Then by enforcing z_1(J_0)=z_2(J_0) it gets:


    C_4m^(1)τ_4m^+/τ_4m^+-1=D_4m^(2)τ_4m^-/τ_4m^--1+
    +C_4m^(2)r^τ_4m^+τ_4m^+/τ_4m^+-1-8/π4m/(4m)^2-1

That reduces to (<ref>) when ñ=n. And by combining eqns. (<ref>), (<ref>), (<ref>) and (<ref>) an equation can be obtained to calculate  D_4k-2^(3). Here is taken the particular case with ñ→∞ and n=1. With that choice, the relationships for the coefficients are finally:


    C_4m^(1)=-2/mπ+D_4m^(2)+r^4mC_4m^(2)



    D_4m^(2)=r^4mC_4m^(2)



    C_4m^(2)=2/mπ-D_4m^(2)r^4m+8/π∑_k=1^∞4mD_4k-2^(3)/(4m)^2-(4k-2)^2



    C_4m^(2)=D_4m^(2)r^4m-8/π∑_k=1^∞D_4k-2^(3) (4k-2) /(4m)^2-(4k-2)^2


Combing the above equations  (<ref>), (<ref>) and (<ref>) it is obtained:



    r^8m-1/4m=∑_k=1^∞D_4k-2^(3)(4m(1-r^8m)+(4k-2)(1+r^8m))/(4m)^2-(4k-2)^2,


that allows calculating D_4m-2^(3) and then D_4m^(2), C_4m^(2) and C_4m^(1). The r parameter is again found numerically by enforcing continuity of the different solutions, and the current streamlines and magnetic field can be seen in FIG.<ref>(d). The d-line for this case can also be found in FIG.<ref>.



§ DEFECT IN INFINITE HALF PLANE

A stream function for the infinite half plane x≥ 0 with a triangular defect that satisfies the boundary conditions is given by:


    ψ_1/J_0=∑_m=1^∞ C_4m^(1)( J/J_0)^τ_4m^+sin (4mθ),


For J<J_0. And for J≥ J_0:


    ψ_2/J_0=∑_m=1^∞ D_4m-2^(2)( J/J_0)^τ_4m-2^-sin ((4m-2)θ+3π/2),



Notice the different units of distance of the C and D coefficients now. The z functions are given by:


    z_1=-1+ia-e^iθ∑_m=1^∞C_4m^(1)/τ_4m^+-1( J/J_0)^τ_4m^+-1[4mcos (4mθ) .
    .-i τ_4m^+sin(4mθ)],




    z_2=-(1-a)+
     -e^iθ∑_m=1^∞D_4m-2^(2)/τ_4m-2^- -1( J/J_0)^τ_4m-2^- -1[(4m-2)cos ((4m-2)θ+3π/2) .
    .-i τ_4m-2^-sin((4m-2)θ+3π/2)],


Applying continuity conditions for ψ and the terms involving the sine functions of z it gets:


    D_4m-2^(2)=8/π∑_k=1^∞4kC_4k^(1)/(4k)^2-(4m-2)^2



    τ_4m-2^-/τ_4m-2^–1 D_4m-2^(2)=8/πa/(4m-1)(4m-3)+
    +8/π∑_k=1^∞4kC_4k^(1)/(4k)^2-(4m-2)^2τ_4k^+/τ_4k^+-1

Combing eqns. (<ref>) and (<ref>):

    ∑_k=1^∞4kC_4k^(1)/(4k)^2-(4m-2)^2(τ_4m-2^-/τ_4m-2^--1-τ_4k^+/τ_4k^+-1)=
    a/(4m-1)(4m-3)

When n→∞ this expression simplifies to:

    ∑_k=1^∞4kC_4k^(1)/(4k)^2-(4m-2)^2(-1/(4k)^2-1)=
    a/(4m-1)(4m-3)

That can be better inverted as:


    C_4k^(1) =-((4k)^2-1)×
    64a/π^2∑_m=1^∞4k/(4k)^2-(4m-2)^21/(4m-1)(4m-3)
    
    =-32a/πk,


by using the fact that 


    a_k=(8/π)^2 ∑_m,l=1^∞4k/(4k)^2-(4m-2)^24l/(4l)^2-(4m-2)^2a_l,


obtained by applying eq. (<ref>) to a expression of the kind ∑_k=1^∞sin(4kθ)a_k=∑_m=1^∞sin((4m-2)θ+3π/2)b_m, and the sum in (<ref>) is calculated explicitly. The current streamlines and magnetic field can be seen in FIG.<ref>(a).


 §.§ Flux flow

Again the case with creep exponent ñ for J<J_0 and n for J≥ J_0 is investigated, when J_0=J_c. Continuity conditions now produce:


    D_4m-2^(2)=8/π∑_k=1^∞4kC_4k^(1)/(4k)^2-(4m-2)^2



    τ_4m-2^-/τ_4m-2^–1 D_4m-2^(2)=8/πa/(4m-1)(4m-3)+
    +8/π∑_k=1^∞4kC_4k^(1)/(4k)^2-(4m-2)^2τ_4k^+/τ_4k^+-1

And the equation for C_4k^(1):

    ∑_k=1^∞4kC_4k^(1)/(4k)^2-(4m-2)^2(τ_4m-2^-/τ_4m-2^--1-τ_4k^+/τ_4k^+-1)=
    a/(4m-1)(4m-3)

And again for the hybrid case, a very high ñ is considered for which the limit ñ→∞ is taken in the above expression, and n=1, which gives:

    ∑_k=1^∞4kC_4k^(1)/(4k)^2-(4m-2)^2(4m-2/4m-1-(4k)^2/(4k)^2-1)=
    a/(4m-1)(4m-3)

The current streamlines and magnetic fields and d-line can be found in FIG.<ref>(b) and FIG.<ref> respectively.


§ DISCUSSION

The triangular defect here studied shares several properties with the planar defect analyzed in <cit.>. In particular, the current density diverges on the tip of the triangle for finite n due to the existence of the hodograph region (3), opposed to the critical state n→∞ in which the current density remains finite in that spot. For the infinite half plane it can be verified that J=J_0 at the triangle tip, by evaluating z_1(θ=3π/4-δ, J=J_0) from equation (<ref>) with the aid of expression (<ref>) for C_4k^(1):


    z_1=-1+ia-e^i(3π/4-δ)8a/π∑_m=1^∞4m^2/4m^2-1(-1)^m cos(4mδ)+
    
     +ie^i(3π/4-δ)8a/π∑_m=1^∞4m^3/4m^2-1(-1)^m sin(4mδ),


And by taking the limit δ→0 it is obtained z_1=-1+a, probing that J=J_0 there. In the case of the critical state in the stripe, it can be equally calculated that J=J_m at the tip of the triangle, the same as on the opposite side of the stripe. In fact J=J_m on all the line connecting both points, since J=J_m on both sides of the constriction, and all the current is passing through a zone where the maximum J is J_m. It also follows that r=(l-a)/l.
 
 In the hybrid case for the stripe, most of the current significantly above J_0 passes near the triangle tip, and the current density J is very close to J_m in most of the constriction. Additionally, J_m is also very close to J_0,  as it can be seen from the r=0.99999301 value. 
 
In the case of the infinite half plane with hybrid creep exponents, it is seen that perturbations on the currents are only important near the indentation, similar to the purely ohmic case, virtually disappearing at distances larger than ∼ 10a, as it can be seen in the current streamlines in FIG.<ref>(d)  and their d-lines in FIG.<ref>(b). The d-lines for the stripe and infinite plane also almost identical in that case. 
 
The d-line in the critical state for the infinite half plane traces a curve quite close to the Bean parabola, but identical actually to the result obtained in <cit.> of x≈ y^2/1.94a for the planar defect. In the appendix it is shown that both the planar defect and triangular indentation have actually the same asymptotic d-line. This is done by replacing the expression (<ref>) for  C_4k^(1), obtained here just from continuity conditions, opposed to the calculation for the planar defect <cit.> where C_4k^(1) is obtained from a expansion of current vortices and anti-vortices. The same expansion procedure is also done in the appendix to obtain approximate explicit expressions for the triangle indentation C_4k^(1) and D_4m-2^(2) coefficients, that have similar functional form<cit.> and limit for high n as the planar defect. A key point of the derivation is that the behaviour of the currents is mainly determined by a vortex solution to a London equation in a transformed hodograph plane of coordinates η=(1/√(n))log(E/E_0) , θ. That vortex is centred in η=0, θ=π/2, and is the main contributor to the Fourier expansions of ψ in both the planar defect and in the triangle indentation. Here both defects share only the height a, giving the same far current behavior, a feature already expected in the Bean model<cit.>, where the asymptotic d-line depends only on the height of the defect.




§ CONCLUSION

The hodograph series method is used to solve the problem of currents circulating in a straight stripe of finite width and an infinite half plane with a triangular border indentation, showing the non-local magnetic fields produced by those currents. Both cases present similar properties as the planar defect already discussed in <cit.>, with diverging currents in the defect tip for finite n, and bounded current for the critical state. An hybrid case of flux-flow and creep was also considered, showing the presence of d-lines, but of short range. For the triangular defect in the infinite half plane, the same asymptotic d-line is obtained as that of the planar defect in the critical state, which is explained here by the similar contribution of vortex currents in a modified hodograph space η, θ. 



§ ACKNOWLEDGEMENTS
 Grants ARC  13/18-08 for Concerted Research Actions, financed by the Wallonia-Brussels Federation, and Batwal D8.31404.060-P. J. I. A. also acknowledges FNRS support reference 2016/V 3/5/013 (W) -IB/LS-6967. 


*


§ ASYMPTOTICS OF CURRENTS FOR INFINITE HALF PLANE

As per <cit.>, eq. (<ref>) is the solution of the equation appearing from the hodograph transformation:

    J^2/n∂ ^2ψ/∂ J^2+J∂ψ/∂ J+∂ ^2ψ/∂θ^2=0

Whereas the transforming expressions for space coordinates (<ref>) are obtained from the equations:

    ∂ z/∂ E=-e^iθ(Eσ/J^2∂ψ/∂θ-iρ∂ψ/∂ E)
    ∂ z/∂θ=-e^iθ(-E^2/J∂ψ/∂ E-iρ∂ψ/∂θ)


Where σ=∂ J/∂ E. The following change of variable can be done for (<ref>):

    η   =(1/√(n))log(E/E_0)=√(n)log(J/J_0)
    ψ   =e^-βηh

With β=n-1/2√(n), and eq.(<ref>) transforms into a London equation:

    ∂ ^2h/∂η^2+∂ ^2 h/∂θ^2-β^2 h=0

Considering η and θ as cartesian coordinates, a change of variable to polar coordinates r, χ is performed, with η=rcosχ, θ=θ_0+sinχ. A vortex solution is obtained in these polar coordinates:

    h=√(2)h_0e^-βη/√(r)cosχ/2

The solution for ψ can be constructed by adding different vortices of the same strength on the line η=0 and alternating signs such that ψ=0 at θ=π/4, 3π/4,.... This is done by adding terms of the form (<ref>) with sign (-1)^m+1 centred at points θ_m=mπ/2, and converting back to the η, θ coordinates:

    ψ=h_0∑_m (-1)^m+1√(r_m+η)/r_me^-β(r_m+η)

Where r_m=√(η^2+(θ-mπ/2)). The coefficients  C_4k^(1) can be obtained then by calculating the Fourier terms:

    C_4k^(1)=8/π∫_π/2^3π/4ψ(J_0,θ)sin(4kθ)dθ

The main contributor to that integral is the term with m=1 of (<ref>), so the other m terms are discarded. For n ≫ 1, ψ(J_0,θ) is very peaked around θ=π/2 so the upper limit can be extended to infinite, and the following expressions are obtained:


    C_4k^(1)=h_08/√(2π)√(√(β^2+(4m)^2)-β)/√(β^2+(4m)^2)




    D_4m-2^(2)=h_08/√(2π)√(√(β^2+(4m-2)^2)+β)/√(β^2+(4m-2)^2)


Where the same principle was used to obtain (<ref>). 
The parameter h_0 is obtained by comparing the behaviour of (<ref>) for high n:

    C_4m^(1)≈ h_0 32m/n^3/4√(2/π),

to the value already calculated here in (<ref>) from continuity conditions for n→∞ and equal to -32ma/π, then:

    h_0=-n^3/4a/√(2π)

The d-line can be obtained from the asymptotic behaviour of the current, that can be calculated by adding up the series (<ref>) when replacing (<ref>), for n→∞:


    z_1=-1+ia+8a/πe^iθ∑_m=1^∞4m/(4m)^2-1( J/J_0)^(4m)^2-1×
    [4mcos (4mθ)  -i (4m)^2sin(4mθ)],


This expression can be written as a function of the series I_0:

    z_1=-1+ia+8a/πe^iθ(I_0+i∂ I_0/∂θ)

With:


    I_0=   ∑_m=1^∞(4m)^2/(4m)^2-1( J/J_0)^(4m)^2-1cos (4mθ) 
    
    =   ∑_m=1^∞(1+1/(4m)^2-1)( J/J_0)^(4m)^2-1cos (4mθ)


Then, for τ=log(J_0/J)>0, z≫ a, J≈ J_0, τ≪ 1, it gets (J/J_0)^(4m)^2-1=e^-τ((4m)^2-1)≈ e^-τ (4m)^2 and:

    I_0≈∑_m=1^∞(1+1/(4m)^2-1)e^-τ (4m)^2cos (4mθ)


    I_0≈∑_m=1^∞ e^-τ (4m)^2cos (4mθ) +∑_m=1^∞e^-τ (4m)^2/(4m)^2-1cos (4mθ)

The exponential term in the second sum can be omitted for τ≪ 1, and replacing the angle γ=θ-π/2:

    I_0≈∑_m=1^∞ e^-τ (4m)^2cos (4mγ) +∑_m=1^∞cos (4mγ) /(4m)^2-1

The first sum can be estimated with the Euler-Maclaurin formula <cit.>, giving exp(-γ^2/4τ)√(π)/8√(τ)-1/2 and the second adds to 1/2-πsin(2γ)/4. Then:

    I_0≈√(π)/8√(τ)e^-γ^2/4τ-π/4sin 2γ


    I_0^'=∂ I_0/∂θ≈-√(π)γ/16τ^3/2e^-γ^2/4τ-π/2cos 2γ

Then, for y≫ a, x≫ 1, from (<ref>):

    x   =b(I_0 sinγ +I_0^'cosγ)
    
    y   =b(-I_0 cosγ +I_0^'sinγ)
 
With b=8a/π. By replacing I_0 and I_0^', naming e_1=exp(-γ^2/4τ)√(π)/8√(τ):

    x   =be_1(sinγ +γ/2τcosγ)
    
    y   =be_1(cosγ -γ/2τsinγ)
 
By neglecting the terms not proportional to e_1. Then, eliminating e_1, and solving for τ:

    τ=-γ/2(xsinγ-ycosγ)/(xcosγ+ysinγ)

Combining (<ref>) and (<ref>), replacing back in (<ref>), and expanding for γ≪ 1:

    32/π b^2x^2=(1+ξγ)/(ξ-γ)^3γe^-γ(1+ξγ)/(ξ-γ)

In which ξ=y/x. As discussed in <cit.>, the d-line can be obtained from the line where the space derivatives of γ(x,y) diverge. This can be calculated from the condition ∂ x/∂γ =0. Then (<ref>) is approximated for ξ≪ 1,  derived with respect to γ while keeping ξ constant and equalled to zero, from which it is obtained γ=ξ/2. Replacing this back in (<ref>) and retaining the lowest order in ξ, the expression  x=y^2/(8a/√(2π e)) for the d-line is obtained. 




§ REFERENCES
 

99


Brandt0 E. H. Brandt , Rep. Prog. Phys. 58, 1465 (1995).

Koshelev A. E. Koshelev, I. A. Sadovskyy, C. L. Phillips, A. Glatz, Phys. Rev. B 93, 060508(R) (2016).

Blatter G. Blatter, V. B. Geshkenbein, J. A. G. Koopmann, Phys. Rev. Lett. 92, 067009 (2004).


Jerem1 J. Brisbois, O.-A. Adami, J. I. Avila, M. Motta, W. A. Ortiz, N. D. Nguyen, P. Vanderbemden, B. Vanderheyden, R. B. G. Kramer, A. V. Silhanek, Phys. Rev. B 93, 054521 (2016).

Mints R. G. Mints, E. H. Brandt, Phys. Rev. B 54, 12421 (1996).

Bean1 C. P. Bean, Phys. Rev. Lett. 8, 250 (1962)

Bean2 C. P. Bean, Rev. Mod. Phys. 36, 31 (1964 ).

Bean3 C. P. Bean, J. Appl. Phys. 41, 2482 (1970).

Schuster1 Th. Schuster, M. V. Indenbom, M. R. Koblischka, H. Kuhn, H. Kronmüller, Phys. Rev. B 49, 3443 (1994).

Schuster2 Th. Schuster, H. Kuhn, M. V. Indenbom, Phys. Rev. B 52, 15621 (1995).

Eisenmenger J. Eisenmenger, P. Leiderer, M. Wallenhorst, H. Dötsch, Phys. Rev. B 64, 104503 (2001).

Vinokur1 V. M. Vinokur , M. V. Feigel'man, V. B. Geshkenbein,  Phys. Rev. Lett. 67, 915 (1991).


Rak A. L. Rakhmanov, D. V. Shantsev, Y. M. Galperin, T. H. Johansen, Phys. Rev. B 70, 224502 (2004).

Deni D. V. Denisov, A. L. Rakhmanov, D. V. Shantsev, Y. M. Galperin, T. H. Johansen, Phys. Rev. B 73, 014512 (2006).

Gur1 A. Gurevich, J. McDonald, Phys. Rev. Lett. 81, 2546 (1998).

Gur2 A. Gurevich, M. Friesen, Phys. Rev. B 62, 4004 (2000).

Gur3 M. Friesen, A. Gurevich, Phys. Rev. B 63, 064521 (2001).

Camp1 A. M. Campbell, J. E. Evetts,  Critical Currents in Superconductors, Monographs on Physics. Taylor & Francis Ltd., London (1972).

Kim1 Y. B. Kim, C. F. Hempstead, A. R. Strnad, Phys. Rev. 139, A1163 (1965).



Vest1 J. I. Vestgården , P. Mikheenko, Y. M. Galperin, T. H. Johansen, New J. Phys. 15, 093001 (2013).

Abramowitz Formula 3.6.28 in Handbook of Mathematical Functions With Formulas, Graphs, and Mathematical Tables,  by M. Abramowitz and I. A. Stegun, National Bureau of Standards, Applied Mathematics Series - 55, 1972. 





