Log In Sign Up

Blow up in a periodic semilinear heat equation

Blow up in a one-dimensional semilinear heat equation is studied using a combination of numerical and analytical tools. The focus is on problems periodic in the space variable and starting out from a nearly flat, positive initial condition. Novel results include various asymptotic approximations that are, in combination, valid over the entire space and time interval right up to and including the blow-up time. Preliminary results on continuing a numerical solution beyond the singularity are also presented.


page 1

page 2

page 3

page 4


On central limit theorems for power variations of the solution to the stochastic heat equation

We consider the stochastic heat equation whose solution is observed disc...

On the numerical approximations of the periodic Schrödinger equation

We consider semidiscrete finite differences schemes for the periodic Scr...

Controllability properties from the exterior under positivity constraints for a 1-D fractional heat equation

We study the controllability to trajectories, under positivity constrain...

A Space-Time DPG Method for the Heat Equation

This paper introduces an ultra-weak space-time DPG method for the heat e...

Semi-analytic integration for a parallel space-time boundary element method modeling the heat equation

The presented paper concentrates on the boundary element method (BEM) fo...

Nonlinear inhomogeneous Fokker-Planck models: energetic-variational structures and long time behavior

Inspired by the modeling of grain growth in polycrystalline materials, w...

Cleaning Schedule Optimization of Heat Exchanger Networks Using Particle Swarm Optimization

Oil refinery is one of industries that require huge energy consumption. ...

1 Introduction

The blow-up phenomenon in nonlinear evolution equations has been studied extensively in the literature. Some studies focus on physical applications, such as singularity formation in fluids Braun and Kluwick (2004, 2005); Hocking et al. (1972); Lushnikov et al. (2021), runaway in thermal processes Dold (1991); Herrero and Velázquez (1993); Lacey (1983), and biological applications Jabbari and King (2013)

. Others deal with numerical aspects such as the computation of blow-up profiles, estimation of blow-up times, and singularity tracking 

Berger and Kohn (1988); Budd et al. (1996); Keller and Lowengrub (1993); Tourigny and Grinfeld (1994); Weideman (2003). For a general review, see Galaktionov and Velázquez (2002).

The present paper is a continuation of the papers Keller and Lowengrub (1993) and Weideman (2003). The equation considered in these papers is the nonlinear heat equation


(although more general nonlinearities and more than one space dimension were also considered in Keller and Lowengrub (1993)). In this paper we consider the initial condition


i.e. periodic, positive, and nearly flat. (A broader class of nearly flat initial data is considered in Appendix A.) Figure 1 shows a typical blow-up scenario for this case.

Figure 1: Typical solution profiles of equation (1), in the case of a nearly flat initial condition such as (2). The solution stays almost flat for a long while until a point blow-up occurs over a relatively short period of time.

In Weideman (2003), it was shown by numerical computation that the approach to blow up is not necessarily uniform. That is, there may be times when the diffusive term dominates, leading to a flattening of the solution profile. At other times, particularly near blow up, the nonlinearity dominates, leading to a steepening of the profile. By numerically continuing the solution into the complex plane, it was shown that this behaviour can be associated with the dynamics of the complex singularities of the solution. When diffusion dominates, they typically move farther from the real axis while the opposite is true when nonlinearity dominates. A point blow up occurs when the singularities reach the real axis. For the nearly flat initial conditions considered in this paper, however we show that there is no simple correspondence between the flatness of the solution profile on the real interval and the proximity of the nearest singularity, except in the blow-up limit, where the steepening solution profile is associated with the impingement of singularities onto the real axis at .

In Keller and Lowengrub (1993), the focus was not on complex singularities but on asymptotic estimation of the blow-up time and the solution profile near blow up. The starting point in that paper was the substitution , which transforms (1) into


This transformation has advantages for both analysis and computation. Advantages for analysis are spelled out in Dold (1991). For numerical computation, it is easier to deal with solutions tending to zero than to infinity. This avoids the need for specialized rescaling algorithms or moving mesh methods Berger and Kohn (1988); Budd et al. (1996). The downside is that the simple polynomial nonlinearity in (1) has been replaced by the more complicated nonlinearity in (3). While this can be ameliorated by multiplication by , it raises the red flag of division by zero as becomes small. However, it was shown in Keller and Lowengrub (1993) and reaffirmed in section 3 of the present paper that the nonlinear term remains bounded as approaches zero. This further suggests the possibility of integration through the zero of , i.e. through the singularity of , but Keller and Lowengrub (1993) reported a failed effort. We continued those investigations and announce preliminary findings here.

