Method for Estimating the Heat Load Imposed on a Cryogenic Refrigerator, Associated Program Product, and Method for Controlling the Refrigerator

Abstract

The invention relates to a method for estimating a heat load imposed on a cryogenic refrigerator, to an associated computer program product, and to a method for controlling the cooling power output by said refrigerator. As the refrigerator (1, 1′) includes a phase separator (40, 40′) comprising a bath (41, 41′) of refrigerant, the method for estimating the heat load imposed on said refrigerator includes a step in which said heat load is estimated using a program executed by a computer, said program being based on a mass balance carried out on the phase separator for expressing variations in the time drift of the height of the bath of refrigerant in the phase separator.

Description

The present invention relates to the field of plant refrigeration.

The present invention in particular relates to the cryogenic refrigeration of plants.

The present invention also relates to the cryogenic refrigeration of plants capable of operating in a variable regime.

A variable operating regime may be encountered in many applications.

This is for example the case in plants using strong magnetic fields.

An example of a plant using strong magnetic fields is a tokamak. A tokamak is a magnetic confinement chamber intended to control a plasma in order to study the possibility of power generation by nuclear fusion.

In tokamaks, a variable, pulsed operating regime may be employed. In this case, although the overall profile of the pulsed loads applied to the refrigerator is known, it is not exactly known when a pulse of load will occur. In addition, other unforeseeable perturbations may occur, which perturbations are related to the operation of the tokamak.

To generate strong magnetic fields without destroying the electromagnets, superconducting electromagnets are used. For an electromagnet to operate as a superconductor, its temperature must be kept below its critical temperature.

This is achieved by virtue of a cryogenic refrigerator, for example cooled by helium.

FIG. 1 is a schematic showing the operation of such a cryogenic refrigerator.

The cryogenic refrigerator 1 comprises a compressor 10 allowing a gas, in this case helium, at room temperature (T₀≈300 K) to be compressed from atmospheric pressure P_a to a pressure P of about 20 bar.

A number of heat exchangers 20, 21, 22, 23, 24 are placed in parallel upstream of a phase separator 40 comprising a bath of liquid helium at a temperature T_f=4.5 K. In this case, five heat exchangers have been provided.

These countercurrent heat exchangers allow the temperature of the helium flowing in the various circuits of these heat exchangers to be gradually decreased. Moreover, for each heat exchanger, the pressure in the two exchange circuits of the exchanger differ, so that the “upstream” circuit is a hot, high-pressure circuit and the “downstream” circuit is a colder, low-pressure circuit. Between the “upstream” circuit of the first heat exchange 20 and the “downstream” circuit of the last heat exchanger 24, the pressure thus increases from P=20 bar at the outlet of the compressor to a pressure slightly higher than atmospheric pressure, to prevent problems with low-pressure cavitation.

The bath 41 of liquid helium is therefore at atmospheric pressure.

The cryogenic refrigerator 1 also comprises a number of means for extracting work.

In the present case, these means consist of two turbines 30, 31.

The first turbine 30 has work done on it at the outlet of the low-pressure circuit of the first heat exchanger 20 and reinjects this work at the low-pressure inlet of the second heat exchanger 21. The second turbine 31 has work done on it at the outlet of the low-pressure circuit of the third heat exchanger 22 and reinjects this work at the low-pressure inlet of the fourth heat exchanger 23.

These turbines 30, 31 are complementary to the heat exchangers and participate, via the work done on them, to the cooling of the helium.

Lastly, the cryogenic refrigerator 1 comprises a Joule-Thomson valve 50, placed between the outlet of the low-pressure circuit of the last heat exchanger 24 and the bath 41 of liquid helium at 4.5 K and atmospheric pressure. This valve 50 liquefies the gaseous helium obtained at the outlet of the low-pressure circuit of the last heat exchanger 24 via an expansion which is accompanied by a drop in the temperature of the helium.

The bath 41 of liquid helium then allows cooling power to be delivered in order to keep the electromagnets of the plant employing strong magnetic fields, for example a tokamak, operating as superconductors, on account of the heating power applied by the plant to this bath 41 of liquid helium.

The heating power applied by the plant to the bath 41 of liquid helium is also called the thermal load.

The design and dimensions of this type of refrigerator are tailored to operating regimes in which the strong magnetic field generated by the plant is stable or varies slowly, i.e. to a regime of permanent or almost permanent thermal operation of the refrigerator.

