Multimodal freezing system for cryogenic 3D printing

Cryogenic 3D printing (Cryo-3DP) technique is used to realize the frozen 3D objects by layered deposition of the liquids inside an insulated chamber. The operating temperature of the process ranges from − 20 to − 25 °C. Cryo-3DP demands high cooling rates initially to reach the process temperature followed by lower cooling rates to sustain it. Multimodal freezing system proposed in this study uses more than one cooling mode that helps in achieving the variable cooling rates. The present system uses two modes of cooling, viz., vapor compression refrigeration (VCR) and CO2 injection. A detailed design of the multimodal freezing system is carried out for the chamber size of 200 mm × 200 mm × 195 mm. The performance of the system is analysed numerically. The results were experimentally validated using the prototype of the multimodal system developed as a part of the present study. The results show that the multimodal system reduces the initial cooling time substantially by providing a high cooling rate that rapidly cools the chamber to initiate the cryo-3DP. VCR system provides a low but sustained cooling rate that maintains the desired temperature. The present study proves the significance of multimodality in cooling system that is fit to deploy in a commercial cryogenic 3D printer.


Introduction
The process of Cryogenic 3D printing (Cryo-3DP) is carried out at sub-zero temperatures, typically around − 20 to − 25 °C [1][2][3][4]. The process starts with the cooling of the work chamber from the room temperature to − 20 to − 25 °C. Once the temperature of the work chamber is achieved, the material deposition subsystem spreads the liquid layer selectively on the work platform. Motion subsystem manoeuvres the deposition head as per the required layer geometry. The deposited layer instantly solidifies. Another layer of the liquid is selectively spread on the previous layer, that undergoes freezing. The process continues till many such layers are deposited onto one another, realizing the three-dimensional object [1,[5][6][7][8]. Figure 1 shows the schematic diagram of the typical cryogenic 3D printer. The deposition and motion systems are enclosed in the work chamber that is insulated to maintain the sub-zero temperatures. The cooling is provided by the refrigeration system.
The process of cryogenic 3D printing finds applications in the field of investment casting (see Fig. 2) [2], photopolymer additive manufacturing for the support structure [3], drug encapsulation in the pharma industry [8,9], micro -manufacturing [4] and architectural models [5].
The proposed system works in two steps. In the first step, a high cooling rate is achieved to cool the work chamber from the room temperature to the desired working temperature of − 20 to − 25 °C. In the second step, a lower cooling rate is maintained to sustain the work chamber's temperature. The temperature vs time characteristics of the system are hypothesized in Fig. 3. Since more than one mode of cooling is being used, the authors term the proposed system as multimodal freezing system (MFS).
To achieve the different cooling rates, it is proposed to use CO 2 injection system to rapidly cool the chamber to − 20 °C, followed by the vapour compression refrigeration (VCR) system for sustained cooling of the chamber at a constant temperature. VCR is a widely available, popular cooling  (a) 3D Printed Ice boundary (b) Block of Ice as a Pattern (c) Cast Part system due to its high efficiency and cost-effectiveness. CO 2 gas is inexpensive and readily available in the bottled form. Figure 4 shows the schematic diagram of the MFS. The evaporator coils and CO 2 inlet are located in the freezing chamber.

The architecture of the work chamber
As shown in Fig. 5a-c, the work chamber walls are made of polyurethane foam (PUF) of thickness 80 mm lined with a stainless steel sheet from the outer surface. The work platform of area 200 mm × 200 mm and thickness 10 mm is attached to the freely movable platform with four support rods. The stepper motor operated platform moves the table up and down inside the freezing chamber. The motion system is properly insulated against the sub-zero temperatures. The chamber is covered with the bellows at the top (see Fig. 5a, b). The printhead is a delicate part with several micro-nozzles. It is exposed to the sub-zero temperatures. It may result in freezing of water inside the nozzles, choking them up. A heater is used to maintain the temperature of the nozzle plate above freezing point. The details of the set up are mentioned in the Sect. 4.

Cooling load estimation
Total cooling load Q Total comprises five different parts, as shown in Table 1. Static cooling loads are single-time cooling loads that arise in step 1 (A and B). Continuous cooling loads arise in step 2, where sustained cooling is required for freezing the water and minimize the heat from the nozzle heater (C, D and E). The ambient temperature is 35 °C, and the process temperature is − 20 °C.

Mathematical formulation
The objective of the mathematical formulation is to determine the time required for cooling the components. (1)

CO 2 system
In the proposed MFS, the CO 2 system is used for the workspace's initial cooling. Initial cooling is required for to A and B part of the cooling load. CO 2 inlet is inside the chamber. The isenthalpic expansion of the CO 2 brings about the Joule-Thomson cooling effect. The drop in the temperature due to pressure drop can be determined with Eq. 2.
(2) ΔT = JT ΔP In Eq. 2, μ JT is the Joule-Thomson coefficient. μ JT is 13 K (MPa) −1 for CO 2 [6]. For the pressure drop of the 5 MPa in the given case, the temperature change is 65 °C, i.e. at the outlet, the temperature is − 30 °C.
The cooling capacity offered by CO 2 can be calculated using Eq. 3. CO 2 gas is converted to solid particles (dry ice) during expansion, but the percentage of the dry ice formation is as low as 2% [7]. It can be neglected as the CO 2 system is used only for a short duration. Specific isobaric heat capacity for CO 2 at − 30 °C and atmospheric pressure is 0.798 kJ kg −1 K −1 ( Table 2). The orifice condition determines the flow rate as per Eq. 4 [7]. The geometrical and flow parameters required to determine mass flow rate are given in Table 2.
It is assumed that CO 2 immediately replaces the air inside the chamber at − 30 °C; thus, cooling of chamber air can be neglected. Heat transfer is predominantly convective. Convective coefficient is found out from the Nusselt number (Nu) correlation (Eq. 5). Equation 5 is a standard Dittus-Boelter correlation for the flow over a flat plate that is chosen since the flow over a horizontal flat work platform is considered here [8].
where Re < 10 5 , 0.6 < Pr < 2000. The Nusselt number correlation gives the heat transfer coefficient of 236.25 W m −2 K −1 . With the given heat transfer coefficient, the heat transfer rate can be determined by Eq. 6.
The heat transfer rate in the CO 2 system is 1.06 kW; thus, it takes 88 s to cool the work platform to − 20 °C.

Vapour compression refrigeration (VCR) system
In the VCR system, the first heat transfer occurs from the chamber medium (air) to evaporator coils and then from the work platform to the chamber medium. Figure 6 shows a typical VCR cycle where process 1-2 is compression, 2-3 Heater for nozzles Q D = qkW Film heater of 10 W 0.010 kW E Air infiltration (Assumed 10% of the continuous load) is isobaric condensation, 3-4 is the isenthalpic expansion, and 4-1 is evaporation. The enthalpy of the refrigerant is the function of the saturation temperature (Eq. 7) [9]. Enthalpy of the refrigerant at various temperatures can be expressed with the empirical relations given by Eqs. 8 and 9.
T s (°C) is saturation temperature, a 1 , a 2 , a 3 , a 4 and b 1 , b 2 , b 3 , b 4 are coefficients, h sl and h sv are the enthalpies of the refrigerants in the saturated liquid and vapor conditions. T sl (°C) and T sv (°C) are the saturated liquid and vapor temperatures of the refrigerant [10].
As per the ASHRAE enthalpy datum, the enthalpy polynomial coefficients for saturated liquid R134a are expressed by Eq. 10. Equation 10, that gives the enthalpy of R134a in liquid state, is used for the present system since it gives the accurate curve-fit relation for R134a [10]. Equation 7 helps to predict the enthalpy data with ± 0.2 kJ kg −1 accuracy for the temperature range − 40 °C to 70 °C, which is within the present range of the study. Here, T l is the liquid temperature of the refrigerant where ∆T b = T sat -T l (°C) ≥ 0.
The saturation liquid temperature for the refrigerant is measured after 32 °C.
Similarly, the ASHRAE enthalpy datum for saturated vapor can be predicted with the Eq. 11 [10], where T is defined in °C.
As the compression is isentropic, the saturation temperature can be assumed evaporation temperature of the refrigerated enclosure, i.e., − 23.3 °C.
For superheated vapor, the enthalpy datum can be calculated with Eq. 12 [10]. Here, T s is the superheat temperature and ΔT s = T s − T sat .
The superheated temperature ( T s ) is 54.4 °C, and the evaporation of the refrigerant occurs at − 23.3 °C. The refrigerant used is R134a. The gas pressure after compression is found to be 1. For small compartment cooling applications, a compressor with a cooling capacity of 281 W is selected. The mass flow rate capacity of the compressor is 36 g s −1 .
If VCR alone were used, the time required to handle the load A and B and cool the system to − 20 °C would be higher than that of the CO 2 system.
where h is the heat transfer coefficient, A is the surface area of the heat exchange and T 2 is evaporator coil temperature and T 1 is the air temperature. The cooling takes place by natural convection. The air passes over the cooling coils of the evaporator. The value of h can be estimated using the Nusselt number correlation for the flow over a horizontal cylindrical tube, i.e. cooling coils of the VCR system (Eq. 16) [11]. The heat transfer coefficient is 3.07 Wm −2 K −1 . Thus, Eq. 15 yields the value of heat transfer rate Q 2 as 0.016 kW. Therefore, if we use VCR for part A and B of the cooling load, it takes 97.23 min under ideal insulation and absence of infiltration for chamber medium (air) to reach the desired temperature. However, there are parts C, D and E of the cooling loads that play a crucial role that extends the cooling time for VCR system, to study their effects in details and to (12) determine the cooling time precisely, numerical simulation is carried out that is discussed subsequently.
The h is calculated based on Eq. (18) for the plate heat transfer, and found to be 0.46 Wm −2 K −1 [12]. The heat transfer rate between the chamber medium and the work platform is given by Eq. 17. The combined effect of Q 2 and Q 3 determines the total time taken for the cooling of the work platform. The time is determined with numerical study subsequently.
The results of the mathematical analysis have been correlated and discussed in the Sect. 5.

Experimental setup and numerical procedure
The time to reach the target temperature is predicted by mathematical analysis as 88 s (work platform) and 97 min (chamber medium) for the CO 2 system and VCR system. The experimental and numerical study is carried out to ascertain the time.

Experimental setup
The experimental setup, as shown in Fig. 7a  head's motion in the X and Y directions. Evaporator coils are positioned around the chamber walls uncovered to increase the heat transfer efficiency as shown in Fig. 7b. Thermocouple measures the temperature variation (Fig. 7c). A 10 W polyimide foil heater is sandwiched between the nozzle plate and the base plate as shown in Fig. 7d. Foil heater is used for maintaining the nozzle plate at 60 °C that prevents the nozzles from freezing. Experiments show that the nozzle plate maintained at 60 °C is effective in keeping the nozzles warm without melting the uppermost frozen layer of ice.
Printhead is the crucial part of the architecture since it dispenses the water on the substrate. Printhead is a delicate part that demands careful handling to keep its nozzles from blockage. Xaar 1201 printhead is used for the present study (Fig. 7d). It is a piezoelectric technology-based printhead with 1280 nozzles arranged in four rows of 320 nozzles each. Each nozzle has a tiny chamber made up of the piezoelectric walls (Fig. 8a).
The piezoelectric walls undergo expansion and contraction based on the voltage signals received from the controlled that ejects the water from the chamber through the nozzles. The water is ejected in the desired pattern by the controller that is synchronized with the X-Y motion system. Figure 8b shows the assembled printhead and Fig. 8c shows the openings of the piezoelectric nozzles with the diameter of 20 µm at 100 × magnification. CO 2 is released in the chamber, and the temperature is recorded using a K-type thermocouple. The change in the temperature is measured every 5 s as the temperature changes rapidly. Similarly, the VCR system is tested, and the temperature is noted after every 2 min as the change in the temperature is slow. The observations of the experiment are mentioned in Fig. 10.

Uncertainty of measurement
The sensitivity of the measuring instruments gives rise to the variation in the measured values. The measurement errors are addressed in Table 3.

Numerical methodology
Three-dimensional numerical simulations are carried out for the MFS using flow simulation CFD package by Solid-works® available commercially. The software solves the Navier-Stokes equations that are formulations of mass, momentum, and energy conservation (Solidworks Flow Simulation Technical Reference, 2012).
Conservation of mass: Conservation of momentum: Conservation of energy: where H = h + u 2 2 and is density of the fluid,u j , u j are the velocity components, x, y and z are the coordinates, ij is the viscous shear stress tensor, P is the pressure, H is the total energy, q i is the diffusive heat flux, Q H is the heat source or sink per unit volume, is the rate of dissipation of the turbulent kinetic energy, S i is the summation of mass distributed external force per unit mass due to porous media, gravity, and rotation where S i = S porous i h is the thermal enthalpy. Subscripts denote the summation over three coordinates. ij is the Kronecker delta function,P B is turbulent generation due to buoyancy, t is turbulent eddy viscosity, is dynamic viscosity, f 1 andf 2 are viscosity factors, C 1 ,C 2 , , k are constants defined empirically. Pressure 60 Bar ± 0.5 Bar Fig. 9 Section view of the geometry for the simulation The geometry is as shown in Fig. 9. An adaptive meshing approach is used. The results of the grid independence study are shown in Table 4. The minimum temperature after 60 s is recorded for the CO 2 system and after 3600 s for the VCR system. After 85,000 cells, the temperature variation is minimal (0.03% and 0.1% for CO 2 and VCR systems, respectively). For CO 2 system simulation, the inlet is set up as a pressure inlet with a total pressure of 60 bar, and the outlet is set up as a pressure outlet with environmental pressure, 1.013 bar. The surface goals for the worktable surface are set as minimum solid temperature. The heat transfer coefficient is assigned as 236 Wm −2 K −1 for the worktable surface. The work chamber walls are considered ideal walls. (see Fig. 9). For the VCR system, the geometrical and flow parameters are mentioned in Table 5. Since it is a natural convection problem, the h values for the cooling coil and the work surface are estimated as 3.07 and 0.48 Wm −2 K −1 , as mentioned in Sect. 3.2. The minimum cooling coil temperature obtained by the VCR system at the evaporator is − 23.3 °C. The work chamber walls are considered ideal walls. (see Fig. 9). The temperature is recorded for the study surface. The results are mentioned in Fig. 10.

Results and discussion
The results show the temperature vs time characteristics of the CO 2 and VCR systems. The temperature profiles of the chamber are also plotted for different time instances. Figure 10 shows the temperature-time graph obtained by the numerical simulation and an experiment of the CO 2 and the VCR system. Mathematical analysis (Sect. 3) gives an estimate for the cooling time. It is observed that the experimental and numerical results agree with each other as estimated by mathematical analysis. It is observed in Fig. 10 that the target temperature − 20 °C is obtained in 103 s with the help of the CO 2 system. Mathematical analysis estimates the time to be 88 s which fairly agrees with the experiment. Mathematical analysis of VCR system yields the cooling time of 97.23 min; however, it is under the condition of perfect insulation and absence of heat sources. In real system, the heat loads C, D and E play a major role since the VCR system acts slowly with a low cooling rate. Numerical analysis is carried out to precisely determine the cooling time and verified with the experiments. Under real conditions with air infiltration and nozzle heater present in the chamber, VCR system takes 3 h to reach from 35 °C to 0 °C and then the cooling rate further slows down that takes another 3 h to reach − 20 °C, as observed experimentally.

Variation of the temperature with time
In the CO 2 system, the CO 2 gas acts as a cooling medium and directly exchanges the heat with the work chamber. However, in the VCR system, the heat is first exchanged with the chamber medium (air) and then the cold air exchanges the heat with the chamber contents, i.e. two convective heat transfer coefficients are involved. It results in slower heat transfer as compared to the CO 2 system.   Figure 11 shows the temperature contours of the solid parts and the chamber medium fluid (CO 2 ). It is observed that the temperatures of the CO 2 gas at 130 s range from − 18.88 to − 78.63 °C and the work platform from − 8 to − 22 °C. Figure 12 shows the work platform temperature is in the range of 22.04 to 28.52 °C at 12,600 s. The experimental results suggest the same.

Demonstration of the ice 3D printing
Ice structures were 3D printed using the present experimental set-up. The model material is demineralized water. Undercut-free prismatic geometries are printed to validate the present system. The 3D printed ice parts are shown in Fig. 13. Figure 13a shows the CAD model of the embossed letters of the institute's acronym. Figure 13b is the 3D printed ice part from the CAD model. Figure 13c shows the CAD model of the ornamental parts such as pendants and Fig. 13d shows the 3D printed ice part.

Applications of Cryo-3DP
Cryo-3DP process has applications in the field of investment casting, polymer 3D printing, pharma etc. AM parts of ice are used as a pattern in investment casting. Figure 14a shows the ice AM pattern of a gear-like shape, Fig. 14b shows the mould and Fig. 14c shows the casting made out of the mould.
Ice is an ideal support material for the polymer jet AM processes since it melts away readily without any residue. This approach is very valuable to obtain precious miniature polymeric parts as the support removal for miniature part is cumbersome with traditional support material. It also reduces the wastage of the precious polymeric material as a support that is usually discarded after use. Figure 15a shows a model of an ant completely printed with ice as a support, Fig. 15b shows the part when the support is melted.
Ice structures produced using cryo-3DP act as templates for microfluidic assay kits and drug delivery devices. Subzero photopolymer liquid poured over these templates is cured by flood light. Later, the ice templates are melted away [14]. Figure 16a depicts the schematic diagram of the process to fabricate pre-sealed drug delivery chips using ice templates. The ice pillars of various heights are shown in Fig. 16b to test the capability of the process. Figure 15c shows the pre-sealed drug delivery chip made by process shown in (a), where the drug is released when the chip is pressed by the finger. Figure 15d shows the pre-sealed mixing reactor, where reagent in the reservoir II flows into the reservoir I after pressed by the finger; and flows back to reservoir II once finger is released.

Conclusion
Cryo-3DP process demands the high cooling rates for the initial cooling of the chamber so that deposition can begin as soon as the machine is started. However, traditional VCR system operates slowly and it takes several hours to reach the desired temperature before the machine could be operational; hence, multimodality using CO 2 gas expansion along with VCR system is proposed in this study.
The numerical as well as experimental results confirm that the proposed system reduces the chamber cooling time significantly. As a standalone VCR system takes up to 6 h to reach the working temperature of − 20 °C due to air infiltration, nozzle heater in the chamber and warm inlet water; multimodal system can reach the working temperature within 103 s. CO 2 is used intermittently for initial cooling and instantaneous cooling when required. Although CO 2 is involved, it used intermittently, hence its consumption is low. The sensors are required to constantly check CO 2 concentration in the air. The safe limit for indoors is below 5000 ppm.
VCR provides the sustained cooling and maintains the constant work chamber temperature throughout the process. Multimodality of the cooling system is beneficial for the commercialization of the cryo-3D printers since it gives the instantaneous cooling and variable cooling rates.
As listed in the Sect. 5.4, there are several important applications where the proposed system is useful. The multimodality helps in rapid cooling that can be useful for the cryo-3DP systems for the commercial use.  Ice as a support material for the photopolymer 3D printing [3] (a) Schematic view of the fabrication process (b) micro-pillars printed on the surface (c) pre-sealed drug delivery chip (d) pre-sealed mixing reactor Fig.16 Microfluidic applications of Cryo-3DP [14]