The Spacetime Metric
STM-D-0908Paper2010Published and peer-reviewed

Numerical Experiments Providing New Insights into Plasma Focus Fusion Devices

Sing Lee · Sor Heoh Saw

Open licence · full text · Creative Commons Attribution 3.0 (CC BY 3.0), as stated at the foot of the paper: ’This article is an open-access article distributed under the terms and conditions of the Creative Commons Attribution license’.

In one page

A plasma focus is a capacitor bank fired into a coaxial electrode: the current sheet sweeps down the tube, collapses onto the axis, and for a few billionths of a second makes a pinch hot and dense enough to fuse deuterium. Sing Lee and Sor Heoh Saw ran thousands of simulations with the Lee model code, which is tuned to each machine by fitting its computed current waveform to the measured one, and report four findings. Lowering a machine's fixed circuit inductance raises the peak current without limit, but the current that actually reaches the pinch peaks and then falls — so peak current is not a proxy for pinch current. Neutron and soft X-ray yields scale as steep powers of the pinch current. And the long-standing worry that neutron output saturates above a few hundred kilojoules has a cause they name: the moving current sheet acts as a fixed resistance of about seven milliohms, which caps the current once the bank impedance drops below it. Their prescription is to raise the voltage instead.

Why it matters hereChapter 12 needs devices that make fusion happen in a compact, repeatable machine rather than a reactor hall, and the plasma focus is the oldest of them; chapter 9 is about self-organised plasma structures, and the focus pinch is the laboratory's most studied one. This paper matters because it replaces a discouraging empirical wall with a circuit quantity you can engineer around.

What it claims

  1. 01Reducing a plasma focus device's static inductance raises the peak circuit current without limit, but the pinch current reaches a broad maximum and then falls: for PF1000 at 35 kV in 3.5 torr deuterium the peak current rose from 1.66 to 4.37 megamperes as static inductance fell from 100 to 5 nanohenries, while the pinch current peaked at 1.05 megamperes near 40 to 30 nanohenries and then declined.Section 3 and Section 3.1, Table 1 and Figure 7

    Published and peer-reviewed
  2. 02Over storage energies from about 8.5 kilojoules to 24.5 megajoules the computed neutron yield scales as the pinch current to the power 4.5 and as the peak current to the power 3.8, and the pinch-current scaling holds with less scatter than the peak-current scaling when machines with mixed parameters are added.Section 4.2, Figures 8 and 9

    Published and peer-reviewed
  3. 03Neon soft X-ray yield over 0.2 kilojoules to 1 megajoule scales as the pinch current to the power 3.6 and as the peak current to the power 3.2, and only the pinch-current index survives unchanged when the circuit and model parameters of real machines are substituted.Section 5.2, Figure 11 and the summarising equations

    Published and peer-reviewed
  4. 04Because the moving current sheet makes the tube inductance change with time, half of that rate of change acts as a genuine resistance — the dynamic resistance — and for deuterium at the universal end-axial speed of about 10 centimetres per microsecond it takes a near-constant value of about 7 milliohms in the axial phase.Section 6.2, Equations 5 to 7; validated numerically in Section 2.3 against Figure 4j

    Published and peer-reviewed
  5. 05That fixed 7 milliohms is the cause of the apparent neutron saturation: as the bank capacitance grows, the surge impedance falls below the dynamic resistance, the peak current tends to the asymptote of charging voltage divided by dynamic resistance, and the steep yield scaling turns a current ceiling into a neutron ceiling.Section 6.2, Table 3 and Figure 14

    Published and peer-reviewed
  6. 06The authors propose escaping that ceiling by raising the operating voltage rather than the stored energy: the saturation level is proportional to voltage, and at 1 megavolt with parallel Blumleins matched to the radial-phase dynamic resistance they estimate peak currents beyond 15 megamperes and pinch currents above 6 megamperes, combined with krypton doping and current-stepping.Section 6.3

    What to watch

Read it

Abstract

Recent extensive and systematic numerical experiments have uncovered new insights into plasma focus fusion devices including the following: (1) a plasma current limitation effect, as device static inductance is reduced towards very small values; (2) scaling laws of neutron yield and soft x-ray yield as functions of storage energies and currents; (3) a global scaling law for neutron yield as a function of storage energy combining experimental and numerical data showing that scaling deterioration has probably been interpreted as neutron 'saturation'; and (4) a fundamental cause of neutron 'saturation'. The ground-breaking insights thus gained may completely change the directions of plasma focus fusion research.

Keywords: plasma focus; numerical experiments; scaling laws; neutron saturation

Authors and affiliations. Sing Lee (Institute for Plasma Focus Studies, Chadstone, Victoria, Australia; INTI University College, Nilai, Malaysia; Nanyang Technological University, National Institute of Education, Singapore) and Sor Heoh Saw (Institute for Plasma Focus Studies; INTI University College). Energies 2010, 3, 711–737. Received 14 February 2010; accepted 3 March 2010; published 12 April 2010. Review article.

1. Introduction

The Lee model code couples the electrical circuit with plasma focus dynamics, thermodynamics, and radiation, enabling a realistic simulation of all gross focus properties. The basic model, described in 1984, was successfully used to assist several projects. Radiation-coupled dynamics was included in the five-phase code, leading to numerical experiments on radiation cooling. The vital role of a finite small disturbance speed discussed by Potter in a Z-pinch situation was incorporated together with real gas thermodynamics and radiation-yield terms. This version of the code assisted other research projects and was web published in 2000 and 2005. Plasma self-absorption was included in 2007, improving the soft X-ray yield simulation. The code has been used extensively in several machines including UNU/ICTP PFF, NX2, and NX1 and has been adapted for the Filippov-type plasma focus DENA. A recent development is the inclusion of the neutron yield using a beam–target mechanism, incorporated in recent versions of the code (versions later than RADPFV5.13), resulting in realistic neutron yield scaling with the pinch current. The versatility and utility of the model are demonstrated in its clear distinction of the pinch current from the peak current and the recent uncovering of a plasma focus pinch current limitation effect as static inductance is reduced towards zero. Extensive numerical experiments had been carried out systematically resulting in the uncovering of neutron and soft X-ray scaling laws over a wider range of energies and currents than attempted before. The numerical experiments also gave insight into the nature and cause of 'neutron saturation'. The description, theory, code, and a broad range of results of this "Universal Plasma Focus Laboratory Facility" are available for download.

A brief description of the 5-phase model is given in the following.

2. The 5-Phase Lee Model Code

The five phases are summarized as follows:

a. Axial Phase (see Figure 1 left part). Described by a snowplow model with an equation of motion which is coupled to a circuit equation. The equation of motion incorporates the axial phase model parameters: mass and current factors fm and fc. The mass swept-up factor fm accounts for not only the porosity of the current sheet but also for the inclination of the moving current sheet shock front structure, boundary layer effects, and all other unspecified effects which have effects equivalent to increasing or reducing the amount of mass in the moving structure, during the axial phase. The current factor fc accounts for the fraction of current effectively flowing in the moving structure (due to all effects such as current shedding at or near the back-wall, and current sheet inclination). This defines the fraction of current effectively driving the structure, during the axial phase.

