Labyrinth ice pattern formation induced by near-infrared irradiation

See allHide authors and affiliations

Science Advances  29 Mar 2019:
Vol. 5, no. 3, eaav1598
DOI: 10.1126/sciadv.aav1598


Patterns are broad phenomena that relate to biology, chemistry, and physics. The dendritic growth of crystals is the most well-known ice pattern formation process. Tyndall figures are water-melting patterns that occur when ice absorbs light and becomes superheated. Here, we report a previously undescribed ice and water pattern formation process induced by near-infrared irradiation that heats one phase more than the other in a two-phase system. The pattern formed during the irradiation of ice crystals tens of micrometers thick in solution near equilibrium. Dynamic holes and a microchannel labyrinth then formed in specific regions and were characterized by a typical distance between melted points. We concluded that the differential absorption of water and ice was the driving force for the pattern formation. Heating ice by laser absorption might be useful in applications such as the cryopreservation of biological samples.


Spontaneous self-organization in a system can lead to pattern formation. Numerous examples are found in nature. One such example is the division of cells in a developing embryo, as described by Turing (13). Turing developed a theory, called the reaction-diffusion (RD) model, by which some systems can spontaneously develop a spatial pattern. He defined two chemicals, an activator and an inhibitor, that react with each other and diffuse differently through a cell tissue. Solving the partial differential equations of such a system, he found six possible solution classes. One of these classes involves the formation of stationary patterns. This model may be modified to address other pattern formation phenomena. Further examples are observed in vegetation patterns in dry environments (4, 5), the skin color of animals (6, 7), wrinkle formation in ferrofluid mixtures under a magnetic field (8), drying paint films (9), the formation of brine channels in sea ice (1012), and the dendritic growth of crystals (13).

Dendritic patterns of ice are the snowflake shapes that result from the dendritic growth of ice crystals in supercooled water or supersaturated vapor. Snowflakes have been widely investigated both experimentally (14, 15) and mathematically (1619).

Another phenomenon that produces dendritic patterns is the internal melting of superheated bulk ice under infrared irradiation, as described by Tyndall (2024). Internal melting creates water cavities with a form that typically resembles dendritic ice crystal growth in supercooled water. At low growth rates, internal melting occurs as a circular disk shape. These water structures are often referred to as “inverse snowflakes,” “Tyndall flowers,” or “Tyndall figures.”

Heating ice with near-infrared (NIR) irradiation has been described in several studies. For example, the laser illumination of ice has been used to create microchannels (25). In that work, a collimated laser beam of 1064 nm was used to create channels with controlled sizes in an ice cube. The channel diameters were 40 to 350 μm, depending on the laser power. In another study of antifreeze proteins, a 980-nm laser was used to melt ice crystals within microfluidic channels (26).

Here, we present a previously unreported ice pattern that develops under uniform NIR irradiation. This pattern is generated when thin ice crystals in solution near equilibrium are exposed to phase-selective irradiation, leading to the spontaneous development of a dynamic labyrinth of holes and microchannels (Fig. 1). Irradiation in the NIR spectrum was selected to exploit the differential absorption spectra of ice and water (27, 28). Illumination at the appropriate wavelength permitted the selective heating of one phase over the other in a two-phase system (ice and solution) (29).

Fig. 1 Labyrinth ice pattern formation.

Labyrinth ice patterns developed under 1540-nm irradiation (14 W/cm2) over a time period of 60 min in a 13-μm-thick 10% sucrose solution layer. T = −2°C.


Labyrinth ice pattern evolved under ice-selective NIR irradiation

We irradiated a thin layer of ice crystals in a sucrose solution using 1540-nm laser illumination at an intensity of 14 W/cm2. At this wavelength, the absorption coefficient of ice is three times than that of water (Fig. 2A and fig. S1A) (27, 28). We prepared a two-phase system by cooling a 13-μm-thick sample of a 10% sucrose solution to −25°C using a custom-made temperature-controlled stage (Fig. 2B). The presence of sucrose slowed the pattern formation dynamics, as discussed below, although the sucrose did not change the water absorption spectrum compared to that of pure water at 1540 nm (fig. S1B). Ice growth and melting changed the solution concentration, which influenced the melting temperature. This process moved the partially frozen system toward global equilibrium. Warming the frozen sample to a temperature of −2°C resulted in the formation of isolated 10- to 100-μm-diameter ice crystals. Irradiating the sample while maintaining the stage at this temperature resulted in the spontaneous development of a spatial pattern (Fig. 1). An example of the pattern evolution is shown in movie S1. Circular holes with a radius of 3.3 ± 0.2 μm formed in the ice crystals. The ice crystals expanded at the expense of the hole areas. The holes dynamically emerged and disappeared periodically as the system was subjected to irradiation (Fig. 3, A and B, and movie S2). We described the hole formation process as follows: The ice became thinner, while a circular boundary emerged (Fig. 3A, a to c), and once the thin layer breached (Fig. 3A, d), the ice melted rapidly to form a hole (Fig. 3A, e). Thereafter, the radius decreased at a rate that increased from 0.021 ± 0.003 μm/s to 0.040 ± 0.004 μm/s (Fig. 3, A, e to a, and C). The holes either remained closed or repeated the opening and closing cycle (Fig. 3A, a to g). A hole that formed (Fig. 3A, e) in a channel of thinner ice could merge with adjacent holes and create an open solution channel (Fig. 3A, h to j). The typical distance between the centers of adjacent holes was 22.4 ± 0.3 μm (fig. S2A and movie S3). This distance indicated repulsion between the holes, probably due to reduced photon absorption at the hole site, which reduced the temperature around the hole.

