Rotationally invariant isospinning baby-skyrmions dressed by fermions

We couple Fermions to the isospinning (2+1) baby-skyrme model. We show that consistent rotationally invariant localised solutions can be found but the Fermionic solutions to the equations of motion are not in general eigenstates of the Hamiltonian. These do however become full eigenstates in a particular limit. We also find the corresponding Fermionic eigenstates of the Hamiltonian, for which our Fermion ansatz is also well suited. For localised Fermionic states, a novel constraint on the internal isospin frequency exists. Throughout our studies, we include the backreaction of the Fermions on the baby-skyrmion.

charges [4]. The model can be thought to arise as an effective field theory from the integration of high energy Fermionic "quark" modes. In this sense the theory does not intrinsically possess any fermionic degrees of freedom. While phenomenologically appealing, the full Skyrme theory is not particularly easy to study. It is a highly non-linear theory in four dimensions and most theoretical progress has to be made numerically. While being a simpler model overall, the baby-Skyrmion model [5] [6] [7] captures the essence and main features of its higher dimensional parent remarkably well. This lower dimensional model is particularly suited for the study of Skyrmion interactions in condensed matter systems [8] [9] [10] and more recently for quantum computation purposes [11] [12] [13]. This model also yields solitonic solutions with Baryon charge but, unlike the parent theory, it is not a-priori built from integrating Fermionic degrees of freedom and therefore can be naturally coupled to them.
The static baby-skyrmion solutions admit an isospinning generalization [14] [15] [16] [17]. This endows the baby-skyrmion with time dependence, in the form of rotations in internal space parametrised by a frequency ω. This frequency is severly restricted: it must obey positivity of the effective kinetic terms in the Lagrangian and be such that perturbations of the effective potential are not tachyonic. The presence of isospin provides the baby-skyrmion with a new instability channel. In particular, for Baryon numbers B ≥ 2, the internal rotation can cause a decay into B = 1 constituent baby-skyrmion solutions. The rotationally symmetric B = 1 solutions are stable under such isospinning rotations.
The topic of Fermion couplings to baby-Skyrmions is a very interesting one, neatly parametrized by the Fermion-baby-Skyrmion coupling constant g. Most of the theoretical work regarding this coupled system is performed in the decoupling limit of g → 0 [18] [19] , where one can safely assume the Fermions do not modify the baby-Skyrmion profile. Fermion scatterings on baby-Skyrmion solutions were previously investigated in [20]. Supersymmetric extensions of this system were also previously considered [21]. Recently however, back reaction of the Fermions was included by raising the value of g and solving the resulting coupled system numerically [22] [23]. This was done for B = 1 baby-skyrmion solutions with no internal isospin, the results show a continuous space of Fermion eigenstates coupled to the baby-skyrmion.
In this paper we wish to extend this analysis to include internal isospin. Restricting to B = 1 rotationally symmetric stable baby-skyrmion solutions, we wish to investigate how one can couple these to Fermions. Unlike the previous case where the Fermions are eigenstates of the resulting Hamiltonian, we will show that the inclusion of isospin means such Fermions are no longer generically eigenstates, becoming so only in a very particular limit. We will therefore interpret the resulting energy of the system as an expectation value of the mean energy of the system, and not as an energy eigenvalue of the Hamiltonian. Throughout our studies, we include the backreaction on the baby-skyrmion.
The paper is structured as follows: in section 2 we introduce the action of the baby-skyrme system coupled with Fermions, section 3 is dedicated to the solutions of this system and finally in section 4 we provide the conclusions of our research and some avenues of future investigation.

