Novel vortex structures in the three-dimensional superconductor under the helical magnetic field from the chiral helimagnet
Saoto Fukui Masaru Kato Yoshihiko Togawa Osamu Sato


APS/123-QED












^1Department of Mathematical Sciences, Osaka Prefecture University, 1-1, Gakuencho, Sakai, Osaka 599-8531, Japan
^2Department of Physics and Electronics, Osaka Prefecture University, 1-1, Gakuencho, Sakai, Osaka 599-8531,Japan
^3Osaka Prefecture University College of Technology, 26-12, Saiwaicho, Neyagawa, Osaka 572-8572, Japan





 








 









 
 





             


We have investigated vortex structures in three-dimensional superconductors under a helical magnetic field from a chiral helimagnet numerically.
In order to obtain vortex structures, we solve three-dimensional Ginzburg-Landau equations with the finite element method.
The distribution of the helical magnetic field is assumed to be proportional to the distribution of the magnetic moments in the chiral helimagnet.
Then, the magnetic field is the same direction in the yz-plane and helical rotation along the helical axis.
Under this helical magnetic field, vortices appear to be perpendicular to the surface of the superconductor.
But we have found that there are tilted vortices toward the helical axis, although there is no component of the magnetic field along the helical axis.
This vortex structure depends on the chirality of the distribution of the helical magnetic field.















74.25.Uv, 74.20.De, 74.81.-g, 75.50.Cc


                             

                              
Novel vortex structures in the three-dimensional superconductor under the helical magnetic field from the chiral helimagnet
    Osamu Sato^3
    December 30, 2023
===========================================================================================================================








§ INTRODUCTION

A chiral helimagnet is one of the attractive magnetic materials.
In this material, magnetic moments form either left- or right-handed helical rotation.

The schematic magnetic structure of the chiral helimagnet is shown in Fig. <ref>(a).
This helical configuration of magnetic moments comes from the competition between two interactions: the ferromagnetic exchange interaction between nearest neighbor magnetic moments and the Dzyaloshinsky-Moriya (DM) interaction.
The ferromagnetic exchange interaction causes nearest neighbor magnetic moments to be parallel, while the DM interaction causes nearest neighbor magnetic moments to be perpendicular to each other<cit.>. 
The DM interaction determines the direction of the rotation, left-handed or right handed.

The magnetic structure in the chiral helimagnet has been observed in CrNb_3S_6<cit.>, CsCuCl_3<cit.>, and Yb(Ni_1-xCu_x)_3Al_9<cit.> experimentally.
Then, properties of chiral helimagnets have been investigated experimentally<cit.> and theoretically<cit.>.

For example, under an external magnetic field perpendicular to the helical axis, the helical configuration of magnetic moments changes into an incommensurate magnetic structure<cit.>.
This magnetic structure is called a chiral soliton lattice, which is shown in Fig. <ref>(b).
The chiral soliton lattice consists of ferromagnetic domains periodically partitioned by 360^∘ domain walls.
The chiral soliton lattice has been also observed experimentally with the Lorentz microscopy and the small-angle electron diffraction<cit.>.
Other interesting phenomena are the giant magnetoresistance<cit.>, the magneto-chiral dichroism<cit.>, the creation of the spin current<cit.>, and the Berezinskii-Kosterlitz-Thouless (BKT) transition<cit.>.


From other point of view, we may expect this peculiar magnetic structure affects other materials.
Therefore, in this paper, we focus on effects on superconductors, in particular, vortex structures in type-II superconductors.


These vortex structures in the type-II superconductor are important for a critical magnetic field and a critical current.
In general, under an uniform magnetic field, vortices appear and form a triangular lattice, which is called the Abrikosov lattice <cit.>.
When the external current flows in the superconductor and vortices move, the electric resistance occurs, which leads to break superconductivity.
So, controlling vortex states is a key factor for applications of superconductors.