In Weideman (2003) periodic boundary conditions were considered while Keller and Lowengrub (1993) looked at the pure initial value problem on the entire real line. In this paper we continue with the periodic case, as this gives us access to highly accurate Fourier spectral methods (which can also be applied to problems on the infinite line, but with a substantial penalty in accuracy). One contribution here is a conversion of the analysis of Keller and Lowengrub (1993) to the periodic situation, which is not just a triviality but contributes significant new results as we shall discuss.

For the numerical computations of this paper, a full spectral method based on a Fourier series


is used. The derivatives on the right-hand side of (3) are computed by analytical differentiation of this series, and the nonlinear terms by convolution and de-convolution. The result is an infinite system of ODEs for the evolution of the Fourier coefficients , which we truncate at . To integrate this system, we use ode45, MATLAB’s standard ODE solver. It has adaptive time stepping to maintain accuracy, and the tolerance parameters for this are set to a stringent . To compute the time at which the numerical solution blows up (i.e. ), we use the fact that the initial conditions considered here lead to blow up at and hence we check when the sum of the equals zero. This can be done by the ‘event’ option in ode45.

We use this numerical solution as reference solution for the purposes of checking the various asymptotic estimates. Because of the relatively smooth nature of , even at the critical time (as will be discussed), we believe it is sufficiently close to the true solution for all verification purposes.

The three main sections of the paper can be summarized in a nutshell as: before blow up, at blow up, and after blow up. More specifically, in section 2 we present a perturbation analysis that approximates the solution to (3) accurately to on the entire periodic space domain and the whole time interval until a time that is close to blow up. Beyond that time up to the blow-up point it has to be modified and this is done by matched asymptotic expansions, the details of which are contained in Appendix A. In this appendix we also analyse the dynamics of the singularities of the solution and in Appendix B the relation between the proximity of the singularities to the real axis and the steepness of the solution profile on is clarified. In sections 2 and 3, the analyses in the appendices are confirmed by numerical experiments. Section 4 is a preliminary report on our efforts in integrating through the singularity at the blow-up time and the subsequent evolution.

2 Two-mode perturbation analysis

The analysis of Keller and Lowengrub (1993) was based on the truncated Taylor expansion


By substituting into (3) and dropping powers of the problem was reduced to a dynamical system in the variables and ; see (18) below. Here we follow an analogous procedure, but consider instead a truncated Fourier expansion


Under the assumption of strictly positive solutions, i.e., , both of these approximations blow up (in the variable ) in finite time at .

Substitution of (6) into (3) and neglecting contributions gives the system


or, in explicit form,


(Note that the assumption precludes the vanishing of the denominators.) The phase plane of this system is shown in Figure 2.

Figure 2: Phase plane of the system (8) in the domain . The analysis of this section approximates solution curves close to . Blow up at in the -equation (1) corresponds to solution trajectories intersecting the dashed line .

Consider solution curves in Figure 2 that originate near , say


As an explicit solution of the system (7) appears not to exist, we settle for a perturbation analysis, by expanding


where , , . Substitution into (7) gives, to zeroth order


At ,


This gives, to , the approximate solution


Plugging this expression into (3) gives


which confirms that the -equation is satisfied to uniformly in , for all bounded away from zero. Figure 3 illustrates that the perturbation approximation (13) and the (numerically computed) solution to the two-mode system (8) are good approximations to the numerical reference solution of (3) on most of the interval , where is the critical time at which blow up occurs. As , the assumption underlying the perturbation analysis and the two-mode approximation, , is no longer valid, however, and the approximations lose accuracy.

Figure 3: The maximum relative error on of the perturbation approximation (13) (blue) and the numerical solution to the two-mode system (8) (red) as a function of . The errors are calculated with reference to the numerical solution of (3) mentioned in the final paragraph of section 1.

Various estimates can be obtained from the perturbation solution (13). By setting at and excluding terms, for example, one obtains the following estimate for the blow-up time


In Appendix A, using the method of matched asymptotic expansions, a higher-order estimate of is derived (see (43) and (67)):


the constants being defined in (66). The accuracy of these estimates is verified in Table 1. (The values of listed in the table were computed by the method described in section 1, and are believed to be correct to all digits shown.)

