Prędkość dzwięku w troposferze

[ELZA]

Prędkość dżwięku jak wiadomo jest zależna od temperatury, w troposferze temperatura spada do około - 60 stopni celsjusza, przy czym jej spadek jest liniowy o okolo 0,6 stopnia na każde 100m w góre...



za "chciałbysiem" wykopałam dane pochodzą z sondażu atmosfery przeprowadzonego przez 26 października 2009r w łebie (tak to dzisiaj)



http://weather.uwyo.edu/cgi-bin/soundin ... STNM=12120



dostępne powyżej



ogólny rozkład temperatury przedstawia się następująco:







punkty powyżej 10500m odrzuciłam ponieważ od tego momentu przypuszczalnie balon wlatuje do strefy tropopauzy gdzie zmiany temperatury są niewielkie.



liniowy spadek temperatury od 0 do około 10500m pokazuje wykres:





wykorzystując wzór na prędkość dźwięku :







gdzie:



k - wykładnik adiabaty- Cp/Cv, który dla powietrza wynosi 1,4

R - indywidualna stała gazowa dla powietrza suchego R=287,05 [J/(kg*K)]

T - temperatura w [K]



udało mi się wykreślić zależność prędkości dźwięku w funkcji wysokości dla troposfery:







W tym zakresie temperatur krzywa pierwiastkowa daje się fajnie przybliżyć prostą, na wykresie dodałam linie trendu z równaniem.



Trzeba jeszcze nad tym dalej popracować, ponieważ wiadomo że powietrze suche nie jest i to wszystko jest jakimś tam przybliżeniem, ale to już coś

pomyślę jeszcze nad tym ale nie dzisiaj



PZDR

Elka

[robercik]

Cześć

Dzięki ELZA

Jak się rozkręcisz to może dasz radę wykonać zestaw zależności na prędkość dźwięku do wysokości całych 100km. Wstawimy do ShizoRocketApplication żeby się ładnie malował wykres w Machach.

Robota jest, trzeba uwzględnić skład powietrza itd...



pzdr

[ELZA]

mozna jeszze pomyslec nad stratosfera, to do okolo 50 km, ale w tym momencie nie mamy dostepnych takich danych meteorologicznych jak tutaj i to bedzie dosc spore przyblizenie... wyzej to juz wlasciwie proznia jest, poszperam popatrze

nad ta stala gazowa tez trzeba jeszcze pomyslec bo ona dotyczy powietrza bez zadnej zawartosci wilgoci...

[Grafit_alpha3]

http://en.wikipedia.org/wiki/Internatio ... Atmosphere

Polecam References. Znajdziesz tam np. U.S. Standard Atmosphere

[ELZA]

dzieki za linka zapoznam sie wieczorem

[kbosak]

Rozwiązanko ultymatywne w C++

Kod: Zaznacz cały

//http://en.wikipedia.org/wiki/Ellipse
double EllipseRadius(double major_radius, double minor_radius, double angle_degrees)
{
    const double angle_radians=angle_degrees*DEG2RAD;
    const double dx=major_radius*cos(angle_radians);
    const double dy=minor_radius*sin(angle_radians);
    return (major_radius*minor_radius)/sqrt((dx*dx) + (dy*dy));
}

//http://en.wikipedia.org/wiki/Earth_radius
double EarthRadius(double latitude_degrees)
{
    return EllipseRadius(6378.1370e3, 6356.7523e3, latitude_degrees);
}

//http://en.wikipedia.org/wiki/Barometric_formula
double TemperatureAtAltitudeCelsius(
    double h, 
    double T0_celsius=15.0, 
    double T0_level=0.0
    )
{
    if(h<11000.0)
    {
	   //return T0_celsius + h*(-6.5/1000.0);//Only if T0_celsius is exactly 15.0C
	   //must meet: T0_celsius at h=T0_level, -56.5 at 11000
	   //http://en.wikipedia.org/wiki/Lapse_rate
	   const double lapse_rate=(-56.5 - T0_celsius)/(11000.0-T0_level);
	   return T0_celsius + (h-T0_level)*lapse_rate;
    }
    if(h<20000.0)
    {
	   return -56.5;
    }
    if(h<32000.0)
    {
	   return -56.5 + (h-20000.0)*(1.0/1000.0);
    }
    if(h<47000.0)
    {
	   return -44.5 + (h-32000.0)*(2.8/1000.0);
    }
    if(h<51000.0)
    {
	   return -2.5;
    }
    if(h<71000.0)
    {
	   return -2.5 + (h-51000.0)*(-2.8/1000.0);
    }
    /*
    if(h>84852.0)
    {
	   h=84852.0;
    }
    if(h<=84852.0)
    */
    {
	   return -58.5 + (h-71000.0)*(-2.0/1000.0);
    }
}

