0

私はモデリカで発射体の軌道をモデル化しようとしていますが、高さに対する air_pressure依存性をモデル化すると、奇妙な結果が得られます. air_pressureのコードだけを別のファイルにコピーするだけで、かなり良い出力が得られます. どうしてそうなの？これは発射体の元のコードです。

class proj
  //position
  Real x(start = 6371000);
  Real y(start = 0);
  Real z(start = 0);
  Real r(start = 6371000);
  Real row(start = 6371000);
  Real h(start = 0);
  //velocity
  Real v_x(start = 2000);
  Real v_y(start = 0);
  Real v_z(start = 0);
  Real v_mag(start = sqrt(4000000));
  Real v_row(start = sqrt(4000000));
  //acceleration
  Real a_x();
  Real a_y();
  Real a_z();
  Real a_mag();
  //total force
  Real ft_x();
  Real ft_y();
  Real ft_z();
  Real ft_mag();
  //Gravitational force
  Real fg_x();
  Real fg_y();
  Real fg_z();
  Real fg_mag();
  //Drag force
  Real fd_x();
  Real fd_y();
  Real fd_z();
  Real fd_mag();
  Real energy();
  //Fictitious force
  Real ff_x();
  Real ff_y();
  Real ff_z();
  Real ff_mag();
  //gravitation variables
  Real g(start = G * earth_mass / earth_radius ^ 2);
  //Drag variables
  Real air_density(start = surface_density);
  Real air_pressure(start = surface_pressure);
  //_________________________________________________________________________
  //Earth parameters
  parameter Real earth_radius = 6371000;
  parameter Real earth_mass = 5.972 * 10 ^ 24;
  parameter Real earth_rot_vel = 7.292115 * 10 ^ (-5);
  //Missile parameters
  parameter Real missile_effective_area = 0.785398;
  parameter Real missile_mass = 200;
  parameter Real Cd = 0.02;
  //Air parameters
  parameter Real temp_lapse_rate = 0.0065;
  parameter Real air_molar_mass = 28.97 / 1000;
  parameter Real surface_temp = 303.15;
  parameter Real surface_pressure = 101325;
  parameter Real surface_density = 1.1644;
  //constants
  parameter Real pi = 3.14159265358979;
  parameter Real G = 6.67384 * 10 ^ (-11);
  parameter Real R = 8.3144621;
  //__________________________________________________________________________
equation
  //position
  der(x) = if h < 0 then 0 else v_x;
  der(y) = if h < 0 then 0 else v_y;
  der(z) = if h < 0 then 0 else v_z;
  r = sqrt(x ^ 2 + y ^ 2 + z ^ 2);
  row = sqrt(x ^ 2 + y ^ 2);
  h = r - earth_radius;
  //velocity
  der(v_x) = a_x;
  der(v_y) = a_y;
  der(v_z) = a_z;
  v_mag = sqrt(v_x ^ 2 + v_y ^ 2 + v_z ^ 2);
  v_row = sqrt(v_x ^ 2 + v_y ^ 2);
  //acceleration
  a_x = ft_x / missile_mass;
  a_y = ft_y / missile_mass;
  a_z = ft_z / missile_mass;
  a_mag = ft_mag / missile_mass;
  //total force
  ft_x = fg_x + fd_x + ff_x;
  ft_y = fg_y + fd_y + ff_y;
  ft_z = fg_z + fd_z + ff_z;
  ft_mag = sqrt(ft_x ^ 2 + ft_y ^ 2 + ft_z ^ 2);
  //Gravitational force
  fg_x = -1 * fg_mag * row / r * x / r;
  fg_y = -1 * fg_mag * row / r * y / r;
  fg_z = -1 * fg_mag * z / r;
  fg_mag = G * earth_mass * missile_mass / r ^ 2;
  //Drag force
  fd_x = if v_mag > 0 then -1 * fd_mag * v_row / v_mag * v_x / v_mag else 0;
  fd_y = if v_mag > 0 then -1 * fd_mag * v_row / v_mag * v_y / v_mag else 0;
  fd_z = if v_mag > 0 then -1 * fd_mag * v_z / v_mag else 0;
  fd_mag = 0.5 * air_density * v_mag ^ 2 * Cd * missile_effective_area;
  //Fictitious force
  ff_x = missile_mass * earth_rot_vel * (2 * v_y + earth_rot_vel * x);
  ff_y = missile_mass * earth_rot_vel * (-2 * v_x + earth_rot_vel * y);
  ff_z = 0;
  ff_mag = sqrt(ff_x ^ 2 + ff_y ^ 2 + ff_z ^ 2);
  energy = 0.5 * missile_mass * v_mag * v_mag - G * earth_mass * missile_mass / r;
  //Gravitation variables
  g = fg_mag / missile_mass;
  //Drag variables
  air_density = air_pressure * air_molar_mass / (R * surface_temp);
  air_pressure = if h > 46600 then 0 else surface_pressure * (1 - temp_lapse_rate * (r - earth_radius) / surface_temp) ^ (g * air_molar_mass / (R * temp_lapse_rate));
end proj;

これは抽出された air_pressure のコードで、単独では正常に動作しますが、発射体コード全体では奇妙な動作をします。

class check
  Real g;
  Real h(start = 0);
  Real pressure;
  //Earth parameters
  parameter Real earth_radius = 6371000;
  parameter Real earth_mass = 5.972 * 10 ^ 24;
  //Air parameters
  parameter Real temp_lapse_rate = 0.0065;
  parameter Real air_molar_mass = 28.97 / 1000;
  parameter Real surface_temp = 303.15;
  parameter Real surface_pressure = 101325;
  //constants
  parameter Real pi = 3.14159265358979;
  parameter Real G = 6.67384 * 10 ^ (-11);
  parameter Real R = 8.3144621;
equation
  g = G * earth_mass / (earth_radius + h) ^ 2;
  h = time;
  pressure = if h > 46600 then 0 else surface_pressure * (1 - temp_lapse_rate * h / surface_temp) ^ (g * air_molar_mass / (R * temp_lapse_rate));
end check;

これは、クラスprojからのair_pressure対時間のグラフです。

ここに画像の説明を入力

これはair_pressureのコードのみを含むclass checkからのair_pressure対time(=height)グラフです

ここに画像の説明を入力

4

2 に答える 2

0

システム内の の値を確認する必要がありhます。で何かおかしなことが起こっているのではないかと思いhます。

于 2013-12-02T16:10:15.160 に答える
0

x が何を指しているのかは明確ではありません。開始値から、x は地球の中心からの距離を参照しているように見えます。これが事実である場合、r は何ですか? 私は、x、y & z が 3D 空間での発射体の位置を地面 (地球の表面) 上の基準点と見なしていると想定しています。これが正しければ、x、y を緯度と地表からの高さとしての経度と z。

これにより、h のあいまいさが修正されます。

于 2013-12-09T22:29:34.023 に答える