b. Radial Inward Shock Phase (see Figure 1 right part, also Figure 2). Described by four coupled equations using an elongating slug model. The first equation computes the radial inward shock speed from the driving magnetic pressure. The second equation computes the axial elongation speed of the column. The third equation computes the speed of the current sheath (magnetic piston), allowing the current sheath to separate from the shock front by applying an adiabatic approximation. The fourth is the circuit equation. Thermodynamic effects due to ionization and excitation are incorporated into these equations, these effects being particularly important for gases other than hydrogen and deuterium. Temperature and number densities are computed during this phase using shock-jump equations. A communication delay between shock front and current sheath due to the finite small disturbance speed is crucially implemented in this phase. The model parameters, radial phase mass swept-up and current factors fmr and fcr, are incorporated in all three radial phases. The mass swept-up factor fmr accounts for all mechanisms which have effects equivalent to increasing or reducing the amount of mass in the moving slug, during the radial phase. The current factor fcr accounts for the fraction of current effectively flowing in the moving piston forming the back of the slug (due to all effects). This defines the fraction of current effectively driving the radial slug.

c. Radial Reflected Shock Phase (see Figure 2). When the shock front hits the axis, because the focus plasma is collisional, a reflected shock develops which moves radially outwards, whilst the radial current sheath piston continues to move inwards. Four coupled equations are also used to describe this phase, these being for the reflected shock moving radially outwards, the piston moving radially inwards, the elongation of the annular column and the circuit. The same model parameters fmr and fcr are used as in the previous radial phase. The plasma temperature behind the reflected shock undergoes a jump by a factor close to 2. Number densities are also computed using the reflected shock jump equations.

d. Slow Compression (Quiescent) or Pinch Phase (see Figure 2). When the out-going reflected shock hits the inward moving piston, the compression enters a radiative phase in which for gases such as neon, radiation emission may actually enhance the compression where we have included energy loss and gain terms from Joule heating and radiation losses into the piston equation of motion. Three coupled equations describe this phase; these being the piston radial motion equation, the pinch column elongation equation and the circuit equation, incorporating the same model parameters as in the previous two phases. The duration of this slow compression phase is set as the time of transit of small disturbances across the pinched plasma column. The computation of this phase is terminated at the end of this duration.

e. Expanded Column Phase. To simulate the current trace beyond this point we allow the column to suddenly attain the radius of the anode, and use the expanded column inductance for further integration. In this final phase the snowplow model is used, and two coupled equations are used similar to the axial phase above. This phase is not considered important as it occurs after the focus pinch.

We note that in radial phases b, c and d, axial acceleration and ejection of mass caused by necking curvatures of the pinching current sheath result in time-dependent strongly center-peaked density distributions. Moreover the transition from phase d to phase e is observed in laboratory measurements to occur in an extremely short time with plasma and current disruptions resulting in localized regions of high densities and temperatures. These center-peaking density effects and localized regions are not modeled in the code, which consequently computes only an average uniform density and an average uniform temperature which are considerably lower than measured peak density and temperature. However, because the four model parameters are obtained by fitting the computed total current waveform to the measured total current waveform, the model incorporates the energy and mass balances equivalent, at least in the gross sense, to all the processes which are not even specifically modeled. Hence the computed gross features such as speeds and trajectories and integrated soft x-ray yields have been extensively tested in numerical experiments for several machines and are found to be comparable with measured values.

Figure 1. Schematic of the axial and radial phases. The left section depicts the axial phase, the right section the radial phase. In the left section, z is the effective position of the current sheath-shock front structure. In the right section rs is the position of the inward moving shock front driven by the piston at position rp. Between rs and rp is the radially imploding slug, elongating with a length zf. The capacitor, static inductance and switch powering the plasma focus is shown for the axial phase schematic only.

Figure 2. Schematic of radius versus time trajectories to illustrate the radial inward shock phase when rs moves radially inwards, the reflected shock phase when the reflected shock moves radially outwards, until it hits the incoming piston rp leading to the start of the pinch phase and finally the expanded column phase.

2.1. From Measured Current Waveform to Modeling for Diagnostics

The Lee model code is configured to work as any plasma focus by inputting: bank parameters, the static inductance L0, the capacitance C0 and stray circuit resistance r0; tube parameters b, a and z0; and operational parameters V0 and P0 and the fill gas.

The computed total current waveform is fitted to the measured waveform by varying model parameters fm, fc, fmr and fcr one by one, until the computed waveform agrees with the measured waveform.

First, the axial model factors fm and fc are adjusted (fitted) until the features (1) computed rising slope of the total current trace and (2) the rounding off of the peak current as well as (3) the peak current itself are in reasonable (typically very good) fit with the measured total current trace (see Figure 3, measured trace fitted with computed trace).

Then we proceed to adjust (fit) the radial phase model factors fmr and fcr until features (4) the computed slope and (5) the depth of the dip agree with the measured. Note that the fitting of the computed trace with the measured current trace is done up to the end of the radial phase which is typically at the bottom of the current dip. Fitting of the computed and measured current traces beyond this point is not done. If there is significant divergence of the computed with the measured trace beyond the end of the radial phase, this divergence is not considered important.

In this case, after fitting the 5 features (1) to (5) above, the following fitted model parameters are obtained: fm = 0.1, fc = 0.7, fmr = 0.12, fcr = 0.68.

Figure 3. The 5-point fitting of computed current trace to the measured (or the reference) current trace. Point 1 is the current rise slope. Point 2 is the topping profile. Point 3 is the peak value of the current. Point 4 is the slope of the current dip. Point 5 is the bottom of the current dip. Fitting is done up to point 5 only. Further agreement or divergence of the computed trace with or from the measured trace is only incidental and not considered to be important.

From experience, it is known that the current trace of the focus is one of the best indicators of gross performance. The axial and radial phase dynamics and the crucial energy transfer into the focus pinch are among the important information that is quickly apparent from the current trace.

The exact time profile of the total current trace is governed by the bank parameters, by the focus tube geometry and the operational parameters. It also depends on the fraction of mass swept-up and the fraction of sheath current and the variation of these fractions through the axial and radial phases. These parameters determine the axial and radial dynamics, specifically the axial and radial speeds which in turn affect the profile and magnitudes of the discharge current. There are many underlying mechanisms in the axial phase such as shock front and current sheet structure, porosity and inclination, boundary layer effects and current shunting and fragmenting which are not simply modeled; likewise in the radial phase mechanisms such as current sheet curvatures and necking leading to axial acceleration and ejection of mass, and plasma and current disruptions. These effects may give rise to localized regions of high density and temperatures. The detailed profile of the discharge current is influenced by these effects and during the pinch phase also reflects the Joule heating and radiative yields. At the end of the pinch phase the total current profile also reflects the sudden transition of the current flow from a constricted pinch to a large column flow. Thus the discharge current powers all dynamic, electrodynamic, thermodynamic and radiation processes in the various phases of the plasma focus. Conversely all the dynamic, electrodynamic, thermodynamic and radiation processes in the various phases of the plasma focus affect the discharge current. It is then no exaggeration to say that the discharge current waveform contains information on all the dynamic, electrodynamic, thermodynamic and radiation processes that occur in the various phases of the plasma focus. This explains the importance attached to matching the computed total current trace to the measured total current trace in the procedure adopted by the Lee model code. Once matched, the fitted model parameters assure that the computation proceeds with all physical mechanisms accounted for, at least in the gross energy and mass balance sense.

2.2. Diagnostics — Time Histories of Dynamics, Energies and Plasma Properties Computed from the Measured Total Current Waveform by the Code