The System
The action we use describing a baby Skyrme system coupled to a Fermion is where g is the coupling constant related to the baby-Skyrmion-Fermion interaction and the a 1 , a 2 are constant coefficients. ⃗ Φ is a vector field of dimension n = 3. The isospin matrices are defined as τ i = I ⊗ σ i , where σ i are the Pauli matrices. The spin matrices are defined as γ µ =γ µ ⊗ I. We also include a potential term for the baby-Skyrmion with its corresponding coupling constant κ 0 . The general literature also considers different kinds of potentials given by higher powers of the one chosen above. In this paper, we restrict our attention to this potential only.
The mass dimensions of the parameters appearing in the Lagrangian are In this work, we use the flat metric η µν = (1, −1, −1) and we choose the spin matricesγ 1 = −iσ 1 , γ 2 = −iσ 2 andγ 0 = σ 3 such that γ µ respect the Clifford algebra For the baby-Skyrmion, we will use the following parametrization where F and G are real-scalar fields. On the other side, for the fermion field we choose where and Ω 4 = Ω 4 (x µ ) are real-scalar field components. The coefficient N i is a normalization coefficient which ensures that Ψ † Ψ d 2 x = 1. In detail, this evaluates to The field Ψ is a spin-isospin spinor, meaning that it describes the wave function of a particle with spin one-half and two possible state of iso-spin: up and down. Without loss of generality, we also take the standard choice for the vacuum ⃗ n = (0, 0, 1) .
Given our choice for the metric, the positivity of the energy of the baby skyrmion sector is guaranteed requiring Finally, the topological charge of the baby Skyrmion is The equations of motion for the baby skyrmion fields are In the above equations c 1 = a 2 a 1 . Similarly, the Dirac equations for the Fermions read

Symmetries
The system admits at least two internal exact symmetries. The first symmetry of the action (1) is an SO(2) ∼ U (1) rotation that acts only on the field components G, Ω 2 and Ω 4 as where α is a constant phase. Equivalently, this rotation acts on the fundamental field Φ and Ψ as where and it corresponds to an iso-spin rotation. The related conserved current J µ and the Fermion current component J µ The previous definition of the isospin operator Q can help us to better clarify the physical meaning of the spin-isospin spinor Ψ in 2 + 1 dimensions. The last ingredient that we need to this end is the definition of the spin operator that, following [22], reads Now, given the basiŝ it is a simple computation to show that the isospin and spin operator Q and S act on these vectors as The spinor Ψ describes therefore a particle state with two possible value of spin and isospin. In this case, the complex functions ] represent the spatial wave-functions associated with each pair of spin-isospin numbers.
Besides the isospin rotation, the system is also invariant under the phase transformation that correspond to Ψ → Ψ ′ = e iβ Ψ. The associated Fermionic current J µ F is The conservation of this current corresponds to the conservation of quantum probability since the time-component of the four-vector (29) represents the probability density (ρ F from now) The divergence of the current (29) gives therefore the continuity equation 3 Combined Skyrmion-Fermion ansatz

The baby-Skyrmion ansatz
In order to solve the equations of motion for the baby Skyrmion-Fermion system, from this section we will use spatial polar coordinates {x, φ}, where x is the radial coordinate and φ is the angular coordinate.
For what concerns the baby-Skyrmion fields, we use the following ansatz Note that this ansatz implies a non trivial periodic dependence on the time coordinate. One of the fundamental aspects of this kind of ansatz is that even though the field profiles depend explicitly on time, the energy density does not, as we will show in section (3.5).
As shown in [14] [17] [16], insertion of this ansatz in the baby-Skyrme lagrangian (ignoring the Fermions for the time being) leads to an effective Lagrangian in which the isospinning contribution causes two main channels of instability. The first is related to the effective kinetic term for the baby-Skyrmion fields This term becomes singular when ω = 1. The second is related to the effective potential for the scalar field, In order to have a real value for the "pion mass", and thus an asymptotically decreasing behaviour for the baby skyrmion profile-function, we require the effective mass of baby-Skyrme excitations to be positive.
We define here an effective angular momentum for the baby-Skyrmion as The topological charge density (11) with the ansatz (32) simplifies to Given the interval for the x−coordinate x ∈ [0, ∞] and the boundary conditions for F the baby Baryon charge B is then B = p. Since in this work we want to consider the axially symmetric ansatz (32) for the baby-Skyrmion, we reduce our analysis to the case Using the baby-Skyrmion ansatz, the Lagrangian for the Fermions can be reduced to where the operator D has components 6