Fig. 2 Microscopy system.

(A) The absorption coefficients of ice and water in the NIR (27, 28). The blue solid line indicates the ice absorption spectrum. The red dashed line indicates the water absorption spectrum. The gray vertical lines indicate the water-selective wavelength (980 nm) and the ice-selective wavelength (1540 nm). (B) The experimental setup. Additional technical information is provided in the Materials and Methods section. vis, visible; CCD, charge-coupled device.

Fig. 3 Analysis of the periodic opening and closing of holes during pattern formation in the 13-μm-thick sample in a 10% sucrose solution under 1540-nm irradiation (14 W/cm2) at a temperature of −2°C.

(A) Single hole dynamics. Left images show the experimental results. Right images show schematic diagrams of the system states. Additional visualizations are provided in movie S2. (B) The radius of the hole described in (A) versus the time. This hole opened four times and then merged with an adjacent hole. The blue rectangle in the top graph indicates a magnified view. Each opening had an initial radius of 3.3 ± 0.2 μm, and the boundaries gradually became clearer. The radius decreased slightly, and then, a hole opening event was observed as a morphology jump. The hole became slightly larger and then gradually decreased. The letters a to j indicate the images at (A). The asterisk in the magnified view indicates the starting point for the calculation of the rate at which the radius decreased. (C) The hole radius decreased at a rate that increased as the radius became smaller. Data were taken from the results obtained from 13-μm-thick samples in a 10% sucrose experiments under 1540-nm irradiation (14 W/cm2). Seventeen holes were analyzed. Error bars represent the SE.

In addition to observing the hole formation, we observed channel formation. We noticed that holes formed preferentially at the channels. Adjacent holes merged and created open channels of solution in the ice (fig. S2B and movie S4). The 1540-nm irradiation induced pattern formation under other conditions, such as in a 5% sodium chloride solution (fig. S3A) or in a 20% sucrose solution at −3° and −5°C (fig. S4) and in pure water (fig. S5 and movie S5). In a 10% sucrose solution at low temperatures (fig. S6), the pattern evolution decreased gradually as the temperature decreased. We tested ice-selective laser illumination at a different wavelength, 1060 nm (128 W/cm2), using a 40-μm-thick sample comprising a 10% sucrose solution at −1.8°C. The labyrinth ice pattern formation was evident under these conditions (fig. S3B).

Higher-temperature regions in the ice crystals corresponded to the topological skeleton

We observed crystals bearing a single hole. These crystals were typically 30 to 50 μm in diameter. This diameter allowed only for a single hole, taking into account the measured distance between holes in larger crystals. Under irradiation, the bulk ice acted as a heat source compared to water; therefore, the higher-temperature regions in the ice appeared to be farthest from the ice-water interface because of the greater distance of heat diffusion. The higher-temperature regions in the ice corresponded, therefore, to the topological skeleton of the ice. The skeleton of a shape is a thinner representation of the shape that follows a specific central path equidistance from the shape boundaries. The skeleton of a shape preserves the complete topological information of the original shape (30, 31). The heat transfer simulations revealed superheating at the center of a circular thin ice crystal (Fig. 4A). Experimentally, we observed the formation of a hole at the center of the ice. As a result, the crystal assumed a toroidal shape (Fig. 4A and movie S6). Hole formation [rupture of thin ice (Fig. 3A, d)] caused a topological change, and a new skeleton emerged, indicating the formation of a new region with a higher temperature (Fig. 4B). Microchannels of thinner ice formed, following the skeleton pattern (Fig. 4B and movie S7). New holes tended to form within these channels. Adjacent holes merged and transformed thinner ice channels into solution channels.

Fig. 4 Heat transfer simulations and skeleton characterization of a single ice crystal under 1540-nm irradiation (14 W/cm2).

(A) Left side of the arrow: The initial ice crystal (top left). The white line indicates the skeleton, i.e., the region at a higher temperature (top right). The temperature distribution in the ice crystal indicated that the temperature at the center was four thousandths of a degree higher than the temperature at the interface (bottom left). Side view of the temperature distribution (bottom right). Right side of the arrow: The ice crystal with a hole formed under irradiation. (B) Left side of the arrow: An ice crystal with a hole formed under irradiation (top left). The white line indicates the skeleton, i.e., the region at a higher temperature (top right). The temperature distribution in the ice crystal indicates that the temperature in the water hole was lower than the temperature around the hole (bottom left). Side view of the temperature distribution (bottom right). Right side of the arrow: The ice melted around the water hole.