During every adjustment of each of the model parameters the code goes through the whole cycle of computation. In the last adjustment, when the computed total current trace is judged to be reasonably well fitted in all five waveform features, computed time histories are presented, in Figure 4a to 4o as an example, as follows: for the NX2 operated at 11 kV, 2.6 Torr neon.

Figure 4. Computed plasma and electrodynamic properties of the NX2, in fifteen panels: (a) fitted computed total current; (b) computed total current and plasma current; (c) computed tube voltage; (d) axial trajectory and speed; (e) radial trajectories; (f) length of elongating structure; (g) radial shock front, piston and elongation speeds; (h) computed tube inductance; (i) inductive energy; (j) piston work and work dissipated by the dynamic resistance; (k) dynamic resistance; (l) ion number density, maximum and averaged uniform; (m) electron number density; (n) plasma temperature, maximum and averaged uniform; (o) neon soft X-ray power.

2.3. Comments on Computed Quantities

The computed total current trace typically agrees very well with the measured because of the fitting. The end of the radial phase is indicated in Figure 4a. Plasma currents are rarely measured. We had done a comparison of the computed plasma current with measured plasma current for the Stuttgart PF78 which shows good agreement of our computed to the measured plasma current. The computed plasma current in this case of the NX2 is shown in Figure 4b. The computed tube voltage is difficult to compare with measured tube voltages in terms of peak values, typically because of poor response time of voltage dividers. However the computed waveform shape in Figure 4c is generally as expected. The computed axial trajectory and speed agree with experimentally obtained time histories. Moreover, the behavior with pressure, running the code again at different pressures, agrees well with experimental results. The radial trajectories and speeds are difficult to measure. The computed trajectories in Figure 4e agree with the scant experimental data available. The length of the radial structure is shown in Figure 4f. Computed radial shock front and piston speeds and speed of the elongation of the structure are shown in Figure 4g. The computed inductance (Figure 4h) shows a steady increase of inductance in the axial phase, followed by a sharp increase, rising by more than a factor of 2 in a radial phase time interval about one tenth the duration of the axial phase for the NX2.

The inductive energy peaks at 70% of initial stored energy, and then drops to 30% during the radial phase, as the sharp drop of current more than offsets the effect of sharply increased inductance (Figure 4i). In Figure 4j is shown the work done by the magnetic piston, computed using the force integrated over distance method. Also shown is the work dissipated by the dynamic resistance, computed using dynamic resistance power integrated over time. We see that the two quantities and profiles agree exactly. This validates the concept of half Ldot as a dynamic resistance (see Section 6.2). The piston work deposited in the plasma increases steadily to some 12% at the end of the axial phase and then rises sharply to just below 30% in the radial phase. Dynamic resistance is shown in Figure 4k. The values of the dynamic resistance in the axial phase, together with the bank surge impedance, are the quantities that determine the peak current. The ion number density has a maximum value derived from shock-jump considerations, and an averaged uniform value derived from overall energy and mass balance considerations. The time profiles of these are shown in Figure 4l. The electron number density (Figure 4m) has similar profiles to the ion density profile, but is modified by the effective charge numbers due to ionization stages reached by the ions. Plasma temperature too has a maximum value and an averaged uniform value derived in the same manner; these are shown in Figure 4n. Computed neon soft x-ray power profile is shown in Figure 4o. The area of the curve is the soft x-ray yield in Joule. Pinch dimensions and lifetime may be estimated from Figure 4e and 4f. The model also computes the neutron yield, for operation in deuterium, using a phenomenological beam-target mechanism. The model does not compute a time history of the neutron emission, only a yield number.

Thus as is demonstrated above, the model code when properly fitted is able to realistically model any plasma focus and act as a guide to diagnostics of plasma dynamics, trajectories, energy distribution and gross plasma properties.

Moreover, using such simulation, series of experiments have been systematically carried out to look for behavior patterns of the plasma focus. Insights uncovered by the series of experiments include: (i) pinch current limitation effect as static inductance is reduced; (ii) neutron and soft X-ray scaling laws; (iii) a global scaling law for neutrons versus storage energy combining experimental and numerical experimental data; and (iv) insight into the nature and a fundamental cause of neutron saturation. These are significant achievements accomplished within a period of twenty months of intensive numerical experimentation.

3. Insight 1 — Pinch Current Limitation Effect as Static Inductance Is Reduced towards Zero

In a recent paper, there was expectation that the large megajoule plasma focus PF1000 in Warsaw could increase its discharge current, and its pinch current, and consequently neutron yield by a reduction of its external or static inductance L0. To investigate this point, experiments were carried out using the Lee Model code. Unexpectedly, the results indicated that whilst the peak current indeed progressively increased with reduction in L0, no improvement may be achieved due to a pinch current limitation effect. Given a fixed C0 powering a plasma focus, there exists an optimum L0 for maximum pinch current. Reducing L0 further will increase neither the pinch current nor the neutron yield. The numerical experiments leading to this unexpected result are described below.

A measured current trace of the PF1000 with C0 = 1332 µF, operated at 27 kV, 3.5 torr deuterium, has been published, with cathode and anode radii b = 16 cm, a = 11.55 cm and anode length z0 = 60 cm. In the numerical experiments we fitted external (or static) inductance L0 = 33.5 nH and stray resistance r0 = 6.1 milliohms, a damping factor RESF of 1.22. The fitted model parameters are: fm = 0.13, fc = 0.7, fmr = 0.35 and fcr = 0.65. The computed current trace agrees very well with the measured trace through all the phases, axial and radial, right down to the bottom of the current dip indicating the end of the pinch phase as shown in Figure 5.

Figure 5. Fitting computed current to measured current traces to obtain fitted parameters fm = 0.13, fc = 0.7, fmr = 0.35 and fcr = 0.65. The measured current trace was for the PF1000 at 27 kV, storage capacity of 1332 µF and fitted static inductance of 33.5 nH.

We carried out numerical experiments for PF1000 using the machine and model parameters determined from Figure 5. Operating the PF1000 at 35 kV and 3.5 Torr, we varied the anode radius a with corresponding adjustment to b to maintain a constant ratio c = b/a = 1.39 and in order to keep the peak axial speed at 10 cm per microsecond. The anode length z0 was also adjusted to maximize the pinch current as L0 was decreased from 100 nH progressively to 5 nH.

As expected, the peak current increased progressively from 1.66 to 4.4 MA. As L0 was reduced from 100 to 35 nH, the pinch current also increased, from 0.96 to 1.05 MA. However, then unexpectedly, on further reduction from 35 to 5 nH, the pinch current stopped increasing, instead decreasing slightly to 1.03 MA at 20 nH, to 1.0 MA at 10 nH, and to 0.97 MA at 5 nH. The neutron yield also had a maximum value of 3.2 times 10¹¹ at 35 nH.

To explain this unexpected result, we examine the energy distribution in the system at the end of the axial phase just before the current drops from its peak value and then again near the bottom of the almost linear drop to the pinch phase indicated by the arrow pointing to 'end of radial phase'.

