API Reference

Primary solvers

Both solvers follow the same lifecycle:

flex.initialize()
flex.run()
w = flex.w          # read deflection before finalize clears it
flex.output()       # optional: save to file or display plots
flex.finalize()     # releases w, qs, and the coefficient matrix

Warning

finalize() deletes w, qs, and the cached coefficient matrix. Read w (and call output() if needed) before calling finalize. Accessing w afterwards raises AttributeError.

class gflex.F2D(filename=None)[source]

Bases: Flexure

Two-dimensional lithospheric flexure solver.

Computes the deflection w(x, y) of a thin elastic plate overlying an inviscid fluid (mantle) given a surface load stress qs. Supports spatially variable elastic thickness Te.

Set instance attributes, then call initialize(), run(), and finalize() in sequence. Read flex.w before calling finalize(); finalize clears all model state including w.

method

Solution method. 'fd' (finite difference, supports variable Te), 'fft' (spectral, requires scalar Te; 2-D only), 'sas' (superposition of analytical solutions, constant Te only), or 'sas_ng' (SAS on an ungridded point cloud).

Type:

str

solver

Linear solver: 'direct' (sparse LU, default).

Type:

str

g

Gravitational acceleration [m s⁻²].

Type:

float

E

Young’s modulus [Pa].

Type:

float

nu

Poisson’s ratio.

Type:

float

rho_m

Mantle density [kg m⁻³].

Type:

float

rho_fill

Infill material density [kg m⁻³] (0 for air, ~1000 for water, ~2700 for rock).

Type:

float

T_e

Elastic thickness [m]. A scalar is broadcast to the full grid.

Type:

float or ndarray of shape (M, N)

qs

Surface load stress [Pa].

Type:

ndarray of shape (M, N)

dx

Grid spacing in the x (column) direction [m].

Type:

float

dy

Grid spacing in the y (row) direction [m].

Type:

float

bc_west, bc_east, bc_north, bc_south

Boundary conditions on the west, east, north, and south edges. FD options: 'zero_displacement_zero_slope' (alias 'clamped'), 'zero_displacement_zero_moment', 'zero_moment_zero_shear' (alias 'free'), 'zero_slope_zero_shear' (alias 'mirror'), 'periodic'. SAS option: 'no_outside_loads' (the default when unset).

Type:

str

sigma_xx

Normal in-plane stress in the x-direction \(\sigma_{xx}\) [Pa]. Supported by fd and fft. Default 0.

Type:

float, optional

sigma_yy

Normal in-plane stress in the y-direction \(\sigma_{yy}\) [Pa]. Supported by fd and fft. Default 0.

Type:

float, optional

sigma_xy

In-plane shear stress \(\sigma_{xy}\) [Pa]. Supported by fd and fft. Default 0.

Type:

float, optional

quiet

Suppress timing output. Default False.

Type:

bool

verbose

Print progress messages. Default True.

Type:

bool

Examples

Minimal finite-difference run:

import numpy as np
from gflex import F2D

flex = F2D()
flex.quiet = True
flex.method = 'fd'
flex.solver = 'direct'
flex.g = 9.8
flex.E = 65e9
flex.nu = 0.25
flex.rho_m = 3300.
flex.rho_fill = 0.
flex.T_e = 30e3 * np.ones((50, 50))
flex.qs = np.zeros((50, 50))
flex.qs[20:30, 20:30] = 1e6   # 10 × 10 cell load
flex.dx = flex.dy = 5000.     # 5 km grid
flex.bc_west = flex.bc_east = flex.bc_south = flex.bc_north = 'zero_moment_zero_shear'
flex.initialize()
flex.run()
flex.finalize()
deflection = flex.w           # (50, 50) array, negative downward
initialize(filename=None)[source]

Validate inputs and prepare the solver.

Must be called once before run(). If a configuration-file path was passed to the constructor (or to this method), parameters are read from that file; otherwise they are taken from the instance attributes set by the caller.