The accelerated radius decrease was consistent with the Gibbs-Thomson effect

A single hole in an ice crystal could be treated mathematically as the opposite of ice crystal growth and melting. The melting temperature of a small ice particle decreased as the ratio of the surface change (dA) to the volume change (dV) increased, as described by the Gibbs-Thomson effect (3236). The melting point depression due to the ice particle surface curvature is given byΔT=Tm()Tm(r)=ΩγTm()ΔHdAdV(1)where Tm(r) is the equilibrium melting temperature of a particle with a finite radius, Tm(∞) is the freezing temperature of bulk ice, Ω is the volume of one mole of water molecules, γ is the surface tension, and ΔH is the molar latent heat of fusion. For a spherical ice crystal, dAdV=2r; thus, ΔT=2ΩγTm()rΔH. For a water disk in ice, as observed in our results, the curvature of the ice was negative, dAdV=1rw. Thus, the melting point wasTm(rw)=Tm()(1+ΩγrwΔH)(2)where rw is the radius of the hole. Equation 2 suggested that Tm(rw) > Tm(∞), i.e., the melting temperature of the hole interface was higher than the melting temperature of the flat ice interface. Therefore, the observed acceleration of the decrease in radius (Fig. 3C), as the radius became smaller, is consistent with the Gibbs-Thomson effect.

In the absence of irradiation, only Ostwald ripening was observed. Larger ice crystals grew at the expense of the smaller crystals, and no spatial pattern developed (fig. S7 and movie S8). Ostwald ripening is driven by decreasing the surface energy and, hence, minimizing the surface area, unlike the tendency of ice to increase the surface area under irradiation for efficient heat release. We suggest that these two opposing tendencies drove the periodic opening and closing of the holes during the pattern formation (Fig. 3, A and B).

As mentioned, ice crystals in pure water (instead of a 10% sucrose solution) under NIR irradiation also developed a spatial pattern (fig. S5 and movie S5). The average radius of the water holes was 3.1 ± 0.2 μm, similar to the hole sizes observed in 10% sucrose solution, but, the rate of decrease in the radius was 1.7 ± 0.1 μm/s, two orders of magnitude faster than the rate observed in the 10% sucrose experiments. We noted that holes were moving slowly within the ice (movie S5). We assumed that the movement was toward the hotter area in the sample, and this is a demonstration of thermodynamic buoyancy described in (37). The patterns in pure water were less stable than those in the 10% sucrose solution.

Suppression of the pattern evolution by a low concentration of black ink

We hypothesized that the pattern developed because of the differential absorption of the ice and solution, which overheated the ice. Thus, we assumed that the pattern would not develop if the absorption of the ice and the liquid were identical. To eliminate the difference in the absorption, i.e., to increase the absorption of the liquid phase, we added black ink particles to the sucrose solution. The ink is a suspension of carbon black particles that absorb visible and infrared (IR) wavelengths (38). We measured the 1540-nm absorption of a series of ink suspensions with a range of ink concentrations (from 0 to 10%) (Fig. 5A) and found that the 1% (v/v) ink solution provided an absorption cross section at 1540 nm equal to that of the ice. In a partially frozen system, the initial solution is concentrated relative to the nonfrozen solution. For this reason, we examined a lower ink suspension concentration, 0.4% ink in 10% sucrose. We found that this ink concentration suppressed pattern development (Fig. 5B). Only a few holes formed in the ice crystals, unlike the labyrinth structures that developed in the absence of the ink.

Fig. 5 The effect of black ink on the water absorption and the pattern formation.

(A) Water absorption at 1540 nm increased with the ink concentration. The black line indicates the fit to the equation A=1eαcz,as αc=αw+CαcC, where αw is the water absorption coefficient and αc is the water absorption coefficient as a function of the ink concentration, C. The width of the absorption cell was measured using a caliper micrometer to be z = 0.017 ± 0.001 cm; therefore, αcC=22±1.4cm1·%ink1. The ink concentration needed to achieve an absorption cross section equal to the ice absorption was calculated using αc = αi. Thus, αw+CαcC=αi. For αw = 12 cm−1 (27) and αi = 34.7 cm−1 (28), C was approximately 1%. The red arrow indicates the calculated absorption value of 0.017 cm thick ice sample. All measured values were calibrated against an empty sample. The error bars represent the SD. (B) Suppression of the ice pattern formation in a 13-μm-thick sample of a solution containing 10% sucrose and 0.4% black ink. T = −2°C. Laser intensity: 14 W/cm2 (1540 nm).

Experiments conducted using a much higher ink concentration, 10% ink in 10% sucrose, in which the solution absorption coefficient was sevenfold higher than that of the ice, revealed the pattern formation (fig. S8 and movie S9). The radius of the holes was smaller, and the rate at which the radius decreased was one order of magnitude larger than the corresponding values observed in experiments conducted without ink.

Temporal oscillations of ice and water under ice-selective NIR irradiation

