Buoyancy¶

The upward acceleration of a block of one density in a surrounding of another, Archimedes' reduced gravity g' = g (rho_s - rho) / rho, gravity positive down: a lighter block rises, a denser one sinks, a neutral one is still.

Source: Kundu & Cohen, Fluid Mechanics, 4.18; Turner, Buoyancy Effects in Fluids, ch. 1. Checked against a reference in numpy, which has its own tests.

Source

Kundu & Cohen, Fluid Mechanics, 4.18 (buoyancy); Turner, Buoyancy Effects in Fluids, ch. 1.

Proof

reference, invariant, limit.

buoyancy.trait

namespace lusine.laws;

float<m/s2> buoyancy(float<kg/m3> density, float<kg/m3> surrounding, float<m/s2> gravity)
{
  return gravity * (surrounding - density) / density;
}

Uses, without proving it here: properties.gas_density.

The upward acceleration of a block of one density in a surrounding of another, Archimedes' reduced gravity g' = g (rho_s - rho) / rho.

In [1]:

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
import matplotlib.pyplot as plt
from lusine import trait_run
EPS = float(np.finfo(np.float32).eps)

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 import matplotlib.pyplot as plt from lusine import trait_run EPS = float(np.finfo(np.float32).eps)

The kernel¶

The law as the library states it, and as this page runs it.

In [2]:

print(trait_run.source('laws/buoyancy.trait'))

print(trait_run.source('laws/buoyancy.trait'))

// lusine.laws.buoyancy: the upward acceleration of a block of one density in a surrounding of another, Archimedes'
// reduced gravity g' = g (rho_s - rho) / rho, gravity positive down: a lighter block rises, a denser one sinks, a
// neutral one is still.
// Source: Kundu & Cohen, Fluid Mechanics, 4.18; Turner, Buoyancy Effects in Fluids, ch. 1.
// Checked against a reference in numpy, which has its own tests.

namespace lusine.laws;

float<m/s2> buoyancy(float<kg/m3> density, float<kg/m3> surrounding, float<m/s2> gravity)
{
  return gravity * (surrounding - density) / density;
}

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.

In [3]:

print(trait_run.law_program('laws/buoyancy.trait', 'buoyancy'))

print(trait_run.law_program('laws/buoyancy.trait', 'buoyancy'))

namespace probe;

lusine.version version = { { number = 1 } };

struct item {
  host float<kg/m3> density;
  host float<kg/m3> surrounding;
  host float<m/s2> gravity;
  float<m/s2> result_0;
};
lattice<item, 1> items(spacing = 1 m, boundary = wall);

void main() { items.result_0 = lusine.laws.buoyancy(items.density, items.surrounding, items.gravity); }

Reference: air a few degrees warmer¶

For air at a fixed pressure g' = g ΔT / T: a block ten degrees warmer than 290 K air rises at about 0.34 m/s². The two densities are the library's lusine.properties.gas_density, for dry air's 28.9645 g/mol.

In [4]:

g = 9.80665
rho = trait_run.law('properties/gas_density.trait', 'gas_density', pressure=101325.0, temperature=[300.0, 290.0],
                    molar_mass=0.0289645)
g1 = trait_run.law('laws/buoyancy.trait', 'buoyancy', density=rho[0], surrounding=rho[1], gravity=g)[0]
print(f"{g1:.4f} m/s^2 against g dT/T = {g * 10 / 290:.4f}")
assert abs(g1 - g * 10 / 290) < 0.002

g = 9.80665 rho = trait_run.law('properties/gas_density.trait', 'gas_density', pressure=101325.0, temperature=[300.0, 290.0], molar_mass=0.0289645) g1 = trait_run.law('laws/buoyancy.trait', 'buoyancy', density=rho[0], surrounding=rho[1], gravity=g)[0] print(f"{g1:.4f} m/s^2 against g dT/T = {g * 10 / 290:.4f}") assert abs(g1 - g * 10 / 290) < 0.002

0.3382 m/s^2 against g dT/T = 0.3382

Invariant: the sign is the density's¶

Lighter rises, denser sinks, and the reduced gravity grows with the difference.

In [5]:

d = np.linspace(0.8, 1.6, 17)
a = trait_run.law('laws/buoyancy.trait', 'buoyancy', density=d, surrounding=1.2, gravity=9.80665)
assert np.all((a > 0) == (np.float32(d) < np.float32(1.2))) and np.all(np.diff(a) < 0)
fig, ax = plt.subplots(figsize=(7, 3))
ax.plot(d, a, 'o-')
ax.axhline(0, lw=0.5)
ax.set_xlabel('block density (kg/m3)')
ax.set_ylabel("g' (m/s^2)")
plt.show()

d = np.linspace(0.8, 1.6, 17) a = trait_run.law('laws/buoyancy.trait', 'buoyancy', density=d, surrounding=1.2, gravity=9.80665) assert np.all((a > 0) == (np.float32(d) < np.float32(1.2))) and np.all(np.diff(a) < 0) fig, ax = plt.subplots(figsize=(7, 3)) ax.plot(d, a, 'o-') ax.axhline(0, lw=0.5) ax.set_xlabel('block density (kg/m3)') ax.set_ylabel("g' (m/s^2)") plt.show()

Limit: neutral is still¶

A block of the surroundings' own density feels nothing.

In [6]:

a = trait_run.law('laws/buoyancy.trait', 'buoyancy', density=[1.2, 1000.0], surrounding=[1.2, 1000.0],
                  gravity=9.80665)
assert np.all(a == 0.0), a
print('neutral:', a)

a = trait_run.law('laws/buoyancy.trait', 'buoyancy', density=[1.2, 1000.0], surrounding=[1.2, 1000.0], gravity=9.80665) assert np.all(a == 0.0), a print('neutral:', a)

neutral: [0. 0.]