earthsciencekit.geophysics#

Domain hub · History · Examples · API reference

Inspectable geophysics models: gravity, the figure of the Earth and geomagnetism; SI units unless stated otherwise.

earthsciencekit.geophysics.block_magnetic_anomaly(x, left, right, top, bottom, magnetization)[source]#

Magnetic anomaly of a vertically magnetized two-dimensional block.

Strips of normally and reversely magnetized crust give the alternating anomalies over mid-ocean ridges.

Parameters:

  • x (float or array_like) – m; constraint: real. Horizontal position of the observer, perpendicular to the strike.

  • left (float or array_like) – m; constraint: left < right. Edges of the block.

  • right (float or array_like) – m; constraint: left < right. Edges of the block.

  • top (float or array_like) – m; constraint: 0 < top < bottom. Depths below the observation level.

  • bottom (float or array_like) – m; constraint: 0 < top < bottom. Depths below the observation level.

  • magnetization (float or array_like) – A/m; constraint: real. Positive downward (normal polarity at the north magnetic pole).

Returns:

Model value [T]: vertical field anomaly, equal to the total-field anomaly reduced to the pole.

Return type:

numpy scalar or ndarray

Notes

Infinite strike length, uniform vertical magnetization and a vertical ambient field. A laterally infinite layer has no external field. Equation: MU0 * magnetization / (2 * np.pi) * (phi(top) - phi(bottom)) with phi(z) = np.arctan((right - x) / z) - np.arctan((left - x) / z). Invalid, empty, nonreal or nonfinite inputs raise ValueError. NumPy broadcasting applies; singular and overflowing outputs raise ValueError.

earthsciencekit.geophysics.bouguer_anomaly(observed_gravity, normal_gravity, height, density, radius)[source]#

Simple Bouguer gravity anomaly.

Removing latitude, height and the rock slab under a station leaves the attraction of buried density contrasts.

Parameters:

  • observed_gravity (float or array_like) – m/s^2; constraint: positive.

  • normal_gravity (float or array_like) – m/s^2; constraint: positive. Normal gravity on the reference surface at the station’s latitude.

  • height (float or array_like) – m; constraint: real. Station height above the reference surface.

  • density (float or array_like) – kg/m^3; constraint: positive. Reduction density of the slab.

  • radius (float or array_like) – m; constraint: positive.

Returns:

Model value [m/s^2]. Scalar input returns a NumPy scalar; arrays broadcast.

Return type:

numpy scalar or ndarray

Notes

Free-air and infinite-slab reductions only, without terrain correction. The free-air anomaly is the same expression without the slab term. Equation: observed_gravity - normal_gravity + 2 * normal_gravity * height / radius - 2 * np.pi * G * density * height. Invalid, empty, nonreal or nonfinite inputs raise ValueError. NumPy broadcasting applies; singular and overflowing outputs raise ValueError.

earthsciencekit.geophysics.bouguer_correction(density, thickness)[source]#

Attraction of an infinite horizontal slab.

Rock between the station and the reference level adds gravity that does not depend on the slab’s extent.

Parameters:

  • density (float or array_like) – kg/m^3; constraint: positive.

  • thickness (float or array_like) – m; constraint: real. Signed station height above the reference level.

Returns:

Model value [m/s^2]. Scalar input returns a NumPy scalar; arrays broadcast.

Return type:

numpy scalar or ndarray

Notes

Flat, uniform, laterally infinite slab; no terrain or curvature correction. About 0.112 mGal per metre at 2670 kg/m^3. Equation: 2 * np.pi * G * density * thickness. Invalid, empty, nonreal or nonfinite inputs raise ValueError. NumPy broadcasting applies; singular and overflowing outputs raise ValueError.

earthsciencekit.geophysics.centrifugal_ratio(angular_speed, radius, gravity)[source]#

Ratio of centrifugal to gravitational acceleration at the equator.

Rotation flattens a planet in proportion to how strongly it competes with gravity.

Parameters:

  • angular_speed (float or array_like) – rad/s; constraint: nonnegative.

  • radius (float or array_like) – m; constraint: positive. Equatorial radius.

  • gravity (float or array_like) – m/s^2; constraint: positive. Equatorial gravity.