In addition to the spatial patterning, we also observed temporal oscillations in the experiments conducted using pure water under the 1540-nm irradiation (fig. S9). Experiments performed using ice formed from pure water resulted in ice that partially melted under illumination. Air bubbles in the ice migrated toward this melted area and combined to form a large bubble. Water drops within this bubble condensed on the cold glass surface. These drops were supercooled. At low temperatures, the water drops froze as ice around the air bubble contacted the supercooled drops. At this point, we observed that the ice crystals inside the air bubble oscillated rapidly (at a frequency of a few hertz) between frozen and partially melted states during irradiation (fig. S9 and movie S10). The oscillations stopped when the laser illumination was turned off. We interpreted these oscillations according to the following mechanism: The ice was heated by irradiation to a temperature above melting, became superheated, and started to melt. As the ice melted to form water, the latent heat cooled the system, and the photon absorption decreased. The cooling then overshot the melting point, and the water became supercooled. At that point, the ice began to grow, and the temperature increased because of the latent heat release and the increase in absorption, until the ice became superheated again and the cycle started over. This experiment did not produce oscillations when conducted under water-selective (980 nm) illumination.

Patterns observed under water-selective conditions

We applied water-selective heating using two methods. First, we applied water-selective irradiation using laser illumination at 980 nm, in which the water absorption coefficient was higher than the ice absorption coefficient (Figs. 2A and 6A) (27, 28). Because of the low absorption at that wavelength compared to the absorption at 1540 nm, we used an 80-μm-thick layer of a 10% sucrose solution. We observed a pattern of ordered straight plates, separated by a constant distance of 15.1 ± 2.4 μm, using this heating method (Fig. 6B and movie S11).

Fig. 6 Pattern formation under water-selective conditions.

(A) Absorption of ice and water at 980 nm versus sample thickness, calculated from the absorption coefficients: α (water, 980 nm) = 0.482 cm−1 (27) and α (ice, 980 nm) = 0.139 cm−1 (28). The absorption equation A = 1 − e−α∙z. The blue solid line refers to the ice. The red dashed line refers to the water. The black dotted line refers to the ice-to-water absorption ratio. The gray vertical line in the magnified plot indicates a sample thickness of 80 μm used in the experiments. (B) Ice pattern obtained under 980-nm irradiation (97 W/cm2) in an 80-μm-thick layer of a 10% sucrose solution at a temperature of −1°C. Additional visualization is provided in movie S11. (C) Ice pattern formed in an 80-μm-thick layer of a 10% sucrose solution containing 10% black ink at a temperature of −3°C under 100-W halogen illumination. Additional visualization is provided in movie S12. (D) Ice pattern of an 80-μm-thick layer of a 10% sucrose solution containing 0.4% black ink at a temperature of −1°C under 100-W halogen illumination. (E) Ice crystals in an 80-μm-thick layer of a 10% sucrose solution at a temperature of −1°C under 100-W halogen illumination. No patterns were observed under these conditions.

A second method of applying water-selective heating involved adding black ink particles to the 10% sucrose solution and subjecting this solution to halogen light illumination (100-W microscope halogen lamp). The halogen illumination applied a broad spectrum with a peak at 1000 nm (39). The black ink particles present in the liquid phase increased the liquid photon absorption compared to that of the ice absorption. To further increase the photon absorption, we used an 80-μm-thick layer. At 10% ink, ice crystals formed a pattern of ordered straight plates separated by a constant distance of 7.8 ± 1.1 μm (Fig. 6C and movie S12). At 0.4% ink, we observed patterns with a typical width of 11.9 ± 3.2 μm (Fig. 6D). The 10% sucrose solution prepared without ink formed no such pattern (Fig. 6E). We noted that the water-selective patterns observed under 980-nm illumination were similar to those observed in 10% ink solution under halogen illumination.


In this study, we described a previously unreported ice and water pattern, achieved by subjecting a two-phase system to an illumination heating source. We demonstrated a process in which one phase was selectively heated more than the other. This heating mode resulted in pattern formation, and we suggested that the driving force involved efficient heat release. As the ice was heated to a greater extent than the liquid phase, a pattern of dynamic holes and microchannel labyrinths was observed (Fig. 1). As the liquid phase was heated to a greater extent than the ice phase, a pattern formation was also observed, albeit with different features. In this case, we observed the assembly of ice stripes and plates (Fig. 6, B to D). When the liquid photon absorption was equal to the ice photon absorption, the spatial pattern development was suppressed (Fig. 5B).

We emphasize that the labyrinth ice pattern formation differed from the Tyndall figures, which were obtained once the bulk single ice crystal was subjected to IR irradiation. These figures resulted from the internal melting of the superheated ice. By contrast, the labyrinth ice was generated when irradiation was applied to a system coupled to a cooling source.

Several mathematical approaches were developed in an effort to understand the different pattern formation processes. For example, the RD model, developed by Turing (2), explains how two chemical substances react and diffuse to create patterns. Turing’s description included differential equations that described the diffusion, creation, and elimination of two materials that influenced one another. These equations yielded a variety of solutions, including spatial pattern formation and temporal oscillations (3).

