Workspace
Jonathan Brooks/

Mass to Drag by Density and Radii (Comparative Calculator)

0
Beta
Spinner

Mass to Drag by Density and Radii (Comparative Calculator)

from math import pi
import numpy as np
import pandas as pd
import seaborn as sns
import matplotlib.pyplot as plt
from scipy.stats import uniform

def main():
    '''This function takes these inputs:
    
        Enter on new lines:

        1. density of sphere,
        2. density of fluid,
        3. velocity of motion of the sphere relative to the fluid,
        4. a number of example spheres or 0 for exact spheres,
        5. a number of exact spheres (if 0 on 4),
        6. length of radius or radii of exact spheres (if 0 on 4 and not 0 on 5,).

        It plots Mass to Drag ratio and Kinetic Energy relative to radius length.'''
    
    sphere_density = float(input('Sphere Density in kg/m^3 (water ice is around 920):'))
    fluid_density = float(input('Fluid Density in kg/m^3 (the atmosphere at sea level is around 1.21):'))
    velocity = float(input('Velocity in m/s (gravity acceleration at sea level is around 9.8 m/s^2):'))
    samp_spheres = int(input('Number of example spheres (0 for exact):'))
    
    if samp_spheres == 0:
        spec_spheres = int(input('Number of exact spheres:'))
        if spec_spheres == 1:
            radii = [float(input(f'Radius of sphere in micrometers:'))]
        elif spec_spheres > 1:
            radii = [float(input(f'Radius of sphere {i+1} in micrometers:')) for i in range(spec_spheres)]  # manual
    elif samp_spheres > 0:
        np.random.seed(samp_spheres)
        radii = uniform.rvs(1, 1e5, size=samp_spheres)  # automatic
    
    radii.sort()
    radii_normz = norm(radii)
    vols = sphere_volume(radii)  # um^3
    masses = mass(sphere_density, vols)  # kg
    kinergy = [(0.5 * mass * velocity**2) for mass in masses]
    ke_normz = norm(kinergy)
    ke_mean = np.mean(kinergy)
    drags = sphere_drag(fluid_density, velocity, radii)  # m*s^3
    mtds = mass_to_drag(masses, drags)  # kg/m^2
    mtds_normz = norm(mtds)
    mini, maxi = min(radii_normz), max(radii_normz)
    mtd_mean = np.mean(mtds)
    alpha = 0.5
    
    # Create a pandas DataFrame
    data = {'Radius': radii, 'Mass_to_Drag_Ratio': mtds, 'Kinergy': kinergy, 'Normed_Radius': radii_normz, 'Normed_Kinergies': ke_normz}
    df = pd.DataFrame(data)
    corcof = df['Radius'].corr(df['Mass_to_Drag_Ratio'])
    print(f'Correlation Coefficient: {corcof}')
    print(f'Average Mass to Drag Ratio (exact if 1 sphere): {mtd_mean:.3g} kg/m*s^3')
    print(f'Average kinetic energy (exact if 1 sphere): {ke_mean:.3g} Joules')


    # Create subplots with shared x-axis
    fig, ax1 = plt.subplots(figsize=(12, 9))

    # Scatter plot for Mass to Drag Ratio
    sns.scatterplot(data=df, x='Radius', y='Mass_to_Drag_Ratio', hue='Normed_Radius', size=radii_normz, sizes=(mini, maxi), alpha=alpha, palette='blend:#7AB,#EDB', legend=False, ax=ax1)

    # Create a second y-axis for Kinergy
    ax2 = ax1.twinx()
    ax2.scatter(df['Radius'], df['Kinergy'], marker='o', color='black', alpha=alpha, label='Kinetic Energy')

    # Set y-axis limits and label for the first axis
    ax1.set_ylim(0, 1.5 * max(mtds_normz))
    ax1.set_ylabel('Mass to Drag Ratio (kg/m*s^3)')

    # Set y-axis limits and label for the second axis
    ax2.set_ylim(0, 1.1 * max(ke_normz))
    ax2.set_ylabel('Kinetic Energy in Joules (kg*(m/s)^2)')

    # Set title
    plt.title('Mass to Drag Ratio and Kinetic Energy vs. Sphere Radius')

    # Display legend
    ax2.legend(loc='upper right')

    plt.xlim(0, 1.5 * max(radii_normz))
    plt.show()


