- Email: [email protected]
Simplified modeling of active magnetic regenerators
i n t e r n a t i o n a l j o u r n a l o f r e f r i g e r a t i o n 3 6 ( 2 0 1 3 ) 9 3 2 e9 4 0
Available online at www.sciencedirect.com
w w w . i i fi i r . o r g
journal homepage: www.elsevier.com/locate/ijrefrig
Simplified modeling of active magnetic regenerators Thomas Burdyny*, Andrew Rowe Department of Mechanical Engineering, Institute for Integrated Energy Systems, University of Victoria, 3800 Finnerty Rd, Victoria, B.C. V8W 3P6, Canada
article info
abstract
Article history:
A simplified model for predicting the general trends of active magnetic regenerator devices
Received 30 July 2012
is presented. To reduce computational demands it is assumed the regenerator has suffi-
Received in revised form
ciently high convective interaction such that a one-phase regenerator approximation is
16 October 2012
sufficient. A corrective term is subsequently added to the thermal conductivity to account
Accepted 21 October 2012
for the convective heat leak. Losses internal and external to the regenerator are also
Available online 31 October 2012
defined for a generalized AMR device. The steady state temperature span across the regenerator is then evaluated using an energy balance and real material properties with
Keywords:
the cooling capacity, work inputs and COP post-calculated using the result. A comparison
Active magnetic refrigerator cycle
of the simplified model with experimental results was performed using single-layer and
Thermodynamics
multilayer regenerators of Gd, Gd0.74Tb0.26 and Gd0.85Er0.15. The model showed good
Magnetocaloric effect
conformity in magnitude and sensitivity in all of the scenarios which varied by heat load,
Magnetic refrigeration
frequency and utilization. ª 2012 Elsevier Ltd and IIR. All rights reserved.
Mode´lisation simplifie´e des re´ge´ne´rateurs magne´tiques actifs Mots cle´s : Cycle du re´frige´rateur magne´tique actif ; Thermodynamique ; Effet magne´tocalorique ; Froid magne´tique
1.
Introduction
An active magnetic refrigerator (AMR) uses magnetically induced temperature changes in solid materials to vary the temperature of a heat transfer fluid. This may be used to pump heat or to generate power. The design of devices using AMR cycles is complicated by the time varying heat transfer interactions between the fluid and solid, non-linear material properties, and the many geometric and operational parameters. Performance is strongly impacted by the field and temperature dependence of magnetic material properties.
One means of predicting the performance of an AMR is a higher-order model that takes solid/fluid interactions into account on a nodal basis. A number of these models have been produced in recent years with a broad range of applications and resolutions (Engelbrecht, 2008; Legait et al., 2009; Risser et al., 2010; Nielsen et al., 2010; Tusek et al., 2010; Liu and Yu, 2011). A recent description of published models has been reported by Nielsen et al. (2011). Detailed physics approaches, while very useful in predicting the temperature span and cooling capacities of a regenerator, are computationally demanding. In addition, research based on these approaches
* Corresponding author. Tel.: þ1 250 661 4400; fax: þ1 250 721 6051. E-mail address: [email protected] (T. Burdyny). 0140-7007/$ e see front matter ª 2012 Elsevier Ltd and IIR. All rights reserved. http://dx.doi.org/10.1016/j.ijrefrig.2012.10.022
i n t e r n a t i o n a l j o u r n a l o f r e f r i g e r a t i o n 3 6 ( 2 0 1 3 ) 9 3 2 e9 4 0
Nomenclature a, b, c, d A B Bi c d DF Fo h H k L m n Pr p Q R Re s t T DT U
Discrete points in a cycle, e Surface area, m2 Magnetic Field, T Biot number, e Specific heat, J kg1 K1 Diameter, m Degradation factor, e Fourier number, e Convection coefficient, enthalpy, Wm2 K1, k J kg1 Enthalpy flux, W Thermal conductivity, Wm1 K1 Length, m Mass, kg MCE scaling exponent, e Prandtl number, e Pressure, Nm2 Heat transfer, net enthalpy flux, W Thermal mass ratio, e Reynold’s number, e Entropy, k J kg1 K1 Non-dimensional time coordinate, s Temperature, K Magnetocaloric Effect, K Utilization, e
focuses primarily on thermodynamic and heat transfer optimization. This is essential for understanding the physics of the AMR but deters optimization of the entire refrigeration process due to the complexity of the regenerator alone. A simpler model based on work by Rowe (2011a, 2011b) was created to aid the general AMR design process. This model is designed to capture the main characteristics and sensitivities of an AMR but simplify the numerical problem. It can be used to analyze systems and focus more on issues such as component sizing, cost analysis and overall system optimization. Previously, the model was used to determine the performance of AMR refrigerators using idealized material properties Rowe (2011b). The current paper uses the same thermodynamic formulation, but real material properties and regenerator characteristics. Model predictions are compared to experimental results from a superconducting AMR test apparatus (SC-AMRTA) (Tura, 2005). Governing equations and underlying assumptions are presented in the next section. Various loss modes in the regenerator and the device are described and included in the energy balance. Finally, numerical results are validated against sets of experimental data which vary by regenerator material, applied heat load, frequency and utilization.
2.
Governing equations
As shown by Rowe (2011a) the following equations describe an AMR when the thermal coupling between the solid and fluid is
W x
933
Power, W Non-dimensional spatial coordinate, m
Greek a b G k r s s F
Thermal diffusivity, m2 s1 Balance, e Porosity, e Geometric form factor, e Non-dimensional conductance, e Density, kg m3 Symmetry, e Period, s Utilization, e
Subscript C c eff f H h m p s
Cold or cooling capacity, e Cycle, e Effective, e Fluid, e Hot or high-field, e Hydraulic, e Magnetic, e Constant pressure, parasitic, e Solid, e
3
Superscript 0 Per unit length, e
sufficiently large to assume a single temperature represents both substances. The spatial position is normalized by the regenerator length, L, so that x varies between 0 and 1. H¼
m0s Lcs ha DT k dT FRH ð1 bÞ þ þ f1 UH FRH dx sc cp RH
(1)
m0s cs DT n 1 dDT DT cRH ðRH 1Þ 1 n þ 1 s dT RH sc T dT þ f2 þ ðRH 1Þf3 UH dx
(2)
0 ¼ Wm
0 Wm þ Qp0 ¼
1 dH L dx
(3)
Eq. (1) describes the net transfer of heat at any location in the AMR, Eq. (2) is the local rate of magnetic work, and Eq. (3) is the periodic steady-state energy balance. Parameters f1ef3 are: UC 1 1 dDT b UC 1þ f1 h RH dT 2 UH UH 2 f2 h1 þ
1 n UC RH n þ 1 UH
dDT dT
UC 1 n 1 UC 1 dDT dDT f3 h 1 1þ RH dT dT UH RH n þ 1 s UH
(4)
(5)
(6)
The governing equations are derived assuming a step-wise variation in field and fluid flow and the temperatures at different points in the cycle are determined in reference to point a e the local temperatures at the beginning of the cold
934
i n t e r n a t i o n a l j o u r n a l o f r e f r i g e r a t i o n 3 6 ( 2 0 1 3 ) 9 3 2 e9 4 0
blow. In the most general case, one can numerically solve Eq. (3) to determine the temperature distribution for specified boundary temperatures and then calculate work and heat transfer.
2.1.
Second-order refrigerants near room temperature
The experiments simulated using the model use second-order refrigerants in room temperature devices. This allows for some simplifications of the governing equations used in the model. Because the thermal capacity of the heat transfer fluid depends on temperature and pressure, the possibility exists for the flow to be thermally imbalanced between the hot and cold blows. For common heat transfer fluids near room temperature however the variation in specific heat and density is small enough that the balance can be assumed to be equal, giving b ¼ 1. It is also assumed that the magnetocaloric effect of second-order materials scales proportional to the field strength, therefore n ¼ 1. With these assumptions Eqs. (1) and (2) are inserted into Eq. (3) which produces the differential equation describing the temperature profile immediately before the fluid flow commences at low field. 2
d T dT sc þ a0 ¼ þ a1 Qp dx2 dx ms Lcs RH 1 1 dDT DT2 1 a0 ¼ cRH RH 2s dT T DT dDT d f1 UH FRH f2 þ ðRH 1Þf3 a1 ¼ T dT dT RH a2
(7)
a2 ¼ k f1 FRH UH d f 1 UH dUC 1 dUH 1 dDT b dUC ¼ 1þ 2 dT RH dT 2 dT dT dT 2 UH 1 d DT þ UC RH dT2 2
(8)
AMR loss mechanisms
Loss mechanisms existing in an AMR refrigerator include parasitic losses within the regenerator and losses external to the regenerator. The general AMR device schematic in Fig. 1 shows the thermodynamic exchanges occurring within a typical device.
3.1.
Parasitic losses
The parasitic losses appearing in Eqs. (3) and (7) are those occurring within the regenerator and across its boundaries. This includes eddy currents in the magnetic material, ambient heat leaks through the regenerator shell and viscous dissipation due to fluid pressure drop in the regenerator matrix. Therefore the regenerator losses per unit length are, 0 0 0 þ Qamb þ Qviscous Qp0 ¼ Qeddy
(10)
These losses are specific to the device design and operating conditions. The viscous dissipation term is calculated using the pressure drop across the regenerator as Qviscous ¼
_ f dp m rf L dx
(11)
As in the previous equations the spatial coordinate is nondimensional. For particles the Ergun equation (Ergun, 1952) can be used to calculate the pressure drop in conjunction with constants found by Kaviany (1995). Alternatively, experimental data can be used. Due to the magnitude of the pressure drop in the regenerator configurations considered here, the heat leak and eddy current terms are negligible in comparison to viscous dissipation losses.
3.2.
Solution of Eq. (7) requires two boundary conditions, Ta0 and Ta1, which results in the temperature distribution of a single state point, Ta(x), across the regenerator. Knowledge of the magnetic cycle allows for the other three state temperatures to be determined in reference to point a (Rowe, 2011a). The average fluid temperature in Eq. (9) is subsequently assumed to be the average of the four states of the cycle. 1 Tf ¼ ðTa þ Tb þ Tc þ Td Þ 4
3.
Regenerator effectiveness
In addition to the aforementioned parasitic losses, performance is also affected by the effectiveness of the regenerator itself. This includes the efficiency of transferring heat between the refrigerant and the working fluid as well as conduction between the ends of the regenerator.
(9)
The fluid temperatures at the ends of the regenerator, Tf (x ¼ 0) and Tf (x ¼ L), are then comparable to the fluctuating temperature measurements taken during experiments. These are then identified as TC and TH, respectively, which define the temperature span of the regenerator. With the temperature span and profile determined, the cooling capacity for the given operating parameters can be evaluated. This requires accounting for losses external to the regenerator due to eddy currents and ambient heat leaks. It is also important to quantify the unavoidable losses occurring within the regenerator itself.
Fig. 1 e Schematic of the general components in an AMR device (right) and thermodynamic fluxes (left) within the system include work inputs, losses and enthalpy fluxes across the regenerator.
i n t e r n a t i o n a l j o u r n a l o f r e f r i g e r a t i o n 3 6 ( 2 0 1 3 ) 9 3 2 e9 4 0
The conductive component is accounted for in the governing expressions through the non-dimensional conductivity term present in Eqs. (1) and (7) and is thus present in both the calculation of enthalpy flux and the differential equation describing the temperature profile. The nondimensional conductivity is calculated using Eq. (12). k¼
sc keff A ms Lcs L
(12)
where keff is the effective thermal conductivity of the regenerator. Because the model assumes perfect heat transfer between the solid and fluid, the effect of convection is removed from the energy balance differential equations. One approach to estimating the impact of finite convection and thermal mass on a regenerator’s effectiveness is to relate it to a heat leak from the hot side to the cold side (Nam and Jeong, 2003). This approach allows for convection losses on the AMR’s cooling power to be post-calculated, without it affecting the temperature profile of the regenerator. Instead, to better approximate imperfect heat transfer in a one-phase passive regenerator model, Vortmeyer and Schaefer (1974) proposed converting convection into an equivalent thermal conductivity. This equivalent term combines with the conductive component and explicitly includes the effects of convection on AMR performance. The effective thermal conductivity can then be written as keff ¼ kcond þ kconv
(13)
The first term, kcond, takes into account the combined conductivity of the solid and fluid phases in the regenerator. This is calculated using the static component of the thermal conductivity used by Dikeos (2006) and Engelbrecht (2008). The dispersive component used in their models however is not considered in the single phase approach as it is internally accounted for within the kconv term proposed by Vortmeyer and Schaefer (1974). The equivalent convective conductivity is, kconv ¼
_ f cf m AR
2
dh 4heff 3
(14)
where dh is the hydraulic diameter, A is the regenerator crosssectional area, R is the thermal mass ratio and heff is the corrected convection coefficient. The following empirical correlation derived by Wakao and Kaguei (1982) describes the convection coefficient for fluids passing through packed beds. 1=3 2 þ 1:1Re0:6 kf f Prf (15) h¼ dpart where Ref is the Reynold’s number based on the particle size of the porous media, Prf is the Prandtl number and dpart is the characteristic particle diameter. The convection coefficient is then corrected in Eqs. (16)e(20) using a degradation factor to account for internal temperature gradients which may exist in the particles (Engelbrecht, 2006), heff ¼ h DF 1 DF ¼ Bi 1 þ fH 5
Bi ¼
fH ¼
hdpart 2ks 1
935
(18) 4 35Fo
as Fo ¼ 2 dpart 2
(19)
(20)
With the non-dimensional thermal conductivity in Eq. (12) fully defined, the system losses due to regenerator effectiveness are accounted for in the model.
3.3.
External losses
The final loss mechanisms in an AMR system are due to effects external to the regenerator. As seen in Fig. 1 these may include eddy currents in the surrounding materials (for example, the cold heat exchanger), heat leaks from the environment to the cold side due to imperfect insulation and heat leaks from the hot side of the system through the structure. Eddy current heating losses can occur within the device due to the presence of the time-varying magnetic field and metallic materials. This is calculated using the following approximation by Kittel (1990) for the electrically conducting components subjected to the time varying field. Qeddy;c ¼
n X GAV 2 B 32r i i¼1
(21)
where G is a geometric form factor, A is the area enclosed by the current loop, V is the material volume, r is the electrical resistivity, B is the magnetic flux and n is the number of parts in the device. Ambient heat leaks occur throughout the entire refrigeration device and across the regenerator shell. This can occur in piping, heat exchangers or generated through bearing friction. This loss varies depending on the overall design of the device and can be estimated through both modeling and experimental methods. Eddy current losses in the AMRTA components which have direct impact on the cooling capacity (i.e. cold heat exchanger) are negligible. Based on work by Tura et al. (2006) an ambient heat leak through the AMRTA structure was estimated as, Qamb;c þ Qleak y 0:08WK1 ðTH TC Þ
(22)
With the external losses defined the total cooling power for the device is defined as QC ¼ Hðx ¼ 0Þ Qeddy;c Qamb;c Qleak
(23)
Using Eq. (23) and the temperature span defined by Eq. (9) the model can be tested against experimental data.
4.
Model description
(16) (17)
The AMR model uses the developed theory and correlations in conjunction with material data, geometric parameters and operating conditions to determine the temperature distribution across the AMR. From this, performance parameters such
936
i n t e r n a t i o n a l j o u r n a l o f r e f r i g e r a t i o n 3 6 ( 2 0 1 3 ) 9 3 2 e9 4 0
as expected cooling capacity and work input can be determined. Required inputs include material properties, effective high and low field strength during a blow, utilization and temperatures at the ends of the regenerator. Since the regenerator end temperatures are an input to the simulation, the temperature span of the regenerator is artificially created while the shape of the distribution is determined by Eq. (7). The cooling capacity and work for a given temperature span is then post-calculated.
4.1.
Fundamental operation
In the numerical model the regenerator’s temperature distribution, Ta(x), is determined by Eq. (7) using a discrete number of points. This equation is solved along the entire domain using MATLAB’s bvp5c ODE solver in conjunction with the user-specified boundary conditions Ta(x ¼ 0) and Ta(x ¼ L). Solution of this equation requires multiple iterations as several of the terms in the ODE are temperature dependent. This includes the adiabatic temperature change, DT, which is simulated by assuming isentropic magnetization of the material such that, DT ¼ Tðs1 ; H2 Þ Tðs1 ; H1 Þ
(24)
Entropy values for Eq. (24) are determined using the integral of the material’s specific heat. The material data was provided by AMES Iowa National Lab for Gd, Gd0.74Tb0.26 (GdTb) and Gd0.85Er0.15 (GdEr) which corresponds to the materials used in the AMRTA experiments. This data is interpolated in MATLAB using the TriScatteredInterp function allowing for the user to select the high and low field experienced by the magnetic material. Furthermore, if the Curie temperature and magnetocaloric effect are known at two or more field strengths for a certain material, this same technique can use the interpolated Gd data to map the material’s data over a wide operating range. This is due to the similarities in the shape of the second-order material data. This technique is used for GdEr and GdTb where data is available only at 0 T and 2 T. By mapping against Gd, approximate data for GdEr or GdTb can be generated for a wide range of field choices. With Ta(x) computed the other three state points are found using knowledge of the magnetic cycle (Rowe, 2011a); this allows for the fluid temperature across the regenerator to be calculated. Cooling capacity is then found using the theory described in the preceding sections. Other useful parameters such as magnetic work, pumping power, efficiency and figures of merit can also be post-calculated using the temperature profile and cooling capacity. Model inputs are selected such that the outputs can be compared to experimental data. The model takes approximately 2 s to solve for a single set of boundary conditions on a laptop computer with a 2.40 GHz Intel Core i5 processor.
4.2.
Experimental data
Experimental data used for validation includes the temperature span across the regenerator, Tspan, at a given hot end temperature, TH, for a specified heat load or cooling capacity, Qc. This allows for two important curves to be created which provide a measure of performance for a given material,
operating parameters and regenerator configuration. These curves are Qc vs Tspan for a fixed TH, and Tspan vs TH for a fixed Qc (typically Qc ¼ 0 W). Furthermore, experimental data sets used to compare the model results vary by material, utilization, heat load and frequency. Fixed and variable parameters used in the experiments are described in Table 1. The regenerators are made of crushed particles which are formed into pucks containing between 40 and 45 g of magnetocaloric material each. Regenerators in the experiments of interest are composed of one, two or three pucks made from Gd, GdEr or GdTb (Tura, 2005). Working fluid properties are found using NIST’s REFPROP program which allows temperature dependent variables to vary between data points. These include the dynamic viscosity, specific heat and density. These properties are relatively constant for the fluid and conditions used in the SC test apparatus for room temperature experiments. The magnetic field used for all of the tests is provided by a 2 T superconducting magnet. The magnetic field seen by the material as a function of time is approximated as a sinusoid. For this reason, the discrete applied field in the model uses RMS values to relate actual operating conditions during a blow period to the step-wise waveforms assumed in deriving Eqs. (1)e(3).
4.3.
Results
A variety of experimental conditions are evaluated to ensure the model is tested in more than a single operating region of the AMR. Once validated, a higher level of confidence then exists in predicting previously untested operating conditions or regions currently unavailable to the experimental device. As such, the following results vary by material type, material mass, utilization, frequency and heat load.
Table 1 e Temperature independent parameters of the superconducting AMRTA used in producing the experimental results. Test parameters Frequency (Hz) Stroke length (m) Effective displacer area (mm2)
0.65e0.8 0.21 2513
Regenerator Material mass (g per puck) Porosity Regenerator length (mm per puck) Regenerator diameter (mm) Characteristic particle diameter (mm) Particle type Number of regenerators
40e45 w55% 25 25 560 Crushed particles 2
Magnetic material Thermal conductivity (W K1) Density (kg m3) Reference specific heat (J kg1 K1)
10.3 7900 381
Working fluid (Helium) Specific heat, Cp (J kg1 K1) Thermal conductivity (W K1) Pressure (atm)
5139 0.153 3e9.7
937
i n t e r n a t i o n a l j o u r n a l o f r e f r i g e r a t i o n 3 6 ( 2 0 1 3 ) 9 3 2 e9 4 0
Fig. 2 compares the model results with all of the aforementioned losses mechanisms included to single pucks of Gd, GdTb and GdEr and under no-load (Qc ¼ 0 W) conditions. In all three data sets the model predictions follow the same general trend as the experimental data. The primary discrepancy exists at high TH values relative to the Curie points where the temperature span begins to drop abruptly (best seen for GdTb). While the general trends are correct for all three materials, the model tends to significantly over predict the performance of the GdTb regenerator. As can be seen from the GdTb plots, decreasing the assumed magnetocaloric effect for this material by a factor 0.77 as compared to the raw data corrects the prediction. As with GdEr, the model results still under-predict the experimental data when operating at temperatures well above the Curie point. It is interesting to note that another more detailed numerical model of the AMR (Dikeos, 2006) produced a similar result for the GdTb pucks. In that case, a better fit was obtained by increasing the static conductivity of the bed. The fact that both assumptions are effective in correcting model predictions is supported by a simplified study by Rowe (2011b) where a non-dimensional parameter determining the cooling power and efficiency of an AMR was found: I¼
F DT k T
(25)
In the work by Dikeos (2006), the impact of increasing the static conductivity decreases I; likewise, in this work, decreasing the adiabatic temperature change leads to a decrease in I. In effect, the adjustments produce the same result in the governing energy balance. Both approaches may be justified given the lack of detailed measurements for the particular batch of GdTb used within experiments as compared to the batch used in measuring the material properties. One possibility is differences in the manufacturing process between batches of GdTb, which would explain why
no adjustments are needed for Gd and GdEr. Furthermore, the irregular particle shapes result in some uncertainty regarding the regenerator thermal effectiveness. The model predictions can also be compared to layered regenerators where experimental data is available for two and three layer regenerators using Gd, GdEr and GdTb. Based on the results for the single GdTb bed, the corrected adiabatic temperature change (0.77 scaling) for this material is used in the model. Fig. 3 demonstrates that the model predicts both the layering aspect and the combination of materials relatively well. The model tends to be less accurate the more the operating temperature is increased above the average Curie temperature of the materials in the bed and as the utilization decreases. Because the gas displaced through the regenerator is constant and adding a material increases the overall regenerator mass, the utilization decreases from w0.3 for the single puck experiments down to w0.1 for the three material tests. The sensitivity of the temperature span to applied load is shown in Fig. 4a for a single Gd puck at TH ¼ 292 K. From the figure the experimental results are in good agreement with the model. The model has been used to also show the effect that system losses have on the predicted cooling capacity. Recalling Eq. (23), the net cooling power is determined by the enthalpy flow at x ¼ 0, H(x ¼ 0), less the parasitic heat leaks from ambient, Qamb,c, and from the hot side of the system, Qleak. The impacts of regenerator thermal ineffectiveness, internal losses and the external losses can be easily shown. In this case internal losses are due to only viscous dissipation. The top curve represents the temperature span for a system with perfect regenerator effectiveness, no viscous dissipation and no ambient cooling losses. Below this curve is the same scenario but with regenerator ineffectiveness included. The next curve then shows the cooling power when only the external losses are excluded while the final solid curve accounts for all system losses and is the predicted cooling capacity for the device. This analysis of system losses is important in optimization studies for determining where
35 55
30
GdTb GdEr
45
20 TSpan [K]
T
Span
[K]
25
50
Gd
15
Gd−GdEr
40 35
10
30
5
25
0 250
Gd−GdTb−GdEr
Gd−GdTb
Gd−Gd
20
260
270
280 290 TH [K]
300
310
Fig. 2 e TSpan vs TH for single pucks of Gd, GdTb and GdEr at Qc [ 0 W, 0.65 Hz and a charge pressure of 9.5 atm. Solid lines represent experimental data while dashed lines are simulations from the model. In the case of GdTb pucks the temperature span is plotted using both original material data and data altered to reduce the MCE.
15 280
285
290
295 300 T [K]
305
310
315
H
Fig. 3 e TSpan vs TH for two and three puck multilayered regenerators at Qc [ 0 W, 0.65 Hz and a charge pressure of 9.5 atm. Solid lines represent experimental data while dashed lines are simulations from the model. Note the vertical axis is offset from 0 K for clarity.
938
i n t e r n a t i o n a l j o u r n a l o f r e f r i g e r a t i o n 3 6 ( 2 0 1 3 ) 9 3 2 e9 4 0
a
30 25
0W
20 TSpan [K]
experimental devices can be improved and potential effects these improvements have on the temperature span and cooling capacity. As is demonstrated in Fig. 4b, for the same operating parameters but at a temperature above the Curie point (TH ¼ 305 K), viscous dissipation losses have a much greater impact on the temperature span than at lower temperatures. In Fig. 5 the temperature span is plotted against three different TH values for two heat load cases, 0 W and 6 W. While both data sets follow the trend of the experimental curve, the model under-predicts experiments at the higher frequency. Potential causes of differences between model predictions and experiments are discussed in the next section. The final experimental result accounts for a variation in the utilization. For the AMRTA the utilization is varied by altering the pressure of the helium in the system thereby changing the density and subsequent thermal mass of the displaced fluid. Increasing the utilization raises the pressure
15 6W
10 5 0
290
295 T [K]
300
305
H
Fig. 5 e TSpan vs TH for single-pucks of Gd at 9.5 atm, 0.8 Hz and varying load. Solid lines represent experimental data while dashed lines are simulations from the model.
20 drop in the system and, in the case of the AMRTA, increases the temperature span achievable by the regenerator up to certain value. Shown in Fig. 6 are the temperature span curves at three different charge pressures, or utilizations. Similar to the heat load sensitivity curves in Fig. 5 the temperature span follows the experimental trends. Furthermore, at a pressure of 3 atm the temperature span is trending lower at TH ¼ 304 K.
Q [W]
15
10
5
Exp. H(x=0)−Q
−Q
amb,c
(all losses)
leak
H (with int. losses, no ext. losses) H (no int. losses, no ext. losses) H (100% effectiveness, no losses)
0 0
5
10
T
Span
b
5. 15
[K]
20
20
In all of the tested scenarios the model results followed the trends of the experimental results. In the case of varying utilization and heat load some differences are seen between the experiments and the model. Differences are most
25
15
9.5 atm
20
10
5
H(x=0)−Q
−Q
amb,c
leak
TSpan [K]
Q [W]
Discussion
(all losses)
10
H (with int. losses, no ext. losses) H (no int. losses, no ext. losses) H (100% effectiveness, no losses)
0 0
5
10 T
Span
[K]
15
6 atm
15
3 atm 5
20
Fig. 4 e Qc vs TSpan for single-pucks of Gd at 0.65 Hz and a utilization corresponding to a charge pressure of 9.5 atm, (a) TH [ 292 K and (b) TH [ 305 K. Results are plotted to show the effect that internal and external losses have the cooling capacity and temperature span. The difference in the impact of viscous dissipation (internal losses) is also shown between the two plots.
0
270
280
290 T [K]
300
310
H
Fig. 6 e TSpan vs TH for single-pucks of Gd at 0.65 Hz, Qc [ 0 W and utilizations corresponding to charge pressures of 3 atm, 6 atm and 9.5 atm. Solid lines represent experimental data while dashed lines are simulations from the model.
i n t e r n a t i o n a l j o u r n a l o f r e f r i g e r a t i o n 3 6 ( 2 0 1 3 ) 9 3 2 e9 4 0
significant when the operating temperatures are increased above the Curie points of the materials in the AMR and when the utilization is decreasing. In general, when these conditions arise the model tends to under-predict performance. These discrepancies are expected due to the simplicity of the model relative to the physics occurring at the time of the experiments. Examples of this are seen in the physical representation of the AMR inside the model including the characteristic particle diameter and simulated RMS field strength. The model also assumes that pure Gd, GdTb and GdEr exist in the regenerator. In the experiments however the regenerator materials are coated in a thin layer of epoxy which will increase interfacial heat transfer resistance, act as a parasitic thermal mass, and may cause flow-channeling. Additionally, since the cooling power of the AMRTA experiments are small (w20 W), any errors in modeling the losses can have a large impact on the expected temperature span as compared to experiments with large cooling capacities. Some errors also exist in the replication of GdEr and GdTb specific heat curves from Gd data as the specific heat curves of second-order materials are similar but not exactly the same. In the operating space defined by the SC AMRTA the model is found to produce repeatable and continuous curves based upon the results in Fig. 2 through Fig. 6. In general, TH values below the Curie point solve faster than those above the Curie point. The model has difficulty solving for the temperature distribution at TH values high above the Curie point when both large pressure drops and low utilizations (<0.15) are present. This is due to the non-dimensional thermal conductivity becoming too low to solve the energy balance differential equation. It was found that setting a minimum value on the conductivity in the differential equation rectified this problem. It was noted by Tura (2005) that obtaining data with TH values high above the Curie points is in general experimentally unstable and difficult to replicate, particularly for varying heat load scenarios at low utilizations. Overall, the results sufficiently follow the tested scenarios enough to provide a basis for proceeding with the model in optimization and design cases where general trends of the performance are of primary importance. The next stage of verification for the model involves systems where the working fluid is a liquid and the expected cooling loads are much higher. Similar to the results used in this paper, an adequate amount of experimental data is available to verify the revised model. In parallel, the model is to be further used in optimization programs and mapping of various operating spaces. Providing quick solutions that accurately predict the trends of the AMR allows for parameters such as the regenerator length to be linked to the length of the permanent or superconducting magnet. The device performance at a given field strength can then be linked with the magnet cost necessary to generate the matching size and magnitude of the given field. This type of design balances the importance of performance and cost requirements on future AMR devices.
6.
Conclusions
Through the use of experimental results it has been shown that a simplified AMR model using real material properties
939
can predict the general trends of a device. This is verified from results using a 0e2 T magnetic field and three different types of materials, Gd, GdTb and GdEr, which were evaluated as single layer and multilayer regenerators. These scenarios showed good sensitivity with experiments both with no-load and an applied load as well as at two different frequencies, 0.65 Hz and 0.8 Hz, and three charge pressures, 9.5 atm, 6 atm and 3 atm. The model also provides flexibility in observing the individual predicted losses of the AMR signifying areas of device improvement. In addition, operating regions which are difficult for the model to predict have been identified and provide direction for further improvements. The model will now be tested against experimental systems using liquid working fluids and will be implemented into costing and optimization programs.
Acknowledgments Support in the form of a Canadian Graduate Scholarship from the Natural Sciences and Engineering Research Council of Canada is greatly appreciated. Partial support for this work was provided by H2Can Strategic Network. The authors also acknowledge the unwavering support of colleagues and friends whose stimulating conversations allowed this paper to come to fruition.
references
Dikeos, J., 2006. Development and Validation of an Active Magnetic Regenerator Refrigerator Cycle Simulation. MASc thesis, University of Victoria, Victoria. Engelbrecht, K., 2006. The effect of internal temperature gradients on regenerator matrix performance. Trans. ASME 128, 1060e1069. Engelbrecht, K., 2008. A Numerical Model of an Active Magnetic Regenerator Refrigerator with Experimental Validation. PhD thesis, University of Wisconsin-Madison. Ergun, S., 1952. Fluid flow through packed column. Chem. Eng. Prog. 48, 125e127. Kaviany, M., 1995. Principles of Heat Transfer in Porous Media. Springer, New York, NY. Kittel, P., 1990. Eddy current heating in magnetic refrigerators. In: Fast, R.W. (Ed.), Adv. Cryog. Eng. Plenum, vol. 35B, pp. 1141e1148. Legait, I., Kedous-Lebouc, A., Rondot, L., 2009. Numerical simulation and analysis of the refrigerant bed behavior using fluent software. In: Egolf, P.W. (Ed.), Third International Conference on Magnetic Refrigeration at Room Temperature. International Institute of Refrigeration. Liu, M., Yu, B., 2011. Numerical investigations on internal temperature distribution and refrigeration performance of reciprocating active magnetic regenerator of room temperature magnetic refrigeration. Int. J. Refrigeration 34, 617e627. Nam, K., Jeong, S., 2003. Measurement of cryogenic regenerator characteristics under oscillating flow and pulsating pressure. Cryogenics 43, 575e581. Nielsen, K.K., Bahl, C.R.H., Smith, A., Pryds, N., Hattel, J., 2010. A comprehensive parameter study of an active magnetic regenerator using a 2D numerical model. Int. J. Refrigeration 33, 753e764.
940
i n t e r n a t i o n a l j o u r n a l o f r e f r i g e r a t i o n 3 6 ( 2 0 1 3 ) 9 3 2 e9 4 0
Nielsen, K.K., Tusek, J., Engelbrecht, K., Schopfer, S., Kitanovski, A., Bahl, C.R.H., Smith, A., Pryds, N., Poredos, A., 2011. Review on numerical modeling of active magnetic regenerators for room temperature applications. Int. J. Refrigeration 34, 603e616. Risser, M., Vasile, C., Engel, T., Keith, B., Muller, C., 2010. Numerical simulation of magnetocaloric system behavior for an industrial application. Int. J. Refrigeration 33 (5), 973e981. Rowe, A., 2011a. Thermodynamics of active magnetic regenerators: part I. Cryogenics 52, 111e118. Rowe, A., 2011b. Thermodynamics active magnetic regenerators: part II. Cryogenics 52, 119e128. Tura, A., 2005. Active Magnetic Regenerator Experimental Optimization. MASc thesis, University of Victoria, Victoria.
Tura, A., Roszmann, J., Dikeos, J., Rowe, A., 2006. Cryogenic active magnetic regenerator test apparatus. Adv. Cryog. Eng. 823, 985e992. Tusek, J., Sarlah, A., Zupan, S., Prebil, I., Kitanovski, A., Poredos, A., 2010. A numerical optimization of a packed bed AMR. In: Egolf, P.W. (Ed.), Fourth International Conference on Magnetic Refrigeration at Room Temperature. International Institute of Refrigeration. Vortmeyer, D., Schaefer, R.J., 1974. Equivalence of one- and twophase models in packed beds: one dimensional theory. Chem. Eng. Sci. 29, 485e491. Wakao, N., Kaguei, S., 1982. Heat and Mass Transfer in Packed Beds. Gordon and Breach Science Publishers, New York, NY.
Available online at www.sciencedirect.com
w w w . i i fi i r . o r g
journal homepage: www.elsevier.com/locate/ijrefrig
Simplified modeling of active magnetic regenerators Thomas Burdyny*, Andrew Rowe Department of Mechanical Engineering, Institute for Integrated Energy Systems, University of Victoria, 3800 Finnerty Rd, Victoria, B.C. V8W 3P6, Canada
article info
abstract
Article history:
A simplified model for predicting the general trends of active magnetic regenerator devices
Received 30 July 2012
is presented. To reduce computational demands it is assumed the regenerator has suffi-
Received in revised form
ciently high convective interaction such that a one-phase regenerator approximation is
16 October 2012
sufficient. A corrective term is subsequently added to the thermal conductivity to account
Accepted 21 October 2012
for the convective heat leak. Losses internal and external to the regenerator are also
Available online 31 October 2012
defined for a generalized AMR device. The steady state temperature span across the regenerator is then evaluated using an energy balance and real material properties with
Keywords:
the cooling capacity, work inputs and COP post-calculated using the result. A comparison
Active magnetic refrigerator cycle
of the simplified model with experimental results was performed using single-layer and
Thermodynamics
multilayer regenerators of Gd, Gd0.74Tb0.26 and Gd0.85Er0.15. The model showed good
Magnetocaloric effect
conformity in magnitude and sensitivity in all of the scenarios which varied by heat load,
Magnetic refrigeration
frequency and utilization. ª 2012 Elsevier Ltd and IIR. All rights reserved.
Mode´lisation simplifie´e des re´ge´ne´rateurs magne´tiques actifs Mots cle´s : Cycle du re´frige´rateur magne´tique actif ; Thermodynamique ; Effet magne´tocalorique ; Froid magne´tique
1.
Introduction
An active magnetic refrigerator (AMR) uses magnetically induced temperature changes in solid materials to vary the temperature of a heat transfer fluid. This may be used to pump heat or to generate power. The design of devices using AMR cycles is complicated by the time varying heat transfer interactions between the fluid and solid, non-linear material properties, and the many geometric and operational parameters. Performance is strongly impacted by the field and temperature dependence of magnetic material properties.
One means of predicting the performance of an AMR is a higher-order model that takes solid/fluid interactions into account on a nodal basis. A number of these models have been produced in recent years with a broad range of applications and resolutions (Engelbrecht, 2008; Legait et al., 2009; Risser et al., 2010; Nielsen et al., 2010; Tusek et al., 2010; Liu and Yu, 2011). A recent description of published models has been reported by Nielsen et al. (2011). Detailed physics approaches, while very useful in predicting the temperature span and cooling capacities of a regenerator, are computationally demanding. In addition, research based on these approaches
* Corresponding author. Tel.: þ1 250 661 4400; fax: þ1 250 721 6051. E-mail address: [email protected] (T. Burdyny). 0140-7007/$ e see front matter ª 2012 Elsevier Ltd and IIR. All rights reserved. http://dx.doi.org/10.1016/j.ijrefrig.2012.10.022
i n t e r n a t i o n a l j o u r n a l o f r e f r i g e r a t i o n 3 6 ( 2 0 1 3 ) 9 3 2 e9 4 0
Nomenclature a, b, c, d A B Bi c d DF Fo h H k L m n Pr p Q R Re s t T DT U
Discrete points in a cycle, e Surface area, m2 Magnetic Field, T Biot number, e Specific heat, J kg1 K1 Diameter, m Degradation factor, e Fourier number, e Convection coefficient, enthalpy, Wm2 K1, k J kg1 Enthalpy flux, W Thermal conductivity, Wm1 K1 Length, m Mass, kg MCE scaling exponent, e Prandtl number, e Pressure, Nm2 Heat transfer, net enthalpy flux, W Thermal mass ratio, e Reynold’s number, e Entropy, k J kg1 K1 Non-dimensional time coordinate, s Temperature, K Magnetocaloric Effect, K Utilization, e
focuses primarily on thermodynamic and heat transfer optimization. This is essential for understanding the physics of the AMR but deters optimization of the entire refrigeration process due to the complexity of the regenerator alone. A simpler model based on work by Rowe (2011a, 2011b) was created to aid the general AMR design process. This model is designed to capture the main characteristics and sensitivities of an AMR but simplify the numerical problem. It can be used to analyze systems and focus more on issues such as component sizing, cost analysis and overall system optimization. Previously, the model was used to determine the performance of AMR refrigerators using idealized material properties Rowe (2011b). The current paper uses the same thermodynamic formulation, but real material properties and regenerator characteristics. Model predictions are compared to experimental results from a superconducting AMR test apparatus (SC-AMRTA) (Tura, 2005). Governing equations and underlying assumptions are presented in the next section. Various loss modes in the regenerator and the device are described and included in the energy balance. Finally, numerical results are validated against sets of experimental data which vary by regenerator material, applied heat load, frequency and utilization.
2.
Governing equations
As shown by Rowe (2011a) the following equations describe an AMR when the thermal coupling between the solid and fluid is
W x
933
Power, W Non-dimensional spatial coordinate, m
Greek a b G k r s s F
Thermal diffusivity, m2 s1 Balance, e Porosity, e Geometric form factor, e Non-dimensional conductance, e Density, kg m3 Symmetry, e Period, s Utilization, e
Subscript C c eff f H h m p s
Cold or cooling capacity, e Cycle, e Effective, e Fluid, e Hot or high-field, e Hydraulic, e Magnetic, e Constant pressure, parasitic, e Solid, e
3
Superscript 0 Per unit length, e
sufficiently large to assume a single temperature represents both substances. The spatial position is normalized by the regenerator length, L, so that x varies between 0 and 1. H¼
m0s Lcs ha DT k dT FRH ð1 bÞ þ þ f1 UH FRH dx sc cp RH
(1)
m0s cs DT n 1 dDT DT cRH ðRH 1Þ 1 n þ 1 s dT RH sc T dT þ f2 þ ðRH 1Þf3 UH dx
(2)
0 ¼ Wm
0 Wm þ Qp0 ¼
1 dH L dx
(3)
Eq. (1) describes the net transfer of heat at any location in the AMR, Eq. (2) is the local rate of magnetic work, and Eq. (3) is the periodic steady-state energy balance. Parameters f1ef3 are: UC 1 1 dDT b UC 1þ f1 h RH dT 2 UH UH 2 f2 h1 þ
1 n UC RH n þ 1 UH
dDT dT
UC 1 n 1 UC 1 dDT dDT f3 h 1 1þ RH dT dT UH RH n þ 1 s UH
(4)
(5)
(6)
The governing equations are derived assuming a step-wise variation in field and fluid flow and the temperatures at different points in the cycle are determined in reference to point a e the local temperatures at the beginning of the cold
934
i n t e r n a t i o n a l j o u r n a l o f r e f r i g e r a t i o n 3 6 ( 2 0 1 3 ) 9 3 2 e9 4 0
blow. In the most general case, one can numerically solve Eq. (3) to determine the temperature distribution for specified boundary temperatures and then calculate work and heat transfer.
2.1.
Second-order refrigerants near room temperature
The experiments simulated using the model use second-order refrigerants in room temperature devices. This allows for some simplifications of the governing equations used in the model. Because the thermal capacity of the heat transfer fluid depends on temperature and pressure, the possibility exists for the flow to be thermally imbalanced between the hot and cold blows. For common heat transfer fluids near room temperature however the variation in specific heat and density is small enough that the balance can be assumed to be equal, giving b ¼ 1. It is also assumed that the magnetocaloric effect of second-order materials scales proportional to the field strength, therefore n ¼ 1. With these assumptions Eqs. (1) and (2) are inserted into Eq. (3) which produces the differential equation describing the temperature profile immediately before the fluid flow commences at low field. 2
d T dT sc þ a0 ¼ þ a1 Qp dx2 dx ms Lcs RH 1 1 dDT DT2 1 a0 ¼ cRH RH 2s dT T DT dDT d f1 UH FRH f2 þ ðRH 1Þf3 a1 ¼ T dT dT RH a2
(7)
a2 ¼ k f1 FRH UH d f 1 UH dUC 1 dUH 1 dDT b dUC ¼ 1þ 2 dT RH dT 2 dT dT dT 2 UH 1 d DT þ UC RH dT2 2
(8)
AMR loss mechanisms
Loss mechanisms existing in an AMR refrigerator include parasitic losses within the regenerator and losses external to the regenerator. The general AMR device schematic in Fig. 1 shows the thermodynamic exchanges occurring within a typical device.
3.1.
Parasitic losses
The parasitic losses appearing in Eqs. (3) and (7) are those occurring within the regenerator and across its boundaries. This includes eddy currents in the magnetic material, ambient heat leaks through the regenerator shell and viscous dissipation due to fluid pressure drop in the regenerator matrix. Therefore the regenerator losses per unit length are, 0 0 0 þ Qamb þ Qviscous Qp0 ¼ Qeddy
(10)
These losses are specific to the device design and operating conditions. The viscous dissipation term is calculated using the pressure drop across the regenerator as Qviscous ¼
_ f dp m rf L dx
(11)
As in the previous equations the spatial coordinate is nondimensional. For particles the Ergun equation (Ergun, 1952) can be used to calculate the pressure drop in conjunction with constants found by Kaviany (1995). Alternatively, experimental data can be used. Due to the magnitude of the pressure drop in the regenerator configurations considered here, the heat leak and eddy current terms are negligible in comparison to viscous dissipation losses.
3.2.
Solution of Eq. (7) requires two boundary conditions, Ta0 and Ta1, which results in the temperature distribution of a single state point, Ta(x), across the regenerator. Knowledge of the magnetic cycle allows for the other three state temperatures to be determined in reference to point a (Rowe, 2011a). The average fluid temperature in Eq. (9) is subsequently assumed to be the average of the four states of the cycle. 1 Tf ¼ ðTa þ Tb þ Tc þ Td Þ 4
3.
Regenerator effectiveness
In addition to the aforementioned parasitic losses, performance is also affected by the effectiveness of the regenerator itself. This includes the efficiency of transferring heat between the refrigerant and the working fluid as well as conduction between the ends of the regenerator.
(9)
The fluid temperatures at the ends of the regenerator, Tf (x ¼ 0) and Tf (x ¼ L), are then comparable to the fluctuating temperature measurements taken during experiments. These are then identified as TC and TH, respectively, which define the temperature span of the regenerator. With the temperature span and profile determined, the cooling capacity for the given operating parameters can be evaluated. This requires accounting for losses external to the regenerator due to eddy currents and ambient heat leaks. It is also important to quantify the unavoidable losses occurring within the regenerator itself.
Fig. 1 e Schematic of the general components in an AMR device (right) and thermodynamic fluxes (left) within the system include work inputs, losses and enthalpy fluxes across the regenerator.
i n t e r n a t i o n a l j o u r n a l o f r e f r i g e r a t i o n 3 6 ( 2 0 1 3 ) 9 3 2 e9 4 0
The conductive component is accounted for in the governing expressions through the non-dimensional conductivity term present in Eqs. (1) and (7) and is thus present in both the calculation of enthalpy flux and the differential equation describing the temperature profile. The nondimensional conductivity is calculated using Eq. (12). k¼
sc keff A ms Lcs L
(12)
where keff is the effective thermal conductivity of the regenerator. Because the model assumes perfect heat transfer between the solid and fluid, the effect of convection is removed from the energy balance differential equations. One approach to estimating the impact of finite convection and thermal mass on a regenerator’s effectiveness is to relate it to a heat leak from the hot side to the cold side (Nam and Jeong, 2003). This approach allows for convection losses on the AMR’s cooling power to be post-calculated, without it affecting the temperature profile of the regenerator. Instead, to better approximate imperfect heat transfer in a one-phase passive regenerator model, Vortmeyer and Schaefer (1974) proposed converting convection into an equivalent thermal conductivity. This equivalent term combines with the conductive component and explicitly includes the effects of convection on AMR performance. The effective thermal conductivity can then be written as keff ¼ kcond þ kconv
(13)
The first term, kcond, takes into account the combined conductivity of the solid and fluid phases in the regenerator. This is calculated using the static component of the thermal conductivity used by Dikeos (2006) and Engelbrecht (2008). The dispersive component used in their models however is not considered in the single phase approach as it is internally accounted for within the kconv term proposed by Vortmeyer and Schaefer (1974). The equivalent convective conductivity is, kconv ¼
_ f cf m AR
2
dh 4heff 3
(14)
where dh is the hydraulic diameter, A is the regenerator crosssectional area, R is the thermal mass ratio and heff is the corrected convection coefficient. The following empirical correlation derived by Wakao and Kaguei (1982) describes the convection coefficient for fluids passing through packed beds. 1=3 2 þ 1:1Re0:6 kf f Prf (15) h¼ dpart where Ref is the Reynold’s number based on the particle size of the porous media, Prf is the Prandtl number and dpart is the characteristic particle diameter. The convection coefficient is then corrected in Eqs. (16)e(20) using a degradation factor to account for internal temperature gradients which may exist in the particles (Engelbrecht, 2006), heff ¼ h DF 1 DF ¼ Bi 1 þ fH 5
Bi ¼
fH ¼
hdpart 2ks 1
935
(18) 4 35Fo
as Fo ¼ 2 dpart 2
(19)
(20)
With the non-dimensional thermal conductivity in Eq. (12) fully defined, the system losses due to regenerator effectiveness are accounted for in the model.
3.3.
External losses
The final loss mechanisms in an AMR system are due to effects external to the regenerator. As seen in Fig. 1 these may include eddy currents in the surrounding materials (for example, the cold heat exchanger), heat leaks from the environment to the cold side due to imperfect insulation and heat leaks from the hot side of the system through the structure. Eddy current heating losses can occur within the device due to the presence of the time-varying magnetic field and metallic materials. This is calculated using the following approximation by Kittel (1990) for the electrically conducting components subjected to the time varying field. Qeddy;c ¼
n X GAV 2 B 32r i i¼1
(21)
where G is a geometric form factor, A is the area enclosed by the current loop, V is the material volume, r is the electrical resistivity, B is the magnetic flux and n is the number of parts in the device. Ambient heat leaks occur throughout the entire refrigeration device and across the regenerator shell. This can occur in piping, heat exchangers or generated through bearing friction. This loss varies depending on the overall design of the device and can be estimated through both modeling and experimental methods. Eddy current losses in the AMRTA components which have direct impact on the cooling capacity (i.e. cold heat exchanger) are negligible. Based on work by Tura et al. (2006) an ambient heat leak through the AMRTA structure was estimated as, Qamb;c þ Qleak y 0:08WK1 ðTH TC Þ
(22)
With the external losses defined the total cooling power for the device is defined as QC ¼ Hðx ¼ 0Þ Qeddy;c Qamb;c Qleak
(23)
Using Eq. (23) and the temperature span defined by Eq. (9) the model can be tested against experimental data.
4.
Model description
(16) (17)
The AMR model uses the developed theory and correlations in conjunction with material data, geometric parameters and operating conditions to determine the temperature distribution across the AMR. From this, performance parameters such
936
i n t e r n a t i o n a l j o u r n a l o f r e f r i g e r a t i o n 3 6 ( 2 0 1 3 ) 9 3 2 e9 4 0
as expected cooling capacity and work input can be determined. Required inputs include material properties, effective high and low field strength during a blow, utilization and temperatures at the ends of the regenerator. Since the regenerator end temperatures are an input to the simulation, the temperature span of the regenerator is artificially created while the shape of the distribution is determined by Eq. (7). The cooling capacity and work for a given temperature span is then post-calculated.
4.1.
Fundamental operation
In the numerical model the regenerator’s temperature distribution, Ta(x), is determined by Eq. (7) using a discrete number of points. This equation is solved along the entire domain using MATLAB’s bvp5c ODE solver in conjunction with the user-specified boundary conditions Ta(x ¼ 0) and Ta(x ¼ L). Solution of this equation requires multiple iterations as several of the terms in the ODE are temperature dependent. This includes the adiabatic temperature change, DT, which is simulated by assuming isentropic magnetization of the material such that, DT ¼ Tðs1 ; H2 Þ Tðs1 ; H1 Þ
(24)
Entropy values for Eq. (24) are determined using the integral of the material’s specific heat. The material data was provided by AMES Iowa National Lab for Gd, Gd0.74Tb0.26 (GdTb) and Gd0.85Er0.15 (GdEr) which corresponds to the materials used in the AMRTA experiments. This data is interpolated in MATLAB using the TriScatteredInterp function allowing for the user to select the high and low field experienced by the magnetic material. Furthermore, if the Curie temperature and magnetocaloric effect are known at two or more field strengths for a certain material, this same technique can use the interpolated Gd data to map the material’s data over a wide operating range. This is due to the similarities in the shape of the second-order material data. This technique is used for GdEr and GdTb where data is available only at 0 T and 2 T. By mapping against Gd, approximate data for GdEr or GdTb can be generated for a wide range of field choices. With Ta(x) computed the other three state points are found using knowledge of the magnetic cycle (Rowe, 2011a); this allows for the fluid temperature across the regenerator to be calculated. Cooling capacity is then found using the theory described in the preceding sections. Other useful parameters such as magnetic work, pumping power, efficiency and figures of merit can also be post-calculated using the temperature profile and cooling capacity. Model inputs are selected such that the outputs can be compared to experimental data. The model takes approximately 2 s to solve for a single set of boundary conditions on a laptop computer with a 2.40 GHz Intel Core i5 processor.
4.2.
Experimental data
Experimental data used for validation includes the temperature span across the regenerator, Tspan, at a given hot end temperature, TH, for a specified heat load or cooling capacity, Qc. This allows for two important curves to be created which provide a measure of performance for a given material,
operating parameters and regenerator configuration. These curves are Qc vs Tspan for a fixed TH, and Tspan vs TH for a fixed Qc (typically Qc ¼ 0 W). Furthermore, experimental data sets used to compare the model results vary by material, utilization, heat load and frequency. Fixed and variable parameters used in the experiments are described in Table 1. The regenerators are made of crushed particles which are formed into pucks containing between 40 and 45 g of magnetocaloric material each. Regenerators in the experiments of interest are composed of one, two or three pucks made from Gd, GdEr or GdTb (Tura, 2005). Working fluid properties are found using NIST’s REFPROP program which allows temperature dependent variables to vary between data points. These include the dynamic viscosity, specific heat and density. These properties are relatively constant for the fluid and conditions used in the SC test apparatus for room temperature experiments. The magnetic field used for all of the tests is provided by a 2 T superconducting magnet. The magnetic field seen by the material as a function of time is approximated as a sinusoid. For this reason, the discrete applied field in the model uses RMS values to relate actual operating conditions during a blow period to the step-wise waveforms assumed in deriving Eqs. (1)e(3).
4.3.
Results
A variety of experimental conditions are evaluated to ensure the model is tested in more than a single operating region of the AMR. Once validated, a higher level of confidence then exists in predicting previously untested operating conditions or regions currently unavailable to the experimental device. As such, the following results vary by material type, material mass, utilization, frequency and heat load.
Table 1 e Temperature independent parameters of the superconducting AMRTA used in producing the experimental results. Test parameters Frequency (Hz) Stroke length (m) Effective displacer area (mm2)
0.65e0.8 0.21 2513
Regenerator Material mass (g per puck) Porosity Regenerator length (mm per puck) Regenerator diameter (mm) Characteristic particle diameter (mm) Particle type Number of regenerators
40e45 w55% 25 25 560 Crushed particles 2
Magnetic material Thermal conductivity (W K1) Density (kg m3) Reference specific heat (J kg1 K1)
10.3 7900 381
Working fluid (Helium) Specific heat, Cp (J kg1 K1) Thermal conductivity (W K1) Pressure (atm)
5139 0.153 3e9.7
937
i n t e r n a t i o n a l j o u r n a l o f r e f r i g e r a t i o n 3 6 ( 2 0 1 3 ) 9 3 2 e9 4 0
Fig. 2 compares the model results with all of the aforementioned losses mechanisms included to single pucks of Gd, GdTb and GdEr and under no-load (Qc ¼ 0 W) conditions. In all three data sets the model predictions follow the same general trend as the experimental data. The primary discrepancy exists at high TH values relative to the Curie points where the temperature span begins to drop abruptly (best seen for GdTb). While the general trends are correct for all three materials, the model tends to significantly over predict the performance of the GdTb regenerator. As can be seen from the GdTb plots, decreasing the assumed magnetocaloric effect for this material by a factor 0.77 as compared to the raw data corrects the prediction. As with GdEr, the model results still under-predict the experimental data when operating at temperatures well above the Curie point. It is interesting to note that another more detailed numerical model of the AMR (Dikeos, 2006) produced a similar result for the GdTb pucks. In that case, a better fit was obtained by increasing the static conductivity of the bed. The fact that both assumptions are effective in correcting model predictions is supported by a simplified study by Rowe (2011b) where a non-dimensional parameter determining the cooling power and efficiency of an AMR was found: I¼
F DT k T
(25)
In the work by Dikeos (2006), the impact of increasing the static conductivity decreases I; likewise, in this work, decreasing the adiabatic temperature change leads to a decrease in I. In effect, the adjustments produce the same result in the governing energy balance. Both approaches may be justified given the lack of detailed measurements for the particular batch of GdTb used within experiments as compared to the batch used in measuring the material properties. One possibility is differences in the manufacturing process between batches of GdTb, which would explain why
no adjustments are needed for Gd and GdEr. Furthermore, the irregular particle shapes result in some uncertainty regarding the regenerator thermal effectiveness. The model predictions can also be compared to layered regenerators where experimental data is available for two and three layer regenerators using Gd, GdEr and GdTb. Based on the results for the single GdTb bed, the corrected adiabatic temperature change (0.77 scaling) for this material is used in the model. Fig. 3 demonstrates that the model predicts both the layering aspect and the combination of materials relatively well. The model tends to be less accurate the more the operating temperature is increased above the average Curie temperature of the materials in the bed and as the utilization decreases. Because the gas displaced through the regenerator is constant and adding a material increases the overall regenerator mass, the utilization decreases from w0.3 for the single puck experiments down to w0.1 for the three material tests. The sensitivity of the temperature span to applied load is shown in Fig. 4a for a single Gd puck at TH ¼ 292 K. From the figure the experimental results are in good agreement with the model. The model has been used to also show the effect that system losses have on the predicted cooling capacity. Recalling Eq. (23), the net cooling power is determined by the enthalpy flow at x ¼ 0, H(x ¼ 0), less the parasitic heat leaks from ambient, Qamb,c, and from the hot side of the system, Qleak. The impacts of regenerator thermal ineffectiveness, internal losses and the external losses can be easily shown. In this case internal losses are due to only viscous dissipation. The top curve represents the temperature span for a system with perfect regenerator effectiveness, no viscous dissipation and no ambient cooling losses. Below this curve is the same scenario but with regenerator ineffectiveness included. The next curve then shows the cooling power when only the external losses are excluded while the final solid curve accounts for all system losses and is the predicted cooling capacity for the device. This analysis of system losses is important in optimization studies for determining where
35 55
30
GdTb GdEr
45
20 TSpan [K]
T
Span
[K]
25
50
Gd
15
Gd−GdEr
40 35
10
30
5
25
0 250
Gd−GdTb−GdEr
Gd−GdTb
Gd−Gd
20
260
270
280 290 TH [K]
300
310
Fig. 2 e TSpan vs TH for single pucks of Gd, GdTb and GdEr at Qc [ 0 W, 0.65 Hz and a charge pressure of 9.5 atm. Solid lines represent experimental data while dashed lines are simulations from the model. In the case of GdTb pucks the temperature span is plotted using both original material data and data altered to reduce the MCE.
15 280
285
290
295 300 T [K]
305
310
315
H
Fig. 3 e TSpan vs TH for two and three puck multilayered regenerators at Qc [ 0 W, 0.65 Hz and a charge pressure of 9.5 atm. Solid lines represent experimental data while dashed lines are simulations from the model. Note the vertical axis is offset from 0 K for clarity.
938
i n t e r n a t i o n a l j o u r n a l o f r e f r i g e r a t i o n 3 6 ( 2 0 1 3 ) 9 3 2 e9 4 0
a
30 25
0W
20 TSpan [K]
experimental devices can be improved and potential effects these improvements have on the temperature span and cooling capacity. As is demonstrated in Fig. 4b, for the same operating parameters but at a temperature above the Curie point (TH ¼ 305 K), viscous dissipation losses have a much greater impact on the temperature span than at lower temperatures. In Fig. 5 the temperature span is plotted against three different TH values for two heat load cases, 0 W and 6 W. While both data sets follow the trend of the experimental curve, the model under-predicts experiments at the higher frequency. Potential causes of differences between model predictions and experiments are discussed in the next section. The final experimental result accounts for a variation in the utilization. For the AMRTA the utilization is varied by altering the pressure of the helium in the system thereby changing the density and subsequent thermal mass of the displaced fluid. Increasing the utilization raises the pressure
15 6W
10 5 0
290
295 T [K]
300
305
H
Fig. 5 e TSpan vs TH for single-pucks of Gd at 9.5 atm, 0.8 Hz and varying load. Solid lines represent experimental data while dashed lines are simulations from the model.
20 drop in the system and, in the case of the AMRTA, increases the temperature span achievable by the regenerator up to certain value. Shown in Fig. 6 are the temperature span curves at three different charge pressures, or utilizations. Similar to the heat load sensitivity curves in Fig. 5 the temperature span follows the experimental trends. Furthermore, at a pressure of 3 atm the temperature span is trending lower at TH ¼ 304 K.
Q [W]
15
10
5
Exp. H(x=0)−Q
−Q
amb,c
(all losses)
leak
H (with int. losses, no ext. losses) H (no int. losses, no ext. losses) H (100% effectiveness, no losses)
0 0
5
10
T
Span
b
5. 15
[K]
20
20
In all of the tested scenarios the model results followed the trends of the experimental results. In the case of varying utilization and heat load some differences are seen between the experiments and the model. Differences are most
25
15
9.5 atm
20
10
5
H(x=0)−Q
−Q
amb,c
leak
TSpan [K]
Q [W]
Discussion
(all losses)
10
H (with int. losses, no ext. losses) H (no int. losses, no ext. losses) H (100% effectiveness, no losses)
0 0
5
10 T
Span
[K]
15
6 atm
15
3 atm 5
20
Fig. 4 e Qc vs TSpan for single-pucks of Gd at 0.65 Hz and a utilization corresponding to a charge pressure of 9.5 atm, (a) TH [ 292 K and (b) TH [ 305 K. Results are plotted to show the effect that internal and external losses have the cooling capacity and temperature span. The difference in the impact of viscous dissipation (internal losses) is also shown between the two plots.
0
270
280
290 T [K]
300
310
H
Fig. 6 e TSpan vs TH for single-pucks of Gd at 0.65 Hz, Qc [ 0 W and utilizations corresponding to charge pressures of 3 atm, 6 atm and 9.5 atm. Solid lines represent experimental data while dashed lines are simulations from the model.
i n t e r n a t i o n a l j o u r n a l o f r e f r i g e r a t i o n 3 6 ( 2 0 1 3 ) 9 3 2 e9 4 0
significant when the operating temperatures are increased above the Curie points of the materials in the AMR and when the utilization is decreasing. In general, when these conditions arise the model tends to under-predict performance. These discrepancies are expected due to the simplicity of the model relative to the physics occurring at the time of the experiments. Examples of this are seen in the physical representation of the AMR inside the model including the characteristic particle diameter and simulated RMS field strength. The model also assumes that pure Gd, GdTb and GdEr exist in the regenerator. In the experiments however the regenerator materials are coated in a thin layer of epoxy which will increase interfacial heat transfer resistance, act as a parasitic thermal mass, and may cause flow-channeling. Additionally, since the cooling power of the AMRTA experiments are small (w20 W), any errors in modeling the losses can have a large impact on the expected temperature span as compared to experiments with large cooling capacities. Some errors also exist in the replication of GdEr and GdTb specific heat curves from Gd data as the specific heat curves of second-order materials are similar but not exactly the same. In the operating space defined by the SC AMRTA the model is found to produce repeatable and continuous curves based upon the results in Fig. 2 through Fig. 6. In general, TH values below the Curie point solve faster than those above the Curie point. The model has difficulty solving for the temperature distribution at TH values high above the Curie point when both large pressure drops and low utilizations (<0.15) are present. This is due to the non-dimensional thermal conductivity becoming too low to solve the energy balance differential equation. It was found that setting a minimum value on the conductivity in the differential equation rectified this problem. It was noted by Tura (2005) that obtaining data with TH values high above the Curie points is in general experimentally unstable and difficult to replicate, particularly for varying heat load scenarios at low utilizations. Overall, the results sufficiently follow the tested scenarios enough to provide a basis for proceeding with the model in optimization and design cases where general trends of the performance are of primary importance. The next stage of verification for the model involves systems where the working fluid is a liquid and the expected cooling loads are much higher. Similar to the results used in this paper, an adequate amount of experimental data is available to verify the revised model. In parallel, the model is to be further used in optimization programs and mapping of various operating spaces. Providing quick solutions that accurately predict the trends of the AMR allows for parameters such as the regenerator length to be linked to the length of the permanent or superconducting magnet. The device performance at a given field strength can then be linked with the magnet cost necessary to generate the matching size and magnitude of the given field. This type of design balances the importance of performance and cost requirements on future AMR devices.
6.
Conclusions
Through the use of experimental results it has been shown that a simplified AMR model using real material properties
939
can predict the general trends of a device. This is verified from results using a 0e2 T magnetic field and three different types of materials, Gd, GdTb and GdEr, which were evaluated as single layer and multilayer regenerators. These scenarios showed good sensitivity with experiments both with no-load and an applied load as well as at two different frequencies, 0.65 Hz and 0.8 Hz, and three charge pressures, 9.5 atm, 6 atm and 3 atm. The model also provides flexibility in observing the individual predicted losses of the AMR signifying areas of device improvement. In addition, operating regions which are difficult for the model to predict have been identified and provide direction for further improvements. The model will now be tested against experimental systems using liquid working fluids and will be implemented into costing and optimization programs.
Acknowledgments Support in the form of a Canadian Graduate Scholarship from the Natural Sciences and Engineering Research Council of Canada is greatly appreciated. Partial support for this work was provided by H2Can Strategic Network. The authors also acknowledge the unwavering support of colleagues and friends whose stimulating conversations allowed this paper to come to fruition.
references
Dikeos, J., 2006. Development and Validation of an Active Magnetic Regenerator Refrigerator Cycle Simulation. MASc thesis, University of Victoria, Victoria. Engelbrecht, K., 2006. The effect of internal temperature gradients on regenerator matrix performance. Trans. ASME 128, 1060e1069. Engelbrecht, K., 2008. A Numerical Model of an Active Magnetic Regenerator Refrigerator with Experimental Validation. PhD thesis, University of Wisconsin-Madison. Ergun, S., 1952. Fluid flow through packed column. Chem. Eng. Prog. 48, 125e127. Kaviany, M., 1995. Principles of Heat Transfer in Porous Media. Springer, New York, NY. Kittel, P., 1990. Eddy current heating in magnetic refrigerators. In: Fast, R.W. (Ed.), Adv. Cryog. Eng. Plenum, vol. 35B, pp. 1141e1148. Legait, I., Kedous-Lebouc, A., Rondot, L., 2009. Numerical simulation and analysis of the refrigerant bed behavior using fluent software. In: Egolf, P.W. (Ed.), Third International Conference on Magnetic Refrigeration at Room Temperature. International Institute of Refrigeration. Liu, M., Yu, B., 2011. Numerical investigations on internal temperature distribution and refrigeration performance of reciprocating active magnetic regenerator of room temperature magnetic refrigeration. Int. J. Refrigeration 34, 617e627. Nam, K., Jeong, S., 2003. Measurement of cryogenic regenerator characteristics under oscillating flow and pulsating pressure. Cryogenics 43, 575e581. Nielsen, K.K., Bahl, C.R.H., Smith, A., Pryds, N., Hattel, J., 2010. A comprehensive parameter study of an active magnetic regenerator using a 2D numerical model. Int. J. Refrigeration 33, 753e764.
940
i n t e r n a t i o n a l j o u r n a l o f r e f r i g e r a t i o n 3 6 ( 2 0 1 3 ) 9 3 2 e9 4 0
Nielsen, K.K., Tusek, J., Engelbrecht, K., Schopfer, S., Kitanovski, A., Bahl, C.R.H., Smith, A., Pryds, N., Poredos, A., 2011. Review on numerical modeling of active magnetic regenerators for room temperature applications. Int. J. Refrigeration 34, 603e616. Risser, M., Vasile, C., Engel, T., Keith, B., Muller, C., 2010. Numerical simulation of magnetocaloric system behavior for an industrial application. Int. J. Refrigeration 33 (5), 973e981. Rowe, A., 2011a. Thermodynamics of active magnetic regenerators: part I. Cryogenics 52, 111e118. Rowe, A., 2011b. Thermodynamics active magnetic regenerators: part II. Cryogenics 52, 119e128. Tura, A., 2005. Active Magnetic Regenerator Experimental Optimization. MASc thesis, University of Victoria, Victoria.
Tura, A., Roszmann, J., Dikeos, J., Rowe, A., 2006. Cryogenic active magnetic regenerator test apparatus. Adv. Cryog. Eng. 823, 985e992. Tusek, J., Sarlah, A., Zupan, S., Prebil, I., Kitanovski, A., Poredos, A., 2010. A numerical optimization of a packed bed AMR. In: Egolf, P.W. (Ed.), Fourth International Conference on Magnetic Refrigeration at Room Temperature. International Institute of Refrigeration. Vortmeyer, D., Schaefer, R.J., 1974. Equivalence of one- and twophase models in packed beds: one dimensional theory. Chem. Eng. Sci. 29, 485e491. Wakao, N., Kaguei, S., 1982. Heat and Mass Transfer in Packed Beds. Gordon and Breach Science Publishers, New York, NY.