To describe our system, a set of differential equations could be defined. The variables should be the thermal energy, the phase of the water, and the solute concentration. The liquid and the solid phases absorb radiation at different levels. These two phases have different thermal conductivities, and the local temperature influences the phase transition. Therefore, the variables are nonlinearly interconnected. For example, the heat diffusion coefficient, D=kρCp, (where k is the thermal conductivity, ρ is the density, and Cp is the heat capacity), of ice is Dice=1.17mm2s, and that of water is 8.4 times lower, Dwater=0.14mm2s.

Pattern evolution can be studied using a phase field model as a tool to describe the phase transitions. In the phase field model, the phase change is given by the derivative of the free energy with respect to the phase change (40). We modeled certain features of ice crystal growth in our system using two-dimensional (2D) axisymmetric simulations (fig. S10A). The phase field and the heat equations are described in detail in Materials and Methods (Eqs. 5 and 10). Heat transfer simulations showed that ice temperature increased with ice crystal size as the crystal absorbed more energy (fig. S10B). Phase field simulations revealed that the irradiation caused the ice crystals to melt at their edges and then to grow toward the cold sapphire at the bottom (fig. S10C). In addition, the phase field revealed a repulsive interaction between the two ice elements (fig. S10D). This demonstrated a steady-state solution in which there is a typical distance between the ice elements. Under experimental conditions, in which ice was formed in an air bubble and thus relatively thermally isolated, temporal oscillations were observed (fig. S9). These oscillations were described as part of the Turing’s solutions.

Although the pattern formation was observed in pure water (fig. S5), it was difficult to maintain the system close to its melting point. The addition of sucrose brought the system to its melting point at any given temperature because of a combination of effects. (i) The melting point decreased with increasing sugar concentration (Blagden’s law). (ii) Solutes were poorly incorporated into the ice fraction. (iii) The total amount of sugar in the solution remained constant as long as no sugar crystallization occurred. As a result of (ii) and (iii), the solute concentration in the liquid increased when the ice fraction increased, and because of (i), the melting point reduced. Therefore, the ice will grow until the solution reaches the concentration in which its melting point matches the set temperature. Thus, the addition of sucrose allowed us to work close to equilibrium at any given temperature and to slow down the dynamics and stabilize the pattern. Several recent studies demonstrated simulations of ice development in a solute using the phase field approach, in which global constraints fixed the amount of ice in the system (40, 41). Insights from these simulations should be incorporated into the full solution of the observed pattern formation.

Several studies have shown that saccharides can affect the ice morphology (42, 43). In our system, the development of the pattern’s morphology was probably not affected by the sugar. The diffusion coefficient of sucrose in an aqueous solution at −2°C is D=3×102μm2s, extrapolated from (44, 45). Thus, a typical diffusion time, t=x22D, over x = 25 μm is approximately 1 s. Features much smaller than 25 μm undergo dynamics with a typical time scale much longer than 1 s; therefore, we expected that the local sucrose gradient would not affect directly the dynamics of hole opening and closing. Another point to consider is the possible change in irradiation absorption as a function of sucrose concentration. Our measurements of NIR absorption of water/sucrose solutions at 1540 nm showed no dependence in the sucrose concentration (fig. S1B). This finding is supported by the work of Jung and Hwang (46). Thus, in our system, at 1540 nm, the photon absorption does not change because of the local sucrose gradient even if these gradients exist.

Other aspects of pattern evolution might be considered, for example, the interaction between light and the ice and water system. The refractive index of ice is 1.31, and that of sucrose solution is 1.33 to 1.36. This difference might account for local light concentration under certain conditions. Another factor potentially driving patterning is the inhomogeneity of the illumination beam itself. The addition of 0.4% black ink, which eliminated absorption differences between the ice and the water but did not change the refractive index difference, substantially suppressed the pattern (Fig. 5B). Thus, we concluded that the refractive index differences and any potential inhomogeneities in the illumination did not drive the pattern formation.

The significance of our findings, in addition to the discovery of the labyrinth ice pattern formation, is relevant to a variety of fields. In food engineering, ingredients, such as solutes and colloids, can be embedded into the ice at the micrometer level. In the frozen food industry, thawing process is needed for food products that preserved at the frozen state but consumed at the nonfrozen state. Ice-selective irradiation can be incorporated into effective thawing methods that avoid overheating and uncontrolled ice growth, i.e., recrystallization, during warming. Microwave ovens are often used to defrost food, but the microwave heats the ice much slower than the water because of the fixed positions of the water molecules in the ice (47). The absorption coefficient of water is 1000 times higher than that of ice in the microwave oven radiation region at 12.23 cm (47). This differential absorption cross section causes the liquid parts to be overheated while the frozen parts continue to defrost. A partial solution to this problem is provided by the defrost action in the microwave oven, which cycles the microwave power on and off to permit heat diffusion. Reversing the absorption ratio to heat the ice more than the liquid water would be an advantage for achieving efficient thawing. In cryopreservation, a critical problem is ice recrystallization during warming from low storage temperatures. Cell membranes can be mechanically injured by the growth of ice crystals (48), or cells can be dehydrated because of osmotic stress that results from extracellular ice (49). Ice-selective NIR irradiation offers a controlled and improved thawing process without overwarming the samples and undesired ice crystals growth.