Equation 1, given here as a named result. The magnetic energy held at peak current in the static inductance plus the axial tube inductance, scaled by the square of the current factor, is redistributed at the end of the radial phase into three parts: the magnetic energy held at pinch current in the static inductance divided by the square of the current factor, plus the axial and pinch inductances; the energy exchanged with the capacitor during the current drop; and the energy imparted to the plasma as the current sheet moves to the pinch position, which is the integral of half the rate of change of inductance multiplied by the square of the current. We approximate the last as half the pinch inductance multiplied by the square of the pinch current, which is an underestimate for this case. If the duration of the radial phase is short compared to the capacitor time constant, the capacitor is effectively decoupled and the capacitor term may be put as zero.

Equation 2, given here as a named result. From that consideration, the square of the pinch current equals the square of the peak current multiplied by the ratio of two inductance sums: the static inductance plus half the axial inductance, divided by twice the static inductance plus the axial inductance plus twice the pinch inductance. Here the current factor has been taken as 0.7 and its square approximated as 0.5.

Generally, as L0 is reduced, the peak current increases; the anode radius a is necessarily increased leading to a longer pinch length, hence a bigger pinch inductance. Lowering L0 also results in a shorter rise time, hence a necessary decrease in the anode length z0, reducing the axial inductance. Thus, from Equation 2, lowering L0 decreases the ratio of pinch current to peak current. Secondly, this situation is compounded by another mechanism. As L0 is reduced, the L-C interaction time of the capacitor bank reduces while the duration of the current drop increases (see Figure 6, discussed in the next section) due to an increasing anode radius. This means that as L0 is reduced, the capacitor bank is more and more coupled to the inductive energy transfer processes with the accompanying induced large voltages that arise from the radial compression. Looking again at the derivation of Equation 2 from Equation 1, a nonzero capacitor term, in this case of positive value, will act to decrease the pinch current further. The lower the L0 the more pronounced is this effect.

Summarizing this discussion, the pinch current limitation is not a simple effect, but is a combination of the two complex effects described above, namely, the interplay of the various inductances involved in the plasma focus processes abetted by the increasing coupling of C0 to the inductive energetic processes, as L0 is reduced.

3.1. Optimum Static Inductance for Maximum Pinch Current and Neutron Yield

From the pinch current limitation effect, it is clear that given a fixed C0 powering a plasma focus, there exists an optimum L0 for maximum pinch current. Reducing L0 further will increase neither the pinch current nor the neutron yield. The results of the numerical experiments carried out are presented in Figure 6 and Table 1.

With large L0 = 100 nH it is seen (Figure 6) that the rising current profile is flattened from what its waveform would be if unloaded; and peaks at around 12 µs (before its unloaded rise time, not shown, of 18 µs) as the current sheet goes into the radial phase. The current drop, less than 25% of peak value, is sharp compared with the current rise profile. At L0 = 30 nH the rising current profile is less flattened, reaching a flat top at around 5 µs, staying practically flat for some 2 µs before the radial phase current drop to 50% of its peak value in a time which is still short compared with the rise time. With L0 of 5 nH, the rise time is now very short, there is hardly any flat top; as soon as the peak is reached, the current waveform droops significantly. There is a small kink on the current waveform of both the L0 = 5 nH, z0 = 20 cm and the L0 = 5 nH, z0 = 40 cm. This kink corresponds to the start of the radial phase which, because of the large anode radius, starts with a relatively low radial speed, causing a momentary reduction in dynamic loading. Looking at the three types of traces it is seen that for L0 = 100 nH to 30 nH, there is a wide range of z0 that may be chosen so that the radial phase may start at peak or near peak current, although the longer values of z0 tend to give better energy transfers into the radial phase.

The optimized situation for each value of L0 is shown in Table 1. The table shows that as L0 is reduced, the peak current rises with each reduction in L0 with no sign of any limitation. However, the pinch current reaches a broad maximum of 1.05 MA around 40 to 30 nH. Neutron yield also shows a similar broad maximum peaking at 3.2 times 10¹¹ neutrons. Figure 7 shows a graphical representation of this pinch current limitation effect. The curve going up to 4 MA at low L0 is the peak current curve. Thus the peak current shows no sign of limitation as L0 is progressively reduced. However the pinch current reaches a broad maximum. From Figure 7 there is a stark and important message: one must distinguish clearly between peak current and pinch current. In general one cannot take the peak current to be representative of the pinch current.

Figure 6. PF1000 current waveforms computed at 35 kV, 3.5 Torr deuterium for a range of static inductance showing the changes in waveforms as the static inductance varies.

Table 1. Effect on currents and ratio of currents as the static inductance is reduced — PF1000 at 35 kV, 3.5 Torr deuterium.

| L0 (nH) | b (cm) | a (cm) | z0 (cm) | Ipeak (MA) | Ipinch (MA) | Yn (10¹¹) | Ipinch / Ipeak | |---|---|---|---|---|---|---|---| | 100 | 15.0 | 10.8 | 80 | 1.66 | 0.96 | 2.44 | 0.58 | | 80 | 16.0 | 11.6 | 80 | 1.81 | 1.00 | 2.71 | 0.55 | | 60 | 18.0 | 13.0 | 70 | 2.02 | 1.03 | 3.01 | 0.51 | | 40 | 21.5 | 15.5 | 55 | 2.36 | 1.05 | 3.20 | 0.44 | | 35 | 22.5 | 16.3 | 53 | 2.47 | 1.05 | 3.20 | 0.43 | | 30 | 23.8 | 17.2 | 50 | 2.61 | 1.05 | 3.10 | 0.40 | | 20 | 28.0 | 21.1 | 32 | 3.13 | 1.03 | 3.00 | 0.33 | | 10 | 33.0 | 23.8 | 28 | 3.65 | 1.00 | 2.45 | 0.27 | | 5 | 40.0 | 28.8 | 20 | 4.37 | 0.97 | 2.00 | 0.22 |

Figure 7. Effect on currents and current ratio (computed) as the static inductance is reduced, PF1000, 35 kV, 3.5 torr deuterium.

We carried out several sets of experiments on the PF1000 for varying static inductance, each set with a different damping factor. In every case, an optimum inductance was found around 30 to 60 nH with the pinch current decreasing as L0 was reduced below the optimum value. The results showed that for PF1000, reducing L0 from its present 20 to 30 nH will increase neither the observed pinch current nor the neutron yield, because of the pinch limitation effect. Indeed, the pinch current decreases very slightly on further reduction to very small values. We would add that we have used a set of model parameters which in our experience is the most reasonable to be used in these numerical experiments. Variations of the model parameters could occur but we are confident that these variations are not likely to occur with such a pattern as to negate the pinch current limitation effect. Nevertheless these variations should be actively monitored and any patterns in the variations should be investigated.

4. Insight 2 — Scaling Laws for Neutrons

4.1. Computation of Neutron Yield — Describing the Beam-Target Mechanism

The neutron yield is computed using a phenomenological beam-target neutron generating mechanism described recently by Gribkov and colleagues and adapted to yield the following equation. A beam of fast deuteron ions is produced by diode action in a thin layer close to the anode, with plasma disruptions generating the necessary high voltages. The beam interacts with the hot dense plasma of the focus pinch column to produce the fusion neutrons.

Equation 3, given here as a named result. The beam-target yield is proportional to the ion density, the square of the pinch current, the square of the pinch length, the natural logarithm of the ratio of the cathode radius to the pinch radius, and the deuteron–deuteron fusion cross section, divided by the square root of the beam energy. The proportionality constant is treated as a calibration constant combining various constants in the derivation process.