Parameters:

filename (str, optional) – Path to a gFlex YAML configuration file. Overrides any filename supplied to the constructor.

run()[source]

Execute the flexural solution.

Selects and runs the method specified by self.method. The deflection array is stored in self.w on return. Call finalize() afterwards to restore any internally modified state.

For repeated solves (e.g. a coupling loop), set cache_factorization = True before initialize() to reuse the LU factorisation when only qs changes; use "no_check" for maximum throughput when the coefficient matrix is guaranteed stable.

finalize()[source]

Release all model state.

Restores self.T_e to its pre-run value if gFlex padded it internally. Then calls the base finalize, which deletes w, qs, and the cached coefficient matrix. Read w before calling this method.

class gflex.F1D(filename=None)[source]

Bases: Flexure

One-dimensional lithospheric flexure solver.

Computes the deflection w(x) of a thin elastic beam overlying an inviscid fluid (mantle) given a surface load stress qs. Supports spatially variable elastic thickness Te.

Set instance attributes, then call initialize(), run(), and finalize() in sequence. Read flex.w before calling finalize(); finalize clears all model state including w.

method

Solution method. 'fd' (finite difference, supports variable Te), 'fft' (spectral, requires scalar Te), 'sas' (superposition of analytical solutions, constant Te only), or 'sas_ng' (SAS on an ungridded point array).

Type:

str

solver

Linear solver: 'direct' (sparse LU, default).

Type:

str

g

Gravitational acceleration [m s⁻²].

Type:

float

E

Young’s modulus [Pa].

Type:

float

nu

Poisson’s ratio.

Type:

float

rho_m

Mantle density [kg m⁻³].

Type:

float

rho_fill

Infill material density [kg m⁻³] (0 for air, ~1000 for water, ~2700 for rock).

Type:

float

T_e

Elastic thickness [m]. A scalar is broadcast to the full grid.

Type:

float or ndarray of shape (N,)

qs

Surface load stress [Pa].

Type:

ndarray of shape (N,)

dx

Grid spacing [m].

Type:

float

bc_west, bc_east

Boundary conditions on the west (left) and east (right) ends. FD options: 'zero_displacement_zero_slope' (alias 'clamped'), 'zero_displacement_zero_moment', 'zero_moment_zero_shear' (alias 'free'), 'zero_slope_zero_shear' (alias 'mirror'), 'periodic', 'sandbox'. SAS option: 'no_outside_loads' (the default when unset).

Type:

str

sigma_xx

Normal stress applied at the plate ends [Pa]. FD only.

Type:

float, optional

quiet

Suppress timing output. Default False.

Type:

bool

verbose

Print progress messages. Default True.

Type:

bool

Examples

Minimal finite-difference run:

import numpy as np
from gflex import F1D

flex = F1D()
flex.quiet = True
flex.method = 'fd'
flex.solver = 'direct'
flex.g = 9.8
flex.E = 65e9
flex.nu = 0.25
flex.rho_m = 3300.
flex.rho_fill = 1000.
flex.T_e = 30e3
flex.qs = np.zeros(300)
flex.qs[100:200] = 1e6      # 100-cell load
flex.dx = 4000.             # 4 km grid
flex.bc_west = 'zero_displacement_zero_slope'
flex.bc_east = 'zero_moment_zero_shear'
flex.initialize()
flex.run()
flex.finalize()
deflection = flex.w         # (300,) array, negative downward
initialize(filename=None)[source]

Validate inputs and prepare the solver.

Must be called once before run(). If a configuration-file path was passed to the constructor (or to this method), parameters are read from that file; otherwise they are taken from the instance attributes set by the caller.

Parameters:

filename (str, optional) – Path to a gFlex YAML configuration file. Overrides any filename supplied to the constructor.

run()[source]

Execute the flexural solution.

Selects and runs the method specified by self.method. The deflection array is stored in self.w on return. Call finalize() afterwards to restore any internally modified state.