To conclude, we described a phenomenon in which ice can form a spatial pattern when one of the two phases (ice or liquid water) is selectively heated more than the other. The heating energy can be applied as a laser illumination or an incandescent light bulb such as halogen illumination. Further investigations of this phenomenon are needed to more fully explain the observed data.


Solutions preparation for the absorption measurements

Sucrose solutions. Sucrose solutions were prepared in different concentrations: 0 (pure water), 10, 20, and 30% (w/v). Twelve grams of sucrose (lot no. BCBL8329V, Sigma Life Science) was dissolved in 30 ml of double-distilled water (DDW). DDW was added to reach a final volume of 40 ml for the 30% sucrose solution, which was then filtered. The 30% sucrose solution was diluted 1.5-fold to a final concentration of 20% sucrose and 3-fold to a final concentration of 10% sucrose.

Black ink suspensions. Carbon black suspensions were prepared in several concentrations: 0 (pure water), 0.1, 0.4, 1, 2.5, 5, and 10% (v/v) by diluting India ink (29.6 ml dispenser; no. 44201, Higgins Pigment-Based Drawing Ink, Waterproof, Black India) in DDW. The 100% ink solution corresponded to the original undiluted suspension.

Solutions preparation for the pattern formation experiments

5% (w/v) NaCl solution. Two grams of sodium chloride (lot no. 1108881, Bio-Lab Chemicals) was dissolved in 30 ml of DDW. DDW was added to reach a final volume of 40 ml, which was then filtered.

10% (w/v) sucrose solution. Eight grams of sucrose (lot no. BCBL8329V, Sigma Life Science) was dissolved in 30 ml of DDW. DDW was added to reach a final volume of 40 ml. The solution was then filtered and diluted twofold.

Black ink suspensions. Carbon black suspensions in two different concentrations, 0.4 and 10% (v/v), were prepared: One microliter of the original ink suspension was mixed with 5 μl of a 20% sucrose solution and 4 μl of DDW to final concentrations of 10% India ink and 10% sucrose. One microliter of the original ink suspension was mixed with 125 μl of 20% sucrose and 124 μl of DDW to final concentrations of 0.4% India ink and 10% sucrose. The 100% ink solution corresponded to the original undiluted suspension. The particle size ranged from 0.1 to 1 μm (50).

Preparation of ice crystals in solution

A 1.5-μl drop of the tested solution was placed on a 25-mm-diameter sapphire glass and covered with a 12-mm-diameter circular glass coverslip, generating a sample thickness of 13 μm. A 2.3-μl drop was used to achieve a sample thickness of 20 μm, a 4.5-μl drop was used to achieve a sample thickness of 40 μm, and a 9 μl-drop was used to achieve a sample thickness of 80 μm. To avoid water evaporation, the gap between the two glasses was sealed with immersion oil (immersion oil type B, no. 16908-16, Electron Microscopy Sciences). This sample was attached using immersion oil to a copper plate in a custom-built temperature-controlled stage (Fig. 2B). The stage contained thermoelectric units that cooled the copper plate and a thermistor that measured the temperature near the sample (not the sample itself). The stage was connected to a temperature controller (model 3040, Newport Corporation) to permit delicate control over the temperature around the sample holder between room temperature and −30°C, with 0.01°C precision. An external water circuit was connected to the stage to remove heat from the cooling units, and dry air flow was applied to avoid humidity around the device. A 3-mm-diameter hole in the middle of the copper plate served as the viewing window. The temperature was controlled using a custom-made LabVIEW interface (National Instruments). The sample was frozen by cooling to −25°C and then warmed near the melting temperature of the solution. During the warming process, separated ice crystals formed in the solution.

The optical setup

A temperature-controlled stage was integrated with an inverted microscope (TE2000-U, Nikon Eclipse), as shown in Fig. 2B. Laser illumination was directed toward the sample using an optical fiber, a set of mirrors, and a filter. Several laser wavelengths were used: single-mode 10-μm fiber, tunable from 1536 to 1569 nm (2 W, 1540 nm; model KPS-BT2-TFL-1550-20-FA, Keopsys); single-mode 10-μm fiber (5 W, 1060 nm with a spectral range of ±0.5 nm; model KPS-BT2-YFL-1060-050-LC-COL, Keopsys); and multimode 100-μm fiber (5 W, 980 nm with a spectral range of ±3 nm; model SDL-822, Spectra Diode Labs).

Laser intensity calculations. For a Gaussian beam with an optical power P and Gaussian beam radius ω, the peak intensity (on the beam axis) is given by IP=P12πω2 (51). In our system, ω was measured to be 1.5, 0.865, and 1.145 mm for the 1540-, 1060-, and 980-nm laser beams.