The deuteron–deuteron cross-section is sensitive to the beam energy in the range 15 to 150 kV; so it is necessary to use the appropriate range of beam energy to compute it. The code computes induced voltages (due to current motion inductive effects) of the order of only 15 to 50 kV. However it is known, from experiments, that the ion energy responsible for the beam-target neutrons is in the range 50 to 150 keV, and for smaller lower-voltage machines the relevant energy could be lower at 30 to 60 keV. Thus in line with experimental observations the cross section is reasonably obtained by using a beam energy of three times the computed maximum induced voltage. This fit was tested by using the beam energy equal to various multiples of that voltage. A reasonably good fit of the computed neutron yields to the measured published neutron yields at energy levels from sub-kilojoule to near megajoule was obtained when the multiple of 3 was used; with poor agreement for most of the data points when for example a multiple of 1 or 2 or 4 or 5 was used. The model uses a calibration constant of 2.7 times 10⁷ obtained by calibrating the yield at an experimental point of 0.5 MA.

The thermonuclear component is also computed in every case and it is found that this component is negligible when compared with the beam-target component. It might be argued that an adjustment to the thermonuclear component could also be attempted in a similar way to the usage of the multiple to the maximum voltage. However, the usage of that multiple has some experimental basis due to ion energy measurements. Moreover the value of the maximum voltage in each numerical experiment is calculated from the slug model leading to the slow compression phase whilst it is known experimentally that after the slow compression phase, instability effects set in which will increase the electric fields operating within the pinch. These are the basic arguments supporting the view that the operational beam energy has a value above the computed maximum induced voltage. For the thermonuclear component a feasible model to adjust the yield upwards has yet to be suggested.

4.2. Scaling Laws for Neutrons from Numerical Experiments over a Range of Energies from 10 kJ to 25 MJ

We apply the Lee model code to the megajoule machine PF1000 over a range of C0 to study the neutrons emitted by PF1000-like bank energies from 10 kJ to 25 MJ.

As shown earlier the PF1000 current trace has been used to fit the model parameters, with very good fitting achieved between the computed and measured current traces (Figure 5). Once the model parameters have been fitted to a machine for a given gas, these model parameters may be used with some degree of confidence when operating parameters such as the voltage are varied. With no measured current waveforms available for the higher megajoule numerical experiments, it is reasonable to keep the model parameters that we have got from the PF1000 fitting.

The optimum pressure for this series of numerical experiments is 10 torr and the ratio c = b/a is retained at 1.39. For each C0, anode length z0 is varied to find the optimum. For each z0, anode radius is varied so that the end axial speed is 10 cm per microsecond. The numerical experiments were carried out for C0 ranging from 14 µF to 39,960 µF corresponding to energies from 8.5 kJ to 24.5 MJ.

For this series of experiments we find that the neutron yield scaling changes from proportional to the square of the storage energy at tens of kJ, to proportional to the storage energy to the power 0.84 at the highest energies (up to 25 MJ) investigated in this series. This is shown in Figure 8.

Figure 8. Neutron yield plotted as a function of storage energy in log-log scale, showing the scaling change from the square of the storage energy at tens of kJ to the storage energy to the power 0.84 at the highest energies (up to 25 MJ).

The scaling of neutron yield with peak current and pinch current over the whole range of energies investigated up to 25 MJ (shown in Figure 9) are as follows: the yield equals 3.2 times 10¹¹ multiplied by the pinch current in megamperes raised to the power 4.5; and the yield equals 1.8 times 10¹⁰ multiplied by the peak current in megamperes raised to the power 3.8. Here the peak current ranges from 0.3 MA to 5.7 MA and the pinch current ranges from 0.2 MA to 2.4 MA.

Figure 9. Log of neutron yield scaling with log of peak current and log of pinch current, for the range of energies investigated, up to 25 MJ.

This compares to an earlier study carried out on several machines with published current traces and neutron yield measurements, operating conditions and machine parameters including the Chilean PF400J, the UNU/ICTP PFF, the NX2 and Poseidon providing slightly higher scaling laws: the yield proportional to the pinch current to the power 4.7 and to the peak current to the power 3.9. The slightly higher value of the scaling is because those machines fitted are of mixed radius ratio, mixed bank parameters, mixed model parameters and currents generally below 1 MA and voltages generally below 35 kV.

Summarising: over wide ranges of energy, optimizing pressure, anode length and radius, the scaling laws for neutron yield are: the yield equals 3.2 times 10¹¹ multiplied by the pinch current to the power 4.5, and 1.8 times 10¹⁰ multiplied by the peak current to the power 3.8, with peak current from 0.3 to 5.7 MA and pinch current from 0.2 to 2.4 MA; and the yield goes as the square of the storage energy at tens of kJ deteriorating to the storage energy to the power 0.84 at megajoule level, up to 25 MJ.

These laws provide useful references and facilitate the understanding of present plasma focus machines. More importantly, these scaling laws are also useful for design considerations of new plasma focus machines particularly if they are intended to operate as optimized neutron sources.

5. Insight 3 — Scaling Laws for Soft X-ray Yield

5.1. Computation of Neon Soft X-ray Yield — the Equations Used in the Computation

We note that the transition from Phase 4 to Phase 5 is observed in laboratory measurements to occur in an extremely short time with plasma and current disruptions resulting in localized regions of high densities and temperatures. These localized regions are not modeled in the code, which consequently computes only an average uniform density, and an average uniform temperature which are considerably lower than measured peak density and temperature. However, because the 4 model parameters are obtained by fitting the computed total current waveform to the measured total current waveform, the model incorporates the energy and mass balances equivalent, at least in the gross sense, to all the processes which are not even specifically modeled. Hence the computed gross features such as speeds and trajectories and integrated soft x-ray yields have been extensively tested in numerical experiments for several machines and are found to be comparable with measured values.

Equation 4, given here as a named result. In the code, the rate of neon line radiation is proportional to 4.6 times 10⁻³¹ multiplied by the square of the ion number density, the effective charge number, the fourth power of the atomic number, the square of the pinch radius and the pinch length, divided by the temperature. For the temperatures of interest in our experiments we take the soft X-ray yield to be that line radiation integrated over the pinch duration.

Hence the soft X-ray energy generated within the plasma pinch depends on the properties: number density, effective charge number, pinch radius, pinch length and temperature. It also depends on the pinch duration since in our code the line radiation is obtained by integrating over the pinch duration. This generated energy is then reduced by the plasma self-absorption which depends primarily on density and temperature; the reduced quantity of energy is then emitted as the soft X-ray yield. These effects are included in the modeling by computing a volumetric plasma self-absorption factor derived from the photonic excitation number, which is a function of the atomic number, the ion density, the effective charge and the temperature. However, in our range of operation, the numerical experiments show that the self absorption is not significant. It was first pointed out by Liu Mahe that a temperature around 300 eV is optimum for soft X-ray production. Shan Bing's subsequent work and our experience through numerical experiments suggest that around 2 times 10⁶ K (below 200 eV) or even a little lower could be better. Hence unlike the case of neutron scaling, for soft X-ray scaling there is an optimum small range of temperatures to operate.

5.2. Scaling Laws for Neon Soft X-rays from Numerical Experiments over a Range of Energies from 0.2 kJ to 1 MJ

We next use the Lee model code to carry out a series of numerical experiments over the energy range 0.2 kJ to 1 MJ. In this case we apply it to a proposed modern fast plasma focus machine with optimized values for the ratio of the outer to inner electrode radius and the static inductance obtained from our numerical experiments.

