Flux Gradient Research Articles

A flux correction for finite-volume discretizations: Achieving second-order accuracy on arbitrary polyhedral grids

In this paper, we propose a flux correction technique generally applicable to practical finite-volume discretizations of a single flux evaluation per face for achieving second-order accuracy on arbitrary polyhedral grids involving non-planar faces. The proposed technique is derived from the k-exact finite-volume discretization approach originally introduced by Brenner. We take it as a general methodology for constructing a second-order discretization on arbitrary polyhedral grids, and identify a term considered as missing in other practical finite-volume discretizations; it is thus considered as a reason for the loss of second-order accuracy in many of such methods: e.g., a cell-centered discretization loses second-order accuracy on grids with non-planar faces as typical in polyhedral grids. In particular, we write the term as a vector involving flux gradients and show that it can be easily implemented as a correction to a numerical flux at each face. Also, we will show that a consistent definition of a control volume is critically important for achieving second-order accuracy on general polyhedral grids. We then discuss its general applicability and impact on accuracy and efficiency with simple but illustrative examples: the edge-based solver is made to achieve faster iterative convergence on irregular tetrahedral grids with adjustable centroids while preserving second-order accuracy; both the edge-based and cell-centered discretizations are made second-order accurate on irregular prismatic grids having non-planar faces.

Modeling Two-Phase Flow Caused by Hydrate Dissociation in a Reservoir on the Concept of Global Pressure

Summary The classic Darcy’s two-phase flow equation has a variety of mathematically equivalent formats, such as pressure/saturation (PS), two-phase pressure (PP), and global pressure/saturation (GP). Based on the concept of global pressure, we derive a new formulation for two-phase flow caused by hydrate dissociation in porous media, coupled with the mass and energy conservation equations, thereby offering a novel theoretical frame for gas hydrate extraction simulation. The new model provides in-depth insights into complex flow fields: The dissociation of hydrates under thermal stimulation may lead to two flow fields of gas and water with different directions in the reservoir. The simulation results are in good agreement with the recorded data set of Masuda’s and Chong’s experiments, which verifies the correctness and applicability of the new model. Based on numerical simulations of the same hydrate dissociation experiments, the efficiency of the GP method was compared with the PP method published in our previous work. The results show that the GP method has more than two times the timestep size of the PP method for unsteady flow and 1.5 times for steady flow. Furthermore, the relative residual of the GP method is even two orders of magnitude lower than that of the PP method for two-phase flows with complex variations of pressure gradient and fluid flux. Therefore, the GP method is significantly more efficient than the PP method for simulating gas hydrate development. The proposed GP method improves the computational efficiency of hydrate extraction simulation at the laboratory scale and understanding the mechanisms of two-phase flow in reservoirs caused by hydrate dissociation. It may have potential advantages for field-scale simulation of hydrate development, which requires further studies to demonstrate.

PCA Benchmark Analysis with ADVANTG3.0.1. and MCNP6.1.1b Codes

The Pool Critical Assembly Pressure Vessel (PCA) benchmark is a well known benchmark in the reactor shielding community which is described in the Shielding Integral Benchmark Archive and Database (SINBAD). It is based on the experiments performed at the PCA facility in the Oak Ridge National Laboratory (ORNL) and it can be used for the qualification of the pressure vessel fluence calculational methodology. The measured quantities to be compared against the calculated values are the equivalent fission fluxes at several experimental access tubes (A1 to A8) in front, behind, and inside the pressure-vessel wall simulator. This benchmark is particularly suitable to test the capabilities of the shielding calculational methodology and cross-section libraries to predict invessel flux gradients because only a few approximations are necessary in the overall analysis. This benchmark was analyzed using a modern hybrid stochastic-deterministic shielding methodology with ADVANTG3.0.1 and MCNP6.1.1b codes. ADVANTG3.0.1 is an automated tool for generating variance reduction (VR) parameters for Monte Carlo (MC) calculations with MCNP5v1.60 code (and higher versions). It is based on the multigroup, discrete ordinates solver Denovo, used for approximating the forward-adjoint transport fluxes to construct VR parameters for the final MC simulation. The VR parameters in form of the weight windows and the source biasing cards can be directly used with unmodified MCNP input. The underlining CADIS methodology in Denovo code was initially developed for biasing local MC results, such as point detector or a limited region detector. The FW-CADIS extension was developed for biasing MC results globally over a mesh tallies or multiple point/region detectors. Both CADIS and FW-CADIS are based on the concept of the neutron importance function, which is a solution of the adjoint Boltzmann transport equation. The equivalent fission fluxes calculated with MCNP are based on several highenergy threshold reactions from international dosimetry libraries IRDF-2002 and IRDFF-2014, distributed by the IAEA Nuclear Data Section. The obtained results show a good agreement with referenced PCA measurements. Visualization of the deterministic solution in 3D was done using the VisIt code from the Lawrence Livermore National Laboratory (LLNL).