The 1540- and 1060-nm lasers were ice selective, meaning that the absorption coefficients of the ice were higher than those of the water: α (water, 1540 nm) = 12 cm−1 (27), α (ice, 1540 nm) = 34.7 cm−1 (28), α (water, 1060 nm) = 0.149 cm−1 (27), and α (ice, 1060 nm) = 0.243 cm−1 (28). The 980-nm laser was water selective: α (water, 980 nm) = 0.482 cm−1 (27) and α (ice, 980 nm) = 0.139 cm−1 (28).

Experiments conducted using halogen illumination used the halogen lamp of the microscope (100-W lamp; model LHS-H100P-1, Osram halogen 12 V). Experimental data were collected using a charge-coupled device camera (DMK 23G274 GigE camera, Imaging Source).

Absorption measurements

Absorption measurements (shown in Fig. 5A and fig. S1B) were conducted using a sample with a thickness of 0.17 ± 0.01 mm, as it provided a measurable absorption value. The absorption ratio of ice to water was 2.4 at that thickness (Fig. 2A). A polydimethylsiloxane (PDMS) spacer with a thickness of 0.17 ± 0.01 mm and a 6-mm-diameter hole was placed on a 25-mm-diameter sapphire glass. A 4.8-μl drop of the tested sample was placed in the PDMS spacer to fill the hole and was covered with a 12-mm-diameter circular glass coverslip. This sample was placed on the temperature-controlled stage in the laser path (Fig. 2B). The incoming and transmitted 1540-nm laser powers were measured using a laser power meter (model FieldMax-TOP, no. 0440L04, COHERENT), which was placed before the sample for the incoming radiation measurement and after the sample for the transmitted radiation measurement. All measurements were calibrated relative to an empty sample holder (PDMS spacer between the sapphire and the glass coverslip). The absorbances of each sample were measured three times at three different laser powers: 100, 500, and 1000 mW. The absorbance was calculated byA=TreferenceTsampleTreference(3)where Treference is the transmitted radiation of an empty sample holder and Tsample is the transmitted radiation of the tested solution. The errors were calculated as the SD.

Image analysis and calculations

Image analysis and calculations were performed using the ImageJ 1.49p software.

Calculation of the radius and rate of water hole radius decrease. Holes were selected manually, and a specific sequence of images that showed the hole formation and closing was cropped. Image sequences were converted to 8 bit, and a dark background threshold was selected to convert the image sequences to a binary sequence. Particle analysis was used to measure the perimeter and the circularity values for each image. We calculated the rate of decrease in the radius from a point that followed the brief expansion and shrinking of the hole. An example for the chosen starting point is marked with an asterisk in Fig. 3B. The rate of decrease in the radius was calculated using the “SLOPE” command in Excel, which uses the equation b=(xx¯)(yy¯)(xx¯)2. The slope was calculated for each radius and included one value before and one value after, i.e., the slope of three points was calculated for each radius. The 10% sucrose experiments were calculated using the average rate for each radius interval (Fig. 3C). Seventeen holes were measured. The 10% ink experiments were conducted by measuring 12 holes, and pure water experiments were conducted by measuring 4 holes. The errors were calculated as the SEs.

Skeleton analysis. Images were converted to 8-bit images, and a threshold was selected to convert the image to a binary image. The relevant crystal in the image was isolated by coloring the background in black. A skeletonized image was obtained using the “skeletonize” command in ImageJ. Image calculations were conducted by adding the skeleton image to the original image.

Distance between the centers of adjacent holes. Ice crystals with holes were selected from three different experiments. Holes that completely or partially opened were chosen manually and labeled. The x and y values of the hole positions were measured. The distance was calculated by dij=(xixj)2+(yiyj)2 and averaged on ij nearest neighbors.

Computational simulations

Heat transfer and phase field simulations were developed using the COMSOL Multiphysics 5.3 software. 3D heat transfer simulations of a single ice crystal were carried out using the “heat transfer in solids” module, as shown in Fig. 4. The parameters were defined as follows: thermal conductivities ofki=2.22WmK (ice), kw=0.58WmK (water), ks=25WmK (sapphire), and kg=1.38WmK (glass); heat capacities of Cpi=2050JkgK (ice), Cpw=4228JkgK (water), Cps=760JkgK (sapphire), and Cpg=703JkgK (glass); and densities of ρi=918kgm3 (ice), ρw=997kgm3 (water), ρs=3980kgm3 (sapphire), and ρg=2203kgm3 (glass). The geometry and materials were as follows: Above a cylindrical sapphire glass with a 3-mm radius and a 0.5-mm height was placed the water with 1.5-mm radius and 0.01-mm height. An ice crystal 0.05 mm in radius and 0.01 mm in height was placed at the center of the water. A cylindrical glass 3 mm in radius and 0.17 mm in height was placed above the water. An ice crystal with a hole was simulated by adding a water hole with a 0.007-mm radius at the center of the ice crystal. Ice was defined as a heat source according to Qi= Ice absorptionLaser powerSample volume(4)where Qi is the power per volume absorbed by the ice. The temperature of the simulation was 0.004°C below the melting point. The mesh was defined to have an element size “normal” and a resolution of 6 × 10−4 for the maximum element size, 1.08 × 10−4 for the minimum element size, 1.5 for the maximum element growth rate, 0.6 for the curvature factor, and 0.5 for the narrow regions.