0.1 0.161963 3.6e-04 1.0e-02 2.6e-02
0.01 0.242093 1.8e-03 1.2e-04 2.8e-04
0.001 0.249220 2.2e-04 1.2e-06 2.8e-06
0.1 0.955542 2.1e-03 7.7e-03 4.5e-03
0.01 0.996241 9.5e-04 8.1e-05 4.9e-05
0.001 0.999631 1.1e-04 8.1e-07 4.9e-07
0.1 3.996685 5.0e-04 1.5e-03 1.4e-05
0.01 3.999802 5.3e-05 1.5e-05 1.2e-07
0.001 3.999982 5.4e-06 1.5e-07 1.2e-09
Table 1: Blow-up times for various parameter choices in the initial condition (2): is the blow-up time as computed from the reference solution, is the blow-up time as estimated from a numerical solution of the two-mode system (7) and and denote the estimates (15) and (16), respectively.

The singularity dynamics mentioned in the introduction can be estimated as follows. The complex singularity in the -equation corresponds to a complex zero of the -equation. Assuming it is located at with real, setting gives


This approximation is valid, however, only for small values of , i.e., near blow up. This follows from the fact that the remainder term on the right-hand side of (14) grows exponentially with . In Figure 4, the asymptotic estimate (17) is compared to a numerical estimate of the singularity location obtained via the method of Sulem et al. (1983), which is essentially an estimation of the width of the strip of analyticity of the solution, by examining the rate of decay of its Fourier coefficients. To apply the method of Sulem et al. (1983), we use the fact that, to leading order and away from the blow-up time, singularities of solutions to the -equation are second-order poles111In Fasondini et al. we show that these singularities are in fact logarithmic branch points, however the branch point singularity appears in the fourth-order term in the local expansion about the singularity. Therefore, for the purpose of estimating the position of the singularity, we only use its leading-order, second-order pole behaviour. We shall find in section 3 that in the limit , the leading order behaviour of the singularity at , which results from the coalescence of two singularities, is of a more complicated form than that of a second-order pole. . Figure 4 also shows asymptotic estimates of the singularity location that are derived in Appendix A.4 by the method of matched asymptotic expansions. The estimates are valid in the limits and but in between, for , the singularity location of the original PDE (3) is described by the singularity location of a more difficult nonlinear initial value problem (68)–(70) which is not analytically solvable. The asymptotics nonetheless correctly indicate the initial movement of the singularity away from the real axis (at a speed that becomes infinite as ) and the final motion shortly before the singularity collides with the real axis at .

Figure 4: The position on the positive imaginary axis of the singularity of the solution to (3) with initial condition . (In the right frame, the approximation (72) is not visible because it is indistinguishable from (74).)

In Appendix B we investigate the relationship between the position of the complex singularity and the height of the peak of the solution profile in the -variable (which is located at , see Figure 1) relative to the solution value at .

Returning to the analysis of Keller and Lowengrub (1993), which was based on the truncated Taylor approximation (5), we note that the focus in that paper was on the behaviour at , not the evolution on as is the focus here. It is therefore instructive to adapt that analysis here and compare results to (13) and (14).

The dynamical system analogous to (7)–(8) is now


This is a much simpler system, and in fact admits a first integral , although we shall not make use of this result and neither was it used in Keller and Lowengrub (1993). (No such first integral could be found for (7)–(8).)

Proceeding with a perturbation analysis based on (9)–(10) give , , and . Therefore, excluding terms, is approximately


in analogy with (13). Plugging into (3) gives


Comparing with (14) shows an advantage of the periodic analysis, namely that its right-hand side is uniformly in , whereas the right-hand side of (20) is only if . That is, with bounded away from zero (13) provides a valid approximation over the entire (periodic) domain, while (19) is valid only near .

On the other hand, the advantage of the analysis of Keller and Lowengrub (1993) is that it gives a valid description of the structure of the solution close to blow up, as discussed in Appendix A.

3 Solution in the blow-up limit

Appendix A

is devoted to an asymptotic analysis of the solution on three time scales that are progressively closer to the blow-up time, namely

, and . The analysis is performed using the method of matched asymptotics and the results on the first time scale are the same as the those obtained in section 2 via a regular perturbation analysis (in particular, recall the approximation (13) and its failure close to the blow up time, as seen in Figure 3).

On the second time scale, , the following asymptotic approximation is derived for :