Specifically, these plant operating regimes lead to a stable thermal load on the cryogenic refrigerator. Operating the cryogenic refrigerator in a permanent or almost permanent thermal regime makes it possible to provide adequate cooling of the electromagnets, stably and reliably.

However, future plants intended for studying the possibility of generating power by nuclear fusion activated by strong magnetic fields plan to employ pulsed magnetic fields.

This is the case for the ITER (France) and JT60_SA (Japan) tokamak projects.

Variations in the magnetic field in the tokamak then cause similar variations in the thermal load applied to the cryogenic refrigerator.

An example of the variation in the thermal load applied to the cryogenic refrigerator of the future Japanese tokamak JT60_SA, intended to operate in a pulsed regime, is shown in FIG. 2. It will be noted that the pulses are a priori periodic, in this instance having a period of 3000 s, the duration of a pulse being about 100 s. Nevertheless, random processes may change this variation, as was mentioned above.

In FIG. 2, it may be seen that operating the tokamak in a pulsed regime leads to the operating regime of the thermal load applied to the cryogenic refrigerator also being pulsed (curve 100). FIG. 2 also shows (curve 101) the average thermal load applied by the tokamak to the refrigerator.

However, cryogenic refrigerators used in existing plants are not designed to provide adequate cooling under such a pulsed regime.

Specifically, increasing the thermal load leads to an increase in the flow rate of the cooled helium, returned to the electromagnets of the plant.

This cools the entire refrigerator, because the cooling engendered by the thermal load unbalances the exchangers 20, 21, 22, 23 and 24 between the high-pressure section and low-pressure section. In addition, the evaporation caused by the load on the helium bath instantaneously increases the return flow rate on the low-pressure side of each heat exchanger, thereby unbalancing the entire cryogenic refrigerator. A substantial increase in the thermal load may even cause the cryogenic refrigerator to shut down.

To overcome this problem, it has already been suggested to smooth the impact of the pulsed magnetic field on the variation of the thermal load applied to the cryogenic refrigerator.

This smoothing consists in limiting variations in the thermal load actually applied to the cryogenic refrigerator, in order to ensure the nominal operation of the cryogenic refrigerator and therefore especially to prevent the refrigerator from shutting down.

For this purpose, it has been suggested to implement mechanical and/or thermal methods, either by installing dedicated means inside the cryogenic refrigerator, or by installing an additional device between the plant employing strong magnetic fields and its cryogenic refrigerator.

For example, the document by Dauguet et al., “Advances in Cryogenic Engineering: Transactions of the Cryogenic Engineering Conference”, CEC Vol. 53, edited by J. G. Weisend II, pp. 564-569 proposes keeping the thermal load on the cryogenic refrigerator constant.

To do this, the cryogenic refrigerator comprises many heat exchangers and cryogenic valves, these valves being activated in order to keep the cooling power delivered by the refrigerator stable at a value corresponding to the average thermal load applied by the tokamak. Thus, the cryogenic refrigerator operates in a permanent or almost permanent regime, even when the tokamak is operating in a pulsed regime.

A similar solution is proposed in document WO 2009/024705.

According to another example, document “Design of the ITER-FEAT cryoplant to achieve stable operation over a wide range of experimental parameters and operation scenarios”, by Claudet et al., Fusion Engineering and Design, 58-59 (2001), pp. 205-209, suggests installing an intermediate device between the tokamak and the cryogenic refrigerator.

This device allows some of the additional helium flow obtained during a peak in the thermal load applied to the cryogenic refrigerator to be diverted upstream of the refrigerator. Thus, the cryogenic refrigerator does not see the increased helium flow related to the pulsed operating regime of the tokamak.

These proposed solutions are based on mechanical and/or thermal devices intended to smooth the variation in the thermal load liable to be applied to the cryogenic refrigerator.

These solutions work correctly.

However, they require additional parts (heat exchangers, valves, etc.) which may quickly prove to be expensive and they are, sometimes, difficult to implement.

One objective of the invention is to solve at least one of the drawbacks of existing cryogenic refrigerators.

Another objective is to solve at least one of the drawbacks of existing cryogenic refrigerators when the plant to be cooled operates in a variable regime.

As was mentioned above, a variable regime may be encountered in many fields. Therefore, use of the invention is not limited to plants generating strong pulsed magnetic fields, such as a tokamak, but extends to any plant requiring cryogenic refrigeration.