Returns:

Model value [1]. Scalar input returns a NumPy scalar; arrays broadcast.

Return type:

numpy scalar or ndarray

Notes

The geodetic m uses GM / (a b) in place of gravity; the two differ at the 0.5% level. Equation: angular_speed**2 * radius / gravity. Invalid, empty, nonreal or nonfinite inputs raise ValueError. NumPy broadcasting applies; singular and overflowing outputs raise ValueError.

earthsciencekit.geophysics.clairaut_gravity(latitude, equatorial_gravity, flattening, centrifugal_ratio)[source]#

Gravity on a rotating spheroid from Clairaut’s theorem.

The increase of gravity toward the poles measures the flattening without any surveying.

Parameters:

  • latitude (float or array_like) – rad; constraint: -pi/2 <= latitude <= pi/2. Geodetic latitude.

  • equatorial_gravity (float or array_like) – m/s^2; constraint: positive.

  • flattening (float or array_like) – 1; constraint: 0 <= value < 1.

  • centrifugal_ratio (float or array_like) – 1; constraint: nonnegative.

Returns:

Model value [m/s^2]. Scalar input returns a NumPy scalar; arrays broadcast.

Return type:

numpy scalar or ndarray

Notes

First order in flattening and rotation for a spheroid in hydrostatic equilibrium; the sin^4 term is omitted. Equation: equatorial_gravity * (1 + (5 * centrifugal_ratio / 2 - flattening) * np.sin(latitude) ** 2). Invalid, empty, nonreal or nonfinite inputs raise ValueError. NumPy broadcasting applies; singular and overflowing outputs raise ValueError.

earthsciencekit.geophysics.degree_geoid_height(anomaly_coefficient, degree, radius, gravity)[source]#

Geoid height of one spherical-harmonic degree of the gravity anomaly.

Stokes’s integral divides each degree of the anomaly by n - 1, so long wavelengths dominate the geoid.

Parameters:

  • anomaly_coefficient (float or array_like) – m/s^2; constraint: real. Coefficient of the degree-n surface harmonic of the gravity anomaly.

  • degree (float or array_like) – 1; constraint: integer >= 2.

  • radius (float or array_like) – m; constraint: positive.

  • gravity (float or array_like) – m/s^2; constraint: positive. Mean normal gravity.

Returns:

Model value [m]: coefficient of the same surface harmonic in the geoid height.

Return type:

numpy scalar or ndarray

Notes

Spherical approximation of Stokes’s formula; degrees 0 and 1 are excluded by construction. Equation: radius * anomaly_coefficient / ((degree - 1) * gravity). Invalid, empty, nonreal or nonfinite inputs raise ValueError. NumPy broadcasting applies; singular and overflowing outputs raise ValueError.

earthsciencekit.geophysics.dipole_field(latitude, distance, moment)[source]#

Field of a geocentric axial dipole.

A uniformly magnetized Earth gives a field that dips more steeply toward the poles.

Parameters:

  • latitude (float or array_like) – rad; constraint: -pi/2 <= latitude <= pi/2. Geocentric latitude.

  • distance (float or array_like) – m; constraint: positive. Distance from the Earth’s centre.

  • moment (float or array_like) – A m^2; constraint: real. Positive for today’s polarity, whose field points down in the northern hemisphere.

Returns:

Model value [T]: the northward and the downward component.

Return type:

ndarray, shape (…, 2)

Notes

Pure axial dipole, outside its sources. The field at the poles is twice that at the equator. Equation: MU0 * moment / (4 * np.pi * distance**3) * [np.cos(latitude), 2 * np.sin(latitude)]. Invalid, empty, nonreal or nonfinite inputs raise ValueError. NumPy broadcasting applies; singular and overflowing outputs raise ValueError.

earthsciencekit.geophysics.dipole_inclination(latitude)[source]#

Inclination of an axial dipole field.

The dip of a compass needle tells a navigator, or a rock, its magnetic latitude.

Parameters:

latitude (float or array_like) – rad; constraint: -pi/2 <= latitude <= pi/2.

Returns:

Model value [rad]: positive downward, for today’s polarity.

Return type:

numpy scalar or ndarray

Notes

Geocentric axial dipole; tan(I) = 2 tan(latitude). Equation: np.arctan2(2 * np.sin(latitude), np.cos(latitude)). Invalid, empty, nonreal or nonfinite inputs raise ValueError. NumPy broadcasting applies; singular and overflowing outputs raise ValueError.

earthsciencekit.geophysics.direction_angle(declination_a, inclination_a, declination_b, inclination_b)[source]#

Angle between two magnetic directions.

A remanence more than 90 degrees from the present field records a reversed polarity.

Parameters:

  • declination_a (float or array_like) – rad; constraint: real.

  • declination_b (float or array_like) – rad; constraint: real.

  • inclination_a (float or array_like) – rad; constraint: -pi/2 <= inclination <= pi/2.

  • inclination_b (float or array_like) – rad; constraint: -pi/2 <= inclination <= pi/2.

Returns:

Model value [rad] in [0, pi]. Scalar input returns a NumPy scalar; arrays broadcast.

Return type:

numpy scalar or ndarray

Notes

Great-circle distance between unit vectors; intensities are ignored. Equation: np.arccos(np.sin(Ia) * np.sin(Ib) + np.cos(Ia) * np.cos(Ib) * np.cos(Da - Db)). Invalid, empty, nonreal or nonfinite inputs raise ValueError. NumPy broadcasting applies; singular and overflowing outputs raise ValueError.

earthsciencekit.geophysics.earth_mean_density(gravity, radius, gravitational_constant=6.6743e-11)[source]#

Mean density of the Earth.

Surface gravity and the gravitational constant together weigh the Earth.

Parameters:

  • gravity (float or array_like) – m/s^2; constraint: positive.

  • radius (float or array_like) – m; constraint: positive.

  • gravitational_constant (float or array_like) – m^3/(kg s^2); constraint: positive. Defaults to earthsciencekit.constants.G.

Returns:

Model value [kg/m^3]. Scalar input returns a NumPy scalar; arrays broadcast.

Return type:

numpy scalar or ndarray

Notes

Spherically symmetric, nonrotating Earth, so gravity = G M / radius**2. Equation: 3 * gravity / (4 * np.pi * gravitational_constant * radius). Invalid, empty, nonreal or nonfinite inputs raise ValueError. NumPy broadcasting applies; singular and overflowing outputs raise ValueError.

earthsciencekit.geophysics.equivalent_water_height(coefficient_change, degree, load_love_number, radius, mean_density, water_density)[source]#

Water layer behind a change of one Stokes coefficient.

Monthly gravity fields from a satellite pair reveal the moving mass of water, ice and groundwater.

Parameters:

  • coefficient_change (float or array_like) – 1; constraint: real. Change of a fully normalized Stokes coefficient of degree n.

  • degree (float or array_like) – 1; constraint: integer >= 0.

  • load_love_number (float or array_like) – 1; constraint: greater than -1. Elastic yielding k'_n of the Earth under the load.

  • radius (float or array_like) – m; constraint: positive.

  • mean_density (float or array_like) – kg/m^3; constraint: positive. Mean density of the Earth.

  • water_density (float or array_like) – kg/m^3; constraint: positive.

Returns:

Model value [m]: coefficient of the same fully normalized harmonic in the water thickness.

Return type:

numpy scalar or ndarray

Notes

Thin surface layer on an elastic Earth (Wahr et al., 1998); valid while the layer is much thinner than the resolution of the truncated expansion. Equation: radius * mean_density * (2 * degree + 1) * coefficient_change / (3 * water_density * (1 + load_love_number)). Invalid, empty, nonreal or nonfinite inputs raise ValueError. NumPy broadcasting applies; singular and overflowing outputs raise ValueError.

earthsciencekit.geophysics.field_elements(north, east, down)[source]#

Declination, inclination and intensity of a field vector.

Magnetic directions in rocks and observatories are reported as two angles and a strength.