Flux–Gradient Relationships Below 2 m Over a Flat Site in Complex Terrain

The surface–atmosphere turbulent exchange fluxes are experimentally determined using the eddy-covariance method. Their estimation using profiles of the variables of interest is a less costly alternative, although restricted to certain ranges of stability and assumed to hold for relatively flat and homogeneous terrain. It relays usually on the prescription of the roughness lengths for momentum, heat and matter, the latter two being adjustable parameters with unclear physical significance. The relations are derived with data from screen level to a few tens of metres upward. The application of these expressions using data only at one level in the surface layer implies assuming zero wind speed and the land surface temperature at their respective roughness lengths. The latter is a quantity that experimentally can only be determined radiatively with a substantial uncertainty. In this work the flux-profile relationships for momentum and sensible heat are assessed over a flat site in moderately inhomogeneous complex terrain in the southern pre-Pyrenees, using data between 2 m and the surface. The main findings are that (i) the classical expressions hold in the daytime for most of the dataset, (ii) the iterative estimations using the Obukhov length and the direct ones using the bulk Richardson number provide very similar results, (iii) using a second observation of temperature avoids a radiometric measure of land surface temperature and the prescription of a thermal roughness length value, (iv) the estimations over wet terrain with high irradiance depart largely from observations.

A nonlinear, non-dispersive energy balance for surfzone waves: infragravity wave dynamics on a sloping beach

A fully nonlinear non-dispersive energy balance for surfzone waves is derived based on the nonlinear shallow water equations to study the nearshore dynamics of infragravity (IG) waves. Based on simulations of waves on a relatively moderate and mild beach slope with a non-hydrostatic wave-flow model (SWASH), the new theory shows that spatial gradients in IG energy flux are nearly completely balanced by the combined effect of bottom stresses and predominantly nonlinear triad interactions. The new balance confirms many features of existing weakly nonlinear theories, and yields an improved description in the inner surfzone where waves become highly nonlinear. A gain of IG energy flux throughout the shoaling and outer surfzones is driven by triad interactions between IG waves and pairs of sea-swell (SS) waves. The IG energy flux decreased in the inner surfzone, primarily through an energy cascade to the swell-band and superharmonic frequencies where wave energy is ultimately dissipated. Dissipation by bottom friction was weak on both slopes. The IG wave breaking, characterized by triads between three IG or two IG waves and one SS wave, was significant only deep inside the surfzone of the mild slope. Even though IG waves broke on the mild slope, nonlinear interactions between IG waves and pairs of SS waves were responsible for at least half of the net IG flux loss.

High-order well-balanced and positivity-preserving finite-difference AWENO scheme with hydrostatic reconstruction for shallow water equations

The shallow water equations (SWEs) admit still water steady-state solutions in which the flux gradients are exactly balanced by the source term. Furthermore, the no-water dry areas make the numerical simulation of the SWEs more difficult since the water height often loses positivity near the dry areas. In this study, we design a high-order well-balanced and positivity-preserving finite-difference AWENO scheme for solving the SWEs with or without dry areas. In this numerical framework, the hydrostatic reconstruction (HR) method is used with two main advantages: 1) The well-balanced property is achieved using the HR method with arbitrary monotone fluxes and the reformulated source terms. 2) The first-order scheme with the Lax-Friedrichs (LF) flux and the HR method maintains the positivity of the water height with an appropriate time step. Therefore, a simple positivity-preserving limiter conjugating the high-order reconstructed flux with the first-order positivity-preserving LF flux is introduced. One- and two-dimensional classical examples with or without dry areas demonstrate that the proposed scheme is high-order, well-balanced, and positivity-preserving.

Influence of Low and Extreme Heat Fluxes on Thermal Degradation of Carbon Fibre-reinforced Polymers