To achieve at least one of these objectives, the invention provides a method for estimating a thermal load applied to a cryogenic refrigerator comprising a phase separator containing a bath of liquid refrigerant, in which this thermal load is estimated using a computer program, said program being based on a mass balance performed on the phase separator allowing the variation in the time derivative of the height of the bath of liquid refrigerant in the phase separator to be expressed.

The method according to the invention will possibly have other technical features, whether alone or in combination:

• • the thermal load applied to the cryogenic refrigerator is a variable thermal load;
  • the variable thermal load applied to the cryogenic refrigerator is pulsed;
  • a step is provided for determining appropriate variables of the cryogenic refrigerator, representing the variation in the time derivative of the height of the bath of liquid refrigerant in the phase separator, before the step of estimating the thermal load using the computer program;
  • the cryogenic refrigerator comprising a valve at the inlet of the phase separator, said appropriate variables comprise at least the degree of opening of the valve, the temperature upstream of the valve and their respective time derivatives, and the thermal load; and
  • the computer program comprises a step in which the data relating to the variation in the time derivative of the height of the bath of liquid refrigerant in the phase separator are filtered.

To achieve at least one of these objectives, the invention also provides a method for regulating the cryogenic refrigerator subjected to a thermal load, in which the thermal load applied to the refrigerator is estimated using the method for estimating a thermal load applied to this refrigerator, according to the invention, and then at least one operating parameter of the refrigerator is regulated depending on the value of the thermal load estimated beforehand.

The regulating method according to the invention will possibly have other technical features, in particular:

• • a step in which the regulation is implemented by modifying the degree of opening of the valve of the refrigerator, said valve being located at the inlet of the phase separator.

To achieve at least one of these objectives, the invention also provides a computer program product comprising programming code instructions for implementing the method for estimating a thermal load applied to a cryogenic refrigerator, according to the invention.

Other features, aims and advantages of the invention will become apparent from the following detailed description, given with reference to the following figures:

FIG. 3 is a schematic showing a mass balance performed on the phase separator of the cryogenic refrigerator, for example as illustrated in FIG. 1 or 4;

FIG. 4 is a schematic of the cryogenic refrigerator named “station 800W@4.5K” and installed at the Commissariat à l'énergie atomique, on which the method according to the invention was tested and validated;

FIG. 5 is a graph showing, on the one hand, the variation in the thermal load actually applied to the cryogenic refrigerator shown in FIG. 4 during a trial, and on the other hand, the variation of the thermal load estimated using the method according to the invention.

The method especially comprises a step of estimating the thermal load w applied to the cryogenic refrigerator 1, 1′ using a computer program.

Details of the development of this program are given below.

This program especially operates on a mass conservation balance performed on the phase separator 40, 40′ of the cryogenic refrigerator 1, 1′.

This mass balance performed on the phase separator 40, 40′ is shown schematically in FIG. 3. Typically, it is a phase separator 40, 40′ such as shown in FIGS. 1 and 4, respectively, which comprises a bath 41, 41′ of liquid refrigerant, for example of helium.

The mass conservation equation for the phase separator 40, 40′ is written:

$\begin{array}{cc}\frac{h}{t}=m_{+}-m_{-}=f\left(u_1,P,T\right)-\frac{w}{L_v}& \left(\mathrm{Eq}.1\right)\end{array}$

where:

h is the height of the bath 41, 41′ of liquid refrigerant in the phase separator 40, 40′ (% of the maximum attainable height h_max of the bath);

m₊ is the flow rate of gas entering into the phase separator (g/s);

m₋ is the flow rate of gas leaving the phase separator (g/s);

w is the thermal load applied to the phase separator, i.e. to the refrigerator, by the plant (W);

u₁ is the degree of opening of the valve 50, 50′ (%);

L_V is the latent vaporisation heat (J/g); and

f(u₁, P, T) is a function depending on the degree of opening u₁ of the value 50, 50′, on the pressure P upstream of the valve and on the temperature T upstream of this value.

The function f(u₁, P, T) may take various developed forms depending on the precision desired for the estimation of the time derivative

$\frac{h}{t}$

of the height h of the liquid refrigerant in the phase separator 40, 40′, which is denoted {dot over (h)} in the following.