Parameters:

  • north (float or array_like) – T; constraint: real, not all zero.

  • east (float or array_like) – T; constraint: real, not all zero.

  • down (float or array_like) – T; constraint: real, not all zero.

Returns:

Model value: declination [rad] east of north in (-pi, pi], inclination [rad] positive down, intensity [T].

Return type:

ndarray, shape (…, 3)

Notes

Local north-east-down frame. A reversed field has declination near pi and inclination of the opposite sign. Equation: [np.arctan2(east, north), np.arctan2(down, np.hypot(north, east)), np.sqrt(north**2 + east**2 + down**2)]. Invalid, empty, nonreal or nonfinite inputs raise ValueError. NumPy broadcasting applies; singular and overflowing outputs raise ValueError.

earthsciencekit.geophysics.flattening_from_j2(j2, centrifugal_ratio)[source]#

Flattening of a hydrostatic Earth from its J2.

The oblateness term of the gravity field and the rotation fix the shape of the level surface.

Parameters:

  • j2 (float or array_like) – 1; constraint: real. Unnormalized second zonal harmonic -C20.

  • centrifugal_ratio (float or array_like) – 1; constraint: nonnegative.

Returns:

Model value [1]. Scalar input returns a NumPy scalar; arrays broadcast.

Return type:

numpy scalar or ndarray

Notes

First order in flattening (Clairaut’s relation for the external field); second-order terms change 1/f by about 0.1%. Equation: 3 * j2 / 2 + centrifugal_ratio / 2. Invalid, empty, nonreal or nonfinite inputs raise ValueError. NumPy broadcasting applies; singular and overflowing outputs raise ValueError.

earthsciencekit.geophysics.free_air_correction(height, gravity, radius)[source]#

Free-air gravity correction.

Gravity falls off with distance from the Earth’s centre, about 0.3086 mGal per metre at the surface.

Parameters:

  • height (float or array_like) – m; constraint: real. Station height above the reference surface.

  • gravity (float or array_like) – m/s^2; constraint: positive. Gravity on the reference surface.

  • radius (float or array_like) – m; constraint: positive.

Returns:

Model value [m/s^2]: the amount to add to gravity observed at height.

Return type:

numpy scalar or ndarray

Notes

First order in height / radius for a spherical point-mass Earth; no rock between station and reference. Equation: 2 * gravity * height / radius. Invalid, empty, nonreal or nonfinite inputs raise ValueError. NumPy broadcasting applies; singular and overflowing outputs raise ValueError.

earthsciencekit.geophysics.gauss_dipole_moment(g10, g11, h11, reference_radius)[source]#

Dipole moment from the degree-1 Gauss coefficients.

Gauss’s analysis weighed the Earth’s magnet from surface measurements alone.

Parameters:

  • g10 (float or array_like) – T; constraint: real. Schmidt semi-normalized degree-1 coefficients.

  • g11 (float or array_like) – T; constraint: real. Schmidt semi-normalized degree-1 coefficients.

  • h11 (float or array_like) – T; constraint: real. Schmidt semi-normalized degree-1 coefficients.

  • reference_radius (float or array_like) – m; constraint: positive.

Returns:

Model value [A m^2]. Scalar input returns a NumPy scalar; arrays broadcast.

Return type:

numpy scalar or ndarray

Notes

Magnitude of the best-fitting geocentric (tilted) dipole; the tilt is arccos(abs(g10) / B0). Equation: 4 * np.pi * reference_radius**3 / MU0 * np.sqrt(g10**2 + g11**2 + h11**2). Invalid, empty, nonreal or nonfinite inputs raise ValueError. NumPy broadcasting applies; singular and overflowing outputs raise ValueError.

earthsciencekit.geophysics.homogeneous_flattening(centrifugal_ratio)[source]#

Flattening of a slowly rotating homogeneous fluid planet.

A uniform self-gravitating fluid in rotation settles into an oblate spheroid.

Parameters:

centrifugal_ratio (float or array_like) – 1; constraint: 0 <= value <= 0.1.

Returns:

Model value [1]: (a - c) / a.