def norm(radii):
    '''Returns normalized radii for proportionality in sns.scatterplot size param.'''
    normz = [radius / max(radii) for radius in radii]    
    radii_normz = np.array([radius * norm for radius, norm in zip(radii, normz)])
    return radii_normz

def sphere_volume(radii):
    '''Returns sphere volume after unit conversion um^3 to m^3'''
    return [(4/3*pi*radius**3 * 1e-18) for radius in radii]

def sphere_drag(fluid_density, velocity, radii):
    '''Returns sphere drag after unit conversion um^2 to m^2'''
    return [(0.45*0.5*fluid_density*(velocity**2)*pi*radius**2 * 1e-12) for radius in radii]
# D = Cd * . 5 * rho * V^2 * A. Where D is equal to the drag (m/s^2), rho is the fluid density (kg/m^3), V is the velocity (m/s), A is a reference area (m^2), and Cd is the drag coefficient (.45 for a sphere).
def mass(sphere_density, vols):
    '''Returns mass as kg'''
    return [(sphere_density * vol) for vol in vols]

def mass_to_drag(masses, drags):
    '''Returns mass to drag ratio as kg/m*s^3'''
    return [mass/drag for mass, drag in zip(masses, drags)]

if __name__ == '__main__':
    main()

This comparative calculator is presently evolving to evaluate the probability that particulate will remain airborne.

Enter on new lines:

  1. density of sphere
  2. density of fluid
  3. velocity of motion of sphere relative to fluid
  4. number of example spheres or 0 for exact spheres
  5. number of exact spheres (if 0 on #4)
  6. length of radius or radii of exact spheres (if 0 on #4 and not 0 on #5,)

Units are kg/m^3 for density and micrometers (millionths of a meter) for radii. I plan on improving the drag calculation to more accurately model reality from sea level, to extraterrestrial atmospheres, to the intersteller medium (derived from the tables here: https://www.engineeringtoolbox.com/standard-atmosphere-d_604.html), and beyond that. Additional ideas for scaling up are elaborated below. TL:DR space odyssey!

In a Physical Geography course, I learned how aerosols and particulates of tiny sizes can stay aloft indefinitely. The Mass to Drag Ratio of spherical objects has linear proportionality to radii length (correlation coefficient of 1); inversely, Air Drag to Mass Ratios decrease as radii increase in spheres of equal density (mass/volume). Larger spheres with an equal density move through through air more easily, and moreso at higher gaseous viscosities (relates to altitude & temperature).

Of course, all of that is very intruiging, and then it scales to planet formation, in the Oort Cloud and similar stellar torii, and protoplanetary discs around stars; at various distances from stars, the solar wind (outward) force equals the solar gravitation (inward) force, acting on spheroids of varying densities and radii lengths. Solar wind has a plasmic (thermoelectrically charged) viscosity, and similar to air density around the Earth, stellar wind density and viscosity increase on a logarithmic scale, inversely proportional to the distance from solar system center, very similar to gravity.

As the pushing power of solar wind diffuses with distance, so does the pulling power of solar gravity. Newtons (gravity units) have inverse proportionality to the square of the distance (d^-2) of separation between centers of gravity, so when the distance doubles (delta d = 2), gravity is reduced by a factor of 1/d^2 or 1/4.

The force of solar wind similarly reduces, by the entropic dispersion of radiant energy; those physics are vastly more complex than Newtonian Gravity, yet I will strive to get an adequate formulaic approximation to build upon.

Asteroids, protoplanets, and planets amass over aeons, from particulate that had blown to the outer solar system in the solar wind (or had otherwise there arrived); as the sizes increase and mass to drag ratios decrease, the balance of forces slides solar orbits of these objects inward. The lower the mass to drag ratio, the lower the orbit gets toward the gravitational center of the solar system.

  • AI Chat
  • Code