In the example given below, the quantity {dot over (h)} depends on the following parameters: u₁, T·u₁, {dot over (u)}₁, {dot over (T)} and w, where {dot over (u)}₁ and {dot over (T)} are the time derivatives of the degree of opening of the valve 50, 50′ and of the temperature upstream of this valve, respectively. These parameters represent the main parameters influencing the mass balance performed on the phase separator 40, 40′.

In this example, the influence of the pressure P in the function f is not taken into consideration.

Different parameters may be chosen, depending on the desired precision of the model, in particular parameters relating to the pressure P could be incorporated.

According to another example, the model takes into account, in order to define the quantity {dot over (h)}, the following parameters: u₁, T, and their respective time derivatives, and the thermal load w.

It is then sought to express the equation (Eq. 1) in a linear form, i.e. in the form of (Eq. 2):

$\stackrel{.}{h}=\left(u_1,\left({\mathrm{Tu}}_1\right),{\stackrel{.}{u}}_1,\stackrel{.}{T},w\right)·\left(\begin{array}{c}{a}_1\\ a_2\\ a_3\\ a_4\\ a_5\end{array}\right)$

where: a₁, a₂, a₃, a₄ and a₅ are coefficients to be determined.

The equation (Eq. 2) may be written in the form of the equation (Eq. 3):

{dot over (h)}=L _q q+L _z ż+a ₅ w

on account of the following notations:

$q:=\left(\begin{array}{c}{u}_1\\ {\mathrm{Tu}}_1\end{array}\right);$ $z=\left(\begin{array}{c}{u}_1\\ T\end{array}\right);$ $L_q\left(a_1,a_2\right);$ $L_z\left(a_3,a_4\right)$

To estimate the thermal load w applied to the cryogenic refrigerator 1, 1′ by the plant, the quantity of interest is the filtered time derivative {dot over (h)} of the height h of liquid refrigerant in the phase separator, which is expressed in a Laplace field.

Filtering allows any “noise” that is capable of influencing the determination of the variable {dot over (h)} to be removed.

Thus, if the time derivative {dot over (h)} of this height h is denoted H(s) in the Laplace field, and the function obtained after first-order filtering of the quantity H(s) is denoted Y(s), then a priori their relationship can be expressed in the form of the equation (Eq. 4):

$Y\left(s\right)=\frac{s}{1+{\tau }_fs}·H\left(s\right)$

where τ_f is the time constant of this first-order filter, which must be defined.

The time constant is defined in the following way. The valve 50, 50′ is opened with a given opening profile, a stepped profile for example. Then, the variation of the time derivative of the height of the bath of liquid is monitored. It is then possible to determine the time constant in a way known per se to those skilled in the art.

Equation (Eq. 4) then corresponds, in the real field, to the following differential equation (Eq. 5):

$\stackrel{.}{y}=\frac{1}{{\tau }_f}\left[-y+\stackrel{.}{h}\right]$

where: {dot over (y)} is the real variable associated with the filtered function Y(s) of equation (Eq. 4) in the Laplace field, and y is the real integrated value of the variable {dot over (y)}.

By inserting the relationship given by equation (Eq. 3) obtained from the mass conservation balance performed on the phase separator 40, 40′ into equation (Eq. 5) above, equation (Eq. 6) is obtained:

$\stackrel{.}{y}=\frac{1}{{\tau }_f}\left[-y+L_qq+L_z\stackrel{.}{z}+a_5w\right]$

It will be understood that this differential equation expresses the mass conservation balance performed on the phase separator after first-order filtering has been carried out on the variable {dot over (h)}.

Next, the state vector x_w is introduced in the form of the equality (Eq. 7):

x _w:=τ_f·y−L_z·z

This notation then allows equation (Eq. 6) to be written in the form of the following equation (Eq. 8):

{dot over (x)} _w =−y+L _q q+a ₅ w

Next, by replacing the variable y with the variable

$\frac{1}{{\tau }_f}\left[x_w+L_zz\right]$

the following equation is obtained:

${\stackrel{.}{x}}_w=-\frac{1}{{\tau }_f}\left[x_w+L_zz\right]+L_qq+a_5w$

which finally gives the following state representation (Eq. 9):