Return type:

numpy scalar or ndarray

Notes

First order in rotation (Maclaurin spheroid). A centrally condensed Earth is less flattened than this. Equation: 5 * centrifugal_ratio / 4. Invalid, empty, nonreal or nonfinite inputs raise ValueError. NumPy broadcasting applies; singular and overflowing outputs raise ValueError.

earthsciencekit.geophysics.love_number_h2(rigidity, density, gravity, radius)[source]#

Degree-2 displacement Love number of a homogeneous elastic Earth.

The solid Earth rises and falls with the tide by a fraction of the equilibrium tide that measures its rigidity.

Parameters:

  • rigidity (float or array_like) – Pa; constraint: nonnegative. Shear modulus; zero gives a fluid Earth.

  • density (float or array_like) – kg/m^3; constraint: positive.

  • gravity (float or array_like) – m/s^2; constraint: positive.

  • radius (float or array_like) – m; constraint: positive.

Returns:

Model value [1]: radial displacement over the equilibrium tide height.

Return type:

numpy scalar or ndarray

Notes

Uniform, incompressible, self-gravitating elastic sphere (Kelvin’s solution); h2 = 5 k2 / 3. Equation: 2.5 / (1 + 19 * rigidity / (2 * density * gravity * radius)). Invalid, empty, nonreal or nonfinite inputs raise ValueError. NumPy broadcasting applies; singular and overflowing outputs raise ValueError.

earthsciencekit.geophysics.love_number_k2(rigidity, density, gravity, radius)[source]#

Degree-2 potential Love number of a homogeneous elastic Earth.

The tidal bulge adds its own gravitational potential, a fraction k2 of the tidal potential.

Parameters:

  • rigidity (float or array_like) – Pa; constraint: nonnegative. Shear modulus; zero gives a fluid Earth.

  • density (float or array_like) – kg/m^3; constraint: positive.

  • gravity (float or array_like) – m/s^2; constraint: positive.

  • radius (float or array_like) – m; constraint: positive.

Returns:

Model value [1]. Scalar input returns a NumPy scalar; arrays broadcast.

Return type:

numpy scalar or ndarray

Notes

Uniform, incompressible, self-gravitating elastic sphere (Kelvin’s solution); the fluid limit is 3/2. Equation: 1.5 / (1 + 19 * rigidity / (2 * density * gravity * radius)). Invalid, empty, nonreal or nonfinite inputs raise ValueError. NumPy broadcasting applies; singular and overflowing outputs raise ValueError.

earthsciencekit.geophysics.magnetic_decay_time(radius, conductivity)[source]#

Free decay time of the slowest magnetic mode of a conducting sphere.

Without a dynamo the Earth’s field would die in tens of thousands of years, far less than its age.

Parameters:

  • radius (float or array_like) – m; constraint: positive.

  • conductivity (float or array_like) – S/m; constraint: positive.

Returns:

Model value [s]: e-folding time of the field of the slowest dipole mode.

Return type:

numpy scalar or ndarray

Notes

Uniform stationary conductor surrounded by an insulator. Equation: MU0 * conductivity * radius**2 / np.pi**2. Invalid, empty, nonreal or nonfinite inputs raise ValueError. NumPy broadcasting applies; singular and overflowing outputs raise ValueError.

earthsciencekit.geophysics.magnetic_reynolds_number(speed, length, conductivity)[source]#

Magnetic Reynolds number.

A conducting fluid can stretch and regenerate a magnetic field only if advection beats diffusion.

Parameters:

  • speed (float or array_like) – m/s; constraint: nonnegative.

  • length (float or array_like) – m; constraint: positive.

  • conductivity (float or array_like) – S/m; constraint: positive.

Returns:

Model value [1]. Scalar input returns a NumPy scalar; arrays broadcast.

Return type:

numpy scalar or ndarray

Notes

Order-of-magnitude scale; dynamos need values above a geometry-dependent threshold of order 10 to 100. Equation: MU0 * conductivity * speed * length. Invalid, empty, nonreal or nonfinite inputs raise ValueError. NumPy broadcasting applies; singular and overflowing outputs raise ValueError.