One of plausible vortex controlling method is using a ferromagnet.
In a ferromagnet / superconductor bilayer system, vortices appear spontaneously <cit.>.
This phenomenon comes from interaction between magnetic domains in the ferromagnet and magnetic fluxes of vortices.
This spontaneous vortex state enhances superconductivity, in particular, its critical current by a pinning of vortices in the superconductor due to magnetic domains in the ferromagnet. 
From this interference between the ferromagnet and the superconductor, we expect novel effects of peculiar magnetic materials; the chiral helimagnet.
Therefore, we investigate effects of the chiral helimagnet on the superconductor theoretically.

In our previous study, we have investigated vortex structures in two-dimensional superconductors under the magnetic field from the chiral helimagnet<cit.>.
Although, magnetic field is created by the helically structured magnetic moments in the chiral helimagnet only the perpendicular component of the magnetic field to the superconducting surface is effective, and other components are neglected in the two-dimensional superconductor, .
In this paper, we consider three-dimensional superconductors, and taking all components of the magnetic field, investigate vortex structures under the helical magnetic field from the chiral helimagnet completely.
In section II, we introduce numerical methods in order to obtain vortex structures in three-dimensional superconductors.
In section III, we show vortex structures under the helical magnetic field from the chiral helimagnet, and discuss origin of these vortex structures.
Finally, in section IV, we summarize our results.




§ METHOD

We consider a three-dimensional superconductor under the helical magnetic field.

The distribution of the helical magnetic field is assumed to be proportional to the distribution of the magnetic moments in the chiral helimagnet.
We obtain distributions of the order parameter in superconductors by solving the Ginzburg-Landau equations.
We start from the Ginzburg-Landau free energy,

    ℱ(ψ,) = ∫_V ( f_n + α(T) |ψ|^2 + β/2 |ψ|^4 ) dV 
            + ∫_V {1/2m_s| ( -iħ∇ - e_s /c) ψ|^2 
                                                                    
                                                                     + ||^2/8π - ·_ ext/4π} dV,

where ψ is a superconducting order parameter and f_n is a free energy density of the normal state.
α(T) is a coefficient which depends on the temperature T, α(T) = α'(T-T_c).
α' and β is a positive constant and T_c is the critical temperature of the superconductor.
m_s is the effective mass of the superconductor and e_s is an effective charge of electrons in the superconductor.
= ∇× is a local magnetic field and  is a magnetic vector potential.
_ ext is an external magnetic field, which is included the magnetic field from the chiral helimagnet.





The order parameter and the vector potential is normalized as,

    ψ̃ = ψ/√(α)/β,    = 2π/Φ_0.

Φ_0 = ch/2e is the quantum flux and e is an elementary charge.
Using the normalized order parameter and the vector potential, we obtain following equations from Eq. (<ref>),

    ∫_V[ ( i∇ψ̃ - ψ̃) ( -i ∇ (δψ̃) -  (δψ̃) ) .  
          . + ( i∇ (δψ̃) -  (δψ̃) ) ( -i∇ψ̃^∗ - ψ̃^∗) . 
          . + 1/ξ^2( |ψ̃|^2 - 1 ) ( ψ̃ (δψ̃^∗) + ψ̃^∗ (δψ̃) ) ] dΩ = 0, 
    
          ∫_V[ κ^2 ξ^2 (  div · div (δ) +  rot · rot (δ) ) . 
          . |ψ|^2 · (δ) - i/2{ψ̃^∗ (∇ψ̃) - ψ̃ (∇ψ̃^∗) }] dΩ
           = κ^2 ξ^2 ∫_V 2π/Φ_0_ ext· rot (δ) dΩ,

where δψ̃ and δ are variations, or test functions of the order parameter and the vector potential, respectively.
κ = λ/ξ is the Ginzburg-Landau parameter, and λ and ξ are the penetration length and the coherence length, respectively.

In order to solve Eqs. (<ref>) and (<ref>), we use the three-dimensional finite element method (FEM).
In the three-dimensional FEM, we divide the system into tetrahedron elements (see Fig. <ref>).

In a tetrahedron, there are four volume coordinates, which is given as,

    N_i^e = ( a_i + b_ix + c_i + d_iz )/6V    (i=1, 2, 3, 4),