which is obtained by combining (44)–(46) and (65)–(67). It also follows that, by setting in (21), we obtain a representation of the solution profile at the blow-up time that is valid as , with :


Figure 5 verifies the accuracy of the asymptotic approximations (13) (away from the blow-up time) and (21)–(22) (close to and at the blow-up time). The left frame of Figure 6 shows the numerical solution at the blow-up time with the blow-up profile (22) superimposed on it. Hence, we have accurate asymptotic expressions for the solution on the entire spatial interval and from all the way up to and including the blow-up time. For at , however, we shall need another asymptotic approximation, namely (24), to be discussed below.

Figure 5: Errors in the asymptotic approximations (13) (blue) and (21) (red) when compared to the reference solution. The numerical integration method returns approximations at times , , where and . Here and (more digits are listed in Table 1). At each , the relative errors of (13) and (21) are calculated on the interval and plotted against the time step index . Because of the adaptive time-stepping, the values are not equidistant, but the spacing is much denser near . As an indication of this, note that the two curves intersect at and , which is already quite close to the critical time even though roughly 800 more time steps are to be taken. From the asymptotic analysis, we expect the approximation on the first time scale (13) (blue) to start to break down and (21) (red) to be valid on the second time scale when . Here we have that at around , , which is indeed when the error curves of the two approximations start to diverge noticeably.

The rate of decay of the Fourier coefficients of will be of relevance in the next section when we consider the possibility of continuing the solution beyond blow-up. From (22) we deduce that the decay as close to the blow-up time because

which implies that


The right frame of Figure 6 confirms the accuracy of this estimate as .

Figure 6: Left: The solution to the -equation (3), displayed on a log-scale, at the blow-up time (solid curve) compared to the asymptotic estimate (22) (dashed curve on top of the solid curve). Here we would like to draw attention to the fact that, at the origin and at the critical time, the values of are computed to levels close to roundoff error (), which translates into huge values in the solution near blow up. (See for example the third frame in Figure 9.) This is a consequence of solving equation (3) instead of (1). Right: The Fourier coefficients of the reference numerical solution at the blow-up time (dots) compared to the coefficients of the global asymptotic blow-up profile (dashed lines), given in (23). The single solid line shows the asymptotic estimate of the decay of the Fourier coefficients of the local blow-up profile, which is given in (25).

At the blow-up time, for exponentially small with respect to , (22) is no longer valid (see the discussion below (61)). Instead,


which follows from (62) and (67). Figure 7 confirms that (24) has better accuracy than (22) for small222Ideally, we would show the accuracy of (24) for even smaller values of than those in Figure 7. However, this would require high-precision computations. In standard double precision (with a machine precision of approximately ) we cannot compute the errors for smaller because round-off errors prevent the accurate evaluation of from its numerically computed Fourier coefficients. at the blow-up time.

Figure 7: The relative error of the asymptotic approximations (22) (blue) and (24) (red) at the blow-up time for small .

Regarding the strength of the singularity at , it follows from (24) that , for but the third derivative, however, blows up as , . A singularity with two bounded derivatives is consistent with the decay of the Fourier coefficients that was derived from (22)333The ‘global’ blow-up profile (22) has a stronger singularity at than the ‘local’ blow-up profile (24) since its second derivative blows up logarithmically at while for the local profile , . The Fourier coefficients of the ‘global’ profile nevertheless decay at the correct rate since its second derivative, though unbounded, is integrable.. In fact, the asymptotic approximation (24) suggests that the Fourier coefficients of the solution precisely at the blow-up time decay slightly faster than the suggested by (23). To show this, symmetry and integration by parts can be used to obtain

From (24) it follows that , , and by making the change of variable with , one obtains,


The solid line in Figure 7 shows the estimate (25), which, unlike the estimate (23), is independent of and . The numerical Fourier coefficients shown in Figure 7 do not decay as fast as (25), which suggests that one would need to perform high-precision computations (to compute the solution closer to the blow-up time or with a larger number of Fourier coefficients) to observe the rate of decay predicted by (25).

For the possible continuation beyond blow-up one needs to confirm that the nonlinear term in the -equation (3), i.e., , remains bounded at at the blow-up time. In fact, from (24) it follows that the nonlinear term vanishes, according to , for .

4 Integrating through the singularity