The following parameters are kept constant: (i) the ratio c = b/a, kept at 1.5, which is practically optimum according to our preliminary numerical trials; (ii) the operating voltage V0, kept at 20 kV; (iii) static inductance L0, kept at 30 nH, which is already low enough to reach the pinch current limitation regime over most of the range of storage energy we are covering; and (iv) the ratio of stray resistance to surge impedance, kept at 0.1, representing a higher performance modern capacitor bank. The model parameters fm, fc, fmr, fcr are also kept at fixed values 0.06, 0.7, 0.16 and 0.7. We choose the model parameters so they represent the average values from the range of machines that we have studied. A typical example of a current trace for these parameters is shown in Figure 10.

The storage energy is varied by changing the capacitance. Parameters that are varied are operating pressure, anode length and anode radius. Parametric variation at each storage energy follows the order pressure, anode length and anode radius until all realistic combinations are investigated. At each storage energy, the optimum combination is found that produces the biggest soft X-ray yield. In other words at each storage energy, a pressure is fixed, an anode length is chosen and the anode radius is varied until the largest yield is found. Then keeping the same values of storage energy and pressure, another anode length is chosen and the radius is varied until the largest yield is found. This procedure is repeated until for that storage energy and pressure, the optimum combination of length and radius is found. Then keeping the same value of storage energy, another pressure is selected. The procedure for parametric variation of length and radius as described above is then carried out for this storage energy and new pressure until the optimum combination is found.

The procedure is then repeated with a new value of storage energy. In this manner after systematically carrying out some 2000 runs, the optimized runs for various energies are tabulated in Table 2.

Figure 10. Computed total current versus time for L0 = 30 nH and V0 = 20 kV, C0 = 30 µF, stray resistance ratio 0.1, radius ratio 1.5 and model parameters fm, fc, fmr, fcr fixed at 0.06, 0.7, 0.16 and 0.7 for optimized anode radius 2.285 cm and anode length 5.2 cm.

Table 2. Optimized configuration found for each storage energy. Optimization carried out with stray resistance ratio 0.1, radius ratio 1.5, L0 = 30 nH and V0 = 20 kV and model parameters fm, fc, fmr, fcr fixed at 0.06, 0.7, 0.16 and 0.7, respectively. The speeds are the peak axial, radial shock and radial piston speeds, in centimetres per microsecond.

| E0 (kJ) | C0 (µF) | a (cm) | z0 (cm) | P0 (Torr) | Ipeak (kA) | Ipinch (kA) | va | vs | vp | Ysxr (J) | |---|---|---|---|---|---|---|---|---|---|---| | 0.2 | 1 | 0.58 | 0.5 | 4.0 | 100 | 68 | 5.6 | 22.5 | 14.9 | 0.44 | | 1 | 5 | 1.18 | 1.5 | 4.0 | 224 | 143 | 6.6 | 23.3 | 15.1 | 7.5 | | 2 | 10 | 1.52 | 2.1 | 4.0 | 300 | 186 | 6.8 | 23.6 | 15.2 | 20 | | 6 | 30 | 2.29 | 5.2 | 4.2 | 512 | 294 | 8.1 | 24.5 | 15.6 | 98 | | 10 | 50 | 2.79 | 7.5 | 4.0 | 642 | 356 | 8.7 | 24.6 | 15.7 | 190 | | 20 | 100 | 3.50 | 13 | 4.0 | 861 | 456 | 9.6 | 24.6 | 16.0 | 470 | | 40 | 200 | 4.55 | 20 | 3.5 | 1109 | 565 | 10.3 | 24.7 | 16.2 | 1000 | | 100 | 500 | 6.21 | 42 | 3.0 | 1477 | 727 | 11.2 | 24.8 | 16.4 | 2700 | | 200 | 1000 | 7.42 | 63 | 3.0 | 1778 | 876 | 11.4 | 24.8 | 16.5 | 5300 | | 400 | 2000 | 8.70 | 98 | 3.0 | 2079 | 1036 | 11.4 | 24.9 | 16.5 | 9400 | | 500 | 2500 | 9.10 | 105 | 2.9 | 2157 | 1086 | 11.5 | 25.1 | 16.7 | 11,000 | | 1000 | 5000 | 10.2 | 160 | 3.0 | 2428 | 1261 | 11.4 | 25.2 | 16.7 | 18,000 |

We then plot the soft X-ray yield against peak current and pinch current and obtain a yield scaling as the pinch current to the power 3.6 and as the peak current to the power 3.2. The pinch current scaling has less scatter than the peak current scaling. We next subject the scaling to further test when the fixed parameters and model parameters are varied. We add in the results of some numerical experiments using the parameters of several existing plasma focus devices including the UNU/ICTP PFF (stray resistance ratio 0.2, radius ratio 3.4, L0 = 110 nH and V0 = 14 kV with fitted model parameters fm = 0.05, fc = 0.7, fmr = 0.2, fcr = 0.8), the NX2 (stray resistance ratio 0.1, radius ratio 2.2, L0 = 20 nH and V0 = 11 kV with fitted model parameters fm = 0.10, fc = 0.7, fmr = 0.12, fcr = 0.68) and PF1000 (stray resistance ratio 0.1, radius ratio 1.39, L0 = 33 nH and V0 = 27 kV with fitted model parameters fm = 0.1, fc = 0.7, fmr = 0.15, fcr = 0.7). These new data points (white data points in Figure 11) contain wide ranges of radius ratio, voltage, static inductance and model parameters. The resulting soft X-ray yield versus pinch current log-log curve remains a straight line, with the scaling index 3.6 unchanged and with no more scatter than before. However the resulting yield versus peak current curve now exhibits considerably larger scatter and the scaling index has changed.

Figure 11. Soft X-ray yield plotted as a function of pinch current and peak current. The parameters kept constant for the black data points are: stray resistance ratio 0.1, radius ratio 1.5, L0 = 30 nH and V0 = 20 kV and model parameters fm, fc, fmr, fcr at 0.06, 0.7, 0.16 and 0.7 respectively. The white data points are for specific machines which have different values for the radius ratio, static inductance and voltage.

We would like to highlight that the consistent behavior of the pinch current in maintaining the scaling of the yield as the pinch current to the power 3.6, with less scatter than the peak current scaling, particularly when mixed-parameters cases are included, strongly supports the conclusion that pinch current scaling is the more universal and robust one. Similar conclusions on the importance of the pinch current in plasma focus performance and scaling laws have been reported.

It may also be worthy to note that our comprehensively surveyed numerical experiments for Mather configurations in the range of energies 0.2 kJ to 1 MJ produce a pinch current scaling rule for soft X-ray yield not compatible with Gates' rule. However it is remarkable that our pinch current scaling index of 3.6, obtained through a set of comprehensive numerical experiments over a range of 0.2 kJ to 1 MJ, on Mather-type devices, is within the range of 3.5 to 4 postulated on the basis of sparse experimental data, basically just two machines one at 5 kJ and the other at 0.9 MJ, by Filippov, for Filippov configurations in the range of energies 5 kJ to 1 MJ.

It must be pointed out that the results represent scaling for comparison with baseline plasma focus devices that have been optimized in terms of electrode dimensions. It must also be emphasized that the scaling with pinch current works well even when there are some variations in the actual device from L0 = 30 nH, V0 = 20 kV and radius ratio 1.5.