where V is a volume of the tetrahedron.
a_i, b_i, c_i, and d_i are given in the Appendix.
The order parameter ψ̃ and the vector potential  are expanded using these volume coordinates,

    ψ̃()    =   ∑_i,e^N_eψ̃_i^e N_i^e, 
     ()    =   ∑_i,e^N_e_i^e N_i^e,
 
where ψ̃_i^e and _i^e are the order parameter and the vector potential at i-th node in the e-th element, respectively.


N_e is the number of finite elements.
We consider following test functions δψ̃ and δ,

    δψ̃_8(e-1) + 2j-1   =    N_j^e    (x ∈ e-th element)  
     0        (otherwise), 
    δψ̃_8(e-1) + 2j   =    iN_j^e   (x ∈ e-th element) 
     0        (otherwise), 
    δ_12(e-1) + 3(j-1)+i   =   _i N_j^e (x ∈ e-th element)  
     0        (otherwise),

where i=1,2,3 and j=1,2,3,4.
_1, _2, and _3 are basis vectors in the three-dimensional space.
We substitute Eqs. (<ref>)-(<ref>) into Eqs. (<ref>) and (<ref>), we obtain these equations,

    ∑_j [ P_ij({Ã}) + P_ij^2R({ψ̃}) ]  Reψ̃_j^e 
            + ∑_j [ Q_ij ({}) + Q_ij^2({ψ̃}) ]  Imψ̃_j^e = V_i^R({ψ̃}),    
          ∑_j [ -Q_ij({}) + Q_ij^2({ψ̃}) ]  Reψ̃_j^e 
            + ∑_j [ P_ij({}) + P_ij^2I({ψ̃}) ]  Imψ̃_j^e = V_i^I({ψ̃}), 
          ∑_j R_ij ({ψ̃}) Ã_jx + ∑_j S_ij^xyÃ_jy + ∑_j S_ij^xzÃ_jz
            = T_i^x - U_i^x, 
          ∑_j R_ij ({ψ̃}) Ã_jy + ∑_j S_ij^yxÃ_jx + ∑_j S_ij^yzÃ_jz
            = T_i^y - U_i^y, 
          ∑_j R_ij ({ψ̃}) Ã_jz + ∑_j S_ij^zxÃ_jx + ∑_j S_ij^zyÃ_jy
            = T_j^z - U_i^z.

Coefficients are given in Appendix and the reference <cit.>.
Solving Eqs. (<ref>)-(<ref>) self consistently, we obtain real and imaginary parts of the order parameter Re ψ̃, Im ψ̃ and three-components of the vector potential A_x, A_y, and A_z.

The magnetic field _ ext includes the helical magnetic field from the chiral helimagnet, 

    _ ext = _ CHM + _ appl,

where _ CHM is the magnetic field from the chiral helimagnet and _ appl is the homogeneous applied magnetic field.
We consider that a distribution of _ CHM is proportional to the configuration of magnetic moments in the chiral helimagnet.
The configuration of magnetic moments is obtained by the Hamiltonian of the chiral helimagnet <cit.>,

    ℋ   =    -2J ∑_n_n ·_n+1 + ·∑_n _n ×_n+1
           - 2μ_B _ appl·∑_n _n,

where _n is the spin at the n-th site, μ_B is the Bohr magneton.
This Hamiltonian consists of three terms.
The first term is the ferromagnetic exchange interaction with magnitude J (>0).
The second term is the DM interaction with the DM vector .
The last term is the Zeeman energy.
We assume that the helical axis is the x-axis.
In the monoaxial chiral helimagnet, the DM vector is parallel to the helical axis, = (D,0,0), and we assume that the direction of spin is perpendicular to the helical axis, so φ = π/2.
We express n-th spin as _n = S(sinθ_ncosφ, sinθ_nsinφ, cosθ_n).
We set _ appl = (0,0,H_ appl).
In the typical chiral helimagnet CrNb_3S_6, the helical period L=48nm is much longer than the lattice constant.
So, we consider the continuum limit.
We minimize Eq. (<ref>) in the continuum limit with respect to θ(x).
We obtain the Sine-Gordon equation,

    d^2θ(x)/dx^2 - H^∗sinθ(x) = 0,