earthsciencekit.geophysics.nodal_precession_rate(semi_major_axis, eccentricity, inclination, j2, radius, gravitational_parameter)[source]#

Secular regression of a satellite orbit’s node.

The equatorial bulge torques an inclined orbit, so the node drifts at a rate that measures J2.

Parameters:

  • semi_major_axis (float or array_like) – m; constraint: positive.

  • eccentricity (float or array_like) – 1; constraint: 0 <= value < 1.

  • inclination (float or array_like) – rad; constraint: 0 <= inclination <= pi.

  • j2 (float or array_like) – 1; constraint: real.

  • radius (float or array_like) – m; constraint: positive. Equatorial reference radius of J2.

  • gravitational_parameter (float or array_like) – m^3/s^2; constraint: positive. G M of the Earth.

Returns:

Model value [rad/s]: negative (westward) for prograde orbits and positive J2.

Return type:

numpy scalar or ndarray

Notes

First-order secular J2 perturbation, without drag or higher harmonics; perigee must clear the reference radius. Equation: -1.5 * n * j2 * (radius / semi_major_axis) ** 2 * np.cos(inclination) / (1 - eccentricity**2) ** 2 with mean motion n = np.sqrt(gravitational_parameter / semi_major_axis**3). Invalid, empty, nonreal or nonfinite inputs raise ValueError. NumPy broadcasting applies; singular and overflowing outputs raise ValueError.

earthsciencekit.geophysics.normal_gravity(latitude)[source]#

Normal gravity on the GRS80 ellipsoid.

The reference field that gravity anomalies are measured against.

Parameters:

latitude (float or array_like) – rad; constraint: -pi/2 <= latitude <= pi/2. Geodetic latitude.

Returns:

Model value [m/s^2]. Scalar input returns a NumPy scalar; arrays broadcast.

Return type:

numpy scalar or ndarray

Notes

Somigliana’s closed form on the ellipsoid surface with the GRS80 constants gamma_e = 9.7803267715, k = 0.001931851353 and e2 = 0.00669438002290; exact for the level ellipsoid. Equation: gamma_e * (1 + k * np.sin(latitude) ** 2) / np.sqrt(1 - e2 * np.sin(latitude) ** 2). Invalid, empty, nonreal or nonfinite inputs raise ValueError. NumPy broadcasting applies; singular and overflowing outputs raise ValueError.

earthsciencekit.geophysics.paleolatitude(inclination)[source]#

Paleolatitude from a remanent inclination.

A rock’s frozen-in dip records how far from the pole it formed.

Parameters:

inclination (float or array_like) – rad; constraint: -pi/2 <= inclination <= pi/2. Positive downward.

Returns:

Model value [rad]: the magnetic latitude of the site when the rock was magnetized.

Return type:

numpy scalar or ndarray

Notes

Assumes a geocentric axial dipole averaged over secular variation; a reversed direction must be inverted first. Equation: np.arctan2(np.sin(inclination), 2 * np.cos(inclination)). Invalid, empty, nonreal or nonfinite inputs raise ValueError. NumPy broadcasting applies; singular and overflowing outputs raise ValueError.

earthsciencekit.geophysics.rikitake_fixed_point(mu, k)[source]#

Steady state of Rikitake’s two-disc dynamo.

The two coupled discs have a pair of steady states of opposite polarity, which the chaotic flow visits in turn.

Parameters:

  • mu (float or array_like) – 1; constraint: positive. Scaled resistance.

  • k (float or array_like) – 1; constraint: positive. Sets A = mu * (k**2 - 1 / k**2).

Returns:

Model value [1]: the scaled currents x, y and the scaled angular speed z of the positive state; the other state is (-x, -y, z).

Return type:

ndarray, shape (…, 3)

Notes

Nondimensional system of simulate_rikitake_dynamo(). Both states are unstable for the usual parameters. Equation: [k, 1 / k, mu * k**2]. Invalid, empty, nonreal or nonfinite inputs raise ValueError. NumPy broadcasting applies; singular and overflowing outputs raise ValueError.