${\stackrel{.}{x}}_w=\left[\frac{-1}{{\tau }_f}\right]x_w+\left(L_q,\frac{-L_z}{{\tau }_f}\right)\left(\begin{array}{c}q\\ z\end{array}\right)+a_5w$ $y=\left[\frac{1}{{\tau }_f}\right]x_w+\left[\left(0,\frac{L_z}{{\tau }_f}\right)\right]\left(\begin{array}{c}q\\ z\end{array}\right)$

where x_w is the state and q, z and w the inputs and y the output.

For the sake of clarity, the state representation (Eq. 9) may also be expressed in the following form (Eq. 10):

${\stackrel{.}{x}}_w={\mathrm{Ax}}_w+B\left(\begin{array}{c}q\\ z\end{array}\right)+\mathrm{Gw}$ $y={\mathrm{Cx}}_w+D\left(\begin{array}{c}q\\ z\end{array}\right)$

on account of the following matrix notations:

$A=\left[\frac{-1}{{\tau }_f}\right];$ $B=\left(L_q,\frac{-L_z}{{\tau }_f}\right);$ $C=\left[\frac{1}{{\tau }_f}\right];$ $D=\left[\left(0,\frac{L_z}{{\tau }_f}\right)\right];$ $G=a_5$

The state representation (Eq. 9/Eq. 10) finally allows the filtered variable y representative of the height h of the liquid refrigerant in the phase separator to be expressed as a function of q, z, w and the coefficients a₁, a₂, a₃, a₄, a₅ and τ_f.

This state representation can be easily implemented in a programmable controller.

This state representation may be extended to estimate the thermal load applied to the cryogenic refrigerator.

To do this, the extended state of the state representation (Eq. 9) is defined by the following matrix, (Eq. 11):

$\xi :=\left(\begin{array}{c}{x}_w\\ s_w\end{array}\right)$

where the second state vector s_w represents the thermal load w.

Then, the state representation (Eq. 10) can be expressed in the following form, denoted (Eq. 12):

$\stackrel{.}{\xi }:=\left(\begin{array}{c}{\stackrel{.}{x}}_w\\ {\stackrel{.}{s}}_w\end{array}\right)=\underset{\underset{\stackrel{\_}{A}}{}}{\left(\begin{array}{cc}A& G\\ 0& 0\end{array}\right)}\left(\begin{array}{c}{x}_w\\ s_w\end{array}\right)+\left(\begin{array}{c}B\\ 0_{1\times 4}\end{array}\right)·\left(\begin{array}{c}q\\ z\end{array}\right)$ $y=\underset{\underset{\stackrel{\_}{C}}{}}{\left(\begin{array}{cc}C& 0\end{array}\right)}\xi +D\left(\begin{array}{c}q\\ z\end{array}\right)$

This representation assumes that {dot over (s)}_w=0, i.e. that the variable s_w representative of the thermal load w applied to the refrigerator has a constant profile. This assumption is correct when the perturbation that arrives at the bath of liquid takes a step form. This is especially the case for a tokamak.

Of course, other assumptions could be made. For example, it may be assumed that the perturbation arrives in the form of a ramp or a sinusoidal variation. This could be the case for cryogenic refrigerators used in plants other than tokamaks.

If (Eq. 13) is written:

$\left(\begin{array}{c}q\\ z\end{array}\right):=\left(\begin{array}{c}{u}_1\\ {\mathrm{Tu}}_1\\ u_1\\ T\end{array}\right)\text{=:}U\left(T,u_1\right)$

then the state representation of (Eq. 12) can be written in the form of (Eq. 14):

{dot over (ξ)}=[Ā]·ξ+[ B]·U(Y,u₁)

y=[ C]·ξ+[D]·U(T,u₁)

Next, if discrete matrices representative of the matrices Ā and B for a sampling period T_e are denoted A_d and B_d respectively, the above relationships can be written:

ξ⁺=[A_d]·ξ+[B_d]·U(T,u₁)

y=[C]·ξ+[D]·U(T,u₁)

where ξ⁺ is the discrete form of the matrix ξ.

It will be recalled that when the continuous matrices A, B, C, D, of a state representation are known, the associated discrete matrices A_d, B_d, C_d and D_d are written:

• A_d=e^A·T; where T_e is the sampling period;

$B_d={\int }_0^{T_e}{}^{A·s}·B·s;$

• C_d=C; and
• D_d=D.

The relationships expressing ξ⁺ and y can be written in shortened form.

To do this, the observer equation will be recalled, namely:

{circumflex over (ξ)}⁺=A_d{circumflex over (ξ)}+B_dU+L(y−ŷ)

where L is a two-row column matrix comprising real values called the observer gain. The observer gain only depends on the response time t_r chosen for the observer.

Next, the relationship expressing y is used for the variable ŷ. In this case ŷ can be written:

ŷ=[ C]·{circumflex over (ξ)}+[D]·U(T,u₁).

By replacing ŷ with the latter expression in the observer equation, the shortened form mentioned above is finally obtained, namely equation (Eq. 15):

{circumflex over (ξ)}⁺=[A_d−L C]·{circumflex over (ξ)}+[B_d−LD]·U(T,u₁)+L·y

in which the observer gain L is chosen so that the eigenvalues of the matrix A_d−L C are in part real and strictly negative. As should be understood with ease by those skilled in the art, the eigenvalues of the matrix A_d−L C correspond to the poles of the corresponding transfer function.

These poles are especially related to the time constant τ_f, insofar as this parameter appears in the expression of this matrix. Thus, these poles are related to the response time t_r chosen for the observer.

However, it must be noted that the variable y is not present in equation (Eq. 15). It must therefore be determined in some other way. To do this, the differential equation (Eq. 5) is used, namely:

$\stackrel{.}{y}=\frac{1}{{\tau }_f}\left[-y+\stackrel{.}{h}\right].$

Next, the vector:

η:=τ_f·y−h

is defined and equation (Eq. 5) is rewritten in the following form:

$\stackrel{.}{\eta }=-y=\frac{-1}{{\tau }_f}\left[\eta +h\right]$ $y=\frac{1}{{\tau }_f}\left[\eta +h\right]$

which may be expressed in the form of the following state representation:

$\stackrel{.}{\eta }=\left[\frac{-1}{{\tau }_f}\right]·\eta +\left[\frac{-1}{{\tau }_f}\right]·h$ $y=\left[\frac{1}{{\tau }_f}\right]·\eta +\left[\frac{1}{{\tau }_f}\right]·h.$

This state representation may then be expressed in discrete form, for a sampling period T_e, in the following way (Eq. 16):

${\eta }^{+}=\left[R_d\right]·\eta +\left[E_d\right]·h$ $y=\left[\frac{1}{{\tau }_f}\right]·\eta +\left[\frac{1}{{\tau }_f}\right]·h$

where η⁺ is the discrete form of the variable {dot over (η)}, the matrices [R_d] and [E_d] being the discrete matrices of the state representation of equation (Eq. 5). These matrices are defined, as a function of the continuous matrices [R]=[−1/τ_f] and E=[−1/τ_f], by the following relationships:

• R_d=e^R·Te where T_e is the sampling period; and

$E_d={\int }_0^{T_e}{}^{R·s}·E·s.$

By then grouping equations (Eq. 15) and (Eq. 16), a set of equations (Eq. 17) may be written:

${\eta }^{+}=\left[R_d\right]·\eta +\left[E_d\right]·h$ ${\hat{\xi }}^{+}=\left[A_d-L\stackrel{\_}{C}\right]·\hat{\xi }+\left[B_d-\mathrm{LD}\right]·U\left(T,u_1\right)+L·\left(\left[\frac{1}{{\tau }_f}\right]·\eta +\left[\frac{1}{{\tau }_f}\right]·h\right).\hat{w}=\left(0,1\right)\hat{\xi }$

This set of equations may be expressed in the following form (Eq. 18):

$X^{+}=\left(\begin{array}{cc}{R}_d& 0_{1\times 2}\\ \frac{1}{{\tau }_f}L& A_d-L\stackrel{\_}{C}\end{array}\right)X+\left(\begin{array}{cc}{0}_{1\times 4}& E_d\\ B_d-\mathrm{LD}& \frac{1}{{\tau }_f}L\end{array}\right)\left(\begin{array}{c}U\left(T,u_1\right)\\ h\end{array}\right)$ $\hat{w}=\left(0,0,1\right)X$

if the following notation

$X:=\left(\begin{array}{c}\eta \\ \hat{\xi }\end{array}\right)$

is used; or, more simply (Eq. 19):

X ⁺ =A _obs X+B _obs U _obs(T,u₁,h)

ŵ=C_obsX

with the following matrices:

$A_{\mathrm{obs}}=\left(\begin{array}{cc}{R}_d& 0_{1\times 2}\\ \frac{1}{{\tau }_f}L& A_d-L\stackrel{\_}{C}\end{array}\right);$ $B_{\mathrm{obs}}=\left(\begin{array}{cc}{0}_{1\times 4}& E_d\\ B_d-\mathrm{LD}& \frac{1}{{\tau }_f}L\end{array}\right);$ $C_{\mathrm{obs}}=\left(0,0,1\right);$ $\mathrm{and}$ $U_{\mathrm{obs}}\left(T,u_1,h\right)=\left(\begin{array}{c}{u}_1\\ {\mathrm{Tu}}_1\\ u_1\\ T\\ h\end{array}\right).$

Equation (Eq. 19) thus allows the thermal load w applied to the refrigerator to be estimated, in a form that can be easily implemented in a programmable controller.

However the thermal load w can be estimated, using (Eq. 19), only once the matrices A_obs and B_obs have been completely defined, which requires that the coefficients a₁ to a₅ be identified, the values of the matrix L be determined and the time constant τ_f of the first-order filter be chosen.

The coefficients a₁ to a₅ are identified in the following way.

The equation (Eq. 2) is integrated between a reference time t₁ and the time t, and then filtered.

The integration returns equation (Eq. 20):

h(t)−h(t₁)=M(t)*a

where M(t) is the line matrix given by:

$M\left(t\right)=\left(\begin{array}{ccccc}{\int }_{t_1}^tu\left(\tau \right)\tau & {\int }_{t_1}^tT\left(\tau \right)u_1\left(\tau \right)\tau & u_1\left(t\right)-u\left(t_1\right)& T\left(t\right)-T\left(t_1\right)& {\int }_{t_1}^tw\left(\tau \right)\tau \end{array}\right),$

and a is the column vector of the coefficients a₁ to a₅.

Equation (Eq. 20) is then filtered with the following filter:

$F\left(s\right)=\frac{s}{1+{\tau }_fs}.$

The coefficients a₁ to a₅ are then determined by minimizing, using the least squares method, the relationships forming a system of equations (Eq. 21), the unknowns of which are the coefficients of the matrix a:

$\left(\begin{array}{c}F\left(M\right)\left(t_1\right)\\ F\left(M\right)\left(t_2\right)\\ ⋮\\ F\left(M\right)\left(t_N\right)\end{array}\right)·a-\left(\begin{array}{c}0\\ F\left(h\right)\left(t_2\right)-h\left(t_1\right)\\ ⋮\\ F\left(h\right)\left(t_N\right)-h\left(t_1\right)\end{array}\right)$

which allows equation (Eq. 20) to be solved in its filtered form between the initial time t₁ and the final time t_N, where N=5 in this instance.

To do this, it is necessary to measure the values of the coefficients of the other matrices.

These measurements were carried out on the refrigerator 1′ shown schematically in FIG. 4.

This refrigerator 1′ is similar to the refrigerator 1 described with reference to FIG. 1, especially as regards the phase separator 40′. This refrigerator 1′ however differs slightly from the refrigerator described with reference to FIG. 1 in various ways.

Specifically, it only has four heat exchangers 20′, 21′, 22′, 23′ for lowering the temperature of the refrigerant, in this case helium, from 300 K to 4.5 K. The pressure delivered by the compressor is 16 bar. Moreover, only one turbine 31′ is provided and the first heat exchanger 20′ contains an additional heat exchanger 60′ (nitrogen, 80 K).

The coefficients of the matrix a were calculated for a compressor outlet pressure P of 16 bar. Moreover, the temperature upstream of this valve 50′ is T=7 K. The time t₁ considered is t₁=1 s, at which point data recording begins. The interval between measurements is 3 s, i.e. dt=t_i+1−t_i=3 s, for i ranging from 1 to N.

Under these conditions, solving the system of equations (Eq. 21) returns the coefficients in table 1 below:

TABLE 1
Coefficient of the matrix a value
a1                          0.003710768747845
a2                          −0.000396506514909
a3                          0.024617305832635
a4                          0.286747643082885
a5                          −0.000240828045189

Once the coefficients a₁ to a₅ have been obtained, they are entered into the program, as are the time constant τ_f and the coefficients of the observer gain L.

A 5% response time, denoted tr5%, was chosen for the observation. In this instance, tr5%=60 s. This indicates that the observer is delayed in order to obtain a 5% estimate of the thermal load w in 60 s.