For repeated solves (e.g. a coupling loop), set cache_factorization = True before initialize() to reuse the LU factorisation when only qs changes; use "no_check" for maximum throughput when the coefficient matrix is guaranteed stable.

finalize()[source]

Release all model state.

Restores self.T_e to its pre-run value if gFlex padded it internally. Then calls the base finalize, which deletes w, qs, and the cached coefficient matrix. Read w before calling this method.

Output

Flexure.output()[source]

Save deflection to file and/or plot, based on optional attributes.

Does nothing if neither w_out_file nor plot_choice has been set. Set w_out_file to a path ending in '.npy' for a binary NumPy array, or any other extension for an ASCII grid. Set plot_choice to 'q', 'w', 'both', or (1D only) 'combo' to display plots.

In-plane stresses

In-plane stresses are set as attributes directly on the solver instance before calling initialize(). They are not available as configuration file keys.

Attribute

Solvers

Description

sigma_xx

FD, FFT (1-D and 2-D)

Normal stress in the x-direction \(\sigma_{xx}\) [Pa]. Default 0.

sigma_yy

FD, FFT (2-D only)

Normal stress in the y-direction \(\sigma_{yy}\) [Pa]. Default 0.

sigma_xy

FD, FFT (2-D only)

Shear stress \(\sigma_{xy}\) [Pa]. Default 0.

All three default to zero if not assigned; setting any of them with SAS or SAS_NG raises a RuntimeWarning and has no effect. See Theory for the governing equations that include these terms.

Domain-padding utilities

These functions help when running F1D or F2D with a spatially variable elastic thickness grid. A smooth padding zone reduces spurious deflections at the domain boundary caused by sharp rigidity gradients, and ensures that the flexural forebulge can develop freely before reaching the boundary.

2-D (F2D)

gflex.pad_domain(Te, qs, dx, dy=None, n_wavelengths=1.0, Te_out=None, E=65000000000.0, nu=0.25, rho_m=3300.0, rho_fill=0.0, g=9.8)[source]

Pad both the elastic thickness and surface load arrays for use with F2D.

Combines recommended_pad_width() and smooth_pad_Te() into a single call, and zero-pads qs to match. The returned pad width p can be used to trim the deflection output after the run:

w_inner = flex.w[p:-p, p:-p]
Parameters:
  • Te ((M, N) array) – Elastic thickness [m] for the inner domain.

  • qs ((M, N) array) – Surface load [Pa] for the inner domain.

  • dx (float) – Grid cell size in the x-direction [m].

  • dy (float, optional) – Grid cell size in the y-direction [m]. Defaults to dx.

  • n_wavelengths (float, optional) – Padding width expressed as a number of flexural wavelengths. Default 1.0; use 0.5 for a less conservative (narrower) padding.

  • Te_out (float, optional) – Te value at the outer edge of the padding ring. Defaults to Te.mean().

  • E (float, optional) – Young’s modulus [Pa]. Default 65 GPa.

  • nu (float, optional) – Poisson’s ratio. Default 0.25.

  • rho_m (float, optional) – Mantle density [kg m^-3]. Default 3300.

  • rho_fill (float, optional) – Infill density [kg m^-3]. Default 0 (air).

  • g (float, optional) – Gravitational acceleration [m s^-2]. Default 9.8.

Returns:

  • Te_padded ((M + 2p, N + 2p) array) – Smoothly tapered elastic thickness.

  • qs_padded ((M + 2p, N + 2p) array) – Surface load zero-padded to match.

  • p (int) – Pad width in grid cells (same on all four sides).

Examples

>>> import numpy as np
>>> from gflex import pad_domain
>>> Te = 10e3 * np.ones((5, 5))
>>> qs = np.zeros((5, 5))
>>> Te_pad, qs_pad, p = pad_domain(Te, qs, dx=10000.,
...     E=65e9, nu=0.25, rho_m=3300., rho_fill=0., g=9.8)
>>> p
13
>>> Te_pad.shape
(31, 31)