The Fermion ansatz
In the previous subsection, we chose the common rotational-invariant ansatz for the baby-Skyrmion, requiring the profile function F to be dependent only on the radial variable x. Moreover, although the explicit time-dependence of the field G, the baby baryon density ρ B with ansatz (32) does not depend on time and thus, in this sense, we can call that configuration a stationary configuration. Considering this choice for baby-Skyrmion sector, in order to solve the full equations of motion that include the coupling between the baby-Skyrmion and Fermion, we propose an ansatz for the fermion field that is consistent with axial-invariance and stationarity. With this aim, we choose with where the real functions v 1 , v 2 , u 1 and u 2 depend only on the radial coordinate x and in which l ∈ Z and k 1 ∈ R. This particular ansatz leads the Fermion density ρ F = Ψ † Ψ to be rotationalinvariant and time-independent, as requested before. In particular, as we will show in the next subsection, it decouples the angular dependence from the equations, leaving radial equations for the field profile functions only. Note, obviously, that the time-dependence of (41) is different with respect to the case of [22], in which zero isospining-angular-velocity ω was considered. This fact leads to a fundamental difference between the two cases: in [22], a fermionic solution of the equations of motion with rotational-invariance symmetry (and stationary fermionic density) is also simultaneously an Eigenstate of the Hamiltonian; contrarily, as we will show in details in section 3.4, in this case the field (40) cannot possess this property. A solution to the equations of motion below is therefore not generically an Eigenstate of the Hamiltonian. It becomes so only in a particular limit. A natural question would therefore be to relax the assumption of time-independence of the the baby-Skyrmion and Fermion profile functions, along with relaxing the need for cylindrical symmetry, and solve the full set of PDE equations. In that case, different considerations about the properties of such a Fermion state would be needed and we leave this for a future work.

Equations of motion
Using the ansatz (32) and (40)-(41), the equations of motion (12)- (13) and (14) reduce to (we will only be considering the case m = 0 in this paper and a 1 = 1, a 2 = 1) As anticipated before, our ansatz leads to a set of radial equations for the profile-functions, decoupling the angular dependence from the system. In the following, we choose to restrict our attention to the state characterized by l = 0 and p = −1. This choice will help the comparison of the previous equations with the ones for the energy Eigenstates in section 3.4, in which l = 0 and p = −1 correspond to the states with zero total angular momentum.
For the Fermionic profile function it's important to check the asymptotic behavior at x → ∞. Then, one can set f = 0 and solve analytically the equations that result, giving a Fermionic behaviour of the form for some constant C i and where K 0 is a Bessel function. This means in particular that exponentially localised Fermionic states can only be found if g 2 > (1 + k 1 ) 2 ω 2 . Note that this imposes a third and novel constraint on the possible values of ω, since Similarly, the relevant boundary conditions at the origin x → 0 are Given the system (42)-(46), once chosen the integer numbers l = 0 and p = −1, and the angular isospin velocity ω, then there exists a value of k 1 for which the system is solvable within the required boundary conditions, giving smooth solutions for the profile functions. Therefore, in this sense, the resolution of the equations of motion is numerically similar to an Eigenstates problem where in our case k 1 plays the role of the unknown Eigenvalue.
Note that while one may think that just taking ω = 0 naively should bring us back to the case without isospin studied in [22], a more subtle limit must be considered. In fact, the right procedure to reduce our system to the one of [22] is to take the limit ω → 0 together with k 1 → ∞, while keeping k 1 ω = const. Then, the energy Eigenvalue ϵ of such a state would be given by ϵ = −k 1 ω. It is clear from the form of the Fermion ansatz (41) that in this limit the Fermion becomes an exact Eigenstate of a time-independent Hamiltonian, taking the common form Ψ = ψ(x, φ)e −iϵt . It was found in [22], in the system with no isospin, that when the mass of the Fermions is zero, the Fermion equations reduce to two equations since u 1 = −v 2 and u 2 = v 1 . However, when ω ̸ = 0 this is not the case. It is only true in the limit discussed above. The difference in shapes between the field profiles is therefore a good measure of how far one is from the limit of exact Eigenstates.
We therefore wish to solve equations (42)-(46) subject to boundary conditions (47)-(50) and the integral constraint (8). These equations have to be solved numerically. For our numerical solver, we compactify the spatial domain using the change of coordinates x = z 1−z , with z ∈ [0, 1]. Our solver discretises the differential equations using a second order central difference scheme. The resulting set of simultaneous equations are solved using a Newton Method. We vary the number of points on the [0, 1] grid we use for our discretization in accordance to the accuracy required by the scale of our solutions, we checked that all our solutions are accurate up to O 10 −6 . As a further check, we always substituted our solutions back into the equations and calculated their norm, finding that the equations are satisfied to order O 10 −12 . For our normalization, we always use an occupation number N i = 1. The code was fully implemented in Wolfram language.
In figures 1 and 2 we show localised Fermionic solutions. As explained above, the solutions respect the condition g 2 > (1 + k 1 ) 2 ω 2 . In order to raise the back-reaction parameter g, we must be careful to respect the additional constraints posed on ω while making sure the numerical procedure converges on a solution. For fixed values of g, κ 0 and ω, we have not been able to find multiple branches of solutions with different k 1 and different number of nodes on the Fermion functions. However, we cannot exclude the existence of this possibility. In figure 3, we solved the system fixing the value of k 1 and finding ω for different g at fixed κ 0 . As said above and as shown in the figure, we found just one value of ω that makes the system solvable. Importantly, the point corresponding to ω = 0 on this plot is a point of exact coincidence between solutions to the equations of motion and Hamiltonian Eigenstates. This is obvious, since ω = 0 in this case coincides with the results of the no isospin study of [22]. It is a good check of our numerical procedure that the point on the g axis for which ω = 0 is closed to the one presented in the Eigenstate analysis in figure 6 in section 3.4, although a different κ 0 has been used.