earthsciencekit.geophysics.simulate_rikitake_dynamo(mu, k, time_scale, time, *, initial_state=(1.0, 0.0, 0.0))[source]#

Integrate Rikitake’s coupled two-disc dynamo.

Parameters:

  • mu (float) – Positive scalar scaled resistance [1].

  • k (float) – Positive scalar setting the torque asymmetry A = mu * (k**2 - 1 / k**2) [1].

  • time_scale (float) – Positive scalar seconds per nondimensional time unit [s].

  • time (array_like) – Nonnegative strictly increasing sampled times [s].

  • initial_state (sequence of float) – Scaled currents x, y and angular speed z at the first time [1].

Returns:

The scaled currents x and y (whose sign is the polarity) and the scaled angular speed z [1].

Return type:

TimeSeries

Notes

dx/dt = -mu * x + y * z, dy/dt = -mu * y + (z - A) * x, dz/dt = 1 - x * y in units of time_scale. The steady states are rikitake_fixed_point(). Trajectories oscillate about one and jump irregularly to the other, reversing both currents: a chaotic analogue of geomagnetic reversals, not a model of the core.

earthsciencekit.geophysics.sphere_gravity(x, depth, radius, density_contrast)[source]#

Gravitational attraction of a uniform sphere.

A compact body with excess density pulls a plumb line sideways and adds to gravity above it.

Parameters:

  • x (float or array_like) – m; constraint: real. Horizontal distance of the observer from the point above the centre.

  • depth (float or array_like) – m; constraint: real. Depth of the centre below the observer, positive downward.

  • radius (float or array_like) – m; constraint: positive.

  • density_contrast (float or array_like) – kg/m^3; constraint: real. Density minus that of the surroundings.

Returns:

Model value [m/s^2]: the horizontal component, positive toward increasing x, and the downward component.

Return type:

ndarray, shape (…, 2)

Notes

The observer must lie outside the sphere, where it attracts as a point mass M = 4/3 pi radius**3 density_contrast. The horizontal component always points toward the centre for positive contrast. Equation: G * M * [-x, depth] / (x**2 + depth**2) ** 1.5. Invalid, empty, nonreal or nonfinite inputs raise ValueError. NumPy broadcasting applies; singular and overflowing outputs raise ValueError.

earthsciencekit.geophysics.stokes_function(spherical_distance)[source]#

Stokes’s kernel for the geoid.

Weighting gravity anomalies over the whole globe by spherical distance gives the geoid height.

Parameters:

spherical_distance (float or array_like) – rad; constraint: 0 < value <= pi.

Returns:

Model value [1]. Scalar input returns a NumPy scalar; arrays broadcast.

Return type:

numpy scalar or ndarray

Notes

Spherical approximation with no masses outside the geoid; singular as 2 / psi at zero distance. Its Legendre expansion is sum((2 n + 1) / (n - 1) P_n(cos psi), n >= 2), so the geoid is radius / (4 pi gravity) * integral(anomaly * S dsigma). Equation: 1 / s - 6 * s + 1 - 5 * np.cos(psi) - 3 * np.cos(psi) * np.log(s + s**2) with s = np.sin(psi / 2). Invalid, empty, nonreal or nonfinite inputs raise ValueError. NumPy broadcasting applies; singular and overflowing outputs raise ValueError.

earthsciencekit.geophysics.tidal_gravimetric_factor(h2, k2)[source]#

Gravimetric factor of the body tide.

A gravimeter on a yielding Earth sees the tide amplified by its uplift and reduced by the bulge’s attraction.

Parameters:

  • h2 (float or array_like) – 1; constraint: nonnegative.

  • k2 (float or array_like) – 1; constraint: nonnegative.

Returns:

Model value [1]: observed over rigid-Earth tidal gravity change.

Return type:

numpy scalar or ndarray

Notes

Degree-2 tide on a spherically symmetric elastic Earth, without ocean loading. About 1.16 for the real Earth. Equation: 1 + h2 - 1.5 * k2. Invalid, empty, nonreal or nonfinite inputs raise ValueError. NumPy broadcasting applies; singular and overflowing outputs raise ValueError.