This study considers the influence of different irradiation scenarios on the thermal degradation of carbon fibre-reinforced polymers (CFRP). Real threats are simulated, such as fires with long-lasting low heat fluxes and nuclear heat flashes with short-lasting high heat fluxes. For this purpose, coated and uncoated quasi-isotropic samples of the commercially available CFRP HexPly® 8552/IM7 are thermally irradiated from one side by an electrical heater of a cone calorimeter and a xenon short-arc lamp of a laboratory heat flash simulator with heat fluxes between 5 and 175 W/cm2 at varying time intervals. The specimens’ temperature is recorded on the front and back side as well as at different laminate depths. The CFRP are analyzed with ultrasonic testing (UT), infrared spectroscopy (ATR-FTIR), scanning electron microscopy (SEM) and micro-focused computed X-Ray tomography (μCT). Destructive tests are performed to determine the mechanical properties in terms of interlaminar shear, compressive and tensile strength. When samples of CFRP are exposed to higher heat flux, high temperatures and temperature gradient values occur along the cross-section. As a result, extreme damage gradients appear in the material, leading to changes in damage behavior and loss of mechanical properties within seconds. However, to ensure the safety of the material in case of thermal exposure, loading limits are introduced, indicating the threshold for strength collapse. In addition, with the application of coatings, thermal degradation of CFRP can be delayed. Finally, the time-heat flux superposition principle is established to predict the residual strength under different loading scenarios.

Multilevel mesh adaptivity for discrete ordinates transport calculation with spatial-moment-ratio indicators

It is difficult to control discretization errors and reduce computational cost at the same time for multi-scale neutron transport problems, and the adaptive mesh refinement technique is one of the powerful and effective methods for high-resolution transport calculation. We propose a multilevel mesh adaptivity algorithm and develop a discrete ordinates calculation framework on 3-D Cartesian meshes. Data management and transport sweep are optimized based on the multilevel pyramid data structure containing pointer-type and array-type variables. A new spatial-moment-ratio error indicator is derived to measure the leading term of distribution functions and to drive the local mesh refinement. The numerical properties of spatial moment factors are better than the normalized gradient of angular flux. Compared with uniform refinement, our adaptivity algorithm reaches the same accuracy level with a computational cost saving of 60%-80% for heterogeneous problems.

Secondary Reflector and Receiver Positions for Uniform Heat Flux Distribution in Parabolic Trough Solar Thermal Collector

Abstract The distributions of heat flux over the circumference of the receiver tubes have an immense influence on the performance and reliability of the parabolic trough solar thermal collectors. The location of the receiver tube and the secondary reflector configuration may largely influence the performance of the system. Therefore, in this study, the effect of receiver tube position and parabolic secondary reflector configuration has been analyzed, and the non-uniformity of solar flux distribution, heat gradient, and power output has been compared. The results of the Monte Carlo ray-tracing analysis to homogenize the receiver tube flux distribution and maximize the output power, making the use of the cutting edge solar optical simulation tool Tonatiuh, has been presented. A parabolic trough collector with a rim angle of 80 deg and aperture area of 40 m2 have been used for the analysis. It has been confirmed that the circumferential heat flux gradient and the local hot spot could be greatly diminished, while the power output tended to reduce slightly due to the shading effect of the secondary reflector. Under the conditions investigated in this work, although the output power decreased by 4.83%, flux gradient reduced significantly, and the non-uniformity of flux distribution has reduced from 0.9757 to 0.5176. A simple design procedure for receiver tube position and secondary reflector configurations to homogenize the receiver tube temperature distribution has also been proposed.

Detection of Uniaxial Fatigue Stress under Magnetic Flux Leakage Signals using Morlet Wavelet

This paper demonstrates the application of continuous wavelet transform technique for magnetic flux leakage signal generated during a uniaxial fatigue test. This is a consideration as the magnetic signal is weak and susceptible to being influenced by an external magnetic field. The magnetic flux leakage signal response of API steel grade X65 is determined using Metal Magnetic Memory under cyclic load conditions ranging from 50% to 85% of the UTS. To facilitate further signal analysis, the magnetic flux gradient, the dH(y)/dx signal were converted from a length base into time series in this study. Magnetic flux leakage readings indicated a maximum UTS load of 56.5 (A/m)/mm at 85%, where a higher load resulted in a higher reading and the signal contained Morlet wavelet coefficient energy of 1.02×106 µe2/Hz. As increasing percentages of UTS loads were applied, the signal analysis revealed an increasing linear trend in the dH(y)/dx and wavelet coefficient energy. The analysis revealed a strong correlation between the wavelet coefficient energy and the dH(y)/dx amplitude, as indicated by the coefficient of determination (R2) value of 0.8572. Hence, this technique can provide critical information about magnetic flux leakage signals that can be used to detect high stress concentration zones.