Energy Eigenstates
In sections 3.3, we showed an example of time-evolution for a baby Skyrmion-Fermion system with non-zero isospin for the baby-Skyrmion. Using the proper ansatz (32) and (40)-(41), we proved that a time-evolution for such a system with rotationally invariant solutions exist. Moreover, these solutions have stationary Baryonic and Fermionic density, being the profile functions f , v 1 , v 2 , u 1 and u 2 timeindependent. However, as we are going to show in this section, differently from the case studied in [22] with no isospin, the stationary and rotational-invariant solutions found here are not Eigenstates of the Hamiltonian. In the case of [22], being the Hamiltonian of the system H time-independent, a simple  ansatz of the type Ψ = ψ(x, φ)e −iϵt reduces the Fermionic equation of motion (14) to the Eigenvalueequation Hψ = ϵψ. In that case, therefore, once the proper ansatz for ψ is chosen, the Fermionic state Ψ in the time-evolution can be simultaneously an energy Eigenstate and a stationary and rotationalinvariant solution of the equations of motion. This procedure, however, cannot be followed in the case of isospinning baby-Skyrmion discussed in this paper. The Fermionic Hamiltonian here is timedependent (see below) and a decomposition of the type Ψ = ψ(x, φ)e −iϵt does not work to solve the equations of motion. A non-trivial time dependence in the ansatz (41) is indeed required to obtain a time-independent system of equations for the profile functions and thus a stationary solution. As a result, the equations (43)-(46) differ from the energy-Eigenvalue equations, as we show below.
Given the Fermion Lagrangian (2), the Fermionic Hamiltonian reads where the time-dependence is explicit. Since in this work we are considering a specific time-evolving system of the type (32) and (40)    Even if the Hamiltonian explicitly depends on time, at every given t it is always possible to define a basis of energy Eigenstates with the corresponding Eigenvalues, i.e.
(52)  Here, Ψ n (x, φ, t) is nth Eigenstate whose function is different at each time and ϵ n (t) is the corresponding Eigenvalue that, in principle, can depend on time too. The complete set of Eigenstates for H(t) of equation (52) has a precise physical meaning. Let us consider the operation of a measurement of energy for the time-evolving Fermionic state Ψ resulting from eqs. (42)-(46) at a given time t 0 . The Fermionic wave-function Ψ at a given time t 0 can be written as a superposition of energy-Eigenstates Ψ n (x, φ, t 0 ) of the Hamiltonian H(t = t 0 ) as Measuring the energy of the fermion at the time t 0 , the wave-function Ψ will collapse to a given energy-Eigenstate Ψ n (x, φ, t 0 ) with probability | ⟨Ψ n (t 0 )|Ψ(t 0 )⟩ | 2 , giving thus a measured energy of value ϵ n (t 0 ).
In order to calculate the Eigenvalues of the Hamiltonian (51), we propose the following ansatz Inserting (54) and (55) into the equation (52), we finally obtain the system: Note that, due to the lack of t-dependence, the Eigenvalues ϵ n do not depend on time, i.e. ϵ n (t) = ϵ n . Once the Eigenvalue equations are written, we finally compare them with the Fermionic equations (43)-(46) of the time-evolution system for Ψ. It is now evident that, being (1 + k 1 )ω ̸ = k 1 ω for every |ω| > 0, there is no possibility for the stationary and rotational-invariant solution Ψ to be simultaneously an Eigenstate of the energy. On the contrary, as already anticipated in the previous section, only in the limit ω → 0, k 1 → ∞ and ωk 1 = const, the time-evolving state Ψ is effectively an energy-Eigenstate with Eigenvalue ϵ n = −ωk 1 .
The ansatz (54) and (55) simplifies the form of eq. (53) and allows us to see explicitly the stationary fermionic solution (40)-(41) written as a superposition of energy-Eigenstates. In fact, at a given time t 0 , the fermionic solution Ψ can be written as where the coefficient c n is given by withũ 1 ,ũ 2 ,ṽ 1 andṽ 2 the (normalized) profile functions of the nth Eigenstates, whereas u 1 , u 2 , v 1 and v 2 are the (normalized) profile functions for the stationary state Ψ. The value |c n | 2 represents the probability of finding the stationary solution Ψ in the nth Eigenstate once the energy is measured. As discussed before, in the limit ω → 0, k 1 → ∞ and ωk 1 = const, the equations for the time-evolving state Ψ become equal to the Eigenstate-equations for a given mth energy-mode. Therefore, in this limit we reduce to c n = 1 for n = m, c n = 0 for n ̸ = m and also Ψ n (x, φ, t) → Ψ n (x, φ) (being ω → 0). The equation (60) becomes then that is the time-evolution of the mth energy-Eigenstate with Eigenvalue ϵ m = −ωk 1 , as said before. We emphasize again how eq. (62) is exactly the time-evolving fermion ansatz proposed in [22]. This limit will become apparent when comparing the mean energies of our solutions to the energy Eigenstates, in section 3.5 below.
We can now explicitly solve the equation (56)-(59) and thus calculate the energy Eigenstates and Eigenvalues of the Hamiltonian. Although the system (56)-(59) appears functionally the same of [22], we note that in our case the profile-function F = F (x) comes from the solution of equations (42)-(46) in which, in addiction to the Fermion backreaction, also the isospin rotation plays a role in deforming the baby-Skyrmion profile. We expect, therefore, for |ω| > 0 to obtain a deformed spectral flow for the Eigenvalues ϵ n compared to [22].
To perform an explicit calculation, it is convenient to define the total angular momentum operator that counts together the spatial angular momentum, the spin and isospin. This operator commutes with the Hamiltonian (51) and thus it is possible to write a set of common Eigenstates. The ansatz (54) with (55) has exactly this property sincê The total-momentum Eigenvalue κ can be used to classify the energy-Eigenstates together with the index n of ϵ n . Note that the operatorK 3 acts in the same way on the time-evolving state Ψ of ansatz (40)-(41). This means that, although Ψ is not an energy-Eigenstate, it still represents an Eigenstate of the total angular momentum.
To make a better comparison with [22], we restrict our analysis about the energy-Eigenstates to the case of κ = 0, and thus l = 0 and p = −1. With this parameters, the four differential equations reduce to two, due to the relationũ 2 =ṽ 1 andũ 1 = −ṽ 2 .
In figure 6, we plot the value of ϵ/|g| as a function of |g| with g = −|g| for four different modes, called A i and B i corresponding to negative or positive energy Eigenvalues, for two values of ω. Note that the mode A 0 for a particular choice of g becomes a zero-mode, reproducing a well-known result from literature [22]. It is clear that the presence of an isospin ω > 0 deforms the spectral flow, as expected. The presence (and number) of the zero modes is not affected by the iso-spin deformation, in general this should be protected by the index theorem, however the precise flow of the zero mode changes with the isospin parameter. In general, it appears switching on an isospin parameter shifts the solutions towards lower values of g (the shift is small given the strict constraints on the possible values of ω). In figures 7, we plot some zoomed versions of figure 6 in order to better appreciate the isospin deformation of the spectral flow.
A physical interpretation of the results of figure 6 turns out to be more subtle then one may expect. Focusing for example on the positive-energy part of the spectrum, figure 7 (a), we recognize that for g ⪆ 1.8 the introduction of the isospin rotation has the effect of increasing the energy Eigenvalue, as we may expect. However, we note the opposite behaviour for low values of g. This result is probably due to a complicated mix of two phenomena: if from one side the isospin rotation generates an increase of the baby-Skyrmion energy and therefore a probable energy-raise of the whole system, it has also the effect of non-trivially deforming the baby-Skyrmion profile and thus the shape of the "potential" in which the Fermion is trapped. Therefore, in order to give a complete explanation of the effect of the isospin rotation on the Fermion Eigenvalues, further specific studies are needed. We leave this for future work.
In figure 8, we show the profile-function forũ 1 andũ 2 (withṽ 1 =ũ 2 andṽ 2 =−ũ 1 ) for the the mode A 0 at two different values of ω. Together with the profile functions, we plot also the baby-skyrmion profile and the fermionic probability density ρ F of eq. (30). In this case, the profile functions of the Eigenstate A 0 have zero nodes and the probability density ρ F is concentrated around the origin of the plane. It is interesting to notice how the presence of the iso-spin rotation deforms the baby-skyrmion profile F (x) (that we remind it is a solution of the equations of motion (42)-(46)). The deformation of the profile function F modifies the eigen-state equations (56)-(59) and, in the end, the profile-functioñ u 1 ,ũ 2 ,ṽ 1 andṽ 2 . In figure 9, we show the same panel for the mode A 1 at two different values of ω. In this case, the fermionic profile functions have one node and again the probability density ρ F is mainly concentrated around the origin.