Summarising: over wide ranges of energy, optimizing pressure, anode length and radius, the scaling laws for neon soft X-rays are found to be a yield in joules of 8.3 times 10³ multiplied by the pinch current in megamperes to the power 3.6, and 600 multiplied by the peak current in megamperes to the power 3.2, where peak current runs from 0.1 to 2.4 MA and pinch current from 0.07 to 1.3 MA; and the yield goes as storage energy to the power 1.6 in the kilojoule range, deteriorating towards storage energy to the power 0.8 towards the megajoule range.

These laws provide useful references and facilitate the understanding of present plasma focus machines. More importantly, these scaling laws are also useful for design considerations of new plasma focus machines particularly if they are intended to operate as neon soft X-ray sources.

6. Insight 4 — Neutron Saturation

Besides being accurately descriptive and related to wide-ranging experimental reality, desirable characteristics of a model include predictive and extrapolative scaling. Moreover a useful model should be accessible, usable and user-friendly and should be capable of providing insights. Insight however cannot be a characteristic of the model in isolation, but is the interactive result of the model with the modeler or model user.

It was observed early in plasma focus research that the neutron yield goes as the square of the capacitor storage energy. Such scaling gave hopes of possible development as a fusion energy source. Devices were scaled up to higher storage energies. It was then observed that the scaling deteriorated, with the yield not increasing as much as suggested by the square-law scaling. In fact some experiments were interpreted as evidence of a neutron saturation effect as the storage energy approached several hundreds of kJ. As recently as 2006 Krauz and November 2007 Scholz have questioned whether the neutron saturation was due to a fundamental cause or to avoidable machine effects such as incorrect formation of plasma current sheath arising from impurities or sheath instabilities. We should note here that the region of discussion, several hundreds of kJ approaching the MJ region, is in contrast to the much higher energy region discussed by Schmidt at which there might be expected to be a decrease in the role of beam target fusion processes.

6.1. The Global Neutron Scaling Law

Recent extensive numerical experiments also showed that whereas at energies up to tens of kJ the square-law scaling held, deterioration of this scaling became apparent above the low hundreds of kJ. This deteriorating trend worsened and tended towards the storage energy to the power 0.8 at tens of MJ. The results of these numerical experiments are summarized with the solid line representing results from numerical experiments. Experimental results from 0.4 kJ to MJ, compiled from several available published sources, are also included as squares in the same figure. The combined experimental and numerical experimental results (Figure 12) appear to have general agreement particularly with regards to the square-law scaling at energies up to 100 kJ, and the deterioration of the scaling from low hundreds of kJ to the 1 MJ level. The global data of Figure 12 suggests that the apparently observed neutron saturation effect is overall not in significant variance with the deterioration of the scaling shown by the numerical experiments.

Figure 12. The global scaling law, combining experimental and numerical data. The global data illustrates neutron yield scaling deterioration observed in numerical experiments from 0.4 kJ to 25 MJ (solid line) using the Lee model code, compared to measurements compiled from publications (squares) of various machines from 0.4 kJ to 1 MJ.

6.2. The Cause of Neutron 'Saturation' is the Dynamic Resistance

A simple yet compelling analysis of the cause of this neutron saturation has been published. In Figure 1, the left side shows a schematic of the plasma dynamics in the axial phase of the Mather-type plasma focus. In that work the simplest representation was used, in which the current sheet is shown to go from the anode to the cathode perpendicularly. Observation shows that there is actually a canting of the current sheet and also that only a fraction, typically 0.7, of the total current participates in driving the current sheet. These points are accounted for in the modeling by model parameters fm and fc. We now represent the plasma focus circuit as shown in Figure 13.

Figure 13. Plasma focus circuit schematic. The capacitor bank with static inductance L0 and stray resistance r0 is switched into the plasma focus tube where a fraction fc of the circuit current effectively drives the plasma creating a time-varying inductance in the focus tube.

We consider only the axial phase. By surveying published results of all Mather-type experiments we find that all deuterium plasma focus devices operate at practically the same speeds and are characterized by a constancy of energy density (per unit mass) over the whole range of devices from the smallest sub-kJ to the largest MJ devices. The time varying tube inductance is the permeability of free space divided by two pi, multiplied by the natural logarithm of the radius ratio, multiplied by the axial position. The rate of change of inductance is therefore 2 times 10⁻⁷ multiplied by the logarithm of the radius ratio multiplied by the axial speed, in SI units. Typically on switching, as the capacitor discharges, the current rises towards its peak value, the current sheet is accelerated, quickly reaching nearly its peak speed and continues accelerating slightly towards its peak speed at the end of the axial phase. Thus for most of its axial distance the current sheet is travelling at a speed close to the end-axial speed. In deuterium the end-axial speed is observed to be about 10 cm per microsecond over the whole range of devices. This fixes the rate of change of inductance as 1.4 times 10⁻² henries per second for all the devices, if we take the radius ratio as 2. This value changes by at most a factor of 2, taking into account the variation of the radius ratio from low values of 1.4, generally for larger machines, to 4, generally for smaller machines. This typical rate of change of inductance may also be expressed as 14 milliohms.

We need now to inquire into the nature of the change in the tube inductance.

Equations 5 to 7, given here as named results. The induced voltage across a time-varying inductance is the current multiplied by the rate of change of inductance, plus the inductance multiplied by the rate of change of current. The instantaneous power delivered into that inductance is therefore the square of the current multiplied by the rate of change of inductance, plus the inductance multiplied by the current multiplied by the rate of change of current. The instantaneous power associated with the stored inductive energy, which is half the inductance multiplied by the square of the current, is only half of the square of the current multiplied by the rate of change of inductance, plus the same second term.

We note that the power associated with the stored inductive energy is not the same as the total power delivered. The difference, half the rate of change of inductance multiplied by the square of the current, is not associated with the inductive energy stored in the tube. We conclude that whenever the inductance changes with time, the instantaneous power delivered to it has a component that is not inductive. Hence this component of power must be resistive in nature; and the quantity of half the rate of change of inductance, also denoted as half Ldot, is identified as a resistance, due to the motion associated with the changing inductance; which we call the dynamic resistance. Note that this is a general result and is independent of the actual processes involved. In the case of the plasma focus axial phase, the motion of the current sheet imparts power to the shock wave structure with consequential shock heating, Joule heating, ionization, radiation and so on. The total power imparted at any instant is just that amount, with this amount powering all consequential processes. We denote the dynamic resistance of the axial phase as DR0.

We have thus identified for the axial phase of the plasma focus a typical dynamic resistance of 7 milliohms due to the motion of the current sheet at 10 cm per microsecond. It should be noted here that similar ideas of the role of the rate of change of inductance as a resistance were discussed by Bernard and colleagues. In that work the effect was discussed only for the radial phase. In our opinion the more important phase for the purpose of neutron saturation is actually the axial phase for the Mather-type plasma focus.

We now resolve the problem into its most basic form as follows. We have a generator, the capacitor charged to 30 kV, with a surge impedance equal to the square root of the static inductance divided by the capacitance, driving a load with a near constant resistance of 7 milliohms. We also assign a value for stray resistance of one tenth of the surge impedance. This situation may be shown in Table 3 where the static inductance is given a typical value of 30 nH. We also include in the last column the results from a circuit computation, discharging the capacitor with initial voltage 30 kV into a fixed resistance load of 7 milliohms simulating the effect of the axial dynamic resistance and a stray resistance of one tenth of the surge impedance.