The fact that the nonlinear term in (3) remains bounded as approaches zero raises the intriguing possibility of numerically integrating through the blow up. This was tried in Keller and Lowengrub (1993), but the authors found “…the calculation actually continues the solution slightly beyond the blow-up time of the original solution . The method becomes unstable a short time after the blow-up happens, however.” We report here on renewed efforts in this direction. (Details of how our numerical strategies differ from those of Keller and Lowengrub (1993) are postponed to the end of the section.)

Figure 8 shows a series of snapshots of the solution at different times. First, observe that there is no sign of instability as the solution passes smoothly through (third frame). What happens next might be unexpected, namely, the solution turns complex. In addition, uniqueness is lost: the complex conjugate of the solution shown is equally likely to appear when parameters (such as the error tolerances in the time integration scheme or the number of terms retained in the Fourier series (4)) are adjusted slightly. Surprising as these results may be, both the non-uniqueness and the fact that the solution turns complex are consistent with theoretical results of Masuda (1984).

Figure 8: Solution to the -equation (3) corresponding to the initial condition (2) with , . (The corresponding solution is shown in Figure 9, and a movie of the dynamics can be seen in the supplementary material that accompanies this paper.) The two noteworthy times are (third frame) and (seventh frame), where , approximately. At there is a zero of at , after which the solution turns complex: the real part is shown in blue and the imaginary part in red. At there are approximate zeros of near

. The solution shown is not unique (its conjugate is equally probable, as discussed in the text.)

Figure 9 shows the -solution that corresponds to the -solution in Figure 8. After the solution turns complex at the critical time (third frame) the modulus of shows wave-like behaviour, with two waves travelling in opposite directions from the origin until they reach the edge of the domain . Because of periodicity they meet up with similar waves from adjacent intervals and a second blow up almost occurs, this time at (seventh frame). The modulus grows considerably and we conjecture that by varying the parameters in the initial condition a proper secondary blow up may be found. Computing over a longer time interval suggests that asymptotically approaches the constant (and real) solution as .

Figure 9: Same as Figure 8 but here the -solution is shown. Note the blow up at and the near blow up at . The dashed curve is the modulus (shown only after the first blow up). Note that the scales on the vertical axes are different in each frame.

We present the results of Figures 8 and 9 knowing full well that we have little theory to draw on as validation. The theoretical results of Masuda (1984)

suggest the possibility of a complex and non-unique solution post blow-up, but does not allow for a quantitative comparison. Nevertheless, the following heuristic observations support the validity of the results shown here.

Firstly, the fact that the numerical solution approaches as provides some confidence, as this solves (1) exactly. Secondly, as noted before, the singularity in the -equation at the critical time is rather weak, leading to Fourier coefficients that decay at the relatively rapid rate ; recall (25). While still much slower than the typical exponential decay rate for analytic periodic functions, we conjecture that this decay is nevertheless sufficiently rapid that the accuracy loss of the spectral method at the critical time is not disastrous. Figure 10 shows the Fourier coefficients near and at the critical time.

As far as we know, the only published results that deal with numerical computations of post-blow-up solutions are Cho et al. (2016); Takayasu et al. (2022). These authors based their computations on complexification of the -variable. Following this idea we integrated along a path that contains a semi-circle in the complex -plane, centred at the estimated singularity. Using this approach we obtained the same results as those shown in Figures 8 and 9.

It should be noted that the methods of Cho et al. (2016); Takayasu et al. (2022) are based on the -equation, not the -equation. Because the solution to the equation grows without bound near the critical time, these methods cannot compute blow-up solutions in the immediate neighbourhood of the critical time accurately. For the same reason, estimates of blow-up times are unreliable. Adapting the methods of Cho et al. (2016); Takayasu et al. (2022) from the -equation to the -equation might be a worthwhile future project.

Figure 10: Modulus of the Fourier coefficients of the solution shown in Figure 8. The times are just prior to , at , and shortly after. (The value of is given to more digits in Table 1.) In the first frame the coefficients decay exponentially, indicative of a function analytic in a neighbourhood of the real axis. At the coefficients decay algebraically as given by (25), because of the influence of the singularity of as it reaches the origin. Immediately afterwards the exponential decay is recovered.

The observation of a complex solution after the critical time has lead us to conjecture that the failure of the finite difference method of Keller and Lowengrub (1993) was because it was almost surely coded as a real system using real arithmetic. One way to allow for possible complex solutions is to use complex arithmetic (default in MATLAB) or simply to split the -equation (3) into its real and imaginary parts. We integrated such a separated system with the Fourier spectral method mentioned in the first section, which is how the results of Figures 810 were computed. The same idea applies, however, to the finite difference method of Keller and Lowengrub (1993).