Energy and currents
The energy density ρ s for the baby-Skyrmion is given by the Skyrmionic component of energy tensor T µν , this is   With the ansatz given by (32), the baby skyrmion energy density reads On the other side, the fermion energy E f can be defined only as the mean value of the Hamiltonian on the time-evolving state Ψ In figure 10, we show the plots for E f , the fermionic energy mean value, as a function of ω for low and large fixed g. We see that the solutions that we found have negative mean-energy-values, meaning they correspond (at least mainly) to superposition of modes of the negative part of the energy spectrum (these were called A modes in [22] and above).
We notice that in 10-right panel all the points belong to a range for which k 1 ≫ 1. This means that, for this range of ω, we are close to the limit in which these solutions are also Eigenstates of the Hamiltonian. More precisely, the different solutions of Ψ consist of superpositions of energy Eigenstates well-peaked on a single mode (A 1 , A 2 , etc.). To emphasize this behaviour, we draw the lines corresponding to the energy-eigenvalues of the modes A 2 and A 3 (the mode A 3 was not considered in the spectral flow in figure 6). As expected, in the regime k 1 ≫ 1, the various values of E f are close to be energy-eigenvalues. This behaviour is even more evident in the limit of ω vanishing (with k 1 becoming very high) that corresponds to the case analysed in [22]. This corresponds to the leftmost points in the right panel of 10, in which we can see that the solution mean energies are coincident with the Eigenvalues. As one increases ω, one proceeds in the direction of solutions no longer being Eigenstates, therefore the solutions energies no longer coincide (this can be seen in the rightmost points of the right panel of Figure 10). The leftmost panel of 10 shows instead solutions with low value of k 1 . For all ranges of ω, and for the reasons discussed above, these solutions are far from being coincident with any Eigenstate of the Hamiltonian. Having finite mean energy means that they can be decomposed in sums of energy Eigenstates, however this decomposition is not a-priori obvious and therefore one cannot assign a specific energy Eigenstate that is close to these solutions.
Another important comment concerning figure 10-right panel regards the jump in the value of E f at a given point of ω/g. Studying the shapes of the solutions, we notice that the sudden change in the value of E f is accompanied by the addition of a node in the profile functions u 1 , u 2 , v 1 and v 2 of Ψ. This non-smooth modification with additional presence of nodes in the shape of our solutions for different values of ω and g has been already shown in figures 1, 2, 4 and 5, where we solved the equations of motion for different parameters ω (or k 1 ) and g. This is expected as different energy Eigenstates of the Hamiltonian can correspond to different numbers of nodes, and thus different superpositions of Eigenstates can possess different numbers of nodes too. The fact that in figure 10-left panel we do not see any jump coincides with all the solutions corresponding to the different points having the same number of nodes. There remains however the question of why this sudden jump happens in the first place. It is perfectly reasonable to expect that as one adiabatically varies the value of ω/g one would stay within the space of solutions with equal numbers of nodes. The sudden jump seems therefore to be related to a global minimization of the energy, as a function of ω/g. We attempted to remain within the same node sector for increasing values of ω/g by varying this vary slowly, however our numerical procedure always jumps to the new set of solutions.    (20) that results is Finally, we write explicitly the form of the Fermion probability density and the probability density current defined in (29)  In figures 11, 12 and 13, we show the plots of the baby-skyrmion energy (67) and the effective J defined in equation (35) for different values of the coupling g. In all these plots, we compare our results with the so-called rigid rotor approximation, i.e. the approximation in which the energy (and angular momentum) of the isospinning baby skyrmion is calculated using the soliton-profile undeformed by the iso-rotation. In that approximation, the energy and the angular momentum would increase quadratically and linearly in ω, as shown in eqs. (67) and (35). As shown in figures 11, 12 and 13, the approximation works only for small values of ω, as expected.