Table 3. Discharge characteristics of equivalent plasma focus circuit, illustrating the 'saturation' of peak current with increase of storage energy to very large values. The last column presents results using a circuit computation, with a fixed resistance load of 7 milliohms, simulating the effect of the axial dynamic resistance and a stray resistance of one tenth of the surge impedance.

| E0 (kJ) | C0 (µF) | Z0 (milliohms) | DR0 (milliohms) | Total Z (milliohms) | Ipeak from V0 over Z (kA) | Ipeak from circuit computation (kA) | |---|---|---|---|---|---|---| | 0.45 | 1 | 173 | 7 | 197 | 152 | 156 | | 4.5 | 10 | 55 | 7 | 67 | 447 | 464 | | 45 | 100 | 17 | 7 | 26 | 1156 | 1234 | | 135 | 300 | 10 | 7 | 18 | 1676 | 1819 | | 450 | 1000 | 5.5 | 7 | 12.9 | 2321 | 2554 | | 1080 | 2400 | 3.5 | 7 | 10.8 | 2781 | 3070 | | 4500 | 10,000 | 1.7 | 7 | 8.8 | 3407 | 3722 | | 45,000 | 100,000 | 0.55 | 7 | 7.6 | 4209 | 4250 |

Plotting the peak current as a function of storage energy we obtain Figure 14, which shows the tendency of the peak current towards saturation as the storage energy reaches large values; the deterioration of the curve becoming apparent at the several hundred kJ level. This is the case for the peak current computed as the voltage divided by the total impedance and also for the circuit discharge with the simulated value of the axial dynamic resistance. In both cases it is seen clearly that a capacitor bank of voltage V0 discharging into a constant resistance such as the dynamic resistance will have a peak current approaching an asymptotic value of the voltage divided by that dynamic resistance when the bank capacitance is increased to such large values that the surge impedance becomes much smaller than the dynamic resistance. Thus the axial dynamic resistance causes current 'saturation'.

Figure 14. Peak current versus storage energy on log-log scale, illustrating peak current 'saturation' at large storage energy.

In Section 4.2 we had shown that the neutron yield goes as the pinch current to the power 4.5 (Equation 8) and as the peak current to the power 3.8 (Equation 9). Hence saturation of the peak current will lead to saturation of the neutron yield.

At this point we note that if we consider that only 0.7 of the total current takes part in driving the current sheet, as typically agreed upon from experimental observations, then there is a correction factor which reduces the axial dynamic resistance by some 40%. That would raise the asymptotic value of the current by some 40%; nevertheless there would still be 'saturation'.

Thus we have shown that current 'saturation' is inevitable as the storage energy is increased to very large values by an increase in capacitance, simply due to the dominance of the axial phase dynamic resistance. This makes the total circuit impedance tend towards an asymptotic value which approaches the dynamic resistance at infinite values of storage energy. The 'saturation' of current inevitably leads to a 'saturation' of neutron yield. Thus the apparently observed neutron 'saturation' which is more accurately represented as a neutron scaling deterioration is inevitable because of the dynamic resistance. In line with current plasma focus terminology we will continue to refer to this scaling deterioration as 'saturation'. The above analysis applies to the Mather-type plasma focus. The Filippov-type plasma focus does not have a clearly defined axial phase. Instead it has a lift-off phase and an extended pre-pinch radial phase which determine the value of the peak current. During these phases the inductance of the Filippov discharge is changing, and the changing inductance will develop a dynamic resistance which will also have the same current 'saturation' effect as the Filippov bank capacitance becomes big enough.

6.3. Beyond Presently Observed Neutron Saturation Regimes

Moreover the 'saturation' as observed in presently available data is due also to the fact that all tabulated machines operate in a narrow range of voltages of 15 to 50 kV. Only the SPEED machines, most notably SPEED II, operated at low hundreds of kV. No extensive data have been published from the SPEED machines. Moreover SPEED II, using Marx technology, has a large bank surge impedance of 50 milliohms which itself would limit the current. If we operate a range of such high voltage machines at a fixed high voltage, say 300 kV, with ever larger storage energy until the surge impedance becomes negligible due to the very large value of capacitance, then the 'saturation' effect would still be there, but the level of 'saturation' would be proportional to the voltage. In this way we can go far above presently observed levels of neutron 'saturation'; moving the research, as it were, into presently beyond-saturation regimes.

Could the technology be extended to 1 MV? That would raise the peak current to beyond 15 MA and the pinch current to over 6 MA. Also multiple Blumleins at 1 MV, in parallel, could provide driver impedance of 100 milliohms, matching the radial phase dynamic resistance and provide fast rise currents peaking at 10 MA with a pinch current value of perhaps 5 MA. Bank energy would be several MJ. The push to higher currents may be combined with proven neutron yield enhancing methods such as doping deuterium with a low percentage of krypton. Further increase in pinch current might be by fast current injection near the start of the radial phase. This could be achieved with charged particle beams or by circuit manipulation such as current-stepping. The Lee model is ideally suited for testing circuit manipulation schemes.

7. Conclusions

This paper has reviewed the extensive and systematic numerical experiments which have been used to uncover new insights into plasma focus fusion devices including the following. A plasma current limitation effect was unexpectedly found, as the static inductance of any focus device is reduced towards very small values. Scaling laws of neutron yield and soft X-ray yield as functions of storage energies, circuit peak current as well as plasma pinch current, were developed over a wider range of parameters than attempted previously. A global scaling law for neutron yield as a function of storage energy was uncovered combining experimental and extensive numerical data, showing that scaling deterioration has been wrongly interpreted as neutron 'saturation'. However in keeping with conventional terminology, the effect of scaling deterioration will continue to be referred to as neutron 'saturation'. The cause of neutron 'saturation' as device storage energy is increased was found to be the axial phase 'dynamic resistance'. With the fundamental cause discovered, it is suggested that beyond 'present saturation' regimes may be reached by going to higher voltages, and using plasma current enhancement techniques such as current-steps. It is expected that numerical experiments will continue to play a major role complementing and even guiding laboratory measurements and practices. The ground-breaking insights thus gained will completely open up the directions of plasma focus fusion research.

(The fifty-seven references are at the source.)

The way in

https://doi.org/10.3390/en3040711TEXT. Full text below, abstract to conclusions, from the publisher PDF. The fourteen figures are not reproduced; their captions are kept because they carry the machine parameters each run used, and the fifteen-panel Figure 4 is summarised by its caption alone. Tables 1, 2 and 3 are reset as Markdown tables from the extracted columns. Equations 1 to 9 did not survive extraction cleanly and are given here as named results in words; the exact expressions are at the source. The micro sign and the ohm symbol were dropped by the extraction in several places and are restored here — speeds are in centimetres per microsecond, capacitances in microfarads, resistances in milliohms. The reference list is not reproduced.

How to cite it

Sing Lee, Sor Heoh Saw (2010) Numerical Experiments Providing New Insights into Plasma Focus Fusion Devices. doi:10.3390/en3040711

Where it sits in the curriculum

Plasmoids, charge clusters and the orbsLattice confinement fusion

Provenance: Retrieved 2026-09-08 · Summary by The Spacetime Metric editorial rail (AI draft from the source text, 2026-09-07)← The library