//http://en.wikipedia.org/wiki/Gravitational_constant
//http://en.wikipedia.org/wiki/Earth_mass
double GLocal(double h, double latitude_degrees)
{
    const double Gconstant=6.67384E-11;
    const double Mearth=5.97219E24;
    const double R=EarthRadius(latitude_degrees);
    const double dist=R+h;
    return Gconstant*Mearth/(dist*dist);
}

int main()
{
    /*
    for(double h=0.0; h<=100000.0; h+=1000.0)
    {
	   printf("%lf %lf\n", h, TemperatureAtAltitudeCelsius(h) );
    }
    */
    
    const double latitude_degrees=53.0;//Torun test site
    
    const double P0=101325.0;//pressure at sea level
    
    //const double T0=15.0;//15[C] is stdatm
    const double T0=22.0;//temperature at T0_level
    const double T0_level=53.0;//Altitude above mean sea level at which temperature is measured
    
    double M=0.0289644;//Molar mass of dry air [kg/mol]
    
    //const double R=8.31432;//stdatm model
    //http://en.wikipedia.org/wiki/Gas_constant
    const double R=8.3144621;//precise ideal gas constant 
    
    //const double g=9.80665;//stdatm model

    //http://en.wikipedia.org/wiki/Barometric_formula
    double integral_value=0;
    const double integral_dh=1;

    const double alti_step=1000.0;
    int alti_i=0;
    double hnext=0;
    const double hmax=(100000.0+(alti_step*0.5));
    double fvalue_previous;//Trapeze and Simpson's rule
    //http://en.wikipedia.org/wiki/Simpson%27s_rule
    double fvalue_middle;//Simpson's rule only
    double P=0.0;
    for(double h=0.0; h<=hmax; h+=integral_dh)
    {
	   const double g=GLocal(h, latitude_degrees);
	   //http://en.wikipedia.org/wiki/Absolute_zero
	   const double T_celsius=TemperatureAtAltitudeCelsius(h, T0, T0_level);
	   const double T_kelvin=T_celsius+273.15;
	   const double air_density=M/(R*T_kelvin);// 1/(R_specific*T_kelvin)
	   const double fvalue=g*air_density;
	   
	   //http://en.wikipedia.org/wiki/Density_of_air
	   if(P>0)
	   {
		  const double Psat=6.1078*pow(10.0, (7.5*T_celsius)/(T_celsius+237.3));//saturation vapor pressure of water
		  const double RH=0;//0.40;//Relative humidity
		  const double Pvapor=RH*Psat;
		  const double Pdry=P-Pvapor;//Since P_total=Pdry+Pvapor
		  const double Mdry=0.0289644;//Molar mass of dry air
		  const double Mvapor=0.018016;//Molar mass of water vapor
		  const double Mwet=( (Pdry*Mdry) + (Pvapor*Mvapor) )/P;
		  M=Mwet;
	   }
	   
	   {//Simpson's rule only
		  const double hmiddle=h-(integral_dh*0.5);
		  const double g=GLocal(hmiddle, latitude_degrees);
		  //http://en.wikipedia.org/wiki/Absolute_zero
		  const double T_celsius=TemperatureAtAltitudeCelsius(h, T0, T0_level);
		  const double T_kelvin=T_celsius+273.15;
		  const double air_density=M/(R*T_kelvin);// 1/(R_specific*T_kelvin)
		  fvalue_middle=g*air_density;
	   }
	   
	   if(h>0)
	   {
		  //integral_value+=fvalue*integral_dh;//Rectangle rule
		  //integral_value+=((fvalue_previous+fvalue)*0.5)*integral_dh;//Trapeze rule
		  integral_value+=( (fvalue_previous + (4*fvalue_middle) + fvalue)/6.0 )*integral_dh;//Simpson's rule
	   }
	   
	   P=P0*exp(-integral_value);
	   if(h>=hnext)
	   {
		  //http://forum.rakiety.org.pl/topic696.html
		  const double gamma=1.4;//Adiabatic index Cp/Cv
		  //http://en.wikipedia.org/wiki/Gas_constant
		  const double R_specific=R/M;//Specific gas constant for dry air, unless you update M
		  //http://en.wikipedia.org/wiki/Speed_of_sound
		  const double mach_speed=sqrt(gamma*(R_specific*T_kelvin));//[m/s]

		  printf("% 14.6lf % 7.0lf %6.2lf\n", P, h, mach_speed);
		  
		  alti_i++;
		  hnext=alti_i*alti_step;
	   }
	   fvalue_previous=fvalue;
    }
    return 0;
}

[kbosak]