Conclusions
In this paper we investigated the coupling between baby-Skyrmions and Fermions in the presence of internal isospin. We restricted our analysis to the case of rotationally symmetric B = 1 solutions and included the Fermion backreaction on the baby-Skyrme solution. We found that rotational symmetry  is consistent with a Fermionic ansatz in which Fermionic solutions to the Fermion-baby-Skyrmion equations of motion are not generically also Eigenstates of the Fermion Hamiltonian. These solutions do become Eigenstates in a particular limit in which ω → 0 but −k 1 ω = ϵ is held fixed. In general, however, we also solved for the Eigenstates of the time-dependent Hamiltonian. Solutions of the equations of motion can be written as superposition of the energy-Eigenstates and thus, once one measures the energy of the system, the fermionic wave-function will collapse to one of these Eigenstates with a given probability. As discussed in this paper, rotationally invariant isospinning baby-Skyrmions are known to suffer from strict stability constraints on the possible values of the isospin parameter ω. We found that in order to localise Fermionic solutions on them, a further constraint on ω must be imposed, coming from the behaviour of the Fermion functions at large radial values. Within these constraints, we found localised solutions for the Fermion coupled to the baby-Skyrmion, for small and large values of the coupling parameter g. As expected from previous studies [22], we find a tower of energy states (both positive and negative) and a single zero mode whose precise location on the g line depends on the isospin parameter ω. We found numerically that increasing ω has the effect of bringing the zero mode towards lower values of g. Furthermore, plots of the mean energy of our Fermionic solutions (see figure  10) show that, in the case of high g, adiabatic variations of the isospin parameter lead to jumps in the energy mean value corresponding to rapid changes in the types of solutions found, characterised by different number of nodes they have in the radial direction. It is well-known that higher winding number rotationally invariant isospinning baby-Skyrmions are in general unstable towards decay into lower winding number non-rotationally invariant configurations [14]. It is an interesting avenue of research to relax the rotationally invariant constraint on the system imposed in this paper and solve the full 2D PDE system of equations to investigate if and how the isospinning localised Fermionic solutions stabilise the Skyrmion solution.