A solid's stress from its strain¶
A solid's stress from its strain, and whether it fails.
Hooke's law in one dimension, sigma = E epsilon, E the Young's modulus and epsilon the strain (the change in length over the length, positive stretched); the solid fails where the stress's magnitude passes its strength (tensile where stretched, compressive where squeezed: a brittle material's, which a game's bonds break by). Source: Timoshenko and Goodier, Theory of Elasticity, 3rd ed., section 3; Ashby and Jones, Engineering Materials 1, 4th ed., chapters 3 and 8. Checked against a reference in numpy, which has its own tests.
- Source
- Timoshenko and Goodier, Theory of Elasticity, 3rd ed., section 3; Ashby and Jones, Engineering Materials 1, 4th ed., chapter 3.
- Proof
- reference, invariant, limit.
elastic.trait
namespace lusine.laws;
float<Pa> elastic(float<Pa> modulus, float strain)
{
return modulus * strain;
}
bool fails(float<Pa> stress, float<Pa> tensile, float<Pa> compressive)
{
return stress > tensile || -stress > compressive;
}
Hooke's law, sigma = E epsilon, with fails beside it: a stress past the tensile strength stretched or the compressive one squeezed.
import os, sys
root = os.getcwd()
while not os.path.isfile(os.path.join(root, 'docsite.json')):
root = os.path.dirname(root)
sys.path.insert(0, os.path.join(root, 'tools'))
import numpy as np
from lusine import trait_run
F32 = np.float32
U = float(np.finfo(F32).eps) / 2 # binary32's unit roundoff, 2^-24
rng = np.random.default_rng(2026)
The kernel¶
The law as the library states it.
print(trait_run.source('laws/elastic.trait'))
// lusine.laws.elastic: a solid's stress from its strain, and whether it fails. Hooke's law in one dimension, sigma =
// E epsilon, E the Young's modulus and epsilon the strain (the
// change in length over the length, positive stretched); the solid fails where the stress's magnitude passes its
// strength (tensile where stretched, compressive where squeezed: a brittle material's, which a game's bonds break by).
// Source: Timoshenko and Goodier, Theory of Elasticity, 3rd ed., section 3; Ashby and Jones, Engineering Materials 1,
// 4th ed., chapters 3 and 8.
// Checked against a reference in numpy, which has its own tests.
namespace lusine.laws;
float<Pa> elastic(float<Pa> modulus, float strain)
{
return modulus * strain;
}
bool fails(float<Pa> stress, float<Pa> tensile, float<Pa> compressive)
{
return stress > tensile || -stress > compressive;
}
The program¶
Each section runs the kernel over a lattice of items, one a case: the host sets every parameter, and main() binds the kernel over the items.
print(trait_run.law_program('laws/elastic.trait', 'elastic'))
from lusine.kernels.laws.elastic import elastic as reference
def law(**kw):
return trait_run.law('laws/elastic.trait', 'elastic', **kw)
namespace probe;
lusine.version version = { { number = 1 } };
struct item {
host float<Pa> modulus;
host float strain;
float<Pa> result_0;
};
lattice<item, 1> items(spacing = 1 m, boundary = wall);
void main() { items.result_0 = lusine.laws.elastic(items.modulus, items.strain); }
Reference: metals and stone¶
Moduli from 10 to 400 GPa and strains to a hundredth either way: the kernel within binary32's roundings of the reference.
e, eps = rng.uniform(1e10, 4e11, 400), rng.uniform(-1e-2, 1e-2, 400)
got, want = law(modulus=e, strain=eps), reference(e, eps)
err = np.max(np.abs(got - want) / np.abs(want))
assert err <= 8 * U, err / U
print(f'within {err / U:.1f} x 2^-24 over {len(e)} solids')
within 2.3 x 2^-24 over 400 solids
Invariant: linear and odd¶
The stress is linear in the strain, and squeezing by as much as stretching gives the opposite stress.
eps = np.linspace(-1e-3, 1e-3, 21)
s = law(modulus=np.full(21, 2e11), strain=eps)
assert np.allclose(s, -s[::-1], rtol=1e-6) and np.allclose(np.diff(s), np.diff(s)[0], rtol=1e-5)
print('linear, and odd in the strain')
linear, and odd in the strain
Limit: failure, brittle and strong¶
A brittle solid, concrete's 3 MPa stretched and 30 squeezed: 4 MPa stretched fails, 4 squeezed holds, 40 squeezed fails.
broke = trait_run.law('laws/elastic.trait', 'fails', stress=np.array([4e6, -4e6, -4e7]), tensile=3e6, compressive=3e7)
assert list(broke.astype(bool)) == [True, False, True]
print('fails stretched at 4 MPa; holds squeezed at 4; fails squeezed at 40')
fails stretched at 4 MPa; holds squeezed at 4; fails squeezed at 40