where H^∗ = 2μ_B H_ appl/(a^2S^2√(J^2 + ||^2)) is a normalized applied magnetic field and a is the lattice constant.
The solution of Eq.(<ref>) is,

    sin( θ - ϕ/2) =  sn( √(H^∗)/kx | k ),

or,

    θ(x) = 2sin^-1[  sn( √(H^∗)/kx | k ) ] + ϕ.

sn(x|k) is the Jacobi's elliptic function, k is the modulus and ϕ is an initial angle at x=0.
k is determined by the relation,

    πtan^-1(||/J)/4√(H^∗) = E(k)/k,

where E(k) is the complete elliptic integral of the second kind.
Using Eq. (<ref>), the external magnetic field in Eq. (<ref>) is,

    (_ ext)_x(x)    =    0, 
    
     (_ ext)_y(x)    =    H_0 sinθ(x), 
    
     (_ ext)_z(x)    =    H_0 cosθ(x) + H_ appl.

H_0 is a magnitude of the magnetic field from the chiral helimagnet.
We solve the Ginzburg-Landau equations (<ref>)-(<ref>) using the magnetic field in Eqs. (<ref>)-(<ref>) numerically.



§ RESULT

Solving the Ginzburg-Landau equations (<ref>)-(<ref>) self-consistently, we obtain distributions of the order parameter in the superconductor under the chiral helimagnet.
We set the Ginzburg-Landau parameter κ = λ/ξ = 5 and the temperature T=0.3T_c, where T_c is the critical temperature of the superconductor.
The ratio between the ferromagnetic exchange interaction and the Dzyaloshinsky-Moriya interaction is taken from the experimental data<cit.>, ||/J = 0.16.
We consider a parallelepiped as our model shown in Fig. <ref>.
The system size is 1.0L'ξ_0 × 15ξ_0 × 13ξ_0.
ξ_0 is a coherence length at T=0 and L'ξ_0 is a helical period of the chiral helimagnet, which is given as,

    L' = 2π/tan^-1(D/J).