Initialising the imaginary part to strictly zero values, the ODE software signals a singularity at the critical time. Initializing it with a small, random perturbation on the order of machine roundoff level (), however, allowed the solution to pass through the critical without any sign of instability, as shown in Figure 8.

One priority for future investigations is a better understanding of the transition from a real to a complex solution and the associated non-uniqueness. In our case the nonzero imaginary part is triggered by noise at the level of roundoff error. The radomness dictates whether the continuation is with one solution or with its complex conjugate.

5 Conclusions

We investigated, asymptotically and numerically, point blow-up solutions to a periodic nonlinear heat equation (1) with nearly flat initial data by considering the solution in the reciprocal variable . We derived asymptotic approximations for the solution on the entire spatial interval and from up to and including at the blow-up time, for which we also derived a second-order approximation. Due to the high accuracy of the Fourier spectral method (including at the blow-up time at which the -solution has a weak singularity, unlike the -solution), we were able to check numerically the validity of the asymptotics. We believe it is unusual for a single numerical method to confirm asymptotics in multiple regimes since typically numerical methods are weak, or highly inefficient, in most asymptotic limits unless they are highly specialised. The key to the success of the numerical method used here is the fact that it approximates the relatively well-behaved -equation (3) rather than the -equation (1) whose solution becomes unbounded.

The investigations in this paper point to a number of topics for future research, not least of which is the validity of the post-blow-up solutions computed in section 4. In addition, the dynamics of complex singularities of blow-up solutions to (1) for a larger class of initial conditions, including initial data leading to non-generic forms of blow up, will be investigated in Fasondini et al. , also via a combination of asymptotic and numerical methods. In Fasondini et al. we shall also explore the singularity structure of these blow-up solutions on their Riemann surfaces in the complex -plane.