The analysis and computation on nonlocal thermoelastic problems of blend composites via enriched second-order multi-scale computational method

This paper proposes an innovative enriched second-order multi-scale (SOMS) computational method to simulate and analyze the nonlocal thermoelastic problems of blend composites with stress and heat flux gradient behaviors. The multiple periodical heterogeneities and periodic configurations of investigated blend composites in different substructures result in a huge computational cost for direct numerical simulations. The significant characteristics of this study are as follows. (1) The nonlocal properties of blend composites in constitutive equations are converted into the source terms of thermoelastic balance equations. The novel macro-micro coupled SOMS computational model for these transformed nonlocal multi-scale problems is derived on the basis of multi-scale asymptotic analysis. The nonlocal thermoelastic behaviors of blend composites can be merely uncovered in the enriched SOMS solutions. (2) The error analysis in the pointwise sense is presented to elucidate the importance and necessity of establishing the enriched SOMS solutions. Furthermore, an explicit error estimate for the SOMS approximate solutions is obtained in the integral sense for these nonlocal multi-scale problems. (3) A multi-scale numerical algorithm is presented to effectively simulate nonlocal thermoelastic problems of blend composites based on finite element method (FEM). Finally, the capability of the proposed enriched SOMS computational method is demonstrated by typical two-dimensional (2D) and three-dimensional (3D) blend composites, presenting not only the excellent numerical accuracy but also the less computational cost. This work proposes a unified multi-scale computational framework for enabling nonlocal thermoelastic behavior analysis of blend composites.

Permeability of Polymer Membranes beyond Linear Response

The permeability of a membrane to solute penetrants is well defined on the linear response level simply as the ratio of penetrants’ flux and concentration gradient at the membrane boundary layers. However, nonlinearities emerge in the flux–force relation j(f) for large driving forces f, in which the definition of permeability becomes ambiguous. Here, we study nonequilibrium membrane permeation orchestrated by a generic driving force using penetrant- and monomer-resolved computer simulations of transport in a polymer network, supported by exact solutions of the Smoluchowski (drift–diffusion) equation in the stationary state. In the simulations, we consider the transport across a finite polymer membrane immersed in a reservoir of penetrants, addressing one- and two-component penetrant systems. We calculate the f-dependent inhomogeneous steady-state density profiles, boundary layer concentrations, and fluxes of the penetrants. The Smoluchowski approach, using solely coarse-grained equilibrium partitioning and diffusion profiles as input, exhibits remarkable qualitative agreement with our nonequilibrium simulations, which serves for rationalization of the observations. We discuss possible definitions of nonequilibrium, f-dependent permeability, distinguishing between “system” and “membrane” permeabilities. In particular, we introduce the concept of differential permeability as a response to f. The latter turns out to be a highly nonmonotonic function of f for low-permeable systems, demonstrating how a differential permselectivity is substantially tunable by the driving force beyond linear response.

A Stable Condition and Adaptive Diffusion Coefficients for the Coarse-Mesh Finite Difference Method

Coarse-mesh finite difference (CMFD) method is a widely used numerical acceleration method. However, the stability of CMFD method is not good for the problems with optically thick regions. In this paper, a stability rule named the “sign preservation rule” in the field of numerical heat transfer is extended to the scheme of CMFD. It is required that the disturbance of neutron current is positively correlated with that of the negative value of flux gradient. A necessary condition for stability of the CMFD method is derived, an adaptive diffusion coefficient equation is proposed to improve the stability of CMFD method, and the corresponding revised CMFD method is called the rCMFD method. With a few modifications of the code, the rCMFD method was implemented in the hexagonal-Z nodal SN (discrete-ordinates) solver in the NECP-SARAX code system. The rCMFD method and other similar acceleration methods were tested by three fast reactor problems which were obtained by modifying the hexagonal pitches of a benchmark problem. The numerical results indicated that the rCMFD method showed better stability than the traditional CMFD method and the artificially diffusive CMFD (adCMFD) method and a better convergence rate than the adCMFD method and the optimally diffusive CMFD (odCMFD) method for these fast reactor problems.