Here, the uniform applied magnetic field, H_ appl/(Φ_0/ξ_0^2) = 0.00.
For ||/J=0.16, L' is approximately 39.2699.
We assume the superconducting region is surrounded by the vacuum region.
The distance between the superconducting region and the vacuum region is 1.5ξ_0, and the size of the superconducting region is (1.0L' - 3.0)ξ_0 × 12ξ_0 × 10ξ_0.
When we calculate the Ginzburg-Landau equations, we set following boundary conditions,

    · = 0,  | ( -iħ∇ + e/c) ψ| · = 0.




The external magnetic field _ ext is given in Eqs. (<ref>)-(<ref>).
We set magnitudes of the helical magnetic field and the applied magnetic field in Eqs. (<ref>)-(<ref>) as H_0/(Φ_0/ξ_0^2) = 0.15 and H_ appl/(Φ_0/ξ_0^2) = 0.00, respectively.
First, we take the angle at x=0 in Eqs. (<ref>) or (<ref>) as ϕ = -π/2.
Then, the distribution of the helical magnetic field is shown in Fig. <ref>.
In Fig. <ref>, we show distribution of the helical magnetic field (Fig. <ref>(a)), each components of the magnetic field (_ ext)_x, (_ ext)_y, and (_ ext)_z (Figs. <ref>(b)-(d)).



Under the magnetic field in Fig. <ref>, we obtain the distribution of the order parameter shown in Fig. <ref>.
In Fig. <ref>, we show distributions of the order parameter in the xy-plane (Fig. <ref>(a)) and the zx-plane (Fig. <ref>(b)).
The cross sections parallel to the xy-plane at z=1.5ξ_0, 11.5ξ_0 and the zx-plane at y=1.5ξ_0, 13.5ξ_0 are interfaces between the superconducting region and the vacuum region.
In Fig. <ref>(a), two vortices appear in the regions around (x/ξ_0, y/ξ_0) ∼ (10, 7.5) and (30, 7.5), where (_ ext)_y/(Φ_0/ξ_0^2) ∼ 0.00 and (_ ext)_z/(Φ_0/ξ_0^2) ∼± 0.15 in Fig. <ref>.
In this regions, the magnetic field is parallel or anti-parallel to z-axis, so these two vortices have quantum fluxes antiparallel to each other.
On the other hand, in Fig. <ref>(b), one vortex appears in the region around (x/ξ_0, y/ξ_0) ∼ (20, 7.5), where (_ ext)_y/(Φ_0/ξ_0^2) ∼ 0.15 and (_ ext)_z/(Φ_0/ξ_0^2) ∼ 0.00.
In this region, the magnetic field is parallel to the direction of y-axis, so the vortex has a quantum flux parallel to the y-axis.
In total, three vortices appear.
They are separated by 0.25ξ_0.
And the angle between nearest neighbor vortices is π/2. 
































The vortex structure is different from our previous work<cit.>, although the model and numerical parameters are the same.
Only difference is the initial random state of our iteration method for solving Eqs. (<ref>)-(<ref>).
To determine more stable state, we should calculate the free energies for both states.
This is the future problem.

Next, we investigate vortex structures with the other distribution of the helical magnetic field with ϕ = π in Eq. (<ref>), which is shown in Fig, <ref>.
Under this magnetic field, we obtain distributions of the order parameter shown in Fig. <ref>.
In Fig. <ref>, we show the distributions of the order parameter and the phases of the order parameter.
In Fig. <ref>(a), two vortices appear in the region of the magnetic field (_ ext)_z/(Φ_0/ξ_0) > 0.
Positions of these vortices in these cross sections from z=1.5ξ_0 to 11.5ξ_0 change along the x-axis.
So, two vortices tilt toward the x-axis, but the x-component magnetic field is zero, (_ ext)_x/(Φ_0/ξ_0^2) = 0.
Vortices parallel to the y-axis does not appear in Fig. <ref>(b). 
This result comes from the screening current of demagnetization factor of the superconductor.
Next, in order to avoid difference between shielding fields parallel to y- and z-axis, we investigate vortex structure in the system with the larger size.
The system size is 1.0L'ξ_0 × 15ξ_0 × 15ξ_0.
In this system, the cross section parallel to the yz-plane is square.
We take the same numerical parameters κ = 5, T=0.3T_c, and ||/ξ_0 = 0.16, and the same distribution of the helical magnetic field, H_0/(Φ_0/ξ_0^2) = 0.15, H_ appl/(Φ_0/ξ_0^2) = 0.0, and ϕ = π in Eq. (<ref>).
Under these numerical parameters, we obtain vortex structures shown in Figs. <ref>(a) and <ref>(b).
We show vortex structures at top and bottom surface and the center cross section in the superconductor.
In Fig. <ref>(a), we see that two vortices around x ∼ 20ξ_0 tilt toward the x-axis around x ∼ 10ξ_0 and 30ξ_0, while, two vortices around x ∼ 10ξ_0 and 30ξ_0, where the magnetic field (_ ext)_y/(Φ_0/ξ_0^2)= ± 0.15, are parallel or antiparallel to y-axis.
So, these vortices have quantum fluxes parallel to the direction of y-axis. 

Compare to the previous model, sheilding the field along the y-axis costs equal energy with shielding the field along the z-axis.
Then, vortices parallel to the direction of y-axis also appear in this system.






































































Next, we examine how the chirality of the helical rotation of magnetic field affects the vortex structures, i.e. we examine difference between vortex structures under right- and left-handed helical magnetic field.
In order to change the direction of the rotation, we take the opposite DM vector, ||/J=0.16 and  is antiparallel to the x-axis.
The helical magnetic field for this DM vector is shown in Fig. <ref>.
The rotation of the helical magnetic field in Fig. <ref> is opposite to that in Fig. <ref>.
We solve the Ginzburg-Landau equations using this helical magnetic field.
We take the numerical parameter, κ = 5, T=0.3T_c, and the magnetic field, H_0/(Φ_0/ξ_0^2) = 0.15 and H_ appl/(Φ_0/ξ_0^2) = 0.0.
The system size is 1.0L'ξ_0 × 15ξ_0 × 15ξ_0, which is the same system size as Fig. <ref>.
Under these conditions, we obtain vortex structures shown in Figs. <ref>(a) and <ref>(b).
Comparing between Figs. <ref>(a) and <ref>(a), directions of vortices are completely opposite.
Then, the direction of tilt of vortices depends on the chirality of the helical magnetic field.


Experimentally, these vortex structures may appear in superconductor / chiral helimagnet hybrid structures.
For example, our model may be equivalent to the system in which a small superconductor is surrounded by a large chiral helimagnet.
On the other hand, in the chiral helimagnet / superconductor bilayer system, when the superconductor is thin, only the perpendicular component of the magnetic structure is effective.
Then, our previous work on the two-dimensional superconductor is applicable to such bilayer systems<cit.>.



Finally, we discuss the movement of vortices under the external current.
When the external current is applied to vortex structures in Figs. <ref>, <ref>, <ref>, and <ref>, along the helical axis, vortices easily move perpendicular to the x-axis.
On the other hand, when the external current flows along the direction of the y-axis, vortices move in the surface parallel to the zx-plane.

If the z-component of vortex direction is not zero, the vortex moves toward the x-axis.
But the distribution of the helical magnetic field varies along the helical axis (x-axis) spatially.
So, the interaction between the vortex and the magnetic field changes.
Then, motion of vortices are obstructed by this interaction, which leads to the increase of the critical current. 




























§ SUMMARY

We have investigated vortex structures in three-dimensional superconductor under the helical magnetic field from the chiral helimagnet numerically.
We have obtained distributions of the order parameter using the three-dimensional Ginzburg-Landau equations.

When two vortices appear in one magnetic field region (_ ext)_z > 0, vortices tilt toward the x-axis, or the helical axis in spite of (_ ext)_x = 0.
This configuration may come from the interaction between vortices and between the vortex and the boundary in the system.
It is confirmed that when the rotation of the helical magnetic field is reversed, directions of tilt of vortices are also reversed.
These vortex structures do not occur under the uniform magnetic field.
Detailed discussion of stability of these vortex states need more simulations.



We only consider superconductors under the helical magnetic field, but without the uniform magnetic field and the external current.
Under the uniform magnetic field, the magnetic structure in the chiral helimagnet changes into the chiral soliton lattice.
In the microscopic system, uniform magnetic field decreases the number of solitons in the chiral soliton lattice discretely<cit.>.
Then, we expect that this discrete change of the soliton affects the structures of vortices.

If the external current flows in this system, two vortices that tilt toward the helical axis move uniquely because of the complicated distributions of currents in the superconductor.
Investigations and discussions about these phenomena are future works.









§ COEFFICIENTS FOR THE FINITE ELEMENT METHOD

In this appendix, we give coefficients in Eqs. (<ref>), (<ref>)-(<ref>) for the finite element method.
First, a_i, b_i, c_i, and d_i in Eqs. (<ref>) are defined as,

    a_i    =   ϵ_i { x_j(y_kz_l - y_lz_k) + x_k(y_lz_j - y_jz_l) . 
          . + x_l(y_jz_k - y_kz_j) }, 
    
     b_i    =   ϵ_i {  y_j(z_l-z_k) + y_k(z_j-z_l) . 
          . + y_l(z_k-z_j) }, 
    
     c_i    =   ϵ_i {  z_j(x_l-x_k) + z_k(x_j-x_l) . 
          . + z_l(x_k-x_j) } ,
    
     d_i    =   ϵ_i {  x_j(y_l-y_k) + x_k(y_j-y_l) . 
          . + x_l(y_k-y_j) },

where x_i, y_i, and z_i  (i=1,2,3,4) are coordinates for nodes in Fig. <ref>(b) and ϵ_i = 1  (i=1,3) or ϵ_i = -1  (i=2,4).
(i,j,k,l) is a cyclic permutation of (1,2,3,4).

Then, we show coefficients in Eqs. (<ref>)-(<ref>).
They are given as,

    P_ij({})    =   ∑_α K_ij^αα + ∑_i_1,i_2 I_i_1i_2ij∑_α A_i_1α^e A_i_2α^e - 1/ξ(T)^2 I_ij, 
    
    
     P_ij^2R({ψ})    =   1/ξ(T)^2∑_i_1i_2 I_i_1i_2ij (3 Reψ_i_1^e  Reψ_i_2^e +  Imψ_i_1^e  Imψ_i_2^e), 
    
    
     P_ij^2I({ψ})     =   1/ξ(T)^2∑_i_1i_2 I_i_1i_2ij ( Reψ_i_1^e  Reψ_i_2^e + 3 Imψ_i_1^e  Imψ_i_2^e), 
    
    
     Q_ij({})    =   ∑_i_1∑_α (J_ji_1i^α - J_ii_1j^α) A_i_1^α, 
    
     Q_ij^2({ψ})        =   2/ξ(T)^2∑_i_1i_2 I_i_1i_2ij Reψ_i_1^e  Imψ_i_2^e, 
    
     R_ij({ψ})            =   κ^2 ξ(T)^2 ∑_α K_ij^αα + ∑_i_1i_2 I_i_1i_2ijψ_i_1^e∗ψ_i_2^e, 
    
     T_i^α({ψ})            =   ∑_i_1i_2 J_i_1i_2i^α Im (ψ_i_1^e∗ψ_i_2^e), 
    
     V_i^R({ψ})                  =   2/ξ(T)^2∑_i_1i_2i_3 I_i_1i_2i_3i Re(ψ_i_1^e ψ_i_2^e∗)  Reψ_i_3^e, 
    
     V_i^I({ψ})                   =   2/ξ(T)^2∑_i_1i_2i_3 I_i_1i_2i_3i Re(ψ_i_1^e ψ_i_2^e∗)  Imψ_i_3^e, 
    
     S_ij^αβ   =   κ^2 ξ(T)^2 (K_ij^αβ - K_ij^βα), 
    
     U_i^α   =   κ^2 ξ(T)^2 2π/Φ_0( J_i^βH_γ - J_i^γH_β).
 
Here, α, and  β=x,y,z.
In Eq. (<ref>), (α, β, γ) is a cyclic permutation of (x,y,z).
Their coefficients are represented by integrals I_ij, I_i_1i_2i_3, I_i_1i_2i_3i_4, J_i^x_i, J_i_1i_2i_3^x_i, and K_ij^x_ix_j.
These integrals are given by,

    I_ij   =   ∫_V_e N_i^e N_j^e dV, 
    
     I_i_1i_2i_3   =   ∫_V_e N_i_1^e N_i_2^e N_i_3^e dV, 
    
     I_i_1i_2i_3i_4   =   ∫_V_e N_i_1^e N_i_2^e N_i_3^e N_i_4^e dV, 
    
     J_i^x_i   =   ∫_V_e∂ N_i^e/∂ x_i dV, 
    
     J_i_1i_2i_3^x_i   =   ∫_V_e∂ N_i_1^e/∂ x_iN_i_2^e N_i_3^e dV, 
    
     K_i_1i_2^x_ix_j   =   ∫_V_e∂ N_i_1^e/∂ x_i∂ N_i_2^e/∂ x_j dV.



99
 DM_1
  I. Dzyaloshinsky, J. Phys. Chem. Solids 4, 241 (1958).
 DM_2
  T. Moriya, Phys. Rev. 120, 91 (1960).
 CHM_Cr_1
  N. J. Ghimire, M. A. McGuire, D. S. Parker, B. Sipos, S. Tang, J.-Q. Yan, B. C. Sales, and D. Mandrus, Phys. Rev. B 87, 104403 (2013).
 CHM_Cr_2
  K. Tsuruta, M. Mito, H. Deguchi, J. Kishine, Y. Kousaka, J. Akimitsu, and K. Inoue, Phys. Rev. B 93, 104402 (2016).
 CHM_Cr_3
  Y. Togawa, Y. Kousaka, K. Inoue, and J. Kishine, J. Phys. Soc. Jpn. 85, 11201 (2016).
 CHM_Cs_1
  H. Ohsumi, A. Tokuda, S. Takeshita, M. Takata, M. Suzuki, N. Kawamura, Y. Kousaka, J. Akimitsu, and T. Arima, Angew. Chem., Int. Ed. 52, 8718 (2013).
 CHM_Cs_2
  Y. Kousaka, T. Koyama, K. Ohishim, K. Kakurai, V. Hutanu, H. Ohsumi, T. Arima, A. Tokura, M. Suzuki, N. Kawamura, A. Nakao, T. Hanashima, J. Suzuki, J. Campo, Y. Miyamoto, A. Sera, K. Inoue, and J. Akimitsu, Phys. Rev. Material 1, 071402(R) (2017).
 CHM_Yb_1
  S. Ohara, S. Fukuta, K. Ohta, H. Kono, T. Yamashita, Y. Matsumoto, and J. Yamaura, JPS Conf. Proc. 3, 017016 (2014).
 CHM_Yb_2
  T. Matsuura, Y. Kita, K. Kubo, Y. Yoshikawa, S. Michimura, T. Inami, Y. Kousaka, K. Inoue, and S. Ohara, J. Phys. Soc. Jpn. 86, 124702 (2017).
 CHM_theory_1
  M. Shinozaki, S. Hoshino, Y. Masaki, J. Kishine, and Y. Kato, J. Phys. Soc. Jpn 85, 074710 (2016).
 Kishine_theory
  J. Kishine and A. S. Ovchinnnikov, Solid State Phys. 66, 1 (2015).
 Kishine_1
  J. Kishine, K. Inoue, and Y. Yoshida, Prog. Theor. Phys. Suppl. 159, 82 (2005).
 Togawa_PRL_1
  Y. Togawa, T. Koyama, K. Takayanagi, S. Mori, Y. Kousaka, J. Akimitsu, S. Nishihara, K. Inoue, A. S. Ovchinnikov, and J. Kishine, Phys. Rev. Lett. 108, 107202 (2012).
 Togawa_PRL_2
  Y. Togawa, Y. Kousaka, S. Nishihara, K. Inoue, J. Akimitsu, A. S. Ovchinnikov, and J. Kishine, Phys. Rev. Lett. 111, 197204 (2013)
 MCh_Dh
  G. L. J. A. Rikken and E. Raupach, Nature 390, 493 (1997)
 spin_current
  K. Tokushuku, J. Kishine, and M. Ogata, J. Phys. Soc. Jpn. 86, 124701 (2017).
 BKT
  I. Proskurin, A. S. Ovchinnikov, and J. Kishine, J. Phys.: Conf. Ser. 903, 012062 (2017).
 Abrikosov 
  A. Abrikosov, Sov. Phys. JETP 5, 1774 (1957).
 Hess
  H. F. Hess, R. B. Robinson, R. C. Dynes, J. M. Valles, Jr., and J. V. Waszczak, Phys. Rev. Lett. 62, 214 (1989).
 FM_SC_Hybrid
  I. F. Lyuksyutov and V. L. Pokrovsky, Adv. Phys. 54, 1 (2005).
 FM_SC_FSB_ex
  M. Lange, M. J. Van Bael, V. V. Moshchalkov, and Y. Bruynseraede, Appl. Phys. Lett., 81, 322 (2002).
 FM_SC_FSB_theory
  M. A. Kayali and V. L. Pokrovsky, Phys. Rev. B 69, 132501 (2004).  
 Fukui_SUST
  S. Fukui, M. Kato and Y. Togawa, Supercond. Sci. Technol. 29, 125008 (2016)
 Fukui_JPSJ
  S. Fukui, M. Kato, Y. Togawa, and O. Sato, submitted to J. Phys. Soc. Jpn. 
 Fukui_3d_proc
  S. Fukui, M. Kato, Y. Togawa, and O. Sato, submitted to J. Phys.: Conf. Ser.
 DM_Cr
  B. J. Chapman, A. C. Bornstein, N. J. Ghimire, D. Mandrus, and M. Lee, Appl. Phys. Lett. 105, 072405 (2014).
 