The first author is grateful to Saleh Tanveer for stimulating and insightful discussions. This research was started while the three authors were in residence at the Isaac Newton Institute for Mathematical Sciences as part of the programme Complex analysis: techniques, applications and computations. This programme was supported by: EPSRC grant number EP/R014604/1. The work of the first author was also supported by the Leverhulme Trust Research Project Grant RPG-2019-144. A grant to the third author from the H.B. Thom foundation of Stellenbosch University is also gratefully acknowledged.


  • [1] M. Berger and R. V. Kohn (1988) A rescaling algorithm for the numerical calculation of blowing-up solutions. Comm. Pure Appl. Math. 41 (6), pp. 841–863. External Links: ISSN 0010-3640, Document, Link, MathReview (Hong Yuan Fu) Cited by: §1, §1.
  • [2] S. Braun and A. Kluwick (2004) Unsteady three-dimensional marginal separation caused by surface-mounted obstacles and/or local suction. J. Fluid Mech. 514, pp. 121–152. Cited by: §1.
  • [3] S. Braun and A. Kluwick (2005) Blow-up and control of marginally separated boundary layers. Philos. Trans. R. Soc. Lond. Ser. A Math. Phys. Eng. Sci. 363 (1830), pp. 1057–1067. Cited by: §1.
  • [4] C. J. Budd, J. Chen, W. Huang, and R. D. Russell (1996) Moving mesh methods with applications to blow-up problems for PDEs. In Numerical analysis 1995 (Dundee, 1995), Pitman Res. Notes Math. Ser., Vol. 344, pp. 1–18. External Links: MathReview Entry Cited by: §1, §1.
  • [5] C.-H. Cho, H. Okamoto, and M. Shoji (2016) A blow-up problem for a nonlinear heat equation in the complex plane of time. Jpn. J. Ind. Appl. Math. 33 (1), pp. 145–166. External Links: ISSN 0916-7005, Document, Link, MathReview Entry Cited by: §4, §4.
  • [6] J. W. Dold (1991) On asymptotic forms of reactive-diffusive runaway. Proc. Roy. Soc. London Ser. A 433 (1889), pp. 521–545. External Links: ISSN 0962-8444, Document, Link, MathReview Entry Cited by: §1, §1.
  • [7] M. Fasondini, J. King, and J.A.C. Weideman Complex-plane singularity dynamics for blow up in a nonlinear heat equation: analysis and computation. In preparation. Cited by: §A.4, §5, footnote 1.
  • [8] V. A. Galaktionov and J. J. L. Velázquez (2002) The problem of blow-up in nonlinear parabolic equations. Disc. Cont. Dyn. Syst. 8 (2), pp. 399–433. Cited by: §A.3, §1.
  • [9] M. A. Herrero and J. J. L. Velázquez (1993) Plane structures in thermal runaway. Israel J. Math. 81 (3), pp. 321–341. External Links: ISSN 0021-2172, Document, Link, MathReview (Victor A. Galaktionov) Cited by: §1.
  • [10] L.M. Hocking, K. Stewartson, J.T. Stuart, and S.N. Brown (1972) A nonlinear instability burst in plane parallel flow. J. Fluid Mech. 51 (4), pp. 705–735. Cited by: §1.
  • [11] S. Jabbari and J.R. King (2013) Discrete and continuum multiscale behaviour in bacterial communication. In Multiscale Computer Modeling in Biomechanics and Biomedical Engineering, pp. 299–320. Cited by: §1.
  • [12] J. B. Keller and J. S. Lowengrub (1993) Asymptotic and numerical results for blowing-up solutions to semilinear heat equations. In Singularities in fluids, plasmas and optics (Heraklion, 1992), NATO Adv. Sci. Inst. Ser. C Math. Phys. Sci., Vol. 404, pp. 111–129. Cited by: §A.2, §1, §1, §1, §1, §2, §2, §2, §2, §4, §4.
  • [13] A. A. Lacey (1983) Mathematical analysis of thermal runaway for spatially inhomogeneous reactions. SIAM J. Appl. Math. 43 (6), pp. 1350–1366. External Links: ISSN 0036-1399, Document, Link, MathReview (J. W. Bebernes) Cited by: §1.
  • [14] P. M. Lushnikov, D. A. Silantyev, and M. Siegel (2021) Collapse versus blow-up and global existence in the generalized Constantin-Lax-Majda equation. J. Nonlinear Sci. 31 (5), pp. Paper No. 82, 56. External Links: ISSN 0938-8974, Document, Link, MathReview Entry Cited by: §1.
  • [15] K. Masuda (1984) Analytic solutions of some nonlinear diffusion equations. Math. Z. 187 (1), pp. 61–73. External Links: ISSN 0025-5874, Document, Link, MathReview (Jan Prüss) Cited by: §4, §4.
  • [16] C. Sulem, P.L. Sulem, and H. Frisch (1983) Tracing complex singularities with spectral methods. J. Comput. Phys. 50 (1), pp. 138–161. Cited by: §2.
  • [17] A. Takayasu, J.P. Lessard, J. Jaquette, and H. Okamoto (2022) Rigorous numerics for nonlinear heat equations in the complex plane of time. Numer. Math. 151 (3), pp. 693–750. External Links: ISSN 0029-599X, Document, Link, MathReview Entry Cited by: §4, §4.
  • [18] Y. Tourigny and M. Grinfeld (1994) Deciphering singularities by discrete methods. Math. Comp. 62 (205), pp. 155–169. Cited by: §1.
  • [19] J.A.C. Weideman (2003) Computing the dynamics of complex singularities of nonlinear PDEs. SIAM J. Appl. Dyn. Syst. 2, pp. 171–186. Cited by: §1, §1, §1, §1.

Appendix A Analysis via the method of matched asymptotic expansions

a.1 Truncated Fourier expansion

In this section we revisit and expand the analysis based on the two-mode Fourier truncation outlined in section 2. Recall that by substituting (6) into the -equation (3) and neglecting the term the system (7) was obtained. From this, one finds by a self consistency argument that near blow-up, generically (i.e. even for )


where and the blow up time are positive constants, so that the blow-up behaviour associated with (6) takes the form


We emphasise that (26)–(27) describe the blow-up behaviour of (7) rather than that of (3): (27) has features in common with, but is not a valid representation of, the blow-up behaviour of the full PDE (3), a key point to which we shall return.

In keeping with our goal of characterising the behaviour of (3) for near-flat initial conditions, we now develop a fully analytic asymptotic description of the behaviour of the two-mode system (7) for initial conditions (9) with . There are two timescales; the first coincides with (10), so that at leading order (11) and (12) follow.

On the second time scale, and near blow up, we set


with , so that


Since each satisfy