2D axisymmetric heat transfer simulations of multiple ice crystals were carried out using the following setup: Above a 3-mm radius, 0.5-mm-thick sapphire glass was placed a 0.001-mm-thick ice and water layer, and above that was placed a 0.17-mm silica glass slide. Four equal-sized ice crystals were designed in the 0.001-mm-thick water sample (figs. S10, A and B). The cold stage was located at the bottom and was set to −2°C. The melting temperature was set to Tm = 0°C. The material properties were defined as described above for the 3D simulations. The mesh was chosen to have an element size of 0.15 μm. The other domains were meshed using triangles without a fixed size. Therefore, the mesh was tighter at the ice-water interface.

The heat equation was defined asρ(φ)Cp(φ)Tt=(k(φ)T)+Q(φ)(5)where φ is the phase (φ = 0 for water and φ = 1 for ice) and ρ is the density defined asρ(φ)=φρice+(1φ)ρwater(6)Cp is the specific heat defined asCp(φ)=φCp,ice+(1φ)Cp,water(7)k is the thermal conductivity defined ask(φ)=φkice+(1φ)kwater(8)Q is the laser heating defined as the volume heated according to Eq. 4Q(φ)=φQice+(1φ)Qwater(9)where Qice and Qwater are defined as described in Eq. 4.

2D axisymmetric phase field simulations (52) were carried to follow the changes in the shapes of the ice elements (figs. S10, C and D). The setup and the mesh were designed as described above for the 2D axisymmetric heat transfer simulation. The phase field equation was defined using the Allen-Cahn equation (53)φt=Mφ[fφεφ22φ]φ[0,1](10)where Mφ is a factor governing the kinetics of the interface. The mobility could be found using the following expression (54)Mφ=μTm6δL(11)where μ is the interface kinetic coefficient that relates the velocity of the surface and the level of supercooling, δ is a measure of the interface thickness, and L is the latent heat. The parameter Mφ determined how quickly the ice front propagated for a given deviation from equilibrium.

In another 2D axisymmetric phase field simulation, we followed the interactions between two ice elements (fig. S10D). The starting point was the endpoint of the previous simulation (fig. S10C), except that two ice elements were brought closer together, with a gap of 2 μm.


Supplementary material for this article is available at

Fig. S1. Water absorption measurements at 1540 nm.

Fig. S2. Features observed during the ice pattern development in a 13-μm-thick 10% sucrose solution under 1540-nm irradiation (14 W/cm2) at T = −2°C.

Fig. S3. Ice pattern formation observed under other conditions.

Fig. S4. Ice pattern formation in 20% sucrose solution.

Fig. S5. Ice pattern formed in a pure water sample.

Fig. S6. Ice pattern formation at lower temperatures.

Fig. S7. A demonstration of ice crystals growth in a sucrose solution not subjected to NIR irradiation.

Fig. S8. Ice pattern formation at a high ink concentration.

Fig. S9. Ice oscillations in an 80-μm-thick pure water sample under 1540-nm irradiation (20 W/cm2) at T = −1.9°C.

Fig. S10. 2D axisymmetric simulations of multiple ice crystals, with an axis of symmetry at r = 0 and a cold stage on the bottom right side.

Movie S1. Labyrinth ice pattern formation.

Movie S2. Single hole periodic dynamics.

Movie S3. An ice crystal with multiple holes.

Movie S4. Microchannel formation in an ice crystal.

Movie S5. Ice pattern in pure water.

Movie S6. Skeleton characterization of a single ice crystal (part A).

Movie S7. Skeleton characterization of a single ice crystal (part B).

Movie S8. A demonstration of Ostwald ripening.

Movie S9. Ice pattern formation at a high ink concentration.

Movie S10. Ice temporal oscillations.

Movie S11. Pattern formation under water-selective irradiation.

Movie S12. Pattern formation under halogen illumination.

This is an open-access article distributed under the terms of the Creative Commons Attribution-NonCommercial license, which permits use, distribution, and reproduction in any medium, so long as the resultant use is not for commercial advantage and provided the original work is properly cited.


Acknowledgments: Funding: This work was supported by the Israel Science Foundation (grant no. 930/16). Author contributions: S.G.P. performed the experiments. A.Kan., A.Kat., and S.G.P. performed the computational simulations. S.G.P., I.B., H.C., and S.U. designed the experiments. S.G.P., V.Y., I.B., H.C., and S.U. analyzed the data. I.B. and V.Y. developed the LabVIEW software. S.G.P. and I.B. wrote the paper. Competing interests: The authors declare that they have no competing interests. Data and materials availability: All data needed to evaluate the conclusions in the paper are present in the paper and/or the Supplementary Materials. Additional data related to this paper may be requested from the authors.

Stay Connected to Science Advances

Navigate This Article