To run F2D with the padded grids:

flex = gflex.F2D()
flex.T_e = Te_pad
flex.qs = qs_pad
# ... set other parameters and run ...
w_inner = flex.w[p:-p, p:-p]   # trim padding from output
gflex.smooth_pad_Te(Te, pad_width, Te_out=None)[source]

Pad a 2-D elastic thickness array with a smooth linear taper.

When a spatially variable Te grid is padded with a constant value before being passed to F2D, the abrupt step in flexural rigidity D at the inner/outer boundary drives spurious deflections via the D-derivative terms in the vWC1994 stencil (issue #45).

This function eliminates that step by linearly blending the inner-domain edge values toward Te_out across the padding ring, reducing the rigidity gradient at the inner/outer boundary by a factor of ~pad_width compared with an abrupt step.

The corresponding surface load array should be padded with zeros, e.g.:

qs_padded = numpy.pad(qs, pad_width, mode='constant')
Parameters:
  • Te ((M, N) array) – Elastic thickness [m] for the inner domain.

  • pad_width (int) – Width of the padding ring in grid cells. Use recommended_pad_width() to obtain a suitable value.

  • Te_out (float, optional) – Te value at the outer edge of the padding ring. Defaults to Te.mean().

Returns:

Te_padded – Padded elastic thickness with a smooth linear taper.

Return type:

(M + 2*pad_width, N + 2*pad_width) array

Examples

>>> import numpy as np
>>> from gflex import smooth_pad_Te
>>> Te = 35e3 * np.ones((10, 10))
>>> Te_pad = smooth_pad_Te(Te, pad_width=4)
>>> Te_pad.shape
(18, 18)

To run F2D with the padded grid:

import numpy as np
from gflex import F2D, smooth_pad_Te, recommended_pad_width
p = recommended_pad_width(Te, dx=5000.)
Te_pad = smooth_pad_Te(Te, p)
qs_pad = np.pad(qs, p, mode='constant')
flex = F2D()
flex.T_e = Te_pad
flex.qs = qs_pad
# ... set other parameters and run ...
w_inner = flex.w[p:-p, p:-p]   # trim padding from output
gflex.recommended_pad_width(Te, dx, E=65000000000.0, nu=0.25, rho_m=3300.0, rho_fill=0.0, g=9.8, n_wavelengths=1.0)[source]

Return the recommended padding width in grid cells for a variable-Te run.

The padded domain boundary should be at least one flexural wavelength from the load so that the plate’s response to the load is negligible at the boundary. Half a wavelength is often sufficient in practice.

The 2-D flexural wavelength is computed from the maximum Te value via flexural_wavelengths(), giving the most conservative (widest) padding estimate.

Parameters:
  • Te (scalar or 2-D array) – Elastic thickness [m]. The maximum value is used.

  • dx (float) – Grid cell size [m]. Use the smaller of dx and dy if they differ.

  • E (float, optional) – Young’s modulus [Pa]. Default 65 GPa.

  • nu (float, optional) – Poisson’s ratio. Default 0.25.

  • rho_m (float, optional) – Mantle density [kg m^-3]. Default 3300.

  • rho_fill (float, optional) – Infill density [kg m^-3]. Default 0 (air).

  • g (float, optional) – Gravitational acceleration [m s^-2]. Default 9.8.

  • n_wavelengths (float, optional) – Number of flexural wavelengths to use as the padding width. Default 1.0. Use 0.5 for a less conservative estimate.

Returns:

Recommended padding width in grid cells.

Return type:

int

1-D (F1D)

gflex.pad_domain_1d(Te, qs, dx, n_wavelengths=1.0, Te_out=None, E=65000000000.0, nu=0.25, rho_m=3300.0, rho_fill=0.0, g=9.8)[source]

Pad both the elastic thickness and surface load arrays for use with F1D.

Combines recommended_pad_width_1d() and smooth_pad_Te_1d() into a single call, and zero-pads qs to match. The returned pad width p can be used to trim the deflection output after the run:

w_inner = flex.w[p:-p]
Parameters:
  • Te (1-D array) – Elastic thickness [m] for the inner domain.

  • qs (1-D array) – Surface load [Pa] for the inner domain.

  • dx (float) – Grid cell size [m].

  • n_wavelengths (float, optional) – Padding width expressed as a number of flexural wavelengths. Default 1.0; use 0.5 for a less conservative (narrower) padding.

  • Te_out (float, optional) – Te value at the outer edge of the padding. Defaults to Te.mean().

  • E (float, optional) – Young’s modulus [Pa]. Default 65 GPa.

  • nu (float, optional) – Poisson’s ratio. Default 0.25.

  • rho_m (float, optional) – Mantle density [kg m^-3]. Default 3300.

  • rho_fill (float, optional) – Infill density [kg m^-3]. Default 0 (air).

  • g (float, optional) – Gravitational acceleration [m s^-2]. Default 9.8.

Returns:

  • Te_padded (1-D array of length len(Te) + 2p) – Smoothly tapered elastic thickness.

  • qs_padded (1-D array of length len(qs) + 2p) – Surface load zero-padded to match.

  • p (int) – Pad width in grid cells (same on both ends).

Examples

>>> import numpy as np
>>> Te = 10e3 * np.ones(5)
>>> qs = np.zeros(5)
>>> Te_pad, qs_pad, p = pad_domain_1d(Te, qs, dx=10000.,
...     E=65e9, nu=0.25, rho_m=3300., rho_fill=0., g=9.8)
>>> p
19
>>> Te_pad.shape
(43,)

To run F1D with the padded arrays:

flex = gflex.F1D()
flex.T_e = Te_pad
flex.qs = qs_pad
# ... set other parameters and run ...
w_inner = flex.w[p:-p]   # trim padding from output
gflex.smooth_pad_Te_1d(Te, pad_width, Te_out=None)[source]

Pad a 1-D elastic thickness array with a smooth linear taper.

When a spatially variable Te array is padded with a constant value before being passed to F1D, the abrupt step in flexural rigidity D at the inner/outer boundary drives spurious deflections via the D-derivative terms in the FD stencil.

This function eliminates that step by linearly blending the inner-domain edge values toward Te_out across the padding region.

The corresponding surface load array should be padded with zeros, e.g.:

qs_padded = numpy.pad(qs, pad_width, mode='constant')
Parameters:
  • Te (1-D array) – Elastic thickness [m] for the inner domain.

  • pad_width (int) – Width of the padding on each end in grid cells. Use recommended_pad_width_1d() to obtain a suitable value.

  • Te_out (float, optional) – Te value at the outer edge of the padding. Defaults to Te.mean().

Returns:

Te_padded – Padded elastic thickness with a smooth linear taper.

Return type:

1-D array of length len(Te) + 2 * pad_width

Examples

>>> import numpy as np
>>> Te = 35e3 * np.ones(10)
>>> Te_pad = smooth_pad_Te_1d(Te, pad_width=4)
>>> Te_pad.shape
(18,)

To run F1D with the padded arrays:

import numpy as np
from gflex import F1D, smooth_pad_Te_1d, recommended_pad_width_1d
p = recommended_pad_width_1d(Te, dx=5000.)
Te_pad = smooth_pad_Te_1d(Te, p)
qs_pad = np.pad(qs, p, mode='constant')
flex = F1D()
flex.T_e = Te_pad
flex.qs = qs_pad
# ... set other parameters and run ...
w_inner = flex.w[p:-p]   # trim padding from output
gflex.recommended_pad_width_1d(Te, dx, E=65000000000.0, nu=0.25, rho_m=3300.0, rho_fill=0.0, g=9.8, n_wavelengths=1.0)[source]

Return the recommended padding width in grid cells for a 1-D variable-Te run.

The padded domain boundary should be at least one 1-D flexural wavelength from the load so that the plate’s response is negligible at the boundary.

The 1-D flexural wavelength is computed from the maximum Te value, giving the most conservative (widest) padding estimate.

Parameters:
  • Te (scalar or 1-D array) – Elastic thickness [m]. The maximum value is used.

  • dx (float) – Grid cell size [m].

  • E (float, optional) – Young’s modulus [Pa]. Default 65 GPa.

  • nu (float, optional) – Poisson’s ratio. Default 0.25.

  • rho_m (float, optional) – Mantle density [kg m^-3]. Default 3300.

  • rho_fill (float, optional) – Infill density [kg m^-3]. Default 0 (air).

  • g (float, optional) – Gravitational acceleration [m s^-2]. Default 9.8.

  • n_wavelengths (float, optional) – Number of flexural wavelengths to use as the padding width. Default 1.0. Use 0.5 for a less conservative estimate.

Returns:

Recommended padding width in grid cells.

Return type:

int

Examples

>>> recommended_pad_width_1d(Te=35e3, dx=5000.)
94

FD boundary-condition warnings

When running F1D or F2D with the finite-difference solver, gFlex issues UserWarning messages for two categories of potentially problematic boundary conditions.

BC-type warnings fire whenever a side carries a BC whose physical interpretation deserves verification:

  • 'zero_moment_zero_shear' (alias 'free') — assumes a free broken plate end (zero moment and shear force). Physically appropriate for rifted or passive continental margins, subduction trenches with an applied edge load (slab pull), and broken-plate flexure (Turcotte & Schubert). Often applied uncritically elsewhere in the literature — verify that one of these settings applies.

Proximity warnings fire for 'zero_displacement_zero_slope' (alias 'clamped') boundaries when the nearest loaded cell is within one flexural wavelength (\(\lambda = 2\pi\alpha\), where \(\alpha = (4D / \Delta\rho g)^{1/4}\)) of the boundary. Within this distance the flexural forebulge — which peaks at \(\approx \pi\alpha\) from the load — will be suppressed by the zero-displacement condition, contaminating the solution. The warning message reports the distance as a fraction of the local flexural wavelength and directs you to the domain-padding utilities.

Warning deduplication in model-coupling loops

Python’s default warning filter shows each unique warning once per call site per interpreter session. In a time-stepping or iterative-coupling loop such as:

for load in loads:
    flex.qs = load
    flex.run()          # warning fires on first iteration, silenced thereafter

the proximity warning fires on the first iteration and is not repeated, even if the load subsequently moves closer to the boundary. To re-enable the warning on every call:

import warnings
warnings.filterwarnings("always", category=UserWarning, module="gflex")

Suppressing warnings you have verified

Once you have confirmed that a boundary condition is appropriate for your setup, suppress the corresponding warning by message text:

import warnings
warnings.filterwarnings("ignore", message=".*zero_moment_zero_shear.*")

To suppress all gFlex warnings:

warnings.filterwarnings("ignore", module="gflex")

Or use a context manager for a single run:

with warnings.catch_warnings():
    warnings.simplefilter("ignore")
    flex.run()

LU factorization cache

For coupling workflows that call run() (or run()) repeatedly with the same grid, elastic thickness, and boundary conditions, the sparse-LU factorization of the coefficient matrix can be cached to avoid re-factorizing on every call. Set the attribute before calling initialize():

Value

Behaviour

False (default)

No caching. The matrix is factorized on every run() call.

True

Cache with hash check. The factorization is reused when a hash of the coefficient matrix matches the stored hash; it is recomputed when any of Te, dx/dy, BCs, or physical parameters change.

"no_check"

Cache without hash check. The stored factorization is reused on every call without computing a hash. Gives the maximum performance benefit. The coefficient matrix is freed from memory immediately after factorization; only the LU factorization is retained. Smart invalidation (see below) still applies: reassigning a matrix-determining input clears the cache and triggers a rebuild on the next call.

Smart invalidation

Reassigning any matrix-determining attribute — T_e, E, nu, g, rho_m, rho_fill, dx, dy, boundary conditions, or in-plane stresses — automatically clears the cached coefficient matrix and LU factorization. No explicit cache management is needed between solves when only qs changes.

Note

Smart invalidation is triggered by assignment (flex.T_e = new_array), not by in-place mutation of a NumPy array (flex.T_e[5] = 40e3). If you mutate an array in place, reassign it afterwards to ensure the cache is correctly invalidated:

flex.T_e[5] = 40e3
flex.T_e = flex.T_e   # trigger invalidation

Example (coupling loop):

flex.cache_factorization = True
flex.initialize()

for load in load_sequence:
    flex.qs = load
    flex.run()
    w = flex.w
    # ... process w ...

flex.finalize()

The cache is cleared by finalize().

Flexural wavelength

gflex.flexural_wavelengths(Te, E, nu, rho_m, rho_fill, g)[source]

Compute flexural parameters and wavelengths for a thin elastic plate.

Parameters:
  • Te (float) – Elastic thickness [m].

  • E (float) – Young’s modulus [Pa].

  • nu (float) – Poisson’s ratio.

  • rho_m (float) – Mantle density [kg m^-3].

  • rho_fill (float) – Infill density [kg m^-3] (e.g. 0 for air, 1000 for water).

  • g (float) – Gravitational acceleration [m s^-2].

Returns:

Keys alpha_1D, lambda_1D, zero_crossing_1D, alpha_2D, lambda_2D, zero_crossing_2D — all in metres.

Return type:

dict

Examples

>>> from gflex import flexural_wavelengths
>>> r = flexural_wavelengths(Te=30e3, E=65e9, nu=0.25,
...                          rho_m=3300., rho_fill=0., g=9.8)
>>> round(r["alpha_2D"] / 1e3, 1)  # km
46.9
>>> round(r["lambda_1D"] / r["lambda_2D"], 6)
1.414214

BMI interface

BmiGflex exposes the CSDMS Basic Model Interface, enabling gFlex to be coupled with other models in the CSDMS framework. It requires the optional bmipy dependency (pip install gflex[bmi]).

class gflex.BmiGflex[source]

Bases: object

BMI wrapper for gFlex lithospheric flexure.

Implements the CSDMS Basic Model Interface v2 specification. Supports 1-D and 2-D gridded flexure solutions (FD, FFT, SAS methods). The SAS_NG point-load method is not suited to the BMI grid model.

Grid

All variables share grid identifier 0, a uniform rectilinear grid. Grid shape follows NumPy / C order: (nrows,) in 1-D or (nrows, ncols) in 2-D, with spacing (dy, dx) and origin at (0, 0).

Time

gFlex solves instantaneous elastic equilibrium. Time is therefore nominal: start=0, step=1, end=inf. Each call to update() applies the current load and computes deflection.

Variables

Input: lithosphere__load_pressure [Pa] — surface normal stress q0 Output: lithosphere__vertical_displacement [m] — downward deflection w

initialize(config_file: str) None[source]

Initialize gFlex from a configuration file.

Parameters:

config_file (str) – Path to a gFlex YAML configuration file.

update() None[source]

Compute flexural deflection for the current load.

Writes the current lithosphere__load_pressure array into the model’s internal qs field, runs the solver, then copies the result into lithosphere__vertical_displacement.

finalize() None[source]

Tear down the model and release resources.

get_value(name: str, dest: ndarray[tuple[Any, ...], dtype[Any]]) ndarray[tuple[Any, ...], dtype[Any]][source]

Copy the flattened values of variable name into dest and return it.

set_value(name: str, src: ndarray[tuple[Any, ...], dtype[Any]]) None[source]

Overwrite the entire array for variable name with values from src.