Numerical solutions of the flexure equation

Abstract. The 4th order differential equation describing elastic flexure of the lithosphere is one of the cornerstones of geodynamics, key to understanding topography, gravity, glacial isostatic rebound, foreland basin evolution and a host of other phenomena. Despite being fully formulated in the 1940’s, a number of significant issues concerning the basic equation have remained overlooked to this day. We first explain the different fundamental forms the equation can take and their difference in meaning and solution procedures. We then show how numerical solutions to flexure problems in general as they are currently formulated, are potentially unreliable in an unpredictable manner for cases where the coefficient of rigidity varies in space due to variations of the elastic thickness parameter. This is due to fundamental issues related to the numerical discretisation scheme employed. We demonstrate an alternative discretisation that is stable and accurate across the broadest conceivable range of conditions and variations of elastic thickness, and show how such a scheme can simulate conditions up to and including a completely broken lithosphere more usually modelled as an end loaded, single, continuous plate. Importantly, our scheme will allow breaks in plate interiors, allowing for instance, the creation of separate blocks of lithosphere which can also share the support of loads. The scheme we use has been known for many years, but remains rarely applied or discussed. We show that it is generally the most suitable finite difference discretisation of fourth order, elliptic equations of the kind describing many phenomena in elasticity, including the problem of bending of elastic beams. We compare the earlier discretisation scheme to the new one in 1 dimensional form, and also give the 2 dimensional discretisation based on the new scheme.We also describe a general issue concerning the numerical stability of any second order finite difference discretisation of a fourth order differential equation like that describing flexure where contrasting magnitudes of coefficients of different summed terms lead to round off problems which in turn destroy matrix positivity. We explain the use of 128 bit, floating point storage for variables to mitigate this issue.

[1]  J. Braun,et al.  Flexure of the lithosphere and the geodynamical evolution of non-cylindrical rifted passive margins: Results from a numerical model incorporating variable elastic thickness, surface processes and 3D thermal subsidence , 2013 .

[2]  R. Walcott,et al.  Flexural rigidity, thickness, and viscosity of the lithosphere , 1970 .

[3]  D. McKenzie,et al.  Estimates of the effective elastic thickness of the continental lithosphere from Bouguer and free air gravity anomalies , 1997 .

[4]  P. Audet,et al.  Variations in elastic thickness in the Canadian Shield , 2003 .

[5]  B. R. Baliga,et al.  SOME EXTENSIONS OF TRIDIAGONAL AND PENTADIAGONAL MATRIX ALGORITHMS , 1995 .

[6]  Walter H. F. Smith,et al.  The Generic Mapping Tools Version 6 , 2019, Geochemistry, Geophysics, Geosystems.

[7]  P. DeCelles,et al.  The modern foreland basin system adjacent to the Central Andes , 1997 .

[8]  Anthony B. Watts.  The effective elastic thickness of the lithosphere and the evolution of foreland basins , 1992 .

[9]  R. Gunn.  A quantitative study of isobaric equilibrium and gravity anomalies in the hawaiian islands , 1943 .

[10]  B. Chen,et al.  Reconsidering Effective Elastic Thickness Estimates by Incorporating the Effect of Sediments: A Case Study for Europe , 2018, Geophysical Research Letters.

[11]  A. B. WATTS,et al.  Isostasy and Flexure of the Lithosphere , 2001 .

[12]  S. Buiter,et al.  An early Pliocene uplift of the central Apenninic foredeep and its geodynamic significance , 2000 .

[13]  W. Krijgsman,et al.  Regional isostatic response to Messinian Salinity Crisis events , 2009 .

[14]  H. Jeffreys.  The Strength of the Earth's Crust , 1914, The Journal of Geology.

[15]  J. Kley,et al.  The Subhercynian Basin: An example of an intraplate foreland basin due to a broken plate , 2020, Solid Earth.

[16]  A. Watts,et al.  The long-term strength of Europe and its implications for plate-forming processes , 2005, Nature.

[17]  A. Watts,et al.  the long-term strength of continental lithosphere: "jelly sandwich" or "crème brûlée"? , 2006 .

[18]  J. Stewart,et al.  Gravity anomalies and spatial variations of flexural rigidity at mountain ranges , 1997 .

[19]  R. Gunn.  Quantitative aspects of juxtaposed ocean deeps, mountain chains, and volcanic ranges , 1947 .

[20]  K. Furlong,et al.  Tectonic loading and subsidence of intermontane basins: Wyoming foreland province , 1985 .

[21]  R. Gunn.  A quantitative evaluation of the influence of the lithosphere on the anomalies of gravity , 1943 .

[22]  P. DeCelles,et al.  Foreland basin systems , 1996 .

[23]  R. Walcott.  Late Quaternary vertical movements in eastern North America: Quantitative evidence of glacio-isostatic rebound , 1972 .

[24]  J. Wees,et al.  A finite-difference technique to incorporate spatial variations in rigidity and planar faults into 3-D models for lithospheric flexure , 1994 .

[25]  M. Steckler,et al.  Observations of flexure and the rheology of the oceanic lithosphere , 1981 .

[26]  S. S. Egan,et al.  The flexural isostatic response of the lithosphere to extensional tectonics , 1992 .

[27]  C. Beaumont,et al.  Mechanical models of tilted block basins , 2013 .

[28]  R. Walcott.  Isostatic response to loading of the crust in Canada , 1970 .