{'fraction_non_alphanumeric': 0.029337486788674417, 'fraction_numerical': 0.006830950659175613, 'mean_word_length': 4.645396306996608, 'pattern_counts': {'":': 0, '<': 2, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 8, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Quantum Field Theory and the Internal States of Elementary Particles

A new application of quantum field theory is developed that gives a description of the internal dynamics of dressed elementary particles and predicts their masses. The fermionic and bosonic quantum fields are treated as interdependent fields satisfying coupled quantum field equations, all expressed at the same space-time coordinate. Quantization is realized by expanding the quantum fields in terms of fermionic creation and annihilation operators. This approach is applied in a QCD description of the light quarks with a zero Higgs field. Originally massless and pointlike, an isolated quark (described in its own center-of-mass) acquires mass and a finite extent when treated as an interacting system of quark and gluon fields. The binding mechanism of this localized system has a topological character, being a consequence of the non-linear nature of QCD, while being insensitive to the magnitude of the coupling constant to lowest order. To prevent this system from collapsing general relativity is introduced. The quark stabilizes at a radius of 8.8 Planck lengths and acquires a mass of 3.2 MeV, in remarkable agreement with accepted phenomenological values. It is suggested that the two higher generations of quarks are associated with the other two real solutions of the Higgs field equations.


I. INTRODUCTION
The standard model has been extremely successful in describing scattering phenomena between elementary particles. Nonetheless, it is clear that this picture of Nature cannot be complete. Most importantly, the standard model contains too many free parameters to be called fundamental, and can thus best be seen as an effective theory. In modern treatments the masses of the standard model fermions are supposed to arise from the couplings of quarks to various Higgs fields, however, so far this has not lead to a reduction of the number of parameters. Of the 19 free parameters of the standard model, 13 belong to the Yukawa/Higgs sector [1]. These parameters are put in by hand, and there is no explanation for the huge hierarchy of masses ranging over 6 orders of magnitude (not including the neutrinos). Further extensions beyond the standard model often lead to a greater -rather than a smaller -set of unknown parameters.
Given this situation one may wonder whether it is not possible to underlie the standard model by a theory at a more fundamental level that evolves into the standard model through the construction of effective solutions of the bare degrees of freedom. Such an underlying theory should contain many of the elements of the standard model but have a simpler structure, i.e. it should embody a certain unification of the diverse elements of the effective theory. The simplicity should imply fewer degrees of freedom and fewer parameters at the fundamental level. The new application of QFT methods, introduced in this paper, enables such a link between fundamental and effective levels by its ability to construct massive dressed elementary particles from more fundamental bare massless particles.
An important element of such a fundamental theory could be the reduction of the three generations in the effective theory to a single one in the fundamental theory.
Since the three generations in the standard model are not associated with any additional symmetry or interaction, three generations may well be the result of multiple solutions at a more fundamental level. This possible scenario would gain in standing if we could explain the emergence of (exactly) three generations from the exact level. We suggest that a single fundamental Higgs field might explain the threefold split. The Higgs field equations would typically have three real solutions, whose classical approximations could be denoted by φ = 0 and φ = ± |µ 2 /λ|. In the quark-Higgs dynamics the positive and negative solutions lead to distinct solutions, so that these three solutions would indeed lead to three distinct fermion sectors. Since the nature of the Higgs Lagrangian at the fundamental level is not (yet) unambiguously established and a non-zero Higgs field would complicate our basic model considerably, we limit ourselves in this paper to the trivial Higgs solution, namely φ = 0. However, this also eliminates the Higgs parameter µ from the model, so that it is unclear how this basic theory acquires a scale. The question of the basic scale parameters in Nature has been considered in the context of cosmology elsewhere [2]. We find that for the light quarks general relativity has to be introduced to ensure the existence of the quarks, so that for this basic generation the cosmological parameters reviewed in [2] also play a role. The presence of different scale parameters at this fundamental level could give a possible explanation of the hierarchy problem in physics, as it would allow for the existence of a range of different scales at a very basic level.
Which further ingredients would survive in a fundamental theory underlying the standard model? It is to be expected that the SU(3) theory of strong interactions will survive, as this theory is unbroken at the standard model level and is not characterized by a complex multiplet structure. The complexities of the multiplet struc-ture and broken symmetries in the standard model originate from the electro-weak interaction and it is not clear which of these complexities emerge in the transition from fundamental to effective theory, or which ones are already present at the most fundamental level. Therefore it is natural to limit oneself originally to QCD. Clearly, so long as one does not consider the electro-weak interactions one cannot describe leptons as effective systems, so we also focus on quarks in this paper. Being a nonlinear theory, one expects that QCD can supply the necessary self-consistent binding mechanism for the effective quark system, recalling the important role non-linear dynamics plays in soliton particle models in hadron physics [3]. By adding the electromagnetic interaction at a later stage one can break the isospin symmetry of QCD and try to explain the mass differences between up and down quarks. Being a linear theory, the electromagnetic force can be treated perturbatively and is not expected to interfere with the expected non-linear binding mechanism, so that it can be ignored initially.
As stated before, the possibility to link the dressed particles and their internal properties and masses to more fundamental bare particles relies on a new methodology which emphasizes the use of the field equations in QFT. The standard applications of quantum field theory (QFT) deal with scattering processes expressed in term of propagators and vertices, yielding the well-known Feynman diagrams. The dressing of fermions is then expressed in a series of time-ordered interaction diagrams modifying the bare propagator. Most of these diagrams are divergent, which has necessitated the renormalization methods to arrive at meaningful physical results, although outside inputs are still required to fix the physical parameters. This dressing process does not yield any internal properties of the particles, due to the fact that the fermions interact with external bosonic fields. The dressing process which is introduced in this paper is of a very different nature as the fermionic and bosonic fields are treated as mutually dependent fields, the dependence being controlled by the coupled field equations. Hence, the boson fields are not treated as external fields and all fields and interactions refer to the same space-time coordinate. It is this treatment that enables one to describe the internal properties of elementary particles.
It is clear that the field equations play a dominant role in the new formulation. This is not unlike most other areas of physics. However, in QFT their use has been very limited to date. Often the field equations in QFT are denoted as classical, and indeed they are often used as such. For example, the scalar Higgs field φ satisfies a classical field equation with various (constant) classical solutions. The "full" solution is then expanded around this classical solution, and only the additional fields are quantized. Despite their classical nature the constant solutions are attributed physical meaning, and their value is often denoted as the vacuum expectation value of the Higgs field [1]. This is a dubious use of the term expectation value as it suggests that one takes the expectation value of an operator, which one does not. It is not sufficient to call something an operator, one also has to treat it as such. The expansion around classical solutions is also common in the soliton and instanton problem [3]. Here the basic solution also satisfies an approximate non-linear classical equation of motion, and quantization is imposed by linearizing the remainder. We feel that these procedures are fundamentally flawed. If classical solutions have any role to play in QFT then they should emerge as an approximation to the quantum operators, not lie at the basis of the quantization procedure. Quantization should be done up front, not as an afterthought. Hence, we promote a much more rigorous quantum-mechanical use of the field equations. Before solving the field equations one should turn them into quantum equations and ensure that the solutions are quantum field operators.
A rigorous quantization procedure can be realized by turning the classical fields into quantum fields by expanding them into creation and annihilation operators, the latter satisfying well-defined (anti-) commutation rules. The insertion of quantum operators in the field equations also means that the order of the fields in the field equations can make an important difference, a phenomenon well-known in most quantization procedures. Hence, before solving the field equations one must establish the correct order of the fields. This is usually possible by demanding consistency, by imposing symmetries or through other physical requirements. At this point it is of interest to review the standard quantization procedure for a linear field theory such as QED, since this may elucidate why classical solutions have played such a prominent role in traditional quantization procedures, a role which should be avoided when one deals with non-linear equations. To derive the electromagnetic field operator in QED one first determines the classical solutions of the classical linear equation (plane waves in the free case). One then quantizes the field by appending creation and annihilation operators to the classical solutions. This well-known procedure leads to the correct quantized electromagnetic field, and therefore may well be responsible for the practice to use classical solutions in the quantization procedure. However, we would have reached the same result by inverting the order, namely, by first expanding the field in terms of creation and annihilation operators of (as yet) unknown functions, and then solve the (quantum) equations of motion. This procedure would yield the same plane wave states. However, by choosing this latter order one guarantees up front that the solutions are quantum solutions. For non-linear field equations it makes a huge difference which order is applied. If one solves the non-linear equations classically, and then appends creation and annihilation operators to these solutions, then the resulting operator is generally no longer a solution of the quantum field equations, as the operators do not behave like c-numbers. Hence, one must first introduce the field operators and then construct operator solutions to the field equations. If the resulting profile functions, which multiply the creation or annihi-lation operators, are identical -or close -to the classical solutions, then this could justify a classical approximation. However, such a situation only arises if the resulting field operator acts like a unit operator in operator space. In general this situation does not apply and classical solutions cannot be expected to play a fundamental role in non-linear field theories.
As indicated above, creation and annihilation operators play a central role in the quantization procedure. The operator nature of the field now even becomes part of the solution process and in general the simple linear expansion of the fermion and boson fields in terms of creation and annihilation operators is no longer a solution of the full set of field equations. One can start with the simple linear expansion of the fermion field, however, after further iterations a much more complex structure will evolve. For example, to lowest order the gluon operator becomes a bilinear operator in terms of the fermionic creation and annihilation operators. This bilinear expansion of the gluon field in terms of fermionic operators, instead of in single bosonic creation and annihilation operators, does not imply that gluon fields have lost their identity. In fact a similar expansion of the gluon field (or the photon field for that matter) in terms of quark operators is possible when it is considered an external field (and thus fully independent field) and satisfies the homogeneous gluon field equations. So such an expansion in fermionic degrees of freedom is perfectly consistent with the idea of independent gluons, soft gluons and glueballs. However, in the self-consistent case, which we are considering here, the gluon and quark fields are interdependent and neither of these fields can be considered as independent, although it still makes sense to talk about the fermion and gluon degrees of freedom. One should also not interpret this formal expansion of the gluon field as a composite model of the gluon. It is known that the quark-antiquark interaction is repulsive in the color octet gluon channel, however, forces play no role in the operator expansion and this fact plays no role in the expansion.
By further iteration one can construct the whole operator structure of the exact quantum fields. Ultimately one arrives at an infinite sum of operator terms, each additional term containing a higher power of creation and annihilation operators. Despite the elaborate operator structure of the fields and the equations it is possible to construct a closed operator solution. In Section II and Appendix A, this complete operator structure of the fields is derived and is expressed in an elegant concise form, involving a vacuum projection operator. By factoring out the quantum operators the quantum field equations can then be reduced to ordinary coupled equations, similar in structure -but not in content -to the classical equations. The elegance of this reduction is further underlined by the fact that the number of profile functions remains finite, despite the infinite nature of the expansion. The dependence on spin, color and isospin can also be factored out, leaving in the end a set of coupled differential equations for c-number profile functions.
The number of profile functions depends on the symmetries present: for our QCD formulation there are only 5 gluon profile functions. But if symmetries are broken (like isospin in QED), more profile functions are needed. One of the great unexpected benefits of this quantum reduction is that the resulting non-linear equations are more amenable to analytic solutions than the original classical field equations.
After further extensive manipulations of the scalar coupled equations, two simple order O(α 0 s ) solutions emerge. The first one is an uninteresting "trivial" solution, in which all effective potentials vanish. The second one features non-zero potentials that become singular at a finite radius r 0 . This absolute confining mechanism of the quark within its self-generated bag is of great beauty and is made possible by a unique QFT mechanism not seen before in particle theory. Since this solution is enabled by the non-linear nature of the field equations and is independent of the coupling constant, it will be denoted as a topological binding mechanism. However, this solution still features some other problems. First, in the absence of a Higgs field (we set φ = 0) this QCD model does not contain a scale, so there is no indication whether this model refers to a real particle. A second problem is that the expectation value of the energy operator is negative, so that this energy cannot be identified with the mass of the quark. In addition, a system with negative energy will collapse to a point as this minimizes its energy. These problems have an amazing solution. By linking QCD to general relativity (GR), we can exploit the fact that GR contains a definite scale and that gravity will act as a repulsive force if the energy of the system is negative. Hence, after the system has contracted to the Planck scale the repulsive gravitational force becomes big enough to halt the collapse. By adding vacuum energy [4] to the system we can compensate for the negative internal energy, and eventually obtain a net positive energy that can be identified with the mass of the system. Amazingly, this mass comes out to be slightly over 3 MeV, which is close to the phenomenological values found for the light quark masses. Hence, this solution can be identified with the light quark doublet. This solution also resolves the apparent conflict between the smallness of the quark mass and its pointlike nature. Normally, a mass of 3 MeV would correspond to a ridiculous size of 65 fm, but our result of 8.8 Planck lengths lies well within observable limits. Now that we have explained the basic elements in our formulation of the internal bound state dynamics of QFT, it is appropriate to mention earlier approaches in QFT that tried to describe localized states. Soon after the development of the Dirac equation, this equation was used to carry out bound-state calculations, for example by deriving relativistic corrections to the non-relativistic hydrogen atom calculations. The success of these calculations already suggested that QFT can be used fruitfully in bound-state problems, however, the nature of these Dirac calculations (the proton is not treated as a quan-tum field) shows that such a calculation cannot really be called a genuine QFT calculation. Historically, the electron also has been subjected to bound-state calculations of quantum nature. Schweber [5] reviews efforts to describe the electron as a bound state in the period 1940-1955. Dirac [6] made an attempt to describe an electron as a particle with a finite, perfectly conducting, surface in 1962. However, these efforts were hampered by the fact that the electroweak and strong gauge theories had not yet been developed at the time, while one can expect that the non-linearity of these non-Abelean gauge theories is essential for the construction of boundstate solutions. Other methods have been developed to handle bound-state problems within standard QFT, for example, the Bethe-Salpeter equations ( [7], [8], [9], [10], [11]), the Blankenbecler-Sugar equations [12], and the Dyson-Schwinger equations [13]. However, none of these theories could describe the internal properties of the elementary particles as the bosonic fields were treated as external fields. More recently, lattice gauge calculations have been applied to nucleons considered as assemblies of quarks [14]. However, these calculations are also not designed to give insights in the quarks and leptons themselves.
Historically, our efforts to solve the non-linear QCD field equations self-consistently started out as an attempt to describe the binding of three quarks in a bag, expecting to develop an alternative to the bag models popular in the eighties [15]. Solutions were found fairly quickly [16], however, their physical interpretation left much to be desired as the mass/energy scale of the light quarks (a few MeV) corresponds to very extensive objects, which are physically unacceptable. Also, only particle states were considered initially, as was common in the field of intermediate energy physics [17] and other papers dealing with pion corrections to quark bag models [15]. It soon became clear that it is difficult to apply such an approach without anti-particles consistently, and this initial approach was replaced by a formulation in which all fields are expanded in particle and anti-particle creation and annihilation operators. When this formulation lead to the exact operator solution it also became clear that this approach cannot be used for composite systems such as three quark bags, as the operator solutions only survive if they operate on one-body state vectors. Hence, the attention switched from describing composite systems to elementary particles and this led to the realization that these new QFT methods can be used to give a description of the internal dynamics of elementary particles, in particular quarks.
The outline of the paper is as follows. In Section II we discuss the equations of motion and derive the form of the field operators in terms of quark creation and annihilation operators. The formal solution of these operator equations is presented, allowing the reduction of the field equations to ordinary differential equations. After defining the c-number wave function components and c-number gluon profile functions in Section III, we de-scribe the form of the ordinary second-order differential equations, controlling the behaviour of the gluon profile functions in Section IV. These equations are highly nonlinear, and at first sight look fairly untractable. However, after introducing a reduction process in Section V, we can cast these equations in an elegant symmetrical form, which allows simple solutions to zeroth order in α s . The linear equations for O(α s )-corrections are also derived. In Section VI we discuss the (Dirac-like) field equation for the quark and the form of the potentials. In Section VII a self-consistent solution of the coupled quark-gluon equations is presented. In Section VIII we discuss the evaluation of expectation values for the energy, spin and colour of the quark system. These allow us to define sufficient boundary conditions to determine the parameters of the solutions uniquely. In Section IX we discuss the necessity to include coupling to GR for the trivial Higgs sector. By incorporating vacuum energy as well, one is able to make precise estimates of the quark mass and radius. In the final section X we summarize the significance of these results and discuss a possible steps for deriving the internal properties and masses of other particles in the standard model.

II. THE QUANTIZATION OF THE EQUATIONS OF MOTION IN QUANTUM CHROMODYNAMICS (QCD)
As stated before we start from a basic QCD Lagrangian [18] without electro-weak and Higgs interactions: where the covariant derivative is defined as follows: Here the λ a are the Gell-Mann SU (3)-matrices. The inproducts run over 8 components: The gluon field tensor is given by the expression: Here, the SU (3) vector product is defined by: where the f abc are the SU (3) structure constants. The classical equations of motion for the gluon fields are given by: while the equation of motion for the quark field is: In the standard applications of QFT, the fermion and boson fields are expanded in the eigenstates of the free Hamiltonian, while the interaction Hamiltonian is used to calculate the scattering diagrams in the perturbative Feynman series. Its main application lies in scattering problems, and the fermion masses are put in by hand. In the current application of QFT we start from a more basic level where the bare fermions have zero masses (i.e. we set M = 0), and we try to solve the full set of field equations self-consistently for a particular (quark) state vector. The quark is described in its own center-of-mass system (the origin of the coordinate system), and characterized by localized profile functions representing either the gluon fields or internal quark wave function. The actual position of the system that describes the effective particle is thus irrelevant and only becomes of relevance when we enter the quark into the multi-particle environment and consider scattering processes.
Quantization of the formulation is accomplished by expanding the fields in terms of creation and annihilation operators, the latter satisfying certain (anti-) commutation relations. Initially the following expansion of the quark fields is used: where b α is the quark annihilation operator, d † α the antiquark creation operator, and the profile functions φ α (x) and φ a α (x) loosely refer to the particle and anti-particle wave functions of the quark in its own center-of mass system (i.e. the internal wave functions). These wave functions are not known beforehand, and even their existence is uncertain, since it is not known whether a localized system description of the quark state b † α |0 will emerge from the formulation. Since the field equations are solved self-consistently, only one space-time coordinate x enters the considerations, so transitions between states are not associated with different times, and the overall state vector b † α |0 stays the same. Hence, the creation and annihilation operators are also time independent and do not evolve over time. The role of these creation and annihilation operators in the current formulation should be distinguished from that in scattering problems, where they are often associated with a particular fixed functional behavior (a plane wave corresponding to a particular positive or negative hyperboloid). When scattering states undergo time evolution through interactions they change their functional behavior, and creation and annihilation operators defined in this way are mixed up, allowing even the mixing of creation and annihilation operators (Bogoliubov mixing [19]). In our formulation the creation and annihilation operators refer to the actual physical state, represented by b † α |0 . Isolated particles do not change over time (unless they are unstable) and their state does not depend on where they are located, so b † α has no space-time dependence. Only, when we use an external representation of the particles and specify their momentum and/or positions with respect to other particles in scattering processes, could we contemplate time dependent creation and annihilation operators to represent the time evolution of the world state vector. However, even in that case a time independent representation of these operators seems more appropriate in QFT, as the continuous evolution of the physical system between transitions is well represented by the continuous profile functions. The evolution of one world state to another through transitions is of a stochastic nature, and should not be represented by changing creation and annihilation operators. Rather transition operators (e.g. Feynman diagrams) should relate initial and final states uniquely defined in terms of fully specified creation and annihilation operators, i.e. time-dependent operators. Clearly, in our formulation no Bogoliubov mixing occurs, and particle and anti-particle components maintain their unique meaning. This is also important for the later developments in the theory, when we make an important distinction in the algebraic properties of particle and anti-particle operators to restore the symmetry between these entities which is broken in the Lagrangian formulation. In our operator formulation the vacuum also maintains its unique identity as the state without any particles present.
We now continue the discussion of Eq. (8). The index α covers the usual quantum numbers of the standard model: spin, color and charge (also called flavor, but we only consider one generation of quarks). Since we limit ourselves to QCD, we can treat the pair of quark states in each generation as a degenerate isospin doublet. The family/generation quantum number is supposed to be explained by the multiplicity of the solutions and is therefore absent. Both isospin and total angular momentum are good quantum numbers in the present formulation, so that we can limit the wave function space accordingly in the description of the quark state. Since the gluon fields carry color and spin, color and spin transitions can take place, and an expansion in terms of a complete set of color and spin states is required. No expansions in other quantum numbers, such as a principal quantum number, is appropriate since if a multiple of such states exist they should follow as a solution from the given equations, rather than appear in the expansion.
If we use the expansion, Eq.(8), in the gluon field equation Eq.(6), then the gluon field operators also become expansions in quark creation and annihilation operators. In first instance this expansion looks like: However, the reinsertion of this expression in the Dirac equation, Eq.(7), leads to higher order terms in the expansion of ψ(x), etc., which in turn leads to higher order terms in the expansion of the gluon field. It would thus seem that we have to revert to perturbative expansions, as is done in standard QFT. However, we will show below that the operator series can formally be summed to infinity, without increasing the number of profile functions.
To complete the quantization process we first recall the usual anti-commutation rules for the creation and annihilation operators: and while all other anti-commutators are zero. For free fields, when the fermion field is expanded in a complete set of plane waves, one can use these anti-commutation rules (i.e. their continuous generalizations) to derive the equaltime anti-commutator relations. For the current selfconsistent treatment of the full set of coupled equations, such relations are absent. This shows that the anticommutation rules Eqs. (10) and (11) are more fundamental than the equal-time anti-commutator relations. Historically one has often assumed the opposite, namely that the equal-time anti-commutator relations are fundamental to QFT (e.g. by stressing the analogy with the classical Poisson brackets), while Eqs. (10) and (11) were seen equivalent, but secondary. We have pointed out before that it is dangerous to rely on QFT on classical analogies, and this is a case in point.
Since the gluon field operator no longer commutes with the quark field operators, it is necessary to specify the order of the quark and gluon field operators in accordance with the symmetries and physical requirements of the theory. The non-commutativity of boson and fermion field operators is not new to QFT, for example Bjorken and Drell [20] discuss it in the context of QED. However, in the current context the consequences are must more farreaching. The anti-symmetric nature of the gluon field implies that we must replace Eq.(4) by: The quantized gluon field equations can be derived along similar lines and read: The quantized form of the Dirac equation is given by: The physical justification for the order of the ψ and A µ a operator will be given in Section VI. Before we can derive the full operator expansion of the quark and gluon fields we need to refine the algebra of creation and annihilation operators. Expressions such as ψλγ ν ψ display a disturbing asymmetry: the right-hand side operator ψ contains quark annihilation operators, but antiquark creation operators, implying that if such an operator is applied to the vacuum, anti-particles can be created, but particles can not. However, we know that there should be a fundamental symmetry between particles and anti-particles and that the decision to call a quark a particle or an anti-particle is largely arbitrary. So this asymmetric outcome is unacceptable. It does not help to invert the role of particles and anti-particles, as then the problem would be that particles can be formed from the vacuum and anti-particles can not. It seems that the symbolic mathematical language in which we formulate QFT is unable to adequately express the natural symmetry between particles and anti-particles. The solution to this problem is that we can maintain the same field theory language, as long as we invert the order of the anti-particle operators in the complete string of operators between the bra and ket states. This rule can also be seen as giving meaning to the physical role of anti-particles as being particles moving backwards in time. This elegant procedure restores the symmetry between particles and anti-particles. The procedure also eliminates the artifact that QFT gives rise to an enormous vacuum energy, which is contradicted by observation. This procedure will be indicated by the symbol R, and will be called the R-product. Hence, the R-product is a way to introduce ordering in expressions that are "timeless" and to enforce the property that particles and anti-particles have opposite ordering. We should caution, however, that this R-product only applies to quantum operators with identical space-time variable x (as is the case in our current application). Operators referring to different variables x and x ′ are not subject to this product, as their time-order can already be specified by the different times t and t ′ . To ensure the correct application of the R-product one could append a variable label to the operators, so that operators referring to different space-time variables are not confused with those having common variables. Since the product does not apply to operators at different points, it has no effect on most scattering series in QFT, which may explain why it has not been discovered before in the application of QFT. The product is however already used implicitly in standard applications of QFT, when the normal product is inserted to remove unphysical vacuum terms. We actually discovered this rule numerically before its theoretical origin was uncovered, as the coupled equations required unexplained minus signs, which could only be justified after the R-product was introduced. The detailed explanation of this product is given in Appendix A.
Using this R-product one can derive a formal solution of the quantized equations of motion (Appendix B). Instead of the expansions, Eqs. (8) and (9), we get the exact expressions: and where Σ ∞ is the result of taking the operator to the limit n → ∞. The operators N andN are given by: The operator Σ ∞ can be recognized as vacuum projection operator, as it only gives non-zero (unity) results if N = N = 0. As far as we know this operator has never before been encountered in QFT. It will play an important role in the solution of the self-consistent field equations. Since this operator is surrounded by creation and annihilation operators, it effectively acts like a one-body, rather than a vacuum, projection operator. The detailed derivation of these relationships is given in Appendix B. Since the space of internal states for the quark is finite (see Section III), this infinite expansion of Σ ∞ is terminated in practice after a finite number of terms. Clearly, this is not true if the state is characterized by a continuous degree of freedom, such as in the case of scattering, when we have to define N as an integral over all momentum states. An important property of this expansion is that the original number of profile functions in the approximate expressions, Eqs. (8) and (9), does not increase in the transition to the exact expressions (15) and (16). This property is largely responsible for the practical feasibility and the exactness of the current theory, and ultimately for the power of QFT.
We now want to illustrate the effective one-body nature of the projection operator. If one operates with ψ(x) on a one-body state one gets: however, for a many-body state one gets a zero result: and similarly for anti-particle states. This shows that the only fermionic solutions of the coupled equations in the current formulation are of one-particle nature (the term particle is used here as a generic term for particle or anti-particle). This clearly demonstrates that we cannot apply this formulation to a composite problem like a proton system, as this involves at least three quarks.
The only way to generalize this approach to composite systems is to collapse the three quarks into a single colorsinglet point-particle. However, this would eliminate the QCD interaction and thus would not yield a proper description of nucleons. For leptons this might be a possible model, however, this would require the introduction of the electro-weak interaction at this basic level, which we have not attempted yet. Obviously, for protons and other hadronic systems the standard lattice gauge calculations remain the best treatment [21]. Now that we have solved the quantum operator problem formally, we can proceed towards the determination of the profile functions. We will see that QFT has more surprises in store, and that the simplest solution for the single quark system is of great beauty, and realizes the localization of the bound system in an amazing way, never seen before in physical theories.

III. EXPRESSIONS FOR THE PROFILE FUNCTIONS
The next step after having removed the (creation-and annihilation-) operator structure in the field equations is to remove the color and spin dependence, by introducing scalar wave function components and scalar gluon profile functions. We start with the parametrization of the wave function of the 1S-quark state: where χ α is the spin wave function, η α is the isospin wave function, and ξ α is the colour wave function. Here we assumed that the radial wave functions f and g are independent of spin, isospin and colour. This charge independence is no longer valid if electromagnetic forces are introduced, and u-and d-quarks become distinct. The corresponding anti-quark states is given by: Notice the absence of phase factors. Usually, one defines the spin function for the anti-particle by expressions such as (−iσ 2 χ α ). Also, the colour and isospin wave functions for the anti-particle belong to adjoint representations of the relevant SU (3) and SU (2) representations. However, the particle and anti-particle creation and annihilation operators automatically take care of these phase factors. Hence, Eq.(23) does not correspond to the normal definition of the anti-particle wave function, but exploits the simplification enabled by the operator formalism to express the anti-particle wave functions in the same wave function basis as the particle states. The time dependence of these wave functions reflects the fact that we describe stationary states. Although, the quarks belonging to the higher generations will decay into the lighter quarks via several decay mechanisms, we expect that these higher generations are also characterized by a real positive eigenvalue E, and that the decay width is determined in the scattering problem and not in the current bound-state formalism.
The wave functions f and g are normalized according to: Using the wave functions Eq. (22) and Eq. (23), we can now evaluate the source terms in Eq. (13), and thereby determine the form of the different gluon field components in Eq. (16). We express the results in terms of five profile functions: A ap a,αβ (x) = The dummy isospin Kronecker δ-function has been suppressed throughout. This simplification becomes clearly invalid as soon as we introduce isospin breaking QED forces. In the next section we will derive equations for the five profile functions F 0 (r), F (r),F 0 (r), F 1 (r) and F 2 (r). The Dirac equation yields equations for the two functions f and g, and are discussed in Section VI.

IV. DERIVATION OF THE SCALAR GLUON FIELD EQUATIONS
After eliminating the common quantum operators one is left with equations in terms of the remaining SU (n) operators and functions. The SU (3) operators can easily be factored out using the relationship: The spin operators satisfy numerous spin identities, which can be used to remove their explicit occurrence as well. We are left with a set of scalar equations of motion for the gluon profile functions. To simplify our notation we introduce some short-hand notations for the source functions. We define: where α s is the strong coupling constant α s = g 2 s /4π. These source terms are related by the identity: We also introduce the auxiliary profile function K , which is going to play a central role in the dynamics: The differential equation for F 0 reads: As expected, quadratic and cubic terms play an important role in these equations. These terms express the non-linear nature of QCD, and the self-coupling of gluons. The solution of the linearized equation would lead to a Coulomb like potential. For the magnetic particle-particle component F we obtain: An earlier version of Eqs.(37) and (38), without particleanti-particle contributions, was presented in Ref. [16]. These equations did not lead to bag-like bound states. However they already displayed the very interesting nonlinear behaviour, suggesting soliton-like solutions. It has since become clear that the particle-anti-particle coupling is essential for the appearance of singular confining potentials.
The profile functionF 0 , corresponding to the gluon field operator A 0,ap a,αβ (x), satisfies the equation: Finally, the two profile functions F 1 and F 2 appearing in A pa a,αβ (x), satisfy the equations: and These 5 equations for the 5 profile functions do not look very tractable, however, after introducing a number of new auxiliary functions they take on a much more elegant form. We introduce the alternative profile functions: and Furthermore, it is convenient to introduce the following expressions: and We then can write Eqs.(37-40) as follows: This has not just simplified the equations enormously, but also has brought out a strong symmetry between F 3 andF 0 on the one hand, and F 1 and F 4 on the other. The equation for F 2 , Eq.(41), can be added to Eq.(40), to yield a very elegant algebraic expression for K: We can derive another important relationship by differentiating the numerator in Eq.(50), and inserting all the expressions for the second derivatives. We obtain: If we combine this equation with the Dirac equation, we can derive a simple relationship for the source functions, which is related to current conservation. In the next section we continue our treatment of the gluon profile functions exploiting the symmetry apparent in the equations (46-49). These steps are necessary to progress towards a physical solution, however, despite their elegance these steps are of necessity rather technical. Hence, those who want to immediately progress to the treatment of the quark field equations can skip to Section VI.

V. REDUCTION OF THE SCALAR GLUON FIELD EQUATIONS
We can exploit the symmetry of the equations (46-49) to partly decouple the equations. To this end we introduce further auxiliary functions: By adding and subtracting the equations (46-49) we find: We have partly decoupled the equations and can introduce further simplifications by introducing the expression: (58) Because of the symmetries of the equations one can also defineÊ with a minus sign. We can invert this relation by expressing K in terms ofÊ: Again this result is insensitive to the sign ofÊ. We now write: The new functionsX,Ŷ ,Û andV satisfy even simpler coupled equations, as the trailing Z dependent terms have vanished: whereŜ We can also express Eq.(50) in terms of these new functions, and obtain: We can easily solve Eqs.(61-64) to O(α 0 s ), by ignoring the source terms. By imposing the parity requirements of the potentials (see Section VI) we are lead to the following solutions: The solutions for α = 1 3 and α = − 1 3 appear to be physically equivalent, the former corresponding to the negative sign ofÊ. Hence, we will only consider the set with α = −1/3. We will see in Section VIII that the α = 0 solutions lead to infinite energy integrals, and must be excluded. For α = −1/3 the |β| = 2/3 solutions are of special interest. For β = 2/3 we obtain a trivial solution where all potentials are zero (F 0 = F =F 0 = F 1 = 0; K = 0;Ê = 1). For β = −2/3, K -and thusÊ -is singular at a finite radius r 0 (the bag solution). As a consequence, all original profile functions F 3 ,F 0 , F 4 and F 1 are also singular at r 0 . However, the reduced functionŝ X,Ŷ ,Û andV remain finite. Hence, our reduction of the equations in this section does not just lead to a simplification of the equations, but it also leads to a systematic elimination of the singularities from the equations.
The solutions summarized in Eq.(67) are valid in the absence of source terms. The existence of solutions for α s ↓ 0 is an important consequence of the non-linear nature of the field equations in the current bound-state approach. It is an illustration of the power of non-linear Gauge theories in the construction of emergent complex structures. In the standard perturbative approach to QFT, a theory with α s ↓ 0 would make no sense, and there certainly would be no binding.
Although terms higher order in α s do not play an instrumental role in the binding mechanism, they will modify the detailed outcome of the theory. In the following we show how these terms can be treated perturbatively, expressing the corrections in terms of linear equations for the remaining profile functions. Combined with the calculation of observables in Section VIII, these linear equations provide a powerful tool in analyzing the physical solutions in more detail. The combination of exact O(α 0 s )-solutions, with perturbative O(α s )corrections might look somewhat heuristic, however, nonlinear equations require novel methods, and this is a very effective approach The natural place to start the perturbative scheme is the set of regularized equations Eqs.(61-64), which are defined in terms of non-singular hatted profile functions. However, the functionsX,Ŷ ,Û and V have non-zero base values. For a perturbative scheme it is preferable to start from zero functions. This can be done by reverting back to functions like F 0 , etc., but now defined in terms of the hatted functions. Unfortunately, this requires a new set of auxiliary functions. However, the symmetry will not escape the reader, and the scheme remains very elegant and natural. The new function are: For β = 2/3 and α = −1/3 one gets K = 0 andÊ = 1, and the old functions are recovered (i.e.F 0 → F 0 . Expressing the differential equations into the new auxiliary functions, we are led back to differential equations of the familiar form (46-49), but now without the K- where the reduced source functions are defined by: In terms of these new functions the O(α 0 s )-solutions Eq.(67) are given by: where α = 0, ± 1 3 as before. We now expand the differential equations to first order, setting α = −1/3. The O(α s ) linearized equations now read: These equations can be solved formally in terms of the source functions. For the two uncoupled equations we obtain:F where we assumed that the physical solutions have a finite spatial extent specified by r 0 . The parameter λ 0 could be combined with β into a single parameter. However, we prefer to specify the O(α 0 s )-solution by a fixed value of β, and treat the O(α s )-corrections via λ 0 . The other two equations can be decoupled. First we multiply Eq.(76) with r, and then differentiate the result twice. After eliminating the functionF 1 , we can express the result as follows: where P (r) is defined by: and Q(r) is defined in terms of source functions by its derivative: One can re-express this in the original source functions, and show that: where ω is the third integration constant. The solution of Eq.(81) can be expressed in terms of the spherical Bessel functions j 1 and n 1 : where a = 3βE and ν is the fourth -and final -free constant of the O(α s )-solutions. For the special case |β| = 2 3 , Q = ω to O(α s ). For small r a non-zero ω generates a singular component ln(r) in P (r) and a component r 2 ln(r) − r 2 inF 0 . The other functionF 1 can be deduced directly from Eq.(76): 2aF ′ 1 = −(rŜ 1 ) ′ − a 2 r 2 P (r) + Q(r) = = −a 2 r 2 P (r) + arS 2 (r) + ω.
We can get one expression without derivatives from Eq.(76): At this stage we have four free parameters: λ 0 , µ, ν and ω, with dimensions E, E 2 , 1 and E 2 , respectively. These parameters should all be of O(α s ). The imposition of such a constraint on "free" parameters seems contradictory, but the non-linear nature of the full equations limits such freedom in the full solution. We can enforce the O(α s ) constraint by absorbing the free parameters into the left-most integrals in Eqs.(79), (80) and (85), by allowing the upper integration boundary, which is now fixed at r 0 , to vary between 0 and r 0 . In this way the corresponding parameters are always of O(α s ), and also have the right order of magnitude. The need for such a prescription is a consequence of the hybrid nature of the solution scheme: one tries to built a perturbative solution on a non-linear O(α 0 s ) basic solution, which itself contains a free parameter β.
This completes the discussion of the equations for the gluon profile functions. Although we seem to have strayed a long way from the original gluon field tensor defined in Eq.(9), expressing this tensor in terms of the reduced profile functions is straightforward. In the next section we discuss the Dirac field equation, which will give further constraints on the profile and source functions. We stress again that all profile functions and quark wave function components have to be determined selfconsistently. Hence, although we made great strides in deriving simple equations for the gluon profile functions in this section, we need the Dirac equation to complete the solution.

VI. THE QUANTIZED DIRAC EQUATION FOR ISOLATED QUARKS
We have to quantize the classical Dirac equation, Eq.(7). Since both ψ and A µ a are operators expressed in quark creation and annihilation operators, quantization amounts to deciding in which order these operators have to appear in the quantized equation. As stated in Section II, the correct choice is ψA µ a , leading to equation (14). This choice guarantees that the quark interacts with the fields it generates itself. The opposite choice, A µ a ψ, does not yield a bound-state solution. For example, it gives rise to the operator b † α b β b γ , which vanishes when operating on a one-particle state b † δ |0 . Since twoand higher-body terms are eliminated by the "projection" operator Σ ∞ , none of the other terms survive either. On the other hand, ψA µ a , leads to the operator b γ b † α b β , which yields a finite result, when operating on a one-particle state b † δ |0 . Hence, the only way we can get binding potentials and bound normalized states is to choose the indicated order. Again we have to extract sets of terms from the equations with similar operator form. We then can eliminate these operators and are left with algebraic equations. In the reduction of the operator equations we have to take into account the R-product when manipulating the anti-particle terms. Ignoring this aspect would lead to erroneous signs in the reduced equations and inconsistencies.
We can express the reduced Dirac equations as a set of scalar equations for the large (f ) and small (g) wave function components. We find: and The bare quark mass M is zero in our application, however, we leave it here as it could be non-zero when a nonzero Higgs field is introduced. Our aim is to derive the masses of dressed quarks, so it would be inappropriate to introduce phenomenological masses at this point. The formulation gives rise to three types of potentials (if we introduce isospin breaking effects we would get a fourth type of potential, as is allowed by the symmetries of the Dirac equation. This fourth class is also of tensor type). The vector potential is given by: while the tensor potential is given by Finally, the scalar potential is exclusively given in terms of K(r) : In these expressions C is an SU (3) constant with the value 8/3. Each of the potentials receives contributions from the particle-anti-particle gluon profile functions. The latter terms are essential for the emerging elegant self-consistent solutions of the equations and thus indispensable for the existence of dressed quarks that are finite in extent. Although the particle-anti-particle gluon amplitudes fluctuate rapidly with time (they are proportional to exp 2iEt or exp −2iEt), this time dependence is cancelled out in the scalar Dirac equations. This is one of the many elegant features of the current self-consistent theory. The required behaviour of the potentials near r = 0 fixes the parametrization in Eq.(67). In lowest order O(α 0 s ) we demand that V and V s (which are even in r) are constant near r = 0, while the tensor potential (which is odd in r) behaves linearly in r for small r. This is consistent with the small r-behaviour of the odd wave function f (r) ∼ r and the even wave function g(r) ∼ r 2 . We will see later that in O(α s ) the lower component g(r) also has small r components r 2 ln(|r|).
By multiplying Eq.(88) with g(r), and Eq.(89) by f (r), and adding them together, we can derive the following relationship: Combining this with Eq.(51) we obtain: and Eq.(95) also follows immediately from current conservation: if we take the particle component on the right, and the anti-particle component on the left. Notice, that Eq.(96) is trivially satisfied for the particle-particle components. Hence, our derivations so far have been completely consistent with current conservation. We can derive another equation from Eqs.(88-89) in terms of source functions: We now derive the basic non-linear self-consistent solution of the QCD field equations by combining the conditions on the gluon and quark profile functions.

VII. THE BASIC SELF-CONSISTENT SOLUTIONS FOR THE DRESSED QUARK
We anticipate that if a bound solution exists for the dressed quark system, that it is bound by strong, possibly singular forces to ensure that its size is small enough to be consistent with the known properties of quarks. The profile function K may well play a prominent role in this regard, as it is given by a ratio (see Eq.(50)) and thus would become singular if the denominator turns zero. The function K stands alone amongst the other profile functions in being given by an algebraic form. This form is a unique consequence of the non-linear nature of the theory and the inclusion of both particle and anti-particle configurations. The reduced profile functionsF 0 ,F ,F 0 andF 1 satisfy regular second-order differential equations (75-78), and are unlikely to have singular or bag-like characteristics (we will show in the following that the reduced source functions in these equations are also devoid of singularities). We should note that singularities in K orÊ would not just be reflected in V s , but also in the other potentials, because these potentials are implicitly defined in terms ofÊ andÊ −1 .
To analyze the behaviour of K orÊ we cannot directly use Eq.(50), as we have no simple equations for the unreduced profile functions. In addition, the unreduced profile functions may well be singular. Rather, we start by casting Eq.(51) in reduced form: The fact that the new equation has the same structure as the unreduced equation illustrates the elegance of the reduction process. The advantage of the reduced form is that it is expressed in the regularized profile functions with their known regular O(α 0 s ) and O(α s )-solutions. By re-expressing Eq.(98) in terms of the original physical source functions, which must be normalizable and display regular properties, we can derive explicit expressions for E andÊ −1 : and where Here the right-most expressions are the O(α 0 s ) solutions. We can use Eq.(95) to rewrite the derivative term.
For β = 2/3 and the top sign in Eqs.(99, 100) we obtain E = 1 and K = 0, corresponding to the trivial solution with zero potentials. The mirror solution, with β = −2/3 and the lower (negative) sign in Eqs.(99, 100), yields the non-trivial bag solution. For both these solutions the potential singularities ofÊ orÊ −1 , which would be due to a zero of 3βEr ± 6α, cancel out. Notice also that the reduced source functions remain finite, even whenÊ or E −1 become singular when S + or S − passes through zero. Hence, the reduced gluon Eqs.(75-78) are devoid of singularities, and form a good starting point for calculating higher order effects.
Taking β = −2/3 we have: HereÊ is singular for f = −g, which serves to define the radius of the system. This condition coincides with the usual demand in the MIT bag model [22] at the bag surface. Using Eq. (59) we obtain: This result can also be obtained by explicit calculation from Eq.(50). Hence, the scalar potential equals: where C = 8/3 as before. Using Eq.(102) we can also derive the other O(α 0 s )-potentials: and If r = r 0 is the point where f = −g, then the quark system is confined to r < r 0 , as the potentials become infinite at r 0 . Hence, if this solution exists, then it is a bag solution with confining potentials Eqs.(104-106).
Can we construct a solution of the Dirac equations for the potentials Eqs.(104-106)? Clearly, this is a highly self-consistent problem, as the potentials are expressed in the wave functions. We consider the Dirac equations Eqs.(88-89) in more detail. First, we combine the vector and tensor terms. Because of the identity, Eq.(94), we can rewrite these equations as: and where Going back to the specific solution, we insert the potentials Eqs.(104-106) into Eq.(109). We then get: We can now eliminate f ′ and g ′ in this expression by means of the effective equations Eqs.
The wave functions are defined in the range: leading to a normalization constant of: These solutions allow us to calculate the scalar, vector and tensor potentials. In Fig. 1 we display these potentials. The origin of the length and energy units used in this figure will be explained in Section IX. Clearly, the vector and tensor potentials are more singular than the scalar one. In Fig. 2 we show the resulting wave functions, Eqs. (111-114), normalized over the interval [0, r 0 ]. These wave functions (in different units) also appear in the phenomenological MIT quark bag model [22], where they describe the wave functions of the three quarks in a nucleon. Hence, the current theory shows how such a bag model could possibly originate from fundamental QFT considerations, although the restriction of our theory to single-particle states forbids a direct application of our theory. Obviously, the context is also entirely different as the size of an MIT bag is about 1 fm, while our single quark bag has Planck length dimensions.
We derived the solutions for a zero bare quark mass. This is an illustration of the importance of self-consistency in these bound-state QFT solutions, and shows that any perturbative considerations may fail, since they even would suggest the wrong sign of the mass. It shows that the calculation of quark solutions for non-zero Higgs fields is feasible, although they have important further consequences, such as additional terms in the total energy expression.
The solutions obtained in this section also have to satisfy additional consistency constraints. If we describe an isolated quark by the state b † α 0 , then its isospin, color and spin are specified by the index α. However, these properties can also be defined as expectation values and expressed as volume integrals over the quantum fields. Naturally, the two possibilities should yield the same answer. For isospin this condition is trivially satisfied, however, for spin and color things are more complicated. These conditions will be discussed in Section VIII. Our O(α 0 s ) solution for β = −2/3 satisfies the resulting conditions exactly. Hence,we have found an exact solution of the equations of motion for α s = 0. From the analysis of scattering problems we know that QCD is an asymptotically free theory. Hence, it may be quite reasonable to set α s = 0 at the Planck scale. Since this is the scale that seems to be relevant for the light quarks (see Section IX) the α s = 0-case might a complete description of dressed quarks in the QCD context. However, in order to describe the mass splitting between the quarks of different charge we need to introduce QED. This theory is not asymptotically free, so that we cannot ignore this coupling and have to apply a perturbative approach in terms of non-zero coupling constants at the Planck scale, anyway. If we assume that α s = 0, then we have to consider higher order corrections, as well and impose the boundary conditions discussed in VIII. The color condition is the most stringent one and leads to three separate conditions in O(α s ). Hence, it effectively fixes the integration constants ν and ω, while putting indirect conditions on λ 0 , which in turn leads to O(α s ) corrections to the value of the radius x 0 . The value of µ is essentially fixed by the spin condition. Values for ν and ω for the bag solution are given in the next sectionVIII. In that section we will also derive the formula for the total (internal) energy of a quark system. The zeroth order result will turn out to be negative. This presents a serious problem if we want to identify this energy with the effective mass of the dressed quark. However, there are various reasons why the model cannot be complete. In the absence of the Higgs field there is no scale parameter in QCD, hence the scale of the quark system is undetermined and it is a mystery how the small quark mass of the light quarks can be reconciled with the effective pointlike nature of the quarks. In Section IX we suggest a novel way to complete the model of single quarks leading to a solution of these interconnected problems.

VIII. OBSERVABLES AND CONSERVED QUANTITIES
The symmetries of QFT determine the conserved quantities and observables. The invariance under phase transformations of the quark field operator leads to current conservation (cf. Eq.(96)). This relation implies the conservation of baryon number B: where the state vector is either a quark: or an anti-quark: The normalization in Eq.(115) is chosen in such a way that the quark has baryon number 1/3, and the antiquark −1/3. Notice that we need the R-product to obtain the latter result.
The invariance under translations leads to the canonical energy-momentum tensor: By adding a complete differential one can re-express this as the Gauge invariant expression ( [23]): where: so that the spatial integrals of Θ µ0 are conserved quantities, corresponding to total energy for µ = 0 and total momentum for µ = i. It is easy to quantize Eq.(119) by taking A µ a outside the quark current, symmetrizing the gluon term, and introducing the R-product.
The spin of the quark is related to the tensor: where We can show that: The spin of the quark can now be expressed by the following matrix element: The calculated expectation value should equal 1/2 in units .
We now express the conserved quantities in terms of the gluon profile functions and the quark wave functions or source functions. From Eq.(119) we get for the total internal energy of the quark system: (125) where r 0 is the radius of the quark system and I(r) is given by: After reduction this expression simplifies: The internal energy of the system is now expressed completely in terms of reduced profile functions and therefore is finite even if the original profile functions are singular. For the O(α 0 s ) bag solution with β = −2/3 we get E int = E − (32/9)E < 0. In Section IX we will address the problem associated with the negative character of this internal energy.
The total momentum of the system is easily calculated and equals zero, as expected: Since we are considering a system in isolation, it would not make sense to assign a momentum to it. Its momentum only becomes relevant when we embed it in ordinary QFT and the quark interacts with other particles. Let us now consider the spin. We obtain: The first term already provides the exact answer, hence the remaining terms should add up to zero. After reduction we can convert this to the condition: or after integration by parts: To O(α 0 s ) this equation is automatically satisfied by the solutionsX(r) =Ŷ (r) = βEr, e.g. by the bag solution with β = −2/3. Its satisfaction for this bag solution can also be related to the identityŜ 0 = S 0 , which guarantees the validity of Eq. (130) to lowest order. The condition (131) can assist in fixing the four free parameters defined in the perturbative scheme in Section V.
There is another conserved quantity for single quarks, namely color. The reason that this is not usually seen as an observable, is because quarks are usually not considered in isolation (physically they only appear in color singlets). However, since we are considering single quark solutions to the field equations, we clearly have to examine this property as well. Color is associated with the invariance of QCD under local gauge transformations. The relevant infinitesimal transformations are [18]: and The current corresponding to this transformation is given by: it follows immediately that: Hence, we can define the conserved color density (µ = 0): The first term in the right-hand side of Eq.(137) equals the bare color, so that we have to demand that the remaining terms yield zero. Since this leads to a constraint on an integral, this condition is not easy to implement. A simpler condition arises if we use the identity, Eq.(135), again. This leads to the condition: Writing this out in profile functions, we obtain: where r 0 is the radius of the quark bag for the bag solution. After reduction this condition can be written as: Notice the similarity with the spin condition Eq. (131).
As for the spin case, this condition is automatically satisfied for the O(α 0 s ) solutions Eq.(67). We can combine the color and spin condition to formulate the following boundary conditions: Since eitherÊ orÊ is singular at r = r 0 , we need to treat these boundary conditions as limiting processes, rather than as point identities. For the bag solution we have: near the surface x = x 0 . After expandingX in a Taylor expansion for r → r 0 , we then arrive at three conditions, which allow us to fix the perturbative parameters ν and ω. Numerically, we find: The third condition is more tricky as it involves an adjustment of the radius of the quark bag. A rough calculation indicates that we have to adjust x 0 from 2.043 to 2.129, however, a more complete numerical study is required to get precise answers.
In this section we considered the observables which play the role of boundary conditions in our formulation. Other observables, such as the magnetic moment, arise from the interactions with external fields, and must be calculated with standard perturbative methods. Clearly, our formulation has little to add to such calculations except perhaps providing physical cut-off parameters in view of the finite size of the dressed quarks.

IX. THE ROLE OF GENERAL RELATIVITY AND THE DETERMINATION OF THE MASS OF THE LIGHT QUARKS
As seen in the previous section, the total QFT energy of the dressed quark system is negative. Within the context of QFT there is only one obvious contribution which could compensate for this negative energy and that is the Higgs field. We expect that the two finite real Higgs solutions, whose classical approximations could be characterized by the values φ = ± √ µ 2 /λ, would have such an effect and lead to an overall positive energy, which can be associated with the mass of the quarks. However, to enable this treatment, one needs to know the exact nature of the Higgs Lagrangian at the fundamental level and preferably its relationship to the massive vector bosons. These Higgs terms would add an additional set of equations and profile functions and therefore add considerably complexity. It therefore seems more natural to look initially at the trivial Higgs sector with a zero Higgs field. Strictly speaking, this "trivial" Higgs sector is not completely zero, as the quark source term would induce a final Higgs field even if we choose the small Higgs solution in the cubic field equation. However, we believe that it would be a fairly good approximation to approximate this lower branch by a zero Higgs field. The consequence of this assumption is that we need to find other solutions for the following two problems (1) The absence of a scale defining parameter; (2) The negative energy of the system. A system with negative energy and unconstrained scale would contract without bound in order to minimize its energy, unless its scale is fixed, which it is not in QCD.
To prevent this collapse we consider the consequences of general relativity (GR). When the energy approaches the Planck scale, the effects of GR become important, and the large magnitude of the negative internal energy will halt the implosion. To treat these GR effects in a QFT context is non-trivial, and brings us in uncharted territory. However, by treating these effects perturbatively, we avoid some of the tricky problems. Since this inclusion of the effects of GR, in combination with the vacuum energy, leads to a spectacular agreement with experiment, we feel that there is good evidence for the validity of this perturbative approach.
The effect of GR can be represented through the metric tensor, whose spatial component 3 g modifies the spatial energy integral: where the energy density is represented by the component Θ 0 0 (depending on the metric convention, we may also need an additional minus sign). The spatial metric 3 g is controlled by the local energy density inside the quark. The simplest approximation is to replace this energy distribution by an effective mass M ef f at the origin, with M ef f equalling the total internal energy −(23/9)E. By expanding the resulting integral in G we obtain a well-defined integral: where γ = 23/9. At the end of this section we will also discuss results where the point mass at the center is replaced by an integrated density distribution based on the internal quark wave functions f and g. Carrying out the integral we obtain approximately: where we replaced 1/r by its expectation value 1/r = E/δx 0 where δ = 0.6019 · · · for the bag wave function. This expression reaches its minimum for The correction 3γGE 2 /δx 0 to the Minkowski metric equals 1/3, which means that the first order (perturbative) usage of GR is consistent. Both the radius r 0 and the internal frequency E are of Planck scale, with r 0 = 8.8 l Planck = 1.4 × 10 −19 fm and E = 0.23 M Planck = 2.8 × 10 21 MeV. The negativity of the internal quark energy was essential for the stabilization of the system. Hence, a property that appeared like a serious problem of the model actually was necessary to stabilize the system and to fix the scale. However, the relationship between the negative internal energy of Planck mass magnitude and the quark mass, which is positive and lies in the MeV range, is still unexplained.
In order to address this problem we introduce the vacuum energy density ǫ, or what is equivalent a finite cosmological constant. Recently, a cosmological model was developed which is based on the presence of this vacuum (dark) energy [4]. The associated non-perturbative metric of the vacuum universe has a very simple time dependence: where the characteristic time t s follows from the postulated vacuum energy density ǫ: t s = (3/8πGǫ) 1/2 . Another consequence of this cosmological model is that t s = 1/H 0 , where H 0 is the Hubble constant, while t s also represents the age of the universe as measured by a comoving observer. Fitting recent supernovae data one arrives at a vacuum energy density of ǫ = 3.97×10 −47 GeV 4 and an age of the universe of t s = 13.8 billion years. The vacuum metric leads to the effective time-dependent vacuum energy density: We now make the novel assumption that the creation of a quark is associated with the creation of a small vacuum universe (a small big bang). The creation of this bag terminates at a time t c , when the size of the bag matches the size r 0 of the QFT quark bag. Since this creation time t c is very small, the effective vacuum energy density is very large (just like it was t c seconds after the big bang), and the resulting energy can potentially cancel the internal energy, leaving a net quark mass of much smaller magnitude. According to the Heisenberg uncertainty relation [24] the total energy during this time interval t c is uncertain by an amount /2t c , so that it is possible during this time interval to create a mass of magnitude /2t c . Hence, we set: and verify whether this assumption for the quark mass is consistent with the energy balance: where E int is the minimum energy determined earlier.
Since the quark mass will turn out to be infinitesimal compared to the other terms, it is natural to ignore the left-hand term: Clearly, the required cancelation between internal and vacuum energy can only happen if the internal QFT energy is negative, so the negativity of the internal energy is again essential for the construction of a realistic quark model. The cancelation condition leads to a creation time t c : This expresses the creation time -and therefore the mass of the quark -exclusively in terms of the cosmological parameters G and ǫ (or H 0 ), as the remaining parameters (x 0 = 2.0478 · · · , γ = 23/9 and δ ≈ . Reference [25] quotes the following lattice values for the average up and down mass: 3.8 ± 0.8 MeV, while the average value excluding the lattice is quoted as: 4.4 ± 1.5 MeV. Very recently Davies et. al. [21] reported on a new determination of the light quark masses of much higher precision. They quote the following average values for the light quark masses: m q (2 GeV) = 3.40 ± .07 MeV and m q (3 GeV) = 3.07 ± .06 MeV, where the energy in brackets is the scale used in their calculation. The agreement between our result 3.17 MeV and these "experimental" values is extremely good, especially if one takes into consideration the enormous differences between the particle physics scale and the cosmological scale, from which our masses have been derived. Naturally, additional effects, like the electro-weak forces (which splits the degeneracy between u and d-quarks) and O(α s )-terms (if α s = 0), have to be included in more detailed calculations. However, the current level of agreement with experiment is clearly a strong endorsement of the theory presented here and confirms that this solution corresponds to the first generation of quarks.

X. SUMMARY AND DISCUSSION
In this paper we formulate a new set of methods in QFT that can be used to describe the internal properties and masses of elementary particles. The exact nonperturbative operator structure of the fermion and boson quantum fields is constructed by demanding the satisfaction of the full set of coupled quantum field equations. Starting from the basic QCD Lagrangian with bare pointlike massless quarks, we derive the internal wave functions of the dressed quarks and the form of the binding potentials. The binding mechanism of the localized quark-gluon state is topological in nature as it enabled by the non-linear structure of the equations and is independent of the magnitude of the strong coupling constant. Since all equations and fields refer to a single space-time coordinate, an additional operator prescription is required to reflect the opposite ordering of particles and anti-particles. This prescription (the so-called R-product) also resolves the so-called cosmological constant problem, as it leads to a more fundamental and less ambiguous definition of the vacuum in QFT.
The calculation of the quark mass is complicated by the fact that the QFT expectation value of the energy operator is negative. It is possible that by including an elementary Higgs-quark interaction, this problem may be resolved, however, this can only account for two of the three quark generations. The third -and probably basic -generation would then correspond to the trivial (zero) Higgs solution. However, in the absence of a Higgs field we need another entity to counter the negative QFT energy and give the system a scale. We propose a novel solution for this puzzle by involving the perturbative use of gravitational forces at the Planck scale. In this way we arrive at a positive mass which is amazingly close to accepted phenomenological values for the light quarks.
The current theory would get further support if we can deduce the existence of the two higher generations of quarks. However, as we indicated above, this may well require the inclusion of a Higgs quark interaction, which will entail more profile functions and a more involved formulation. But it also offers a chance to verify various models of this force and the associated electro-weak interaction at a more fundamental level than the standard model, by comparing the predictions to experiment. Affirmation of the proposed Lagrangians would undoubtedly lead to a more unified picture of Nature. Another, desirable extension is the perturbative inclusion of QED to explain the mass differences between quarks within one generation. Clearly, QCD cannot describe the dressing of leptons, so that for this sector the non-linear electro-weak interaction must induce the binding mechanism. The application of our formulation to unbroken SU(2) is fairly straightforward (for example we already found that the factor 32/9 is replaced by 3/2), but the chiral nature of the electro-weak interaction is more difficult to include. Again, this application to leptons offers a good opportunity to test various models of the electro-weak interac-tions at the fundamental level, as it is possible that the complicated multiplet structure in the standard model is absent at a more fundamental level. The smallness of the neutrino masses would put extremely strong constraints on the possible theories, and if a model was successful in predicting the lepton masses then it would mark a historic breakthrough in our understanding of Nature.
The application of the new methodology is not completely limited to the description of the internal dynamics of elementary particles. It can also be used to construct free boson propagators using the new quantization method for solving the (homogeneous) equations of motion. The R-product also could play a role in scattering diagrams if there are non-linear vertices with common space-time coordinates (as pointed out repeatedly, the R-product does not affect the algebra of operators defined for different space-time coordinates). In scattering calculations the finite size of dressed fermions (emerging from the current theory) could also play a role in defining physical cut-offs, and thereby ensure the finiteness of otherwise infinite diagrams. However, the size of the light quarks (about 8 Planck lengths) is beyond measurement. If the higher generations of quarks are indeed due to coupling to an elementary Higgs field, then their sizes would be substantially larger, and of the order of 0.004 fm on average, so they might be measurable. Outside of these specific effects, the standard scattering formulation would remain unaltered, and so all the successful scattering applications of the standard model seem fully consistent with the current generalizations of QFT.
Appendix A: The R-product The formulation of QFT contains an inherent asymmetry between particles and anti-particles, which reflects the inadequacy of the mathematical language in which we express this formulation. Since we are constrained by this language we have to add another rule to QFT, which together with the usual rules constitutes a symmetric representation of the world of particles and anti-particles. In this appendix we define this rule. It should form an inherent part of the foundations of QFT. Its absence in the current practice of QFT has led to many misconceptions and is also responsible for the so-called cosmological constant problem, which constitutes the biggest discrepancy ever between theory and experiment in the history of physics.
Consider a typical operator expression in QFT:ψOψ. If we evaluate the vacuum matrix element of this operator and concentrate on the particle contributions we get: as expected. If we look at the expectation value for a particle state we get: expressing the desired link between the state considered and the matrix element. However, when we consider the anti-particle contribution of the operator we have: which is in general non-zero, and can even be infinite, as the sum is unrestricted. Clearly, this asymmetry is unacceptable and illustrates the incompleteness of the standard QFT description (notice that there is no problem with mixed terms such as b † α d † β , which are responsible for the quantum fluctuations). The origin of this asymmetry between particles and anti-particles is that the right-hand operator ψ contains particle annihilation operators and anti-particle creation operators. This particular combination is required by baryon (or lepton) number conservation. However, it leads to the possibility of creating anti-particles -and the impossibility to create particles -when operating with ψ on the ket-vector, even if it is vacuum state. Similarly, the operator on the left,ψ, contains a combination of particle creation and anti-particle annihilation operators, again required by baryon number conservation. This combination displays a similar unacceptable asymmetry when operating to the left on a bra-vector. Let us see how we can restore the symmetry between particles and anti-particles.
To eliminate vacuum terms, like those in Eq.(A3), one usually applies the normal product to the bilinear expression: Clearly, the anti-particle contributions to the vacuum matrix element now also vanish. The usual justification for this heuristic procedure in QFT is that in QFT one is not interested in absolute energies, but rather in relative energies, such as excitation energies. By consistently replacing bilinear expressions by the normal product, one then eliminates such sources of infinities from the start. However, within GR energy leads to a modification of the metric, so that one cannot arbitrarily rescale the energy. Thus the heuristic replacement of bilinear expressions by their normal product is no longer seen as acceptable. This has led to the consideration of QFT vacuum contributions to the energy as a real phenomenon (especially in cosmology), although this is in total disagreement with experimental observations. As we saw above, the normal product is needed to restore the symmetry between particles and anti-particles, so we cannot simply dismiss it because it is a heuristic trick. Rather, we have to explain the origin of this heuristic prescription and should definitely not accept the reality of the vacuum contributions as a fait accompli, especially since these contributions are associated with an unacceptable asymmetry between particles and anti-particles and are in strong disagreement with observation.
The key towards the solution of this problem is to examine the normal product from a different perspective. The normal product is defined as a rearrangement of the operators, so that the creation operators are shifted to the left of the annihilation operators. The normal product is also widely used in the derivations of the scattering expansion in QFT, although there it is used as a convenient mathematical construct and not employed heuristically. In its current heuristic application, Eq.(A4), it is equivalent to the inversion of the order of the antiparticle operators, as the normal product leaves the particle operators unaffected in the QFT expressionψOψ as Hence, an equivalent definition of the normal product for the relevant QFT expressions is that the order of the anti-particle operators must be reversed. We will refer to this re-ordering as the R-product. We will argue that there is a fundamental theoretical justification for this product, which is lacking in the heuristic application of the normal product. Hence, in the bilinear expressions considered one can justify the use of the normal product indirectly via the justification for the R-product. For longer expressions the usage of the R-product, rather than the normal product, leads to different results. However, most of the long expressions involve different space-time coordinates, where neither the R-product, nor the heuristic normal product applies. If we analyze the usage of the heuristic application of the normal product in the existing QFT literature, we see that it is customary to define the (interaction) Hamiltonian heuristically as normal ordered (see e.g. [26]), before it is inserted in the matrix elements of the S-matrix expansion. This means that there is no heuristic normal ordering between operators corresponding to different coordinates x and x ′ (it makes no difference if we let x → x ′ , for the current practice it is only important that the coordinates were different originally). In fact the heuristic usage of normal ordering between different spatial integrals in the S-matrix expansion would destroy the whole structure and meaning of the S-matrix expansion. Hence, the heuristic use of the normal product must be limited to strings of operators that all belong to a single coordinate x (we will later discuss the fundamental origin of this rule). But for such short expressions, corresponding to a single x, the normal product and the R-product are often equivalent, so that either one can be used. This equivalence for the simplest QFT expressions is the likely reason why the R-product has not been discovered before in QFT, as the usage of the normal product was often adequate.
For longer expressions belonging to a single space-time coordinate the usage of the R-product leads to very different results than the normal product, so here it becomes of importance to use the R-product. In the current boundstate application of QFT we encounter many long chains of creation and annihilation operators belonging to a single coordinate x, so our transition to the R-product has major consequences for the current formulation (in standard scattering applications of non-Abelian QFT it might also have consequences, as these also involve strings of operators at the same space-time point). Consider a longer local QFT expression, like d α b β b † γ d † δ and b † γ d † δ d α b β (it is essential to use expressions that actually occur in the usage of QFT, since we are trying to correct a defect of QFT that occurs in its actual representation). Using the normal product we get: whereas the R-product leads to: Under the normal product the distinction between the two initial expressions is lost, while the R-product maintains the integrity of the expression, while still ensuring that the vacuum matrix element vanishes. If we apply the resulting expressions to a one-particle state then we get zero under the normal product, while the R-product gives a finite outcome. Hence, the R-product leads to the survival of non-linear terms in the field equations which are essential for the existence of bound states, such as considered in this paper. The physical justification of the R-product is that it embodies Feynman's interpretation of anti-particles as particles moving backward in time. In the S-matrix expansion the ordering of operators defined at different times is taken care off by the time-ordering procedures of standard QFT. However, when anti-particles all belong to the same time t, this ordering cannot be implemented by explicit time ordering. That is why one has to impose the R-product on local expressions to ensure the correct ordering of interactions and vertices. The terms forward and backward in time may seem a bit inappropriate for local operators defined at a single t, so let us analyze this concept in the current context. In the matrix elements the ket-vector is usually seen as the initial state, while the bra-vector is seen as the final state. From this perspective, the operator right-most in a chain first operates on the ket-vector, and is "earlier" in time than operators further to the left. For a chain of anti-particle operator the left-most operator must first act on the initial state (i.e. the ket vector), etc. This can be accomplished by applying the R-product to expressions dependent on one space-time coordinate.
One consequence of the R-product is that one cannot insert a complete set of states somewhere inside the operator expression. Hence, the manipulation and reduction of an R-product requires special care, as we will also see in the following example. Consider the following reduction: We see that if we commute the anti-particle operators under the R-product, we get an extra minus sign in the anticommutator, compared to the original anti-commutator, Eq. (11). This is conveniently summarized by the following anti-commutation rule: R · · · d † α , d β · · · = R [· · · (−δ αβ ) · · · ] , x = x ′ , (A8) where x and x ′ are the coordinates belonging to the two anti-particle operators. Since, all operators have to be expressed in R-product form, we have to apply this extra minus sign whenever we contract anti-particle terms corresponding to the same space-time point in the equations of motion or in the expressions for the observables. These minus signs are absolutely necessary for deriving consistent equations, in fact we established the need for these signs, before we had discovered their origin in the Rproduct. We will see in Appendix B, that the R-product also enables us to derive the general operator solution of the equations of motion by including operator chains of any length in the expressions.