High order well-balanced conservative finite difference AWENO scheme with hydrostatic reconstruction for the Euler equations under gravitational fields

The system of compressible Euler equations under gravitational fields drew much attention during the past few decades, due to its broad applications in astrophysics and atmospheric sciences. A critical feature of the system is that it admits a hydrostatic equilibrium state when the gradient of fluid flux is exactly balanced by the source term brought by gravitational fields, which is referred to as well-balancedness. In this paper, a well-balanced high order conservative finite difference alternative WENO (AWENO) scheme is proposed for the system. The finite difference AWENO scheme allows us to use an arbitrary consistent Riemann solver, such as LF, HLLC, etc. To preserve hydrostatic equilibrium, the hydrostatic reconstruction technique is adopted. We slightly modify the original hydrostatic reconstruction procedure to achieve the fifth order of accuracy in the framework of finite difference methods, without introducing spurious oscillations. Besides, the discretization is fully conservative without the source term, which is not the case for the methods with original hydrostatic reconstructions. Moreover, we will show that the well-balancedness is exactly guaranteed by our approach at the discrete level. Finally, numerical examples illustrate the effectiveness and robustness of the designed approach.

Thermoelastic response of skin using time-fractional dual-phase-lag bioheat heat transfer equation

In this paper, a time-fractional dual-phase-lag (DPL) bioheat heat conduction model is presented to treat thermoelastic response of skin subjected to sudden temperature shock. The skin tissue is treated as a monolayered thermoelastic material. A modified Pennes bioheat transfer equation with non-Fourier effect is introduced to capture the influence of blood perfusion rate in the skin tissue on temperature distribution. The transient response of the temperature field of skin tissue is derived analytically and calculated numerically. The mechanical response is considered and the elastic displacement induced by temperature change is evaluated. The effects of two relaxation times and fractional orders of heat flux and temperature gradient together with blood flow on the temperature distribution in the skin tissue are analyzed.

Controlling conservation laws II: Compressible Navier–Stokes equations

We propose, study, and compute solutions to a class of optimal control problems for hyperbolic systems of conservation laws and their viscous regularization [17]. We take barotropic compressible Navier–Stokes equations (BNS) as a canonical example. We first apply the entropy–entropy flux–metric condition for BNS. We select an entropy function and rewrite BNS to a summation of flux and metric gradient of entropy. We then develop a metric variational problem for BNS, whose critical points form a primal-dual BNS system. We design a finite difference scheme for the variational system. The numerical approximations of conservation laws are implicit in time. We solve the variational problem with an algorithm inspired by the primal–dual hybrid gradient method. This includes a new method for solving implicit time approximations for conservation laws, which seems to be unconditionally stable. Several numerical examples are presented to demonstrate the effectiveness of the proposed algorithm.

On the use of sacrificial layer in ECDM process for form accuracy

During the electrochemical discharge machining (ECDM) process, major imperfections on work material occurs due to stray cutting caused by thermal energy. Some of these are re-solidification, melt flow over, thermal splashes, micro-cracks, work material breakage, hole over cut (HOC), and sudden taper. These thermal deteriorations are inherent in thermal-based processes like the ECDM process. Post-processing of these machined features at micron-level on glass material is a challenging task. To overcome this issue, in the present article, the depth up to which thermal damages occurred was predicted at given input energy conditions and the same thickness of glass wafer was bonded over the work material during machining. In order to find out the thickness, a pressurized tool feeding and dial gauge attachment was designed. Using this attachment, instantaneous depth was measured at every second of machining time and effective energy supplied to each plane of the work material was calculated analytically. The same was imported in to COMSOL Multiphysics software for simulating the heat flux gradient and temperature distribution on each plane of the work material. Using the simulation data, the depths from the hole entrance that are more prone to thermal deteriorations was identified. This depth was considered as critical depth for respective input energy conditions. Applied voltage, pulse on time, and glass wafer thickness were considered as input parametric conditions, and machined hole depth, HOC, and hole quality on work material were considered as output characteristics. Experiments were performed with various parametric conditions, and it was observed that the glass wafer carried away the major thermal deteriorations and leaving minimal damages on the work material. At the optimal parametric conditions, HOC was reduced by 71.5% and hole quality was improved significantly. Application of glass wafer compromises around 39% hole depth.

