Dissociation of methane from a layer of methane-hydrate particles: A new simple model

A simple model of the dissociation of methane from a layer of methane-hydrate particles is suggested. The model is based on the assumption that methane-hydrate particles form a flat porous film lying on an adiabatic wall. This film is heated by ambient air convection. It is assumed that when the temperature inside a certain region within the methane-hydrate reaches the critical temperature for the release of methane from the methane-hydrate, all methane is released from the methane-hydrate and the latter turns into ice. The contribution of heat required for the release of methane is considered. The model is based on the analytical solution to the one-dimensional heat transfer equation in the two-layer (methane-hydrate/ice) system. This solution is incorporated into the numerical code and used at each time step of the calculations. Model predictions are compared with experimental data for the heating of a layer of methane-hydrate of porosity 0.3 and thickness 5 mm. It is shown that good agreement between the predictions of the model and experimental data is observed at the initial stage of the process. At longer times, the model predicts faster methane release than observed experimentally unless the effect of self-preservation can be ignored. This deviation between modelling results and experimental data is attributed to the main assumption of the model that all methane is instantaneously released from the methane-hydrate when the methane-hydrate temperature reaches the above-mentioned critical temperature. © 2023 The Author(s). Published by Elsevier Ltd. This is an open access article under the CC BY license ( http://creativecommons.org/licenses/by/4.0/ )


Introduction
Methane-hydrates are solid crystalline compounds. Water molecules form the outer framework of the crystalline structure that is made of hydrogen bonds and a cage-like crystal lattice. Guest gas molecules are bound to water molecules by van der Waals forces. Stabilisation of the crystal lattice is realised at a certain temperature and pressure (cf. equilibrium curve linking pressure and temperature) [1][2][3] . To date, there are about 60 known types of gas molecules that form gas-hydrates with various crystal structures. The most common structures are cubic (sI and sII) and hexagonal (sH) [1,2] . Methane-hydrates have the ability to contain a large number of gas molecules per unit volume of clathrate hy- The dissociation of methane-hydrates into ice and methane depends on a large number of key parameters including: pressure, temperature, degree of deviation of temperature (pressure) from the equilibrium curve, the size of methane-hydrate particles and specific surface area, the type of elementary crystalline cell, the peculiarities of the interaction of gas molecules with water molecules, the porosity of both the powder layer and of the ice shell that occurs during the dissociation of the gas hydrate, the methane-hydrate structural characteristics, defects [1,2,23] , the external heat flux [24] , and self-preservation of methane-hydrates during dissociation at low temperatures (below the melting point of ice) [1,2] . Thus, when burning a combustible gas over the surface of the methane-hydrate powder layer, a high heat flux is realised. The latter leads to an increase in the dissociation rate by more than ten times (compared to the rate of dissociation without combustion).
Dissociation of methane-hydrates has been described by a semi-empirical expression obtained on the basis of the Arrhenius equation, when the dissociation rate is proportional to the deviation of the pressure from the equilibrium value [25] . At temperatures above the critical temperature for methane-hydrate dissociation, the activation energy of methane hydrate is known to increase from 78.3 to 81 kJ/mol with an increase in temperature [25,26] . During the dissociation of methane hydrate at a temperature below the melting point of ice, water is not formed (methanehydrate dissociates into methane and ice) and the activation energy becomes significantly lower (33.5 kJ/mol) [27] .
When modelling the processes in methane-hydrates, dissociation is considered to be the process of movement of the dissociation front in the presence of a constant heat flow [28] . The kinetic mechanism of dissociation is discussed in [25] , where an internal kinetic constant and activation energy were obtained based on a summarisation of experimental data. The model described in that paper does not take into account the limitations associated with heat and mass transfer processes. The need to simultaneously account for the kinetics of dissociation and heat transfer is considered in [29] . The results of simulation of methane-hydrate dissociation in a porous medium under thermal stimulation are presented in [30] . The results of modelling methane-hydrate dissociation taking into account mass and momentum processes in different phases (gas, water and methane-hydrate), as well as the kinetics of methane-hydrate dissociation in a porous medium, are presented in [31] . The results of modelling methane-hydrate decomposition using CFD software tools are presented in [32] .
Unlike the case for temperatures above the ice melting point, modelling the dissociation of methane-hydrates at temperatures below this point is complicated by the phenomenon of 'selfpreservation'. In the temperature range 233 K -267 K, the rate of dissociation of methane-hydrate decreases by 2-4 orders of magnitude [33] . Abnormally low rates of dissociation of methanehydrates during annealing are referred to as the phenomenon of self-preservation. The influence of key factors on the rate of dissociation at these temperatures was investigated in [33][34][35] .
The motion of the dissociation front of a spherical particle for a specified diffusion coefficient and pressure drop was investigated by the authors of [ 36 ]. In the temperature range 230 K to 267 K, the dissociation rate decreases by several orders of magnitude. In this range, the diffusion coefficient strongly depends on temperature and cannot be considered to be constant 1 . This limits the applicability of the results presented in [ 36 ]. The results of the modelling of the phenomenon, taking into account the diffusion process, are presented in [ 37,38 ]. The dissociation of methane-hydrate depends on the porosity of the ice shell, which is controlled by kinetic and filtration coefficients [ 39,40 ].
As follows from this brief analysis of the modelling of the dissociation of methane-hydrates, this problem cannot be considered as fully solved. The model described in this paper is not intended to replace the previously developed models but rather to supplement them, focusing on the processes which have been largely overlooked so far. Our analysis shifts the focus from the analysis of spherical methane-hydrate particles to the analysis of a porous methane-hydrate layer formed of pressed methane-hydrate particles. This approach will allow us to directly compare the modelling results with experimental data. In contrast to most previously developed numerical models of the process [41][42][43][44], the new model is based on the incorporation of the newly obtained analytical solution to the heat transfer equation into the numerical code. This approach is similar to the one widely used for modelling of the processes in individual droplets [ 45 ]. The model described in this paper is designed not to replace the above-mentioned purely numerical models but to supplement them.
Basic equations and approximations of the model are described in Section 2 . An analysis of the predictions of the model using values of parameters that are typical of engineering applications is presented in Section 3 . Section 4 focuses on the validation of the model using the available experimental data. The main results of the paper are summarised in Section 5 .

Basic equations and approximations
It is assumed that methane-hydrate particles are pressed together to form a flat film lying on an adiabatic wall. This film is considered to be porous with porosity (volume fraction of pores) of less than or equal to 60%. The upper surface of the film is heated up by ambient air at temperature T g with convection heat transfer coefficient h assumed to be constant (see [ 45,46 ] for the discussion of this assumption). At the start of the process, the temperature of the methane-hydrate is below the critical temperature above which methane can be released from it ( T cr = 193 K), while the ambient temperature T g is above T cr . Self-preservation of methanehydrate and its effect on the temperature distribution was not considered. The heating of a methane-hydrate film is described by the following equation: with the initial condition: where t is time, z is the distance from the wall, κ h is the thermal diffusivity of the methane-hydrate, P (z, t ) describes the volumetric heat source or sink (e.g. due to radiative heating or chemical reactions). The model is developed in such a way that it can potentially use any functional form of P (z, t ) . At this stage of the model's development, however, these processes are not considered, and this function will be assumed to be equal to zero. The effect of the chemical reaction is considered in the conditions at the methanehydrate/ice interface (see Eq. (5) ). The boundary condition for Eq. (1) at the adiabatic wall is presented as: where shows the derivative with respect to z. The boundary condition for Eq. (1) at the methane-hydrate surface is presented as: where z h shows the surface of the methane-hydrate film, k h is the thermal conductivity of the methane-hydrate, h is the convective heat transfer coefficient. After a certain time step t 1 the temperature of the methanehydrate in a thin layer close to the surface of the film is expected to raise above T cr . The thickness of this layer z 1 ≥ 0 is estimated as the distance of the point where T (z h − z 1 ) = T cr from the surface of the film. It is assumed that methane is instantaneously released from this layer and that the methane-hydrate turns into ice within it. This release of methane is accompanied by the absorption of part of the heat supplied by the ambient gas. The heat flux required for this release of methane is estimated as: where ρ h is the density of the methane-hydrate, α m is the mass fraction of methane in the methane-hydrate in layer z 1 , Q m is the specific heat required for this release of methane (in J/kg). α m is constant during each time step but is allowed to change between time steps. The mass flux of methane from the layer surface is estimated as: It is assumed that the density of ice (with cracks forming after the release of methane) is the same as that of the methanehydrate. The dependence of both densities on temperature is ignored. The differences between the thermal conductivity and specific heat capacities of methane-hydrate and ice and their dependence on temperature are considered.
Once methane has been released from the methane-hydrate and a layer of ice has been formed above the layer of methanehydrate, the heat transfer process in a composite methanehydrate/ice film is described by Eq. (1) with κ defined as: where κ h and κ i are the thermal diffusivities of methane-hydrate and ice, respectively; z h and z i are the heights of the upper surfaces of the methane-hydrate and ice layers, respectively. In this case Eq. (1) needs to be solved subject to the following initial condition: The boundary conditions for this equation at the methanehydrate/ice interface are presented as: where | d z h d t | is estimated at the previous time step. | d z h d t | remains constant during each time step. z i remains constant during the whole process as the densities of methane-hydrate and ice are assumed to be constant and equal.
The boundary condition at the surface of the ice is given by Eq. (4) with z h and k h replaced by z i and k i , respectively (see Fig. 1 ).
The analytical solution to Eq. (1) with the above-mentioned boundary and initial conditions is presented as (see Appendix A for the details): where v n and n are given by Expressions (A31) and (A50) , respectively; φ is given by Formula (A11) .
Expression (10) is used at each time step of the calculations, taking into account the changes of | d z h d t | with time.

Analysis
Expression (10) was incorporated into Matlab R2020a and used at each time step of the calculations. These involved 50 terms in the series for expression (10) , time steps of 0.01-1 s, and 20,0 0 0 cells across the methane-hydrate/ice layer. The roots of Equation (A27) for eigenvalues were found using the bisection method with absolute accuracy of 10 −15 (cf. Appendix B). The thermodynamic and transport properties of nitrogen (approximation of air), methane, methane-hydrate and ice, used in calculations, are summarised in Appendix E.
The methane-hydrate initial temperature and ambient gas temperature were assumed equal to 173.15 K and 300 K, respectively. The thickness of the methane-hydrate layer, the critical temperature of methane-hydrate dissociation, porosity and convection heat transfer coefficient were taken equal to 10 mm, 193.15 K, 0.6, and 100 W/(m 2 K), respectively. Ambient pressure was atmospheric (101.325 kPa) [ 47,48 ]. Two cases of the distribution of α m along z were considered. This is the linear distribution given by the equation where α m 0 is assumed equal to 0.1 in this section, and We start our analysis with the verification of the prediction of the model described in Section 2 . This verification was performed based on the comparison between the predictions of the new model, hereafter referred to as the analytical model (AM), and the model based on the purely numerical solution of the heat transfer equation in a methane-hydrate layer using the finite element numerical tools for the simulation of heat conduction and convection in COMSOL Multiphysics, hereafter referred to as the numerical model (NM). The results of this comparison for the time evolution of the temperature at z = 0 ( T c ) and z = z i ( T s ) predicted by the analytical and numerical models are shown in  predicted for distribution (11) are slightly lower than those predicted for distribution (12) .
Plots of temperature versus ξ at various time instants are shown in Fig. 4 for the same input parameters as in Fig. 3 . These plots are consistent with the results presented in Fig. 3 . As follows from Fig. 4 , the level at which this temperature reaches the critical temperature (193.15 K) for methane release (methane-hydrate/ice interface) moves towards the wall with time. At t = 20 s the thickness of the methane-hydrate layer had reduced by about 3.4 mm, at t = 40 s, by about 5.3 mm, at t = 60 s, by about 6.8 mm, and at t = 80 s, by about 8.5 mm.
Our next focus is on the relative reduction of the mass of methane ( m ) in the methane-hydrate predicted by the model for distribution (11) . The initial mass of methane is estimated as: where ρ m is methane density, S h the area of the layer, α m 0 is the maximal mass fraction of methane in the methane-hydrate ( α m 0 = 0 . 1 for distribution (11) ).
During the methane release process, when the thickness of the methane layer reduces from z i to z h , the mass of methane in the  (6) show the methane-hydrate/ice interfaces. Distribution (11) was used for these plots.

methane-hydrate reduces to
where α m is the mass fraction of methane at z = z h .
From (13) and (14) it follows that Plots of m/m 0 versus time for h = 10 W/(m 2 K) and five initial thicknesses of methane-hydrate layer ( z h 0 = z i ), and for the initial thicknesses of methane-hydrate layer 3 mm and five convection heat transfer coefficients h are presented in Plots of the layer surface temperature T s versus time for h = 10 W/(m 2 K) and five initial thicknesses of methane-hydrate layer ( z h 0 = z i ) are shown in Fig. 6 (a). All plots are presented until the time instant when all methane has been released from the methane-hydrate. As can be seen from this figure, in all cases T s increases with time. The rate of this increase is rather large initially, which is attributed to the relatively low thermal conductivity of the methane-hydrate. At the later stage, this rate decreases which is attributed to the conversion of methane-hydrate into ice near the surface of the layer (the thermal conductivity of ice is higher than that of methane-hydrate).
Plots of T s versus time for an initial thickness of the methanehydrate layer equal to 3 mm and five convection heat transfer coefficients h are shown in Fig. 6 (b). As follows from this figure, the rate of increase of T s increases with increasing h as expected.

Modelling results versus experimental data
Experimental studies of the dissociation of methane from methane-hydrates with porositiy ε = 0 . 3 were carried out using laboratory scales. The distribution of the mass fraction of methane in the methane-hydrates was approximated by Expression (11) with α m 0 = 0 . 1 . Methane-hydrate samples in the form of tablets of diameter 25 mm and thickness 5 mm were used in the experiments.
The working area in which the methane hydrate was produced was a cylinder of 40X14 steel, which had a wall thickness of 6 mm and a bottom thickness of 12 mm. This design allowed it to withstand pressures up to 150 bar. A WIKA EN 837-1 316L pressure gauge was used to monitor the pressure. The working space of the stand had dimensions: height 70 mm, diameter 30 mm. Water was loaded into the working space with the obturator open. The apparatus was filled with gas through a shut-off valve. The output part of the valve was connected directly to the gas cylinder by means of an adapter. All parts of the device were made from either heat-treated steel 40X13 (body, obturator, ring), steel 12X18N9T (valve, pressure gauge) or BrZH9-4 (power nut). In the first stage, distilled water and a stand were placed in the refrigerator. A day later, ice was taken out of the refrigerator, and was then placed in a chopper. Next, the resulting ice chips were stirred in the body of the stand. At the next stage, the working space of the stand was sealed. The raw materials were mixed in a chamber at an ambient temperature of about 273 K. Next, gas was pumped into the chamber. The procedure was repeated three times to obtain methane-hydrate with a methane concentration of 13-14%.
Before the start of the experiments, the samples were stored in liquid nitrogen in a Dewar vessel at a temperature of 83 K. During the experiment, a sample was removed from the Dewar vessel and placed on the weighing platform of the OHAUS Pioneer PX125D laboratory scales (with the resolution of 0.0 0 0 01 g, and nonlinearity 0.0 0 01 g). The initial mass of the sample was 1 . 3 ± 0 . 1 g.
During a series of 3 to 5 experiments for each sample, continuous recording of the mass of the methane-hydrate sample was carried out. The data was digitised in increments of 1 s.  (1) and predicted by the model assuming that pores in the methane-hydrate are filled with CH 4 (2) and air approximated by N 2 (3). The distribution of the mass fraction of methane was assumed to be given by Expression (11) ; the initial and ambient temperatures were 83.15 K and 273.15 K, respectively; the critical temperature for methane release was 193.15 K, and pressure was atmospheric (101.325 kPa). It was assumed that α m 0 = 0 . 1 , porosity of methane-hydrate and ice ε were equal to 0.3, and h = 5 W/(m 2 K). The thickness of the layer was 5 mm. The input parameters in simulations and experiments were the same.
Two cases were considered. Firstly, it was assumed that pores in both methane-hydrate and ice are filled with nitrogen (approximating air), which is justified at the initial stage of the process. Secondly, it was assumed that pores in both methane-hydrate and ice are filled with methane, which is justified during most of the process, except its initial stage. The transition between these stages was not considered in our analysis.
Plots of m/m 0 versus time observed experimentally and predicted by the model are shown in Fig. 7 . The value of h = 5 W/(m 2 K) used in the model was a fitting parameter (it cannot be inferred from the experimental data available to us). As follows from Fig. 7 , in both cases good agreement between the predictions of the model and experimental data is observed only at the initial stage of the process (up to about 10 s assuming that pores in the methane-hydrate are filled with CH 4 and up to about 6 s assuming that pores in the methane-hydrate are filled with N 2 ). The agreement between the modelling and experimental results is visibly better for Case 1 than for Case 2 at the initial stage of the process (less than about 3 s), while this agreement is better for Case 2 than for Case 1 at times longer than about 3 s. This result is expected as the pores are filled with air at the initial stage of the process and with methane at longer times. Note that the delay in the onset of dissociation was excluded from the results of calculations presented in Fig. 7 to match the conditions during the recording of methane-hydrate dissociation as a dependence of m (t) /m 0 on time in the experiment.
At longer times, in both cases the model predicts faster methane release than observed experimentally. This deviation between modelling results and experimental data is attributed to the main assumption of the model that all methane is instantaneously released from the methane-hydrate when the methanehydrate temperature reaches the critical temperature 193.15 K. The finite rate of methane release would slow down the process which is consistent with the experimental observations, the results of which are shown in Fig. 7 . In other words, the model described in the paper is expected to predict the maximal possible rate of release of methane rather than its actual release rate except at the very beginning of the process.  Fig. 7 (1) and predicted by the model assuming that pores in the methane-hydrate are filled with CH 4 (2) and air approximated by N 2 (3). The distribution of the mass fraction of methane was assumed to be given by Expression (11) ; the same input parameters as in Fig. 7 were used in calculations.
Note that it is not easy to track the dissociation front experimentally. The methane-hydrate turns into ice without contrasting changes in colour or volume. This tracking would require in situ layer-by-layer analysis of the crystal structure or tomographic study. Neither approach was performed in our experiments. Instead the location of this front was inferred from the mass balance (gas/ice/methane-hydrate). The time evolution of a normalised thickness of the ice crust inferred from experimental data and predicted by our model, assuming that pores in the methane-hydrate are filled with CH 4 and N 2 , is shown in Fig. 8 . As follows from this figure, in both cases good agreement between the predictions of the model and experimental data is observed only at the initial stage of the process (up to about 12 s assuming that pores in the methane-hydrate are filled with CH 4 and up to about 8 s assuming that pores in the methane-hydrate are filled with N 2 ), in agreement with the results presented in Fig. 7 .
In the second experiment, the same setup as described above was used but the ambient temperature was taken equal to 223.15 K (-50 • C). At this temperature the contribution of the effects of self-preservation, which are expected to take place in the temperature range 233-268 K, can be safely ignored.
Plots of m/m 0 and ξ i = (z i − z h ) /z i versus time for h = 5 W/(m 2 K) and initial thickness of the hydrate layer 15 mm are presented in Figs. 9 and 10 , respectively. In both cases, two values of porosity ( ε) were considered (0.7 and 0.6). It should also be noted that during dissociation the height of the powder layer is reduced by about 10-15%, which leads to a decrease in ε from the initial 0.65-0.7 to 0.5-0.6 at the final stage of dissociation. With a decrease in porosity, the heat of dissociation increases, which leads to more intensive cooling of the powder and, accordingly, to a slowdown in the rate of dissociation at the end of the process.
This allows us to consider ε = 0 . 7 to be typical for the beginning of the process, while ε = 0 . 6 is typical for the final stage. As follows from Figs. 9 and 10 , the agreement between the modelling and experimental results is much better for the second experiment than for the first, the results of which are shown in Figs. 7 and 8 . This is attributed to the fact that our model does not consider the effect of self-preservation which is not expected to take place in the second experiment, but can play an important role in the first experiment except at the initial stage of the process. Note that the decrease in the ambient temperature led to a substantial increase in the duration of the process from about 20 s to about 10 minutes.  (1), predicted by the model assuming that porosity ε = 0 . 7 (2), predicted by the model assuming that ε = 0 . 6 (3). The distribution of the mass fraction of methane was assumed to be given by Expression (11) ; the initial and ambient temperatures were 192.15 K and 223.15 K, respectively; the critical temperature for methane release was 193.15 K, and pressure was atmospheric (101.325 kPa). Pores in the methane-hydrate were assumed to be filled with CH 4 .  (1), predicted by the model assuming that ε = 0 . 7 (2), predicted by the model assuming that ε = 0 . 6 (3). The distribution of the mass fraction of methane was assumed to be given by Expression (11) ; the initial and ambient temperatures were 192.15 K and 223.15 K, respectively; the critical temperature for methane release was 193.15 K, and pressure was atmospheric (101.325 kPa). Pores in the methane-hydrate were assumed to be filled with CH 4 .
As follows from Figs. 9 and 10 , the agreement between modelling and experimental results is better for ε = 0 . 7 than for ε = 0 . 6 at the beginning of the process than at its end. The opposite is true at the end of the process. This agrees with the above-mentioned observation that ε = 0 . 7 is typical for the beginning of the process, while ε = 0 . 6 is typical for the final stage. The increased deviation between modelling and experimental results at times longer than 500 s is attributed to weighting errors when the mass of the remaining methane in the powder is low (about 0.01 g).

Conclusions
A simple model for the dissociation of methane from a layer of methane-hydrate particles is suggested. The model is based on the assumption that methane-hydrate particles form a flat porous film lying on an adiabatic wall. This film is heated by ambient air convection. It was assumed that when the temperature inside a certain region within the methane-hydrate reaches the critical temperature for the release of methane from the methane-hydrate ( T cr ), all methane is released from the methane-hydrate and the latter turns into ice. The contribution of heat required for the release of methane is considered in the model.
Two cases of the initial distribution of the mass fraction of methane in the methane-hydrate layer ( α m ) were considered. Firstly, it was assumed that α m linearly decreases from its maximal value at the wall to zero at the surface of the film. Secondly, the distribution was assumed to be constant. The first distribution is closer to that observed experimentally and our analysis was mainly focused on this case. The model is based on the analytical solution to the onedimensional heat transfer equation in a two-layer (methanehydrate/ice) system. This solution was incorporated into the numerical code and used at each time step of the calculations.
It was seen that at the very initial stage of methane-hydrate layer heating no conversion of methane-hydrate into ice was predicted until the methane-hydrate surface temperature reached T cr . Once the hydrate surface temperature reached T cr the speed of the recession of the methane-hydrate/ice interface was predicted to increase rapidly until it reached the maximal value. After that it began to decrease.
The speeds of the predicted release of methane from the methane-hydrate and rise of the layer surface temperature were shown to increase with decreasing initial thickness of the layer and increasing values of the convection heat transfer coefficient.
Model predictions were compared with experimental data relevant to the heating of a layer of methane-hydrate with porosity 0.3 and thickness 5 mm, and where the mass fraction of the methane decreased from 0.1 (at the wall) to zero (at the layer surface). The initial and ambient temperatures were 83.15 K and 273.15 K, respectively. Two cases were considered. Firstly, it was assumed that the pores in the methane-hydrate were filled with air (approximated by nitrogen). Secondly, it was assumed that the pores in the methane-hydrate were filled with methane. The first case refers to the very initial stage of the process. The second case refers to most of the process except its initial stage.
It was shown that in both cases good agreement between the predictions of the model and experimental data is observed at the initial stage of the process. At longer times, in both cases the model predicted faster methane release than observed experimentally unless the effect of self-preservation can be ignored. This deviation between modelling results and experimental data was attributed to the main assumption of the model that all methane is instantaneously released from the methane-hydrate when the methane-hydrate temperature reaches the critical temperature for methane release (193.15 K). Thus, the model is expected to predict the maximal possible rate of release of methane rather than its actual release rate except at the beginning of the process.

Declaration of Competing Interest
The work described has not been published previously and it is not under consideration for publication elsewhere. Its publication is approved by all authors and tacitly or explicitly by the responsible authorities where the work was carried out. If accepted, it will not be published elsewhere in the same form, in English or in any other language, including electronically without the written consent of the copyright-holder.

Data availability
Data will be made available on request.

Appendix A. Derivation of Formula (10)
Let us consider an infinite composite layer of porous methanehydrate and ice placed on an adiabatic wall. The layer of ice, of thickness z i − z h , is placed above the layer of methane-hydrate of thickness z h . The surface of the layer of ice is heated by ambient air at temperature T g with a convection heat transfer coefficient of h .
These assumptions lead us to the following transient heat transfer equation for the composite film with the following initial condition: where z h and z i are the heights of the upper surfaces of the methane-hydrate and ice layers, respectively; κ is thermal diffusivity. Note that The boundary condition for Eq. (A.1) at the adiabatic wall is presented as: The boundary conditions for Eq. (A.1) at the methane-hydrate/ice interface are presented as: where ρ h is the density of the hydrate, α m is the mass fraction of the methane in the methane-hydrate, Q m is the specific heat required for the release of methane (in J/kg), d z h d t ≤ 0 is the speed of the movement of the methane-hydrate/ice interface, estimated at the previous time step. Note that z h and d z h d t remain constant during each time step. z i remains constant during the whole process as the densities of methane-hydrate and ice are assumed to be constant and equal; their temperature dependence is ignored.
The boundary condition at the surface of the ice is presented as Having introduced the new variable and functions and remembering that z i is constant, we re-write (A.1) as The initial and boundary conditions for Eq. (A.5) are presented as: To convert the boundary condition at the ice surface ( Eq. (A.9) ) into a homogeneous one, the new function V (t, ξ ) is introduced via the equation: Expression (A.11) allows us to re-write Eq. (A.5) and initial and boundary conditions (A .6) -(A .9) as We look for the following solution to Eq. (A.13) : where v n (ξ ) form the full set of non-trivial solutions to the eigenvalue problem: subject to boundary conditions: From the last condition in System (A.22) we obtain B 1 = 0 . Three other conditions (second, third, and first in System (A.22) ) lead to the following system of equations To have a nontrivial solution to System (A.25) , the determinant of this system should be equal to zero: Remembering that B 1 = 0 , Expressions (A.24) can be re-written as: (A.29) The following formula was used when deriving (A.28) The explicit expressions for v n for λ n inferred from (A.
Using the following formulae (A. 40) The fact that functions v n are orthogonal with weight b allows us to perform the following expansions: to discrete values of ( n ) subject to the initial condition: Assuming that p n is constant during the time step, the solution to (A.48) is straightforward: This expression is identical to Expression (10) .

Appendix B. Derivation and analysis of Equation (A.27)
The expansion of (A.26) leads to the following expression for : Expressions (B.2) and (B.3) can be presented as: (B.5) Remembering that: Expressions (B.4) and (B.5) can be simplified to: Dividing both sides of (B.9) by −λ cos (λa h ξ h ) cos (λa i (1 − ξ h )) we obtain Eq. (A.27) . A graphical solution of Eq. (A.27) for a typical set of input parameters (see [ 47,48 ]) is illustrated in Fig. 11 . The methane-hydrate initial temperature and ambient gas temperature were taken as 173.15 K and 300 K, respectively. The convection heat transfer coefficient was assumed equal to 100 W/(m 2 K). The thickness of the methane-hydrate layer was taken as 10 mm, ξ h = 1 , the ambient pressure was set to 101.325 kPa. The porosity was assumed equal to 0.6. The mass fraction of methane in the methane-hydrate was inferred from (11) Table 1 Density, thermal conductivity, and specific heat capacity of nitrogen, methane, ice and methane-hydrate at temperatures T = 193 . 15 K for ice and methane-hydrate, and T = 300 K for N 2 and CH 4 [1] . In our analysis the difference between the densities of methane-hydrate and ice was ignored, and the density of methane-hydrate was assumed to be the same as that of ice.