The time constant τ_f is determined using the relationship tr5%=3τ_f. Thus, the time constant τ_f has a value τ_f=20 s.

Moreover, to determine the eigenvalues of the matrix A_d−L C, the poles z of the corresponding transfer function are determined, i.e.

$z={}^{-\frac{3·t}{\mathrm{tr}5\%}}=\left(\begin{array}{c}0.8607\\ 0.8607\end{array}\right)$

insofar as dt=3 s, as mentioned above.

These data especially allow the matrices A_obs and B_obs to be defined.

The thermal load w can then be estimated with this program using the system of equations (Eq. 19).

FIG. 5 shows the estimation of the thermal load obtained by the program described above and the actual variation in the thermal load applied to the refrigerator shown in FIG. 4.

As may be seen in said FIG. 5, although the thermal load applied to the refrigerator 1′ is highly variable, and even random, the estimation is of a very high quality. The estimation follows the actual variation in the thermal load during the peaks in power.

This estimation is also very good in an equilibrium regime, for example between the time t=0 s and the time t=3×10⁴ s in FIG. 5.

Finally, it is possible to correctly estimate, in real time, the variation in the thermal load applied to the cryogenic refrigerator.

This estimation is very good whether the operating regime is a variable or equilibrium regime.

This is particularly advantageous because it is then possible to regulate one or more operating parameters of the cryogenic refrigerator.

It may for example be envisioned to regulate the temperature at the outlet of the turbine 31′, the height of the bath of liquid in the separator, or other parameters.

Considering the above example, it will then be possible to physically adjust variables such as the degree of opening of the valve, the temperature and/or the pressure upstream of this valve, etc. in order to regulate the or each operating parameter considered.

In particular, in the case of pulsed tokamak regimes, it is possible to estimate, in real time, the value of the thermal load applied to the cryogenic refrigerator, without knowledge of events outside of the refrigerator.

This regulation then allows adequate cooling of the plant to be ensured, avoiding the risk of the refrigerator shutting down and consequently the plant itself shutting down.

This regulation is inexpensive because it requires no major hardware.

It should be noted that the coefficients identified for the matrix a are valid for the refrigerator 1′ shown in FIG. 4, taking account of the parameters chosen for equation (Eq. 2).

Different values would have to be identified for the matrix a if other parameters were to be included in equations (Eq. 2). For example, if it were desired to have an even more precise model, it could be envisioned to furthermore take into account the pressure P upstream of the valve 50′, and the time derivative of this pressure P.

Different values would also have to be identified for the matrix a for the same refrigerator operating under different conditions, for example if the outlet pressure of the compressor were different.

Claims

1. A method for estimating a thermal load applied to a cryogenic refrigerator comprising a phase separator containing a bath of liquid refrigerant, in which this thermal load is estimated using a computer program, said program being based on a mass balance performed on the phase separator allowing the variation in the time derivative of the height of the bath of liquid refrigerant in the phase separator to be expressed.

2. The method as claimed in claim 1, in which the thermal load applied to the cryogenic refrigerator is a variable thermal load.

3. The method as claimed in claim 1, in which the variable thermal load applied to the cryogenic refrigerator is pulsed.

4. The method as claimed in claim 1, in which a step is provided for determining appropriate variables of the cryogenic refrigerator, representing the variation in the time derivative of the height of the bath of liquid refrigerant in the phase separator, before the step of estimating the thermal load using the computer program.

5. The method as claimed in claim 1, in which, the cryogenic refrigerator comprising a valve at the inlet of the phase separator, said appropriate variables comprise at least the degree of opening of the valve, the temperature upstream of the valve and their respective time derivatives, and the thermal load.

6. The method as claimed in claim 4, in which the computer program comprises a step in which the data relating to the variation in the time derivative of the height of the bath of liquid refrigerant in the phase separator are filtered.

7. A method for regulating the cryogenic refrigerator subjected to a thermal load, in which the thermal load applied to the refrigerator is estimated using the method as claimed in claim 1, and then at least one operating parameter of the refrigerator is regulated depending on the value of the thermal load estimated beforehand.

8. The method as claimed in claim 7, in which the regulating step is implemented by modifying the degree of opening of the valve of the refrigerator, said valve being located at the inlet of the phase separator.

9. A computer program product comprising programming code instructions for implementing a method as claimed in claim 1.