we find on matching with (11)–(12) that


(a result that relies on the observation that has no contribution for , see (12)). Thus

in (26)–(27). Moreover, while the rescaling (28) leads to a modified balance, i.e. (29) in place of (11) and (12), the result (30) implies that the solution passes unscathed through , with (11)–(12) being recovered for . Thus with the two-mode approximation (6) and (8) continuation through blow up seems to be straightforward, in contrast to that of the -equation (3), as we shall subsequently demonstrate.

a.2 Truncated Taylor expansion

Much of what follows in this subsection revisits results from [12], which we derive using a different approach (namely, matched asymptotic expansions).

As in section 2, we also consider the truncated Taylor approximation (5), which gives rise to the system (18) whose blow-up behaviour takes the form (again, by a self-consistency argument)


for constants and . In contrast to the results of the previous subsection, (18) does capture the blow-up behaviour of the full PDE (3), a point to which we shall also come back.

We now return to initial conditions of the form (9), with the scalings (10) applying for , so expanding in the form

(18) implies


Under the scalings (28) the leading-order balances do not change and the results of relevance below read


the matching into in (32) leading to the contributions. The final scale in this case is then more subtle than those above: we set




The introduction of in (34) is associated with the expansion for in (33) disordering and (18) becomes


so that, on matching with (33),


so that in (31).

For reasons that will become apparent below, it is instructive to record the implications of (37) under the scaling



this will reappear in the analysis below in describing the local blow-up behaviour, but also confirms that (5) cannot describe the spatial profile at blow up (i.e. at ).

a.3 More general near-flat initial data

We now turn to the derivation of the blow-up behaviour, subjecting (3), the full PDE, to the initial data

with, again, . It is striking that this limit allows a near complete analytical description of the transition to blow up through the fully nonlinear regime. Three time scales are required. On the first time scale:


All three terms in this expansion are required for what follows. An immediate result is that

leading to


We shall consider general , but in the case of real-valued -periodic initial conditions (or zero Neumann boundary conditions on a finite domain), we can Fourier decompose in the usual way, so that

We note this special case for two reasons – firstly for its relevance to section A.1 and secondly for the obvious observation that the high-frequency modes are rapidly decaying, providing additional motivation for the analysis of sections 2 and A.1 but being also in some respects counter-intuitive, given that point blow up subsequently ensues.

We define

so that


At next order we have


so that


where (42) serves to define .

On the second time scale:


Here , with and with correction terms to identified below444The notation here differs somewhat from that of the previous subsections.. We shall also need to consider the rescaling (38). We have

so, matching with (40), (42), and defining and via


it follows for that


and hence



We take blow up to occur at , requiring that and take


for constants , with 555The case corresponds to non-generic forms of blow up – we shall not pursue such matters here. and where the requirement that at , implies that in (43) are specified by (47), with being determined by the linear problems (39) and (41).

For , we first generate the required matching conditions from (40), (42) and (47), whereby


wherein we have retained only the required terms in and – in general will also lead to contributions of the form

in but these can be neglected for our purposes, being sub-dominant as with .

We have


and from this we find that the relevant terms for simply reproduce the matching condition (48) obtained from the expansion for . Importantly, the expansion (48) disorders for , as in section A.2, leading us on to our final time scale, as follows:




We emphasise that much of what follows reconstructs known blow-up behaviour (see the review of [8], for example.) The novelty here lies in the focus on the near-flat initial data, which allows more detailed characterization of the full spatial behaviour.

The near self-similar solution ansatz (50) transforms (49) to

so that


implying that




but satisfies (following cancellation of a number of terms using (51))


Since we need to preclude exponential growth of (i.e. to exclude contributions with as ), (54) both requires that


and generates the solvability condition


that the term in is only fully determined via (54) necessitates that the expansion be taken up to . Matching with (48) then requires that


simultaneously confirming matching with the second, fourth and final terms in (48). The calculation of and requires solvability conditions at yet higher orders (which we will not pursue), with (48) requiring that

By analogy with the subdivision into and in (II), we need also to consider , where

though a third scale (namely ) will also require consideration here. Since


gives on matching into



so that, again matching into ,


were not a solution to (56), a term in would also be present here, further clarifying the status of (56) as a solvability condition.

Finally, we can simply set in (44)–(46), noting that corresponds to exponentially small , to obtain the profile at blow up for . Thus as , with we have


The expression (60) has small- behaviour