Optimization of a vibrating MEMS electromagnetic energy harvester using simulations

Nowadays, wireless sensor networks (WSN) are becoming essential in our daily life. However, a major constraint concerns the energy power supply. Indeed, batteries need to be recharged or replaced often which implies a limited lifetime for WSN nodes. One alternative consists in harvesting mechanical energy from surrounding vibrations of the environment. Using finite element simulations, we report here a complete guideline to optimize a MEMS electromagnetic energy harvester consisting of an in-plane vibrating silicon frame supporting an array of micromagnets that faces a static 2D micro-coil. The dimensioning of the magnet array and the specific design of the coil are addressed, considering patterned 50 \(\upmu \)m thick NdFeB films with out of plane magnetization. The optimization of the electromechanical coupling which allows to efficiently convert the energy results from a trade-off between the high magnetic flux gradients produced by the micromagnets and the maximum number of turns of the facing coil.

Trapped magnetic field distribution above a superconducting linear Halbach array

In applications requiring a large magnetic force, permanent magnets with non-parallel magnetization directions can be assembled in a Halbach array to generate a large gradient of magnetic flux density. The saturation magnetization of permanent magnets, however, brings a fundamental limit on the performance of this configuration. In the present work, we investigate experimentally the assembly of cuboid bulk, large grain melt-textured YBa2Cu3O superconductors ( mm3) with orthogonal c-axes so as to form a basic unit of Halbach array. The experiments are carried out at 77 K. The experimental distribution of the magnetic flux density above the array of trapped-field superconductors is compared to a similar array made of permanent magnets. A simple analytical model is developed and is shown to accurately reproduce the main experimental observations. The results suggest that a redistribution occurs in the current flowing in the central sample when the distance between the superconductors is reduced, whereas the neighbouring superconductors are unaffected. It is shown that this current redistribution yields a reduced contribution of the central sample to the magnetic flux density above the centre of the array and a new negative contribution associated with stray fields to the magnetic flux density at this location. This interpretation is confirmed by modelling of the distribution of transport currents in the superconductor using a 3D finite element model.

The direction of landscape erosion

Abstract. The rate of erosion of a landscape depends largely on local gradient and material fluxes. Since both quantities are functions of the shape of the catchment surface, this dependence constitutes a mathematical straitjacket, in the sense that – subject to simplifying assumptions about the erosion process, and absent variations in external forcing and erodibility – the rate of change of surface geometry is solely a function of surface geometry. Here we demonstrate how to use this geometric self-constraint to convert a gradient-dependent erosion model into its equivalent Hamiltonian, and explore the implications of having a Hamiltonian description of the erosion process. To achieve this conversion, we recognize that the rate of erosion defines the velocity of surface motion in its orthogonal direction, and we express this rate in its reciprocal form as the surface-normal slowness. By rewriting surface tilt in terms of normal slowness components and deploying a substitution developed in geometric mechanics, we extract what is known as the fundamental metric function of the model phase space; its square is the Hamiltonian. Such a Hamiltonian provides several new ways to solve for the evolution of an erosion surface: here we use it to derive Hamilton's ray-tracing equations, which describe both the velocity of a surface point and the rate of change of the surface-normal slowness at that point. In this context, gradient-dependent erosion involves two distinct directions: (i) the surface-normal direction, which points subvertically downwards, and (ii) the erosion ray direction, which points upstream at a generally small angle to horizontal with a sign controlled by the scaling of erosion with slope. If the model erosion rate scales faster than linearly with gradient, the rays point obliquely upwards, but if erosion scales sublinearly with gradient, the rays point obliquely downwards. This dependence of erosional anisotropy on gradient scaling explains why, as previous studies have shown, model knickpoints behave in two distinct ways depending on the gradient exponent. Analysis of the Hamiltonian shows that the erosion rays carry boundary-condition information upstream, and that they are geodesics, meaning that surface evolution takes the path of least erosion time. Correspondingly, the time it takes for external changes to propagate into and change a landscape is set by the velocity of these rays. The Hamiltonian also reveals that gradient-dependent erosion surfaces have a critical tilt, given by a simple function of the gradient scaling exponent, at which ray-propagation behaviour changes. Channel profiles generated from the non-dimensionalized Hamiltonian have a shape entirely determined by the scaling exponents and by a dimensionless erosion rate expressed as the surface tilt at the downstream boundary.