earthsciencekit.geophysics.torsion_balance_gravitational_constant(arm_length, separation, deflection, attracting_mass, period)[source]#

Gravitational constant from a torsion balance.

The twist of a suspended rod between large masses, timed by its free swing, measures G.

Parameters:

  • arm_length (float or array_like) – m; constraint: positive. Distance between the two small balls.

  • separation (float or array_like) – m; constraint: positive. Centre distance between each small ball and its large ball.

  • deflection (float or array_like) – rad; constraint: positive. Equilibrium twist of the rod.

  • attracting_mass (float or array_like) – kg; constraint: positive. Mass of each large ball.

  • period (float or array_like) – s; constraint: positive. Free oscillation period of the rod.

Returns:

Model value [m^3/(kg s^2)]. Scalar input returns a NumPy scalar; arrays broadcast.

Return type:

numpy scalar or ndarray

Notes

Point masses, a massless rod and attraction by the near large ball only; the small mass cancels. Equation: 2 * np.pi**2 * arm_length * separation**2 * deflection / (attracting_mass * period**2). Invalid, empty, nonreal or nonfinite inputs raise ValueError. NumPy broadcasting applies; singular and overflowing outputs raise ValueError.

earthsciencekit.geophysics.virtual_geomagnetic_pole(site_latitude, site_longitude, declination, inclination)[source]#

Virtual geomagnetic pole of a site direction.

Poles from rocks of one age on one continent should agree; when they wander, the continent moved.

Parameters:

  • site_latitude (float or array_like) – rad; constraint: -pi/2 <= latitude <= pi/2.

  • site_longitude (float or array_like) – rad; constraint: real.

  • declination (float or array_like) – rad; constraint: real.

  • inclination (float or array_like) – rad; constraint: -pi/2 <= inclination <= pi/2.

Returns:

Model value [rad]: pole latitude and pole longitude in [0, 2 pi).

Return type:

ndarray, shape (…, 2)

Notes

The pole of the geocentric dipole that gives the observed direction at the site; a reversed direction gives a southern pole. Its angular distance from the site is p = np.arctan2(2 * np.cos(I), np.sin(I)), along the great circle leaving the site at azimuth D. Equation: pole = site * np.cos(p) + (north * np.cos(D) + east * np.sin(D)) * np.sin(p) with the site’s unit position vector and local north and east unit vectors; latitude np.arctan2(z, np.hypot(x, y)). Invalid, empty, nonreal or nonfinite inputs raise ValueError. NumPy broadcasting applies; singular and overflowing outputs raise ValueError.

earthsciencekit.geophysics.zonal_field(latitude, distance, degree, coefficient, reference_radius)[source]#

Field of one zonal Gauss coefficient.

Expanding the potential in spherical harmonics separates the dipole from smaller, faster-decaying parts.

Parameters:

  • latitude (float or array_like) – rad; constraint: -pi/2 <= latitude <= pi/2. Geocentric latitude.

  • distance (float or array_like) – m; constraint: positive. Distance from the Earth’s centre.

  • degree (float or array_like) – 1; constraint: integer >= 1.

  • coefficient (float or array_like) – T; constraint: real. Gauss coefficient g_n^0; for zonal terms Schmidt and unnormalized Legendre functions agree.

  • reference_radius (float or array_like) – m; constraint: positive. Usually 6371.2 km.

Returns:

Model value [T]: the northward and the downward component.

Return type:

ndarray, shape (…, 2)

Notes

Internal sources only, potential V = reference_radius * sum((a / r) ** (n + 1) g_n^0 P_n(sin latitude)). Sum several calls for a truncated zonal model; non-zonal terms (m > 0) are not included. Equation: (a / r) ** (n + 2) * coefficient * [lpmv(1, n, x), -(n + 1) * lpmv(0, n, x)] with x = np.sin(latitude), a = reference_radius, r = distance and n = degree. Invalid, empty, nonreal or nonfinite inputs raise ValueError. NumPy broadcasting applies; singular and overflowing outputs raise ValueError.