This calculates frequencies of beams as part of the "Wind Sounds: Which component sings?" post on ianohara.com

natural_frequencies.py

#!/usr/bin/env python3

import dataclasses
from dataclasses import dataclass
import math

@dataclass
class Beam:
    name: str
    width_x: float
    width_y: float
    thickness: float
    length: float
    youngs_modulus: float
    density: float

    def _I(self):
        """
        From wikipedia: https://en.wikipedia.org/wiki/List_of_second_moments_of_area
        X is the b dimension
        Y is the h dimension
        """
        b = self.width_x
        b1 = self.width_x - 2*self.thickness
        h = self.width_y
        h1 = self.width_y - 2*self.thickness

        I_xx = (b * h**3 - b1 * h1**3) / 12
        I_yy = (b**3 * h - b1**3 * h1) / 12

        return (I_xx, I_yy)

    def I_xx(self):
        return self._I()[0]

    def I_yy(self):
        return self._I()[1]

    def area(self):
        return self.width_x * self.width_y

def mm_to_m(mm):
    return mm/1000.0

def calc_w(length, area, I, E, rho, n):
    in_rad_per_s = n**2 * math.pi**2 * ((E * I) / (rho * area * length**4))**(0.5)

    return (180.0/math.pi)*in_rad_per_s

def calc_w_in_x(beam, mode):
    return calc_w(beam.length, beam.area(), beam.I_xx(), beam.youngs_modulus, beam.density, mode)

def calc_w_in_y(beam, mode):
    return calc_w(beam.length, beam.area(), beam.I_yy(), beam.youngs_modulus, beam.density, mode)

ASTM_A847_YOUNGS_MODULUS = 190e9 # 190 GPa
ASTM_A847_DENSITY = 7900.0 # kg / m^3

def nudge_length(beam, frac_change):
    return dataclasses.replace(beam, length=beam.length*(1.0 + frac_change))

def nudge_thickness(beam, frac_change):
    return dataclasses.replace(beam, thickness=beam.thickness*(1.0 + frac_change))

def nudge_width_x(beam, frac_change):
    return dataclasses.replace(beam, width_x=beam.width_x*(1.0 + frac_change))

def nudge_width_y(beam, frac_change):
    return dataclasses.replace(beam, width_y=beam.width_y*(1.0 * frac_change))

def print_beam_freqs(beam):
    print(f"Beam {beam.name}:")
    print(f"  l={beam.length} [m], x={beam.width_x} [m], y={beam.width_y} [m], t={beam.thickness} [m], E={beam.youngs_modulus} [Pa], rho={beam.density} [kg/m^3]")
    print(f"  w_1 in xx={calc_w_in_x(beam, 1):.0f} [Hz]")
    print(f"  w_1 in yy={calc_w_in_y(beam, 1):.0f} [Hz]")

def print_beams_table(beams):
    print("NORMAL")
    for beam in beams:
        print_beam_freqs(beam)

if __name__ == "__main__":
    beams = [
        Beam(
            name="B",
            width_x=mm_to_m(152),
            width_y=mm_to_m(102),
            thickness=mm_to_m(6.4),
            length=10.2,
            youngs_modulus=ASTM_A847_YOUNGS_MODULUS,
            density=ASTM_A847_DENSITY
        ),
        Beam(
            name="C",
            width_x=mm_to_m(254),
            width_y=mm_to_m(152),
            thickness=mm_to_m(9.5),
            length=18.2,
            youngs_modulus=ASTM_A847_YOUNGS_MODULUS,
            density=ASTM_A847_DENSITY
        ),
        Beam(
            name="D",
            width_x=mm_to_m(102),
            width_y=mm_to_m(51),
            thickness=mm_to_m(6.4),
            length=4.9,
            youngs_modulus=ASTM_A847_YOUNGS_MODULUS,
            density=ASTM_A847_DENSITY
        ),
        Beam(
            name="E",
            width_x=mm_to_m(102),
            width_y=mm_to_m(51),
            thickness=mm_to_m(6.4),
            length=9.1,
            youngs_modulus=ASTM_A847_YOUNGS_MODULUS,
            density=ASTM_A847_DENSITY
        )
    ]